Cómo formular lineal modelo mixto para averiguar los efectos de las variables continuas?

Question

Tengo un conjunto de datos con tasa de crecimiento como variable de respuesta (resp en el ejemplo) y de la temperatura, la disponibilidad de alimentos, y de la salinidad como variables predictoras (pred1 través pred3 en el ejemplo). Las variables predictoras son "continua" con intervalos semanales y tienen diferentes unidades. Las mediciones span semanal (con valores perdidos para algunos ejemplos) a lo largo de un año (week en el ejemplo; definido desde el principio del experimento). Tengo varias muestras y quiero cuantificar (por encima de todas las muestras):

1. Cuánto de cada variable predictora explica la variación en la tasa de crecimiento
2. El efecto relativo de cada variable de predicción en la tasa de crecimiento

Entiendo que los lineales de los modelos mixtos podría ser una solución para este problema ya que tengo varias muestras y dependiente de las mediciones a lo largo del tiempo. Mi pregunta es: ¿Cuál sería el modelo óptimo de formulaciones que utilizan lme4 paquete de R?

Ejemplo de datos está disponible aquí. Y aquí está un resumen de la misma:

library(ggplot2)
tmp <- melt(X, id = c("Sample", "weeks"))
ggplot(tmp, aes(x = weeks, y = value)) + geom_line() + facet_wrap(Sample ~ variable, scales = "free_y")

He intentado siguiente:

Como una solución para el punto 1:

library("lme4")
library("MuMIn")

p1 <- lmer(resp ~ pred1 + (1|Sample) + (1|weeks), data = X)
p2 <- lmer(resp ~ pred2 + (1|Sample) + (1|weeks), data = X)
p3 <- lmer(resp ~ pred3 + (1|Sample) + (1|weeks), data = X)

margr2 <- data.frame(Pred = c("pred1", "pred2", "pred3"), marginal.R2 = c(r.squaredGLMM(p1)[[1]], r.squaredGLMM(p2)[[1]], r.squaredGLMM(p3)[[1]]))

ggplot(margr2, aes(x = Pred, y = marginal.R2)) + geom_bar(stat = "identity")

Marginal $R^2$ calculado por el método publicado aquí debe indicar el total de varianza explicada por cada variable predictora hasta donde tengo entendido, y suponiendo que mi modelo formulaciones son correctos.

Para el efecto relativo (punto 2), creo que primero hay que tener las variables predictoras en una misma escala. Sólo entonces puedo comparar por tener todos ellos en el modelo y la eliminación de los interceptos:

Xs <- X
Xs[4:6] <- scale(Xs[4:6])

mod <- lmer(resp ~ pred1 + pred2 + pred3 - 1 + (1|weeks) + (1|Sample), data = Xs)
cis <- confint(mod)[4:6,]

releff <- data.frame(par = rownames(cis), lower = cis[,1], est = fixef(mod), upper = cis[,2])

Con el fin de hacer la interpretación más intuitiva, yo escala de los efectos de máximo valor absoluto de a través de los intervalos de confianza (yo sólo estoy interesado en efecto relativo):

tmp <- c(releff$lower,releff$upper)

add <- 100*releff[c("lower", "est", "upper")]/max(abs(tmp))
colnames(add) <- paste0("rel.", colnames(add))

releff <- cbind(releff, add)

ggplot(releff, aes(x = par, y = rel.est, ymin = rel.lower, ymax = rel.upper)) + geom_pointrange() + geom_hline(yintercept = 0)

Las variables predictoras son "significativos", donde el CIs no cruce la línea horizontal (a mi entender). No estoy seguro de si estos enfoques tienen mucho sentido, y por eso estoy pidiendo ayuda.

Answer 1

Una idea para la mejora de la marginal $R^2$ cálculo utilizada será la de evaluar este con el resto de los predictores incluidos en el modelo. Como se destaca aquí, el marginal $R^2$ cálculo sólo tiene en cuenta un predictor en un momento.

Una alternativa es colocar dos modelos. Un modelo contiene todos los predictores, el otro tiene una predictor caído. Los modelos pueden ser comparadas para ver la disminución en el marginal $R^2$ que es debido a la eliminación de la predictor. Por ejemplo:

m1 <- lmer(resp ~ pred1 + pred2 + pred3 + (1|weeks) + (1|Sample), data = Xs)
m2 <- lmer(resp ~ pred2 + pred3 + (1|weeks) + (1|Sample), data = Xs)

r.squaredGLMM(m1)[[1]]-r.squaredGLMM(m2)[[1]]

Esto indica que el marginal $R^2$ cae un poco, simplemente quitando el primer predictor. Este ecos de su enfoque, pero tiene el beneficio añadido de todos los predictores en el modelo que se utiliza para calcular la bondad de ajuste.

Con respecto a la construcción de un modelo adecuado, ¿por qué han quitado la intersección? Esta es una pieza clave de información. Cuando lo que se está forzando el modelo a pasar por el origen. Específicamente, va a aplicar que cuando los predictores tomar valores de 0 de la predicción de la respuesta debe ser 0. Estoy sospechando de que esto probablemente no es lo que usted desea.

Ya que usted dijo que usted está interesado en los efectos relativos de los factores, la estandarización de sus predictores, como lo han hecho es una buena idea.

Una alternativa es ajustar un modelo con predictores escalados como este:

m3 <- lmer(resp ~ pred1 + pred2 + pred3 + (1|weeks) + (1|Sample), data = X)

Desde que se ha estandarizado el predictores de la estimación de $\beta$s representan el efecto relativo de los predictores del resultado $resp$.

Para comprobar si estas relaciones son propensos a ser cierto, no sólo en la muestra, sino también en la población, un enfoque sensato es llevar a cabo comparaciones de modelos como el cociente de probabilidad de las pruebas, AIC o BIC.

La manera de hacer esto es para quitar los predictores en una forma gradual y comparar los dos modelos con el método de comparación de elección. Si tal comparación revela que un predictor no contribuyan significativamente a una mejora en el ajuste global, a continuación, puede quitar este predictor de su modelo y también considerar la posibilidad de informar de que no parece haber ninguna relación entre ese predictor y su resultado. Mucha información en este sitio para hacer la comparación de modelos.