Loading web-font TeX/Math/Italic

Search

Friday, March 15, 2013

Correlations of fixed effects in linear mixed models

Ever wondered what those correlations are in a linear mixed model? For example:

> (lm.full<-lmer(wear~material-1+(1|Subject), data = BHHshoes))
Linear mixed model fit by REML
Formula: wear ~ material - 1 + (1 | Subject)
Data: BHHshoes
AIC BIC logLik deviance REMLdev
62.9 66.9 -27.5 53.8 54.9
Random effects:
Groups Name Variance Std.Dev.
Subject (Intercept) 6.1009 2.470
Residual 0.0749 0.274
Number of obs: 20, groups: Subject, 10
Fixed effects:
Estimate Std. Error t value
materialA 10.630 0.786 13.5
materialB 11.040 0.786 14.1
Correlation of Fixed Effects:
matrlA
materialB 0.988

The estimated correlation between \hat{\beta}_1 and \hat{\beta}_2 is 0.988.  Note that

\hat{\beta}_1 = (Y_{1,1} + Y_{2,1} + \dots + Y_{10,1})/10=10.360

and 

\hat{\beta}_2 = (Y_{1,2} + Y_{2,2} + \dots + Y_{10,2})/10 = 11.040

From this we can recover the correlation 0.988 as follows:

> b1.vals<-subset(BHHshoes,material=="A")$wear
> b2.vals<-subset(BHHshoes,material=="B")$wear
>
> vcovmatrix<-var(cbind(b1.vals,b2.vals))
>
> covar<-vcovmatrix[1,2]
> sds<-sqrt(diag(vcovmatrix))
> covar/(sds[1]*sds[2])
b1.vals
0.98823

By comparison, in the linear model version of the above:

> summary(lm<-lm(wear~material-1,BHHshoes))
> X<-model.matrix(lm)
> 2.49^2*solve(t(X)%*%X)
materialA materialB
materialA 0.62001 0.00000
materialB 0.00000 0.62001

because Var(\hat{\beta}) = \hat{\sigma}^2 (X^T X)^{-1}.

1 comment:

Thomas said...

Cool!