Why lme4 doesn't throw an error for illegitimate an model
Very interesting Phillip! I'll think about your final question! But, in
your terminology, what is the relationship between "adjustments", "random
effects", and "shrinkage" (of cluster estimates towards the total cluster
mean)?
On the other hand, let's now use an individual-level predictor (i.e.,
"ses": varies within & across clusters) in the random part of the model
WHILE excluding that predictor from the fixed part.
what does `ses` in the output of `VarCorr(m43)` below represent?
hsb <- read.csv('
https://raw.githubusercontent.com/rnorouzian/e/master/hsb.csv')
m43 <- lmer(math ~ female*minority + (ses | sch.id), data = hsb)
VarCorr(m43)
Groups Name Std.Dev. Corr
sch.id (Intercept) 2.1961
ses 2.0963 -0.285
Residual 5.9732
On Wed, Oct 7, 2020 at 6:29 PM Phillip Alday <phillip.alday at mpi.nl> wrote:
Potentially something that is uninteresting or nonsensical or not physically real ... the math for a Klein bottle is well defined, but I doubt you'll find one in the three spatial dimensions we experience. ;) Here it is actually capturing some aspects of the different sectors:
summary(coef(lmer(math ~ (sector | sch.id), data =
hsb))$sch.id[,"sector"]) Min. 1st Qu. Median Mean 3rd Qu. Max. -2.0087 -0.1117 0.4981 0.4908 1.0710 2.8777
summary(lmer(math ~ sector + (1|sch.id), data=hsb))
Linear mixed model fit by REML ['lmerMod']
Formula: math ~ sector + (1 | sch.id)
Data: hsb
REML criterion at convergence: 47080.1
Scaled residuals:
Min 1Q Median 3Q Max
-3.0130 -0.7523 0.0253 0.7602 2.7472
Random effects:
Groups Name Variance Std.Dev.
sch.id (Intercept) 6.677 2.584
Residual 39.151 6.257
Number of obs: 7185, groups: sch.id, 160
Fixed effects:
Estimate Std. Error t value
(Intercept) 11.3930 0.2928 38.907
sector 2.8049 0.4391 6.388
Basically, it's the (shrunken) sector by-school adjustment
from-the-zeroed-out-fixed/population-effect . Since you've omitted
sector from the fixed effects, then the adjustments correspond to
(shrunken) estimates (well, predictions, see previous fine print) of
what the total population + sch.id-level effect would be. Because the
sch.id-level effect is in reality ideally/approximately(*) zero, these
work out to be approximations to the population level effects.
(*) "ideally/approximately" because the variability between schools may
differ between sectors. In other words, the public schools may have
adjustments distributed as N(public_mean, public_sd) and the private
schools may have adjustments distributed as N(private_mean, private_sd),
where the two SDs aren't equal. The means are handled by the fixed
effects, the SDs by the random effects. This is what you get from this
model:
summary(lmer(math ~ sector + (1+sector|sch.id), data=hsb))
Linear mixed model fit by REML ['lmerMod']
Formula: math ~ sector + (1 + sector | sch.id)
Data: hsb
REML criterion at convergence: 47080.1
Scaled residuals:
Min 1Q Median 3Q Max
-3.01392 -0.75219 0.02518 0.76045 2.74806
Random effects:
Groups Name Variance Std.Dev. Corr
sch.id (Intercept) 6.7346 2.5951
sector 0.5322 0.7295 -0.17
Residual 39.1513 6.2571
Number of obs: 7185, groups: sch.id, 160
Fixed effects:
Estimate Std. Error t value
(Intercept) 11.3930 0.2939 38.762
sector 2.8048 0.4387 6.394
summary(coef(lmer(math ~ sector + (sector | sch.id), data =
hsb))$sch.id[,"sector"]) Min. 1st Qu. Median Mean 3rd Qu. Max. 2.531 2.741 2.806 2.805 2.860 3.104 I'll end with a question which I won't answer but it will help you to think about why fitting these models might be useful: how is this related to heteroskedacity? Phillip On 7/10/20 11:59 pm, Simon Harmel wrote:
Also Phillip, what does `sector` in the output of `VarCorr(mn)` below denote, now that you say this model is mathematically defined? mn <- lmer(math ~ ses + (sector | sch.id), data = hsb)
VarCorr(mn)
Groups Name Std.Dev. Corr sch.id <http://sch.id> (Intercept) 2.0256 sector 1.3717 -0.071 Residual 6.0858