¿Cómo interpretar los coeficientes de un modelo mixto multivariante en lme4 sin intercepto global?

Question

Estoy intentando ajustar un modelo mixto multivariante (es decir, de respuesta múltiple) en R . Aparte del ASReml-r y SabreR paquetes (que requieren software externo), parece que esto sólo es posible en MCMCglmm . En el papel que acompaña al MCMCglmm (pág. 6), Jarrod Hadfield describe el proceso de ajuste de un modelo de este tipo como la remodelación de múltiples variables de respuesta en una variable de formato largo y la posterior supresión del intercepto global. Entiendo que la supresión del intercepto cambia la interpretación del coeficiente de cada nivel de la variable de respuesta para que sea la media de ese nivel. Teniendo en cuenta lo anterior, ¿es posible ajustar un modelo mixto multivariante utilizando lme4 ? Por ejemplo:

data(mtcars)
library(reshape2)
mtcars <- melt(mtcars, measure.vars = c("drat", "mpg", "hp"))
library(lme4)
m1 <- lmer(value ~ -1 + variable:gear + variable:carb + (1 | factor(carb)),
    data = mtcars)
summary(m1)
#  Linear mixed model fit by REML 
#  Formula: value ~ -1 + variable:gear + variable:carb + (1 | factor(carb)) 
#     Data: mtcars 
#   AIC   BIC logLik deviance REMLdev
#   913 933.5 -448.5    920.2     897
#  Random effects:
#   Groups       Name        Variance Std.Dev.
#   factor(carb) (Intercept) 509.89   22.581  
#   Residual                 796.21   28.217  
#  Number of obs: 96, groups: factor(carb), 6
#  
#  Fixed effects:
#                    Estimate Std. Error t value
#  variabledrat:gear  -7.6411     4.4054  -1.734
#  variablempg:gear   -1.2401     4.4054  -0.281
#  variablehp:gear     0.7485     4.4054   0.170
#  variabledrat:carb   5.9783     4.7333   1.263
#  variablempg:carb    3.3779     4.7333   0.714
#  variablehp:carb    43.6594     4.7333   9.224

¿Cómo interpretar los coeficientes de este modelo? ¿Funcionaría también este método para los modelos lineales mixtos generalizados?

Answer 1

Tu idea es buena, pero en tu ejemplo, olvidaste modelar diferentes interceptos y diferentes varianzas aleatorias para cada rasgo, por lo que tu salida no es interpretable tal cual. Un modelo correcto sería:

m1 <- lmer(value ~ -1 + variable + variable:gear + variable:carb + (0 + variable | factor(carb))

En ese caso, obtendría las estimaciones de los efectos fijos de cada variable (por ejemplo, variabledrat:gear es el efecto del predictor gear sobre la respuesta drat ), pero también obtendría los interceptos de cada variable (por ejemplo variabledrat para la intercepción de la respuesta drat ) y la varianza aleatoria de cada variable y las correlaciones entre variables:

Groups       Name         Std.Dev. Corr     
 factor(carb) variabledrat 23.80             
              variablempg  24.27    0.20     
              variablehp   23.80    0.00 0.00
 Residual                  23.80       

Una descripción más detallada de estos métodos ha sido redactada por Ben Bolker así como el uso de MCMCglmm en un marco bayesiano. Otro paquete nuevo, mcglm también puede manejar modelos multivariantes, incluso con respuestas no normales, pero tiene que codificar sus matrices de diseño aleatorio. Pronto estará disponible un tutorial (consulte la página de ayuda de R).