khoảng tin cậy hồ sơ trong R: mle2

Question

Tôi đang cố sử dụng lệnh mle2, trong gói bbmle. Tôi đang xem p2 của "Ước tính khả năng tối đa và phân tích với gói bbmle" của Bolker. Bằng cách nào đó tôi không nhập được các giá trị bắt đầu đúng. Dưới đây là mã tái sản xuất:

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)) 

Và đây là lỗi tôi nhận được:

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 

Bất cứ ý tưởng tại sao? Cảm ơn bạn đã giúp đỡ!

Answer 1

Bạn đã bỏ lỡ một vài điều, nhưng điều quan trọng nhất là theo mặc định mle2 có một danh sách của thông số; bạn có thể làm cho nó tham số vector thay vào đó, nhưng bạn phải làm việc chăm chỉ hơn một chút.

Tôi đã chỉnh sửa mã hơi ở những nơi. (Tôi đã thay đổi chức năng log-likelihood thành một hàm log-khả năng tiêu cực, mà không có cái này sẽ không bao giờ hoạt động!)

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) 

Bây giờ chúng ta làm vừa vặn. Điều quan trọng là phải xác định vecpar=TRUE và sử dụng parnames để cho mle2 biết tên của các yếu tố trong vector tham số ...

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)) 

Như đã chỉ ra ở trên cho ví dụ cụ thể này, bạn vừa tái triển khai các probit hồi quy (mặc dù tôi hiểu rằng bây giờ bạn muốn mở rộng này cho phép các biến ngẫu nhiên một cách nào đó ...)

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

là một lưu ý cuối cùng, tôi sẽ nói rằng bạn nên kiểm tra đối số parameters, cho phép bạn để chỉ định mô hình tuyến tính phụ cho bất kỳ tham số nào, và giao diện 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) dường như làm việc tốt (cho TCTD hồ sơ theo yêu cầu) với điều này thiết lập.

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