Como interpretar coeficientes de um modelo misto multivariado no lme4 sem interceptação geral?

10

Estou tentando ajustar um modelo misto multivariado (isto é, resposta múltipla) R. Além dos pacotes ASReml-re SabreR(que requerem software externo), parece que isso só é possível no MCMCglmm. No artigo que acompanha o MCMCglmmpacote (pp.6), Jarrod Hadfield descreve o processo de ajustar um modelo como remodelar várias variáveis ​​de resposta em uma variável de formato longo e suprimir a interceptação geral. Meu entendimento é que suprimir a interceptação altera a interpretação do coeficiente para cada nível da variável de resposta como a média para esse nível. Dado o exposto, é, portanto, possível ajustar um modelo misto multivariado usando lme4? Por exemplo:

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

Como alguém interpretaria os coeficientes desse modelo? Esse método também funcionaria para modelos mistos lineares generalizados?

Chris
fonte

Respostas:

2

Sua ideia é boa, mas no seu exemplo, você esqueceu de modelar diferentes interceptações e diferentes variações aleatórias para cada característica, para que sua saída não seja interpretável como é. Um modelo correto seria:

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

Nesse caso, você obteria as estimativas de efeitos fixos em cada variável (por exemplo, variabledrat:gearé o efeito do preditor gearna resposta drat), mas também receberia as interceptações para cada variável (por exemplo, variabledratpara a interceptação de resposta drat) e a aleatória variância de cada variável e as correlações entre variáveis:

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       

Uma descrição mais detalhada desses métodos foi escrita por Ben Bolker , bem como o uso MCMCglmmem uma estrutura bayesiana. Outro novo pacote mcglmtambém pode lidar com modelos multivariados, mesmo com respostas não normais, mas você precisa codificar suas matrizes de design aleatórias. Um tutorial deve estar disponível em breve (consulte a página de ajuda do R).

Charlotte R
fonte