---
options: 
output: slidy_presentation
---
Multi-level Models (ELM Section 8.8)
========================================================
author: Nancy Reid

date: March 4

The junior schools project data
========================================================
```{r}
library(faraway)
data(jsp)
head(jsp)
jspr = jsp[jsp$year == 2,]
```

Plotting
========================================================

```{r, echo=FALSE}
plot(jitter(math) ~ jitter(raven), data = jspr, xlab = "Raven score", ylab = "Math score")
```

=========================================================
```{r, echo=F}
boxplot(math ~ social, data = jspr, xlab = "Social class", ylab = "Math Score")
```
The data
=========================================================
```{r}
names(jspr)
```
"Raven"  is test score on entering (measure of ability?)

"Social" is measure of socio-economic status

 50 schools; up to 4 classes in each school; 
 
linear model
=========================================================
```{r, tidy = TRUE, width = 1000}
glin <- lm(math ~ raven*gender*social, data = jspr)
anova(glin)
```

simpler model
=========================================================
```{r}
summary(lm(math ~ raven + social, data = jspr))
```

random effects model
=========================================================
```{r}
library(lme4)
schoolmod <- lmer(math ~ raven*social*gender + 
                    (1 | school) + (1 | school:class),
                  data = jspr)
anova(schoolmod)
```
```{r, echo=T, eval = T}
jspr$craven <- jspr$raven - mean(jspr$raven)
schoolmod2 <- lmer(math ~ craven*social + (1|school) + (1 | school:class), data = jspr)
```
output
=========================================================
```{r}
summary(schoolmod2)
```

plots
=========================================================
```{r}
par(mfrow=c(1,2))
qqnorm(resid(schoolmod2), main = "")
plot(fitted(schoolmod2), resid(schoolmod2),xlab = "fitted", ylab = "residuals")
```