No lineal de efectos mixtos de regresión de R

Question

Sorprendentemente, era incapaz de encontrar una respuesta a la siguiente pregunta con Google:

Tengo algunos datos biológicos de varios individuos que muestran una más o menos sigmoidea de crecimiento de comportamiento en el tiempo. Por tanto, quiero el modelo es el uso de un estándar de crecimiento logístico

P(t) = k*p0*exp(r*t) / (k+p0*(exp(r*t)-1))

con p0 siendo el valor inicial en t=0, k siendo el límite asintótico en t->infinito y r es la velocidad de crecimiento. Tan lejos como puedo ver, puedo fácilmente modelo esta usando nls (falta de entendimiento de mi parte: ¿por qué no puedo modelo algo similar utilizando el estándar de la regresión logit mediante la escala de tiempo y de datos? EDIT: Gracias Nick, al parecer la gente hace que, por ejemplo, para proporciones, pero rara vez http://www.stata-journal.com/article.html?article=st0147 . La siguiente pregunta en esta tangente sería si el modelo, posiblemente, puede manejar los datos atípicos >1).

Ahora quiero permitir que algunos fijos (principalmente categórica) y algunos al azar (un individuo de IDENTIFICACIÓN y, posiblemente, también un estudio de IDENTIFICACIÓN) efectos sobre los tres parámetros k, p0 y r. Es nlme la mejor manera de hacer esto? El SSlogis modelo parece sensato pensar que lo que estoy tratando de hacer, es eso correcto? Es cualquiera de las siguientes sensato modelo para empezar? No puedo parecer para obtener los valores de partida de la derecha y update() sólo parece funcionar para los efectos aleatorios, no se fija - en cualquier insinuación?

nlme(y ~ k*p0*exp(r*t) / (k+p0*(exp(r*t)-1)), ## not working at all (bad numerical properties?)
            data = data,
            fixed = k + p0 + r ~ var1 + var2,
            random = k + p0 + r ~ 1|UID,
            start = c(p0=1, k=100, r=1))

nlme(y ~ SSlogis(t, Asym, xmid, scal), ## not working, as start= is inappropriate
            data = data,
            fixed = Asym + xmid + scal ~ var1 + var2, ## works fine with ~ 1
            random = Asym + xmid + scal ~ 1|UID,
            start = getInitial(y ~ SSlogis(Dauer, Asym, xmid, scal), data = data))

Como soy nueva no-lineal de los modelos mixtos, en particular, y los modelos no lineales en general, me gustaría apreciar algunas recomendaciones bibliográficas o enlaces a tutoriales / Faq con preguntas de novato.

Answer 1

Quería compartir algunas de las cosas que he aprendido desde que se hace esta pregunta. nlme parece razonablemente una manera de modelo no lineal de efectos mixtos en R. Comenzar con un simple modelo base:

library(nlme)
data <- groupedData(y ~ t | UID, data=data) ## not strictly necessary
initVals <- getInitial(y ~ SSlogis(t, Asym, xmid, scal), data = data)
baseModel<- nlme(y ~ SSlogis(t, Asym, xmid, scal),
    data = data,
    fixed = list(Asym ~ 1, xmid ~ 1, scal ~ 1),
    random = Asym + xmid + scal ~ 1|UID,
    start = initVals
)

A continuación, utilice actualización para aumentar la complejidad del modelo. El parámetro de inicio es un poco difícil de trabajar, lo que puede tardar algunos retoques para averiguar el orden. Nota cómo el nuevo efecto fijo para var1 en Asym de la siguiente manera regular de efectos fijos para Asim.

 nestedModel <- update(baseModel, fixed=list(Asym ~ var1, xmid ~ 1, scal ~ 1), start = c(fixef(baseModel)[1], 0, fixef(baseModel)[2], fixef(baseModel)[3]))

lme4 parecía más robusto frente a los valores atípicos en mis datos y parecía ofrecer más fiable de la convergencia para los modelos más complejos. Sin embargo, parece que el inconveniente es que la probabilidad de funciones deben ser especificados de forma manual. El siguiente es el modelo de crecimiento logístico con un efecto fijo de var1 (binario) en Asim. Usted puede agregar efectos fijos en xmid y scal en una manera similar. Nota la extraña manera de especificar el modelo utilizando una doble fórmula como resultado ~ de efectos fijos ~ de efectos aleatorios.

library(lme4) ## careful loading nlme and lme4 concurrently
customLogitModel <- function(t, Asym, AsymVar1, xmid, scal) {
    (Asym+AsymVar1*var1)/(1+exp((xmid-t)/scal))
}

customLogitModelGradient <- deriv(
    body(customLogitModel)[[2]], 
    namevec = c("Asym", "AsymVar1", "xmid", "scal"), 
    function.arg=customLogitModel
)

## find starting parameters
initVals <- getInitial(y ~ SSlogis(t, Asym, xmid, scal), data = data)

# Fit the model
model <- nlmer(
    y ~ customLogitModelGradient(t=t, Asym, AsymVar1, xmid, scal, var1=var) ~ 
    # Random effects with a second ~
    (Asym | UID) + (xmid | UID) + (scal | UID), 
    data = data, 
    start = c(Asym=initVals[1], AsymVar1=0, xmid=initVals[2], scal=initVals[3])
)