[R-sig-ME] Why lme4 doesn't throw an error for illegitimate an model

[Simon Harmel]

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 using 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
> >
>

	[[alternative HTML version deleted]]