Gibt es eine Möglichkeit, eine Beziehung zwischen Koeffizienten in der logistischen Regression zu erzwingen?

Ich möchte ein logistisches Regressionsmodell angeben, bei dem ich die folgende Beziehung habe:

$E[Y_i|X_i] = f(\beta x_{i1} + \beta^2x_{i2})$ wobei die inverse Logit-Funktion ist.$f$

Gibt es eine "schnelle" Möglichkeit, dies mit bereits vorhandenen R-Funktionen zu tun, oder gibt es einen Namen für ein Modell wie dieses? Mir ist klar, dass ich den Newton-Raphson-Algorithmus, der für die logistische Regression verwendet wird, modifizieren kann, aber dies ist eine Menge theoretischer und codierender Arbeit, und ich suche nach einer Abkürzung.

BEARBEITEN: Das Abrufen von Punktschätzungen für ist mit optim () oder einem anderen Optimierer in R ziemlich einfach, um die Wahrscheinlichkeit zu maximieren. Aber ich brauche Standardfehler bei diesen Jungs.$\beta$

Dies ist mit der Optimierungsfunktion in R ziemlich einfach. Mein Verständnis ist, dass Sie eine logistische Regression ausführen möchten, bei der y binär ist. Sie schreiben einfach die Funktion und stecken sie dann in optim. Unten ist ein Code, den ich nicht ausgeführt habe (Pseudocode).

#d is your data frame and y is normalized to 0,1
your.fun=function(b)
{

    EXP=exp(d$x1*b +d$x2*b^2)

    VALS=( EXP/(1+EXP) )^(d$y)*( 1/(1+EXP) )^(1-d$y) 
    return(-sum(log(VALS)))
}

result=optim(0,your.fun,method="BFGS",hessian=TRUE)
# estimates 
 result$par
    #standard errors
    sqrt(diag(inv(result$hessian)))
# maximum log likelihood
-result$value

Beachten Sie, dass your.fun das Negativ einer Log-Likelihood-Funktion ist. Optim maximiert also die Protokollwahrscheinlichkeit (standardmäßig minimiert optim alles, weshalb ich die Funktion negativ gemacht habe). Wenn Y nicht binär ist, finden Sie unter http://fisher.osu.edu/~schroeder.9/AMIS900/ch5.pdf multinomiale und bedingte Funktionsformulare in Logit-Modellen.

Die obige Antwort ist richtig. Als Referenz finden Sie hier einen ausgearbeiteten funktionierenden R-Code, um ihn zu berechnen. Ich habe mir erlaubt, einen Abschnitt hinzuzufügen, weil Sie wahrscheinlich einen davon wollen.

## make some data
set.seed(1234)
N <- 2000
x1 <- rnorm(N)
x2 <- rnorm(N)
## create linear predictor
lpred <- 0.5 + 0.5 * x1 + 0.25 * x2
## apply inverse link function
ey <- 1/(1 + exp(-lpred))
## sample some dependent variable
y <- rbinom(N, prob=ey, size=rep(1,N))

dat <- matrix(c(x1, x2, y), nrow=N, ncol=3)
colnames(dat) <- c('x1', 'x2', 'y')

Erstellen Sie nun eine Protokollwahrscheinlichkeitsfunktion zur Maximierung, verwenden Sie diese, dbinomweil sie vorhanden ist, und summieren Sie die Ergebnisse

## the log likelihood function
log.like <- function(beta, dat){
  lpred <- beta[1] + dat[,'x1'] * beta[2] + dat[,'x2'] * beta[2]**2
  ey <- 1/(1 + exp(-lpred))
  sum(dbinom(dat[,'y'], prob=ey, size=rep(1,nrow(dat)), log=TRUE))
}

und passen Sie das Modell mit maximaler Wahrscheinlichkeit an. Ich habe mir nicht die Mühe gemacht, einen Farbverlauf anzubieten oder eine Optimierungsmethode zu wählen, aber vielleicht möchten Sie beides tun.

## fit
res <- optim(par=c(1,1), ## starting values 
             fn=log.like,
             control=list(fnscale=-1), ## maximise not minimise
             hessian=TRUE, ## for SEs
             dat=dat)

Schauen Sie sich jetzt die Ergebnisse an. Die ML-Parameterschätzungen und asymptotischen SEs sind:

## results
data.frame(coef=res$par,
           SE=sqrt(diag(solve(-res$hessian))))

was sein sollte

##        coef         SE
## 1 0.4731680 0.04828779
## 2 0.5799311 0.03363505

oder es gibt einen Fehler (der immer möglich ist).

Es gelten die üblichen Vorbehalte gegen von Hessen abgeleitete Standardfehler.