Por que a introdução de um efeito aleatório de inclinação aumentou o SE da inclinação?

Estou tentando analisar o efeito de Year na variável logInd para um grupo específico de indivíduos (eu tenho 3 grupos). O modelo mais simples:

> fix1 = lm(logInd ~ 0 + Group + Year:Group, data = mydata)
> summary(fix1)

Call:
lm(formula = logInd ~ 0 + Group + Year:Group, data = mydata)

Residuals:
    Min      1Q  Median      3Q     Max 
-5.5835 -0.3543 -0.0024  0.3944  4.7294 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
Group1       4.6395740  0.0466217  99.515  < 2e-16 ***
Group2       4.8094268  0.0534118  90.044  < 2e-16 ***
Group3       4.5607287  0.0561066  81.287  < 2e-16 ***
Group1:Year -0.0084165  0.0027144  -3.101  0.00195 ** 
Group2:Year  0.0032369  0.0031098   1.041  0.29802    
Group3:Year  0.0006081  0.0032666   0.186  0.85235    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.7926 on 2981 degrees of freedom
Multiple R-squared: 0.9717,     Adjusted R-squared: 0.9716 
F-statistic: 1.705e+04 on 6 and 2981 DF,  p-value: < 2.2e-16

Podemos ver que o Grupo1 está diminuindo significativamente, os Grupos2 e 3 estão aumentando, mas não significativamente.

Claramente, o indivíduo deve ter efeito aleatório, por isso introduzo o efeito de interceptação aleatória para cada indivíduo:

> mix1a = lmer(logInd ~ 0 + Group + Year:Group + (1|Individual), data = mydata)
> summary(mix1a)
Linear mixed model fit by REML 
Formula: logInd ~ 0 + Group + Year:Group + (1 | Individual) 
   Data: mydata 
  AIC  BIC logLik deviance REMLdev
 4727 4775  -2356     4671    4711
Random effects:
 Groups     Name        Variance Std.Dev.
 Individual (Intercept) 0.39357  0.62735 
 Residual               0.24532  0.49530 
Number of obs: 2987, groups: Individual, 103

Fixed effects:
              Estimate Std. Error t value
Group1       4.6395740  0.1010868   45.90
Group2       4.8094268  0.1158095   41.53
Group3       4.5607287  0.1216522   37.49
Group1:Year -0.0084165  0.0016963   -4.96
Group2:Year  0.0032369  0.0019433    1.67
Group3:Year  0.0006081  0.0020414    0.30

Correlation of Fixed Effects:
            Group1 Group2 Group3 Grp1:Y Grp2:Y
Group2       0.000                            
Group3       0.000  0.000                     
Group1:Year -0.252  0.000  0.000              
Group2:Year  0.000 -0.252  0.000  0.000       
Group3:Year  0.000  0.000 -0.252  0.000  0.000

Teve um efeito esperado - o SE das inclinações (coeficientes Grupo 1-3 - Ano) agora é menor e o SE residual também é menor.

Os indivíduos também são diferentes em inclinação, então eu também introduzi o efeito de inclinação aleatória:

> mix1c = lmer(logInd ~ 0 + Group + Year:Group + (1 + Year|Individual), data = mydata)
> summary(mix1c)
Linear mixed model fit by REML 
Formula: logInd ~ 0 + Group + Year:Group + (1 + Year | Individual) 
   Data: mydata 
  AIC  BIC logLik deviance REMLdev
 2941 3001  -1461     2885    2921
Random effects:
 Groups     Name        Variance  Std.Dev. Corr   
 Individual (Intercept) 0.1054790 0.324775        
            Year        0.0017447 0.041769 -0.246 
 Residual               0.1223920 0.349846        
Number of obs: 2987, groups: Individual, 103

Fixed effects:
              Estimate Std. Error t value
Group1       4.6395740  0.0541746   85.64
Group2       4.8094268  0.0620648   77.49
Group3       4.5607287  0.0651960   69.95
Group1:Year -0.0084165  0.0065557   -1.28
Group2:Year  0.0032369  0.0075105    0.43
Group3:Year  0.0006081  0.0078894    0.08

Correlation of Fixed Effects:
            Group1 Group2 Group3 Grp1:Y Grp2:Y
Group2       0.000                            
Group3       0.000  0.000                     
Group1:Year -0.285  0.000  0.000              
Group2:Year  0.000 -0.285  0.000  0.000       
Group3:Year  0.000  0.000 -0.285  0.000  0.000

Mas agora, ao contrário do esperado, o SE das inclinações (coeficientes Grupo1-3: Ano) agora é muito maior, ainda mais alto do que sem efeito aleatório!

Como isso é possível? Eu esperaria que o efeito aleatório "coma" a variabilidade inexplicável e aumente a "segurança" da estimativa!

No entanto, o SE residual se comporta conforme o esperado - é menor do que no modelo de interceptação aleatória.

Aqui estão os dados, se necessário.

Editar

Agora eu percebi um fato surpreendente. Se eu fizer a regressão linear para cada indivíduo separadamente e executar a ANOVA nas inclinações resultantes, obtenho exatamente o mesmo resultado que o modelo de inclinação aleatória! Você saberia o porquê?

indivSlope = c()
for (indiv in 1:103) {
    mod1 = lm(logInd ~ Year, data = mydata[mydata$Individual == indiv,])
    indivSlope[indiv] = coef(mod1)['Year']
}

indivGroup = unique(mydata[,c("Individual", "Group")])[,"Group"]


anova1 = lm(indivSlope ~ 0 + indivGroup)
summary(anova1)

Call:
lm(formula = indivSlope ~ 0 + indivGroup)

Residuals:
      Min        1Q    Median        3Q       Max 
-0.176288 -0.016502  0.004692  0.020316  0.153086 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)
indivGroup1 -0.0084165  0.0065555  -1.284    0.202
indivGroup2  0.0032369  0.0075103   0.431    0.667
indivGroup3  0.0006081  0.0078892   0.077    0.939

Residual standard error: 0.04248 on 100 degrees of freedom
Multiple R-squared: 0.01807,    Adjusted R-squared: -0.01139 
F-statistic: 0.6133 on 3 and 100 DF,  p-value: 0.6079

Aqui estão os dados, se necessário.

r mixed-model lme4-nlme random-effects-model Curioso
fonte

Você precisa de um efeito fixo de ano para ter um efeito fixo ano: interação em grupo. Em geral, você não pode incluir um termo de interação sem incluir também os efeitos principais. Você realmente acha que não há componente fixo para o efeito ano? E, se sim, como poderia haver um ano fixo: interação em grupo?

John

E, por que não há interceptação fixa? Você pode ter ambos, fixo e aleatório.

John

@ John, este modelo é completamente válido. Este é apenas um problema da codificação desejada da variável categórica. Dessa forma, é a interceptação do Grupo , e é a inclinação do Grupo . Se o principal efeito de Ano e a interceptação forem incluídos, as estimativas serão as diferenças da interceptação do Grupo Grupo 1, e da mesma forma com as inclinações. Group

i

$i$

i

$i$ Group

i

$i$ :Year

i

$i$

i

$i$

Aniko

@ John, isso é fora de tópico para a minha pergunta, no entanto: acredite em mim, tudo bem, fiz muitas experiências com isso. Meu primeiro modelo de filme é totalmente equivalente a logInd ~ Year*Group, apenas os coeficientes têm formas diferentes, nada mais. Depende do seu gosto e da forma dos coeficientes que você gosta, nada mais. Não há exclusão de "Efeito principal do ano" no meu primeiro modelo enquanto você escreve ... logInd ~ Year*Groupfaz exatamente o mesmo, o Yearcoeficiente não é o principal efeito, mas o Grupo1: Ano.

Curioso

OK, puro, não tinha considerado a interceptação 0 e o Grupo sendo categóricos.

John

Acho que o problema está nas suas expectativas :) Observe que quando você adiciona uma interceptação aleatória para cada indivíduo, o erro padrão das interceptações aumenta. Como cada indivíduo pode ter sua própria interceptação, a média do grupo é menos certa. O mesmo aconteceu com a inclinação aleatória: você não está mais estimando uma inclinação comum (dentro do grupo), mas a média de inclinações variáveis.

EDIT: Por que um modelo melhor não fornece uma estimativa mais precisa?

Vamos pensar sobre o contrário: por que o modelo inicial subestima o erro padrão? Pressupõe independência de observações que não são independentes. O segundo modelo relaxa essa suposição (de uma maneira que afeta as interceptações) e o terceiro relaxa ainda mais.

EDIT 2: relacionamento com muitos modelos específicos de pacientes

Sua observação é uma propriedade conhecida (e se você tivesse apenas dois anos, o modelo de efeitos aleatórios seria equivalente a um teste t emparelhado). Acho que não consigo gerenciar uma prova real, mas talvez a redação dos dois modelos torne o relacionamento mais claro. Vamos ignorar a variável de agrupamento, pois isso apenas complicaria a notação. Usarei letras gregas para efeitos aleatórios e letras latinas para efeitos fixos.

$i$ $j$

Y_{i j} = a + α_{i} + (b + β_{i}) x_{i j} + ϵ_{i j},

$Y_{ij} = a + \alpha_i + (b+\beta_i)x_{ij} + \epsilon_{ij},$

(α_{i}, β_{i})^{'} \sim N (0, Σ)

$(\alpha_i,\beta_i)'\sim N(0,\Sigma)$

ϵ_{i j} \sim N (0, σ^{2})

$\epsilon_{ij}\sim N(0,\sigma^2)$

Y_{i j} = a_{i} + b_{i} x_{i j} + ϵ_{i j},

$Y_{ij} = a_i + b_i x_{ij}+ \epsilon_{ij},$

ϵ_{i j} \sim N (0, σ_{i}^{2})

$\epsilon_{ij}\sim N(0,\sigma_i^2)$

[Nota: o seguinte é realmente apenas uma mão:]

$a_i$ $a+\alpha_i$ $b_i$ $b+\beta_i$ $b_i$ $b$ $\sigma$ $\sigma_i$ $\alpha_i$

Aniko
fonte

Obrigado Aniko. Você está certo, meus cálculos confirmam isso, mas eu gostaria de ver o porquê ... Parece contra-intuitivo. Melhorei os modelos - introduzindo efeitos aleatórios, descrevi melhor a estrutura do erro. Erro residual confirma - é cada vez mais baixo. Então, com esses modelos melhores e mais precisos, eu esperaria ter uma inclinação mais precisa ... Eu sei que estou errado em algum lugar, por favor me ajude a vê-lo.

Curioso

Obrigado Aniko, esse é um ponto de vista interessante! Estou interessado apenas em declives (Grupo *: Ano), não em interceptar aqui .. então, meu primeiro passo para introduzir o efeito aleatório de itcept relaxou essa suposição de independência e levou à SE mais baixa (da inclinação ..) e depois ao próximo passo provavelmente foi demais (??) e fez o contrário (pior ainda, SE ..) .. talvez eu precise pensar sobre isso, obrigado.

Curioso

Agora também estou surpreso com um fato muito interessante - veja minha edição. Você saberia por que isso?

Curioso

Não acho que a suposição de independência tenha sido muito relaxada! Era errado para começar.

Aniko

Tomas, um modelo "preciso" não significa que as estimativas serão mais precisas. Como exemplo extremo, escolha qualquer modelo sem dados que você desejar, como aquele que prevê que todas as respostas são zero. Este modelo é absolutamente certo em sua estimativa de zero. É, portanto, o mais preciso possível - mas provavelmente também o mais errado possível. Dando um modelo a um escopo maior para ajustar parâmetros, geralmente significa que esses parâmetros são ajustados com menos precisão, não mais. Um modelo melhor, porque pode quantificar a incerteza não capturada por um modelo pior, geralmente apresenta erros padrão maiores.

whuber

Por que a introdução de um efeito aleatório de inclinação aumentou o SE da inclinação?

Mas agora, ao contrário do esperado, o SE das inclinações (coeficientes Grupo1-3: Ano) agora é muito maior, ainda mais alto do que sem efeito aleatório!

Editar

Respostas: