R：錯誤與內爾德 - 米德優化選擇初始值

Question

的數據可以發現here

library(nlme) 
library(dfoptim) 
dat0 <- read.table("aids.dat2",head=T) 
dat1 <- dat0[dat0$day<=90, ] # use only first 90-day data 
dat2 <- dat1[!apply(is.na(dat1),1,any),] # remove missing data 
aids.dat <- groupedData(lgcopy ~ day | patid, data=dat2) 
aids.dat$log10copy = log10(aids.dat$lgcopy) 

myfun2 <- function(arg){ 
    s.p1 <- arg[1] 
    s.b1 <- arg[2] 
    s.p2 <- arg[3] 
    s.b2 <- arg[4] 
    model = nlme(log10copy ~ exp(p1-b1*day) + exp(p2-b2*day), 
       fixed = list(p1 ~ 1, b1 ~ 1, p2 ~ 1, b2 ~ 1), 
       random = list(patid = pdDiag(list(p1 ~ 1, b1 ~ 1, p2 ~ 1, b2 ~ 1))), 
       start = list(fixed = c(p1 = s.p1, b1 = s.b1, p2 = s.p2, b2 = s.b2)), 
       data =aids.dat) 
    return(model$logLik) 
} 

nmkb(fn = myfun2, par = c(10,0.5,6,0.005), lower = c(5, 0.1, 5, 0.001), upper = c(15, 1, 10, 0.1)) 

運行上面的代碼中，我碰到幾個錯誤：

Error in nlme.formula(log10copy ~ exp(p1 - b1 * day) + exp(p2 - b2 * day), : 
    step halving factor reduced below minimum in PNLS step 
In addition: Warning message: 
In nlme.formula(log10copy ~ exp(p1 - b1 * day) + exp(p2 - b2 * day), : 
    Singular precision matrix in level -1, block 1 

的模型擬合精與從par = c(10,0.5,6,0.005)盯着值。不過，我認爲隨着隨機算法開始使用lower = c(5, 0.1, 5, 0.001), upper = c(15, 1, 10, 0.1)中的其他起始值，nlme調用會遇到上述問題，因爲它對起始值非常敏感。結果，nmkb調用從不算什麼。

有沒有辦法規避這種情況？

Answer 1

模型日誌似乎需要最大化，但R中的許多優化過程給出了最小化結果。所以要優化的函數必須是負對數似然。因此，它應該是這樣的：

myfun2 <- function(arg){ 
    s.p1 <- arg[1] 
    s.b1 <- arg[2] 
    s.p2 <- arg[3] 
    s.b2 <- arg[4] 
    model = nlme(log10copy ~ exp(p1-b1*day) + exp(p2-b2*day), 
       fixed = list(p1 ~ 1, b1 ~ 1, p2 ~ 1, b2 ~ 1), 
       random = list(patid = pdDiag(list(p1 ~ 1, b1 ~ 1, p2 ~ 1, b2 ~ 1))), 
       start = list(fixed = c(p1 = s.p1, b1 = s.b1, p2 = s.p2, b2 = s.b2)), 
       data =aids.dat) 
    return(-model$logLik) 
} 

雖然還是有很多的警告，有我的機器上沒有更多的錯誤，並且算法收斂成功。

$par 
[1] 13.460199068 0.848526199 7.764024099 0.001513636 

$value 
[1] -358.6631 

$feval 
[1] 197 

$restarts 
[1] 0 

$convergence 
[1] 0 

$message 
[1] "Successful convergence" 

Warning messages: 
1: In nlme.formula(log10copy ~ exp(p1 - b1 * day) + exp(p2 - b2 * day), : 
    Singular precision matrix in level -1, block 1