Os valores que você obtém dos BLUPs não são estimados da mesma maneira que as estimativas BLUE de efeitos fixos; por convenção, os BLUPs são referidos como previsões . Quando você ajusta um modelo de efeitos mistos, o que é estimado inicialmente são a média e a variação (e possivelmente a covariância) dos efeitos aleatórios. O efeito aleatório de uma determinada unidade de estudo (digamos, um aluno) é subsequentemente calculado a partir da média e variância estimadas e dos dados. Em um modelo linear simples, a média é estimada (assim como a variância residual), mas as pontuações observadas são consideradas compostas por isso e pelo erro, que é uma variável aleatória. Em um modelo de efeitos mistos, o efeito para uma determinada unidade também é uma variável aleatória (embora, em certo sentido, ela já tenha sido realizada).
Você também pode tratar essas unidades como efeitos fixos, se quiser. Nesse caso, os parâmetros para essa unidade são estimados como de costume. Nesse caso, no entanto, a média (por exemplo) da população da qual as unidades foram extraídas não é estimada.
Além disso, a suposição por trás dos efeitos aleatórios é que eles foram amostrados aleatoriamente em alguma população, e é a população com a qual você se importa. A suposição subjacente aos efeitos fixos é que você selecionou essas unidades intencionalmente porque essas são as únicas unidades com as quais você se importa.
Se você se virar e ajustar um modelo de efeitos mistos e prever esses mesmos efeitos, eles tendem a ser "encolhidos" em relação à média da população em relação às estimativas de efeitos fixos. Você pode pensar nisso como análogo a uma análise bayesiana em que a média estimada e a variância especificam um prior normal e o BLUP é como a média do posterior resultante da combinação otimizada dos dados com o prior.
A quantidade de encolhimento varia com base em vários fatores. Uma determinação importante de quão longe as previsões de efeitos aleatórios estarão das estimativas de efeitos fixos é a razão entre a variação dos efeitos aleatórios e a variação do erro. Aqui está uma R
demonstração rápida do caso mais simples, com unidades de 5 'nível 2', apenas com ajuste de meios (interceptações). (Você pode pensar nisso como pontuações de teste para alunos das turmas.)
library(lme4) # we'll need to use this package
set.seed(1673) # this makes the example exactly reproducible
nj = 5; ni = 5; g = as.factor(rep(c(1:nj), each=ni))
##### model 1
pop.mean = 16; sigma.g = 1; sigma.e = 5
r.eff1 = rnorm(nj, mean=0, sd=sigma.g)
error = rnorm(nj*ni, mean=0, sd=sigma.e)
y = pop.mean + rep(r.eff1, each=ni) + error
re.mod1 = lmer(y~(1|g))
fe.mod1 = lm(y~0+g)
df1 = data.frame(fe1=coef(fe.mod1), re1=coef(re.mod1)$g)
##### model 2
pop.mean = 16; sigma.g = 5; sigma.e = 5
r.eff2 = rnorm(nj, mean=0, sd=sigma.g)
error = rnorm(nj*ni, mean=0, sd=sigma.e)
y = pop.mean + rep(r.eff2, each=ni) + error
re.mod2 = lmer(y~(1|g))
fe.mod2 = lm(y~0+g)
df2 = data.frame(fe2=coef(fe.mod2), re2=coef(re.mod2)$g)
##### model 3
pop.mean = 16; sigma.g = 5; sigma.e = 1
r.eff3 = rnorm(nj, mean=0, sd=sigma.g)
error = rnorm(nj*ni, mean=0, sd=sigma.e)
y = pop.mean + rep(r.eff3, each=ni) + error
re.mod3 = lmer(y~(1|g))
fe.mod3 = lm(y~0+g)
df3 = data.frame(fe3=coef(fe.mod3), re3=coef(re.mod3)$g)
Portanto, as proporções da variação dos efeitos aleatórios para a variação do erro são 1/5 para model 1
, 5/5 para model 2
e 5/1 para model 3
. Note que usei level significa codificar para os modelos de efeitos fixos. Agora podemos examinar como os efeitos fixos estimados e os efeitos aleatórios previstos se comparam para esses três cenários.
df1
# fe1 re1
# g1 17.88528 15.9897
# g2 18.38737 15.9897
# g3 14.85108 15.9897
# g4 14.92801 15.9897
# g5 13.89675 15.9897
df2
# fe2 re2
# g1 10.979130 11.32997
# g2 13.002723 13.14321
# g3 26.118189 24.89537
# g4 12.109896 12.34319
# g5 9.561495 10.05969
df3
# fe3 re3
# g1 13.08629 13.19965
# g2 16.36932 16.31164
# g3 17.60149 17.47962
# g4 15.51098 15.49802
# g5 13.74309 13.82224
Outra maneira de terminar com previsões de efeitos aleatórios mais próximas das estimativas de efeitos fixos é quando você tem mais dados. Podemos comparar a model 1
partir de cima, com sua baixa proporção de efeitos aleatórios e variação de erro, com uma versão ( model 1b
) com a mesma proporção, mas com muito mais dados (observe isso em ni = 500
vez de ni = 5
).
##### model 1b
nj = 5; ni = 500; g = as.factor(rep(c(1:nj), each=ni))
pop.mean = 16; sigma.g = 1; sigma.e = 5
r.eff1b = rnorm(nj, mean=0, sd=sigma.g)
error = rnorm(nj*ni, mean=0, sd=sigma.e)
y = pop.mean + rep(r.eff1b, each=ni) + error
re.mod1b = lmer(y~(1|g))
fe.mod1b = lm(y~0+g)
df1b = data.frame(fe1b=coef(fe.mod1b), re1b=coef(re.mod1b)$g)
Aqui estão os efeitos:
df1
# fe1 re1
# g1 17.88528 15.9897
# g2 18.38737 15.9897
# g3 14.85108 15.9897
# g4 14.92801 15.9897
# g5 13.89675 15.9897
df1b
# fe1b re1b
# g1 15.29064 15.29543
# g2 14.05557 14.08403
# g3 13.97053 14.00061
# g4 16.94697 16.92004
# g5 17.44085 17.40445
Em uma nota algo relacionado, Doug Bates (o autor do lme4 pacote R) não gosta do termo "BLUP" e usa "modo condicional" em vez (ver pp. 22-23 do seu projecto de lme4 livro pdf ). Em particular, ele aponta na seção 1.6 que "BLUP" pode ser usado apenas de maneira significativa para modelos lineares de efeitos mistos.