Проблемы с достижением сходимости в нелинейной модели смешанных эффектов в книге Пиньеро и Бейтса - PullRequest
Новая
       25

Проблемы с достижением сходимости в нелинейной модели смешанных эффектов в книге Пиньеро и Бейтса

0 голосов
/ 26 декабря 2018

Я работаю над книгой Пиньеро и Бейтса Модели со смешанными эффектами в S и S-Plus в R. У меня проблемы с получением сходящейся в главе 8 модели (стр. 387). 

library(nlme)

fm1Wafer.nlmeR <- nlme(current ~ A + B * cos(4.5679 * voltage) + C * sin(4.5679 * voltage),
                       data = Wafer,
                       fixed = list(A ~ voltage + I(voltage^2), B + C ~ 1),
                       random = list(Wafer = A ~ voltage + I(voltage^2),
                                     Site = pdBlocked(list(A ~ 1, A ~ voltage + I(voltage^2)-1))),
                       start = c(-4.26, 5.62, 1.26, -0.10, 0.10), # starting values taken from fixed effects of another model earlier in the book
                       method = "REML",
                       control = nlmeControl(opt = "nlm"))

Как видите, я попробовал оптимизатор nlm.Оптимизатор nlminb по умолчанию также не работает.Оба выдают это сообщение об ошибке 

Error in nlme.formula(current ~ A + B * cos(4.5679 * voltage) + C * sin(4.5679 *  : 
  maximum number of iterations (maxIter = 50) reached without convergence
In addition: Warning messages:
1: In logLik.reStruct(object, conLin) :
  Singular precision matrix in level -2, block 1
2: In logLik.reStruct(object, conLin) :
  Singular precision matrix in level -2, block 1

Есть предложения?Есть несколько моделей, основанных на этом ниже по течению в книге, поэтому было бы неплохо заставить их сходиться.

1 Ответ

0 голосов
/ 26 декабря 2018

Ошибка предполагает, что установка чего-то вроде nlmeControl(opt = "nlm", maxIter = 2000) поможет, но на самом деле это не так.Я попробовал 4000, но это заняло вечность ...

То, что кажется истинным виновником, это method = "REML".Оставив значение по умолчанию "ML", вы получите ожидаемый результат. 

fm1Wafer.nlmeR <- nlme(current ~ A + B * cos(4.5679 * voltage) + C * sin(4.5679 * voltage),
                       data = Wafer,
                       fixed = list(A ~ voltage + I(voltage^2), B + C ~ 1),
                       random = list(Wafer = A ~ voltage + I(voltage^2),
                                     Site = pdBlocked(list(A ~ 1, A ~ voltage + I(voltage^2) - 1))),
                       start = c(-4.26, 5.62, 1.26, -0.10, 0.10))
...