encuentro el otro, la respuesta más viejo aquí insatisfactoria. Yo distingo entre los casos donde, estadísticamente, la edad no tiene impacto y, a la inversa, encontramos un error de cálculo. Personalmente, he cometido errores de carrera al combinar estos dos casos. R ha señalado esto último y me gustaría profundizar en por qué es así.
El modelo que OP ha especificado es un modelo de crecimiento, con pendientes e intersecciones aleatorias. Se incluye una gran intersección, pero no una pendiente de Grand Age. Una restricción desagradable que se impone ajustando una pendiente aleatoria sin la adición de su término "grandioso" es que está forzando que la pendiente aleatoria tenga 0 significa, lo cual es muy difícil de optimizar. Los modelos marginales indican que la edad no tiene un valor diferente estadísticamente significativo de 0 en el modelo. Además, agregar edad como un efecto fijo no soluciona el problema.
> lme(conc~ age, random=~age|Lot, data=IGF)
Error in lme.formula(conc ~ age, random = ~age | Lot, data = IGF) :
nlminb problem, convergence error code = 1
message = iteration limit reached without convergence (10)
Aquí el error es obvio. Puede ser tentador establecer el número de iteraciones. lmeControl
tiene muchas estimaciones iterativas. Pero incluso eso no funciona:
> fit <- lme(conc~ 1, random=~age|Lot, data=IGF,
control = lmeControl(maxIter = 1e8, msMaxIter = 1e8))
Error in lme.formula(conc ~ 1, random = ~age | Lot,
data = IGF, control = lmeControl(maxIter = 1e+08, :
nlminb problem, convergence error code = 1
message = singular convergence (7)
Así que no es una cuestión de precisión, el optimizador se está quedando fuera de los límites.
Debe haber diferencias clave entre los dos modelos que ha propuesto y una forma de diagnosticar el error que ha encontrado. Un enfoque simple es especificar un ajuste "prolija" para el modelo problemática:
> lme(conc~ 1, random=~age|Lot, data=IGF, control = lmeControl(msVerbose = TRUE))
0: 602.96050: 2.63471 4.78706 141.598
1: 602.85855: 3.09182 4.81754 141.597
2: 602.85312: 3.12199 4.97587 141.598
3: 602.83803: 3.23502 4.93514 141.598
(truncated)
48: 602.76219: 6.22172 4.81029 4211.89
49: 602.76217: 6.26814 4.81000 4425.23
50: 602.76216: 6.31630 4.80997 4638.57
50: 602.76216: 6.31630 4.80997 4638.57
El primer término es la REML (creo). Los términos segundo a cuarto son los parámetros de un objeto llamado lmeSt
de la clase lmeStructInt
, lmeStruct
y modelStruct
. Si usa el depurador de Rstudio para inspeccionar los atributos de este objeto (el eje del problema), verá que es el componente de efectos aleatorios que explota aquí.coef(lmeSt)
después de 50 iteraciones produce reStruct.Lot1 reStruct.Lot2 reStruct.Lot3 6.316295 4.809975 4638.570586
como se ve arriba y produce
> coef(lmeSt, unconstrained = FALSE)
reStruct.Lot.var((Intercept)) reStruct.Lot.cov(age,(Intercept))
306382.7 2567534.6
reStruct.Lot.var(age)
21531399.4
que es el mismo que el
Browse[1]> lmeSt$reStruct$Lot
Positive definite matrix structure of class pdLogChol representing
(Intercept) age
(Intercept) 306382.7 2567535
age 2567534.6 21531399
Por lo tanto, está claro que la covarianza de los efectos aleatorios es algo que está explotando aquí para este optimizador particular. Las rutinas PORT en nlminb
han sido criticadas por sus errores no informativos. El texto de David Gay (Bell Labs) está aquí http://ms.mcmaster.ca/~bolker/misc/port.pdf La documentación de PORT sugiere que nuestro error 7 al usar un billon iter max "x puede tener demasiados componentes libres. Ver §5". En lugar de corregir el algoritmo, nos corresponde preguntar si hay resultados aproximados que deberían generar resultados similares. Es, por ejemplo, son fáciles de instalar un objeto lmList
para llegar a la intersección aleatoria y variación de pendiente aleatoria:
> fit <- lmList(conc ~ age | Lot, data=IGF)
> cov(coef(fit))
(Intercept) age
(Intercept) 0.13763699 -0.018609973
age -0.01860997 0.003435819
aunque idealmente deberían ser ponderados por sus respectivos pesos de precisión:
Para utilizar el nlme
paquete observo que la optimización no restringida utilizando BFGS no produce tal error y da resultados similares:
> lme(conc ~ 1, data=IGF, random=~age|Lot, control = lmeControl(opt = 'optim'))
Linear mixed-effects model fit by REML
Data: IGF
Log-restricted-likelihood: -292.9675
Fixed: conc ~ 1
(Intercept)
5.333577
Random effects:
Formula: ~age | Lot
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
(Intercept) 0.9976 (Intr)
age 0.005647296 -0.698
Residual 0.820819785
Number of Observations: 237
Number of Groups: 10
una declaración sintáctica alternativa de un modelo de este tipo se puede hacer con t él mucho más fácil lme4
paquete:
library(lme4)
lmer(conc ~ 1 + (age | Lot), data=IGF)
que produce:
> lmer(conc ~ 1 + (age | Lot), data=IGF)
Linear mixed model fit by REML ['lmerMod']
Formula: conc ~ 1 + (age | Lot)
Data: IGF
REML criterion at convergence: 585.7987
Random effects:
Groups Name Std.Dev. Corr
Lot (Intercept) 0.056254
age 0.006687 -1.00
Residual 0.820609
Number of obs: 237, groups: Lot, 10
Fixed Effects:
(Intercept)
5.331
Un atributo de lmer
y su optimizador es que las correlaciones de efectos aleatorios que están muy cerca de 1, 0 o -1 se fijan simplemente a esos valores, ya que simplifica sustancialmente la optimización (y la eficiencia estadística de la estimación).
Tomados en conjunto, esto hace no sugieren que la edad no tiene un efecto, como se dijo anteriormente, y este argumento puede ser respaldado por los resultados numéricos.
Eché un vistazo rápido a esto y nada salta a la vista. Puede que tengas mejor suerte en la lista de correo 'r-sig-mixed-models', que tiene una concentración mucho mayor de personas familiarizadas con este paquete ... –
¿Has intentado aumentar los límites de iteración en el primer ejemplo? Ver '? LmeControl'. –
Ver respuesta y comentarios a continuación. Su primer modelo no tiene la edad como un efecto fijo, ni las restricciones de efectos aleatorios que tiene el segundo modelo. –