Playground: Chapter 10

Mixed-Effects Models

This page practises the ideas of Chapter 10 in four parts. Practise the chapter works through the book’s exercises, with starter code, hints, and solutions. Go further goes beyond the book, with new questions and some of R’s built-in datasets. Check your understanding asks short questions whose answers open with a click. Do it yourself gives open tasks with an empty code space and no answers: you choose and check the models yourself, as you will in your own thesis.

The exercises use the lme4 package, which is large, so the first run on this page can take a minute while your browser downloads it. The data frame panel is prepared as in the chapter: one row per student per semester, with the background variables and time counted from 0 (semester 1) to 3 (semester 4). In the exercises, replace each ______ with your own code and press Run Code.

Practise the chapter

These are the exercises at the end of Chapter 10, with the same numbers.

Exercise 1: Sleep over the semesters

Fit a random intercept model for sleep_hours over time. Calculate the ICC, and describe whether sleep changes over the semesters.

Solution

sleep_ri <- lmer(sleep_hours ~ time + (1 | student_id), data = panel)
vc <- as.data.frame(VarCorr(sleep_ri))
vc$vcov[1] / sum(vc$vcov)
confint(sleep_ri, parm = "beta_", method = "Wald")

The ICC is 0.83: 83% of the variation in sleep is between students, so students differ a great deal in how much they sleep, but each sleeps about the same from one semester to the next. The slope of time is almost exactly zero (0.0005 hours per semester), and its confidence interval, from −0.015 to 0.016, is narrow and centred on zero: sleep does not change over the two years, and the interval rules out any change larger than about a minute per semester.

Exercise 2: Is a random slope worth keeping?

Add a random slope of time to the model in Exercise 1, and compare the two models with anova(). Decide whether the random slope is worth keeping.

Solution

sleep_rs <- lmer(sleep_hours ~ time + (time | student_id), data = panel)
anova(sleep_ri, sleep_rs)
VarCorr(sleep_rs)

The likelihood ratio test favours the random slope (\(\chi^2\) = 7.8, 2 df, p = 0.02), and its AIC is slightly lower. But VarCorr() shows that the slopes hardly vary: their standard deviation is 0.05 hours (3 minutes) per semester, against 0.95 hours between students’ overall levels. Statistically detectable, practically unimportant. Either model is defensible; for a thesis, the simpler random intercept model with a sentence noting the test would be a reasonable choice.

Exercise 3: Do part-time students change differently?

Fit a model of wellbeing over time with study mode and their interaction, with random intercepts and slopes, and describe whether part-time students’ wellbeing changes at a different rate.

Solution

wb_mode <- lmer(wellbeing ~ time * study_mode + (time | student_id), data = panel)

Full-time students start at 62.8 and decline by 0.47 points per semester. Part-time students start 3.9 points lower (95% CI −5.9 to −1.9), but their slope differs by only 0.14 points per semester, with an interval from −0.34 to 0.61 that includes zero. Part-time students’ wellbeing is lower from the start, and the gap stays about the same over the two years: there is no evidence that it changes at a different rate.

Exercise 4: Supervisors and study hours

Fit a model with random intercepts for supervisors and students for study_hours, and calculate how much of its variation lies between supervisors.

Solution

study_levels <- lmer(study_hours ~ time + (1 | supervisor_id) + (1 | student_id),
                     data = panel)

The share between supervisors is 0, and R reports a singular fit: the model found no differences between supervisors at all, beyond what differences between their students explain. About 91% of the variation lies between students and 9% within students from semester to semester. Supervisors matter for wellbeing, as the chapter showed, but not for how much students study. The singular fit is the model’s way of saying that the supervisor level can be dropped.

Exercise 5: Supervisors who differ more

Change the simulation at the start of the chapter so that supervisors differ more (sd = 10 for style), and describe what happens to the false alarm rate of the ordinary regression. The page repeats the study 200 times rather than 500, to keep the time in the browser reasonable; it still takes a minute.

Solution

set.seed(12)
rowMeans(replicate(200, simulate_clustered(style_sd = 10)))

The ordinary regression now declares the useless course significant in about 46% of the studies, against about 27% in the chapter; the mixed model stays near 8%. When supervisors differ more, students of the same supervisor resemble each other more (the ICC rises), so the 200 students carry even less independent information: the study is really about 20 supervisors. The ordinary regression still counts 200 independent students, so its standard error is too small by more, and the difference between two sets of 10 supervisors looks like an effect of the course nearly half the time.

Exercise 6: Explaining it to a fellow student

In your own words, explain to a fellow student why an ordinary regression of all the semester records would be wrong. Write your answer first, then open the model answer.

Model answer

“Each student appears up to four times, and a student’s four records are much more alike than records of four different students: someone who sleeps little in semester 1 usually sleeps little in semester 4 too. An ordinary regression treats the 2,326 records as 2,326 independent pieces of evidence, but they are really about 600 people, and together they carry only about as much information as 700 or so independent records. So the ordinary regression thinks it knows more than it does: its standard errors are too small, and its p-values too small, and it will find ‘effects’ that are just differences between a few students. A mixed model knows which records belong to the same student and gets the uncertainty right.”

Go further

These exercises go beyond the book.

Exercise 7: How much information is in the records?

The chapter calculated the design effect for wellbeing. Calculate it for GPA and for sleep, and the number of independent records each set of records is worth.

Solution

design_effect <- 1 + (records_per_student - 1) * icc

With 3.9 records per student, GPA has an ICC of 0.71 and a design effect of 3.0: its 2,325 records are worth about 767 independent observations. Sleep has an ICC of 0.83 and a design effect of 3.4, so its records are worth only about 679. The more stable a variable is within each student, the more its repeated records repeat each other, and the less each new semester adds. For sleep, which barely changes, measuring four times gives little more information than measuring once.

Exercise 8: Averaging first, or a mixed model

One simple way to avoid repeated records is to average each student’s records and run an ordinary regression on the 600 averages. Compare this with a mixed model, and with the wrong analysis of all records, for the difference in wellbeing between part-time and full-time students.

Solution

averages <- panel |>
  summarise(wellbeing = mean(wellbeing), .by = c(student_id, study_mode))

The averages and the mixed model agree almost perfectly: part-time students are 3.7 points lower, with a standard error of 0.98 in both. The ordinary regression of all 2,326 records gives a similar estimate (3.8) but a standard error of only 0.54, far too small, because it counts every record as independent. For a variable that does not change within a student, such as study mode, averaging is a perfectly good solution. It cannot answer questions about change over time, though, and it gives students with one record the same weight as students with four; that is where the mixed model is needed.

Exercise 9: How long a study needs to run

lme4’s sleepstudy data follows 18 people over 10 days of sleep restriction. Fit the random slope model to the first five days only (days 0 to 4) and compare it with the model for all ten days.

Solution

early <- lmer(Reaction ~ Days + (Days | Subject),
              data = sleepstudy |> filter(Days <= 4))

With five days, the average slope is 8.2 ms per day (SE 2.4), against 10.5 (SE 1.5) with ten days: a smaller and much less precise estimate from half the data. The residual standard deviation also shrinks, from 25.6 to 14.9 ms, which suggests that reaction times become more erratic in the later days. A short study would have described a milder effect and missed how people deteriorate as sleep loss accumulates. A slope describes only the period observed, and extrapolating it beyond that period is risky.

Exercise 10: Growth curves for chicks

R’s built-in ChickWeight data records the weight of 50 chicks, measured every few days from birth to day 21. Fit a model in which every chick has its own starting weight and growth rate, and describe how much the chicks differ in growth.

Solution

chick_model <- lmer(weight ~ Time + (Time | Chick), data = ChickWeight)

The average chick gains 8.5 grams a day, and the standard deviation of the growth rates is 3.8 grams a day: a typical chick grows about 45% faster or slower than average, so chicks differ enormously in growth. The model’s average intercept, 29 grams, is below the real birth weight of about 41 grams, and the intercepts and slopes correlate at −0.95. Both are signs that a straight line is too simple: chicks grow faster as they get bigger, so the line starts too low. A curved growth model, with a squared term for time, would describe them better.

Check your understanding

Answer each question in your own words first, then click to see a model answer.

1. Why are repeated records of the same student not independent?

Answer

Because each student brings stable characteristics to every record: their habits, circumstances, personality, and the way they answer questions. Two records of the same student are therefore more alike than records of two different students, and the second record partly repeats the information in the first.

2. What does the intraclass correlation (ICC) measure?

Answer

The share of all the variation that lies between groups, such as students or supervisors, rather than within them. An ICC of 0 means that the grouping does not matter; an ICC near 1 means that members of a group are almost identical. The higher the ICC, the more wrong an analysis that ignores the grouping would be.

3. What is the difference between a fixed effect and a random effect?

Answer

A fixed effect is a quantity of direct interest, estimated as one number: the average effect of time, or the difference between study modes. A random effect describes how groups vary around the average, estimated as a standard deviation: how much students’ starting levels or growth rates differ. The individual students are not of interest in themselves; they are a sample from a population whose variation the random effect describes.

4. What does a singular fit mean, and what should you do about it?

Answer

That the model’s random effects are more complex than the data can support: typically, a variance estimated as zero, or a correlation of exactly ±1. It does not mean that the model is wrong in its fixed effects. The usual response is to simplify the random effects, for example by dropping a random slope or a level that shows no variation, and to report the simplification.

5. In the chapter’s simulation, why did the ordinary regression produce so many false alarms?

Answer

The course was given to supervisors, and there were only 20 of them. Students of the same supervisor shared the supervisor’s style, so the 200 students were far from 200 independent pieces of evidence. The ordinary regression counted them as independent, underestimated the standard error, and mistook chance differences between two small groups of supervisors for an effect of the course.

Do it yourself

These tasks have no starter code and no answers. Each code space below is empty and ready to run: write your own code, as you would for a thesis. The data (students, semesters, and panel) and the dplyr and lme4 packages are already loaded, and R’s built-in datasets are always available.

1. Model the change in exercise_days over the four semesters. Choose the random effects with a likelihood ratio test, calculate the ICC, and write two sentences reporting the result.

2. Test whether the effect of the workshop invitation on wellbeing over time differs between Master’s and PhD students, with a three-way interaction of time, workshop, and programme. Describe the result with predicted wellbeing for the four groups in semester 2.

3. Write the paragraph reporting a mixed model for GPA over the semesters, like the chapter’s reporting paragraph for wellbeing, including the random effects and the number of observations and students, with every number calculated by code.

4. R’s built-in Orange data records the circumference of five orange trees at seven ages. Fit a model in which every tree has its own growth rate, and explain why estimating the spread of the trees’ intercepts and slopes from only five trees is uncertain.

Work on your own computer

Download the Chapter 10 project

The project contains the data and all four parts of this page as an R script: the exercises with blanks, the questions to check your understanding (answers in solutions.R), and the open tasks, each with space to write your code. On your own computer, the simulation in Exercise 5 can use the chapter’s 500 repetitions.

Download chapter10.zip, unzip it, and double-click chapter10.Rproj.
Or type this one line in RStudio’s Console:

usethis::use_course("https://polla-fattah.github.io/data2thesis_r/playground/chapter10.zip")