有没有办法让已经定义的参数在 R 中使用 optim() 函数时显示为丢失?

Is there a way so that an already defined argument can appear as missing when using optim() function in R?

我正在尝试获取 Gumbel 分布的对数似然的最大似然估计量以进行生存分析(我这么说是为了让您不会对对数似然函数感到陌生,我认为它是正确的)。为了做到这一点,我必须通过使用 optim 函数来最大化负对数似然,我尝试这样做但控制台在 fn(par, ...) 中给我一个错误:参数 "b" 是缺失,没有默认值。

我也尝试过以与此 link 的答案类似的方式进行操作:Solve for maximum likelihood with two parameters under constraints,但控制台游戏如下: optim(c(1, 1), function(x) log_lhood(x[1], x[2], d = lung$status, 错误: objective optim 中的函数求值为长度 0 而不是 1。

log_lhood <- function(m,b,d,t){                           
  sum<-0
  for (i in 1:length(lung)){ 
    if (d[i] == 1){
      sum<- sum - log(1-exp(-exp(-(t[i]-m)/b)))
    } else {
      sum<- sum - log((1/b)*exp(-(t[i]-m/b+exp(-(t[i]-m/b)))))
    }
  }
}
#a,b parameter optimization

optim(c(0,0), fn = log_lhood, d = lung$status, t = KM_fit$time) #1st way

optim(c(1, 1), function(x) log_lhood(x[1], x[2],d=lung$status,t=KM_fit_test$time)) #2nd way as in the link

这里有一些问题。该函数的第一个参数应该是一个参数向量。另外你需要 nrow(lung) 而不是 length(lung),最好使用 length(d)。另外你不应该在这里使用循环,它非常低效,使用 ifelse() (在 R 中我们总是尝试向量化一切)。您还需要检查是否可以为所有参数值计算对数似然(例如 b=0)。你也忘了return(sum)。另外 sum 是一个你不应该屏蔽的有用函数。

这会运行。

library(reprex)

lung <- data.frame(status=c(0,0,1,1))
KM_fit <- data.frame(time=c(0,1,2,3))

log_lhood <- function(x,d,t){
    m <- x[1]
    b <- x[2]
    sum <-0
    for (i in 1:nrow(lung)){
        if (d[i] == 1){
            sum <- sum - log(1-exp(-exp(-(t[i]-m)/b)))
        } else {
            sum <- sum - log((1/b)*exp(-(t[i]-m/b+exp(-(t[i]-m/b)))))
        }
    }
    return(sum)
}
#a,b parameter optimization

optim(par=c(0,1), fn = log_lhood, d = lung$status, t = KM_fit$time)
$par
[1] 1.661373 1.811780

$value
[1] 5.318068

$counts
function gradient 
      63       NA 

$convergence
[1] 0

$message
NULL

您可以像这样重写您的函数。

log_lhood <- function(x,d,t){
    m <- x[1]
    b <- x[2]
    s <- ifelse(d==1,
                  -log(1-exp(-exp(-(t-m)/b))),
                  -log((1/b)*exp(-(t-m/b+exp(-(t-m/b)))))
                  )
    return(sum(s, na.rm=TRUE))
}