Mixed models dragon long

example
Dragons
Interested in how smart dragons are, and the relationship with body length.

Data in data folder on Vula.

Exploratory data analysis
In practice
Should standardise variables for LMM

- Helps with model convergence

- Can help with interpretation (it means estimated coefficients all on same scale
- can compare effect sizes)

R: scale()

Stata: ? do by hand for each variable

Fit linear regression (ignore grouping)
Big dragons are smarter.

Ok (?)
Hmmm
Hmmmm...
More EDA (should have done this first..)
Stratify?
Treat group as discrete covariate?
Body length not significant.

Good (?)
Interaction model?
Mixed effects model
- Account for correlation (grouping) with fewer parameters
- Avoid problems with multiple comparisons
- Uses all the data available (no reduction in power due to stratification)
Data structure
Mountain range (we’ve seen) and Site (sites are within mountain ranges)
Fixed effects
We are interested in the smartness ~ length association, so we treat body length
as a fixed effect. Fixed because we want to estimate that association and make
inference on it.

Random effects (always categorical) are grouping factors which we are trying to
adjust for. We want to adjust for mountain range (dragons on same mountain
might be exposed to same ‘stuff’ and hence more alike than dragons on another
mountain) and also might want to account for site.
Random effect of mountain range
R:
R code
mixed.lmer <- lmer(testScore ~ bodyLength2 + (1|mountainRange), data = dragons)

plot(mixed.lmer)

qqnorm(resid(mixed.lmer))
qqline(resid(mixed.lmer))
Stata code
xtmixed testScore bodyLength2 || mountainRange

predict Y, fitted

predict R, rstudent

twoway connected R Y

qnorm R
Variance components
Mountain:

339.7 / (339.7 + 223.8)

= 0.60

60% of variance in score (after

accounting for fixed effects) is
attributable to mountain range

Every unit increase in (scaled) body

length accounts for 0.5 unit
increase in test score.
Nested, crossed, and other effects
Random effects can be crossed or nested.

Nested means what is sounds like (and what we have in this case) factors nested
without other factors. Eg. teachers nested within schools

“Crossed” random effects in this example would be the case where dragons were
not so territorial, and so they could also be measured on different mountain
ranges.

Can differentiate between fully crossed and partially crossed

 Site is nested within mountain range
Which means site and mountain
range are not ‘independent’ effects

This model is WRONG

It is fitting crossed effects when we

want to fit nested effects ***

*** coding matters

Code: multiple random effects
mixed.WRONG <- lmer(testScore ~ bodyLength2 + (1|mountainRange) + (1|site), data = dragons)

summary(mixed.WRONG)

xtmixed testScore bodyLength2 || mountainRange : site

This model is correct

- New variable
- Implicit vs explicit
nesting
- Same code (!?!)
Cross table of mountain x site
‘a’ , ‘b’, and ‘c’ are meaningless on own - they only uniquely
identify a site within a mountain range

This is implicit nesting [coding/data set structure]

New variable ‘sample’

It uniquely identifies site within mountain range

This is explicit nesting [coding/data structure]

This is recommended, helps avoid confusion. Especially with R. Because R

functions allow a lot of flexibility with model specification.
Variance components
Mountain range:

327.56 /
(23.09+327.56+208.58)

= 0.586

Site (within range):

23.09 / (23.09+327.56+208.58)

= 0.041
Body length not significant (?).

Mountain range determines cleverness.

More code
mixed.lmer2 <- lmer(testScore ~ bodyLength2 + (1|mountainRange) + (1|sample), data =
dragons) summary(mixed.lmer2)

R would also allow:

(1|mountainRange/site) or even (1|mountainRange) + (1|mountainRange:site)

xtmixed testScore bodyLength2 || mountainRange : sample

Generic strategy for fitting mixed effects models
Fit the full model (most complex, most complete)

Sort out the random effects structure (do you need nested/crossed/other random
effects - ie if a random effect is accounting for very little variance - why keep it)?

Sort out the fixed effects (are all the necessary confounders and effects accounted
for and in correct form) [keep random effects constant]

Finished.

Yes, look at diagnostics (for the residuals), plot your predictions, plot your random
effects. Is your model useful?
Tutorial seized from: https://gkhajduk.github.io/2017-03-09-mixed-models/