R : mle2의 프로필 신뢰 구간

Question

bbmle 패키지에 mle2 명령을 사용하려고합니다. Bolker는 "bbmle 패키지로 최대 가능성 추정 및 분석"의 p2를보고 있습니다. 어떻게 든 올바른 시작 값을 입력하지 못했습니다. 여기에 재현 코드입니다 :

l.lik.probit <-function(par, ivs, dv){ 
Y <- as.matrix(dv) 
X <- as.matrix(ivs) 
K <-ncol(X) 
b <- as.matrix(par[1:K]) 
phi <- pnorm(X %*% b) 
sum(Y * log(phi) + (1 - Y) * log(1 - phi)) 
} 

n=200 

set.seed(1000) 

x1 <- rnorm(n) 
x2 <- rnorm(n) 
x3 <- rnorm(n) 
x4 <- rnorm(n) 

latentz<- 1 + 2.0 * x1 + 3.0 * x2 + 5.0 * x3 + 8.0 * x4 + rnorm(n,0,5) 

y <- latentz 
y[latentz < 1] <- 0 
y[latentz >=1] <- 1 
x <- cbind(1,x1,x2,x3,x4) 
values.start <-c(1,1,1,1,1) 

foo2<-mle2(l.lik.probit, start=list(dv=0,ivs=values.start),method="BFGS",optimizer="optim", data=list(Y=y,X=x)) 

는 그리고이 오류가 나는 얻을 수있다 :

Error in mle2(l.lik.probit, start = list(Y = 0, X = values.start), method = "BFGS", : 
    some named arguments in 'start' are not arguments to the specified log-likelihood function 

어떤 생각을 왜? 당신의 도움을 주셔서 감사합니다!

Answer 1

당신은 몇 가지를 놓쳤지만, 가장 중요한 것은 기본적으로 mle2이 목록의 매개 변수를 취합니다. 대신 벡터 매개 변수를 사용할 수 있지만 조금 더 세게 작업해야합니다.

코드를 약간 수정했습니다. (나는이 작동하지 않을 것 없이는 음의 로그 - 우도 함수에 로그 - 우도 함수를 변경!)

l.lik.probit <-function(par, ivs, dv){ 
    K <- ncol(ivs) 
    b <- as.matrix(par[1:K]) 
    phi <- pnorm(ivs %*% b) 
    -sum(dv * log(phi) + (1 - dv) * log(1 - phi)) 
} 

n <- 200 

set.seed(1000) 

dat <- data.frame(x1=rnorm(n), 
        x2=rnorm(n), 
        x3=rnorm(n), 
        x4=rnorm(n)) 

beta <- c(1,2,3,5,8) 
mm <- model.matrix(~x1+x2+x3+x4,data=dat) 
latentz<- rnorm(n,mean=mm%*%beta,sd=5) 

y <- latentz 
y[latentz < 1] <- 0 
y[latentz >=1] <- 1 
x <- mm 
values.start <- rep(1,5) 

지금 우리가 적합을한다. 중요한 것은

library("bbmle") 
names(values.start) <- parnames(l.lik.probit) <- paste0("b",0:4) 
m1 <- mle2(l.lik.probit, start=values.start, 
      vecpar=TRUE, 
      method="BFGS",optimizer="optim", 
      data=list(dv=y,ivs=x)) 

으로

이 특정 예를 들어 그냥 프로 빗를 다시 구현 한 위의 지적 ... vecpar=TRUE를 지정하고 mle2 매개 변수 벡터의 요소의 이름을 알려 parnames을 사용하는 것입니다 회귀 (나는 ... 당신은 지금 어떤 방법이 분산을 허용하려면이 옵션을 확장하려는 것으로 알고 있지만) 최종 참고로

dat2 <- data.frame(dat,y) 
m2 <- glm(y~x1+x2+x3+x4,family=binomial(link="probit"), 
    data=dat2) 

, 나는 당신을 허용하는의 parameters 인수를 확인해야한다고 말할 것 임의의 하나의 파라미터에 대한 서브 - 선형 모델을 특정하기 위해, 상기 formula 인터페이스 :

m3 <- mle2(y~dbinom(prob=pnorm(eta),size=1), 
      parameters=list(eta~x1+x2+x3+x4), 
      start=list(eta=0), 
      data=dat2) 

PS confint(foo2)이 셋업과 (요청한 프로파일 CI를주는) 잘 작동 할 것으로 보인다.

ae <- function(x,y) all.equal(unname(coef(x)),unname(coef(y)),tol=5e-5) 
ae(m1,m2) && ae(m2,m3)