Multilevel models

Interpretation + Estimation

Prof. Maria Tackett

Feb 28, 2024

Announcements

  • HW 03 due TODAY at 11:59pm
  • DataFest: March 1 - 3 in Penn Pavilion
    • See Slack for more information
  • Please fill out Campus Culture Survey by March 1
    • You should have received an email with a personalized link

Topics

  • Interpret and inference for multilevel model coefficients

  • Calculate and interpret intraclass correlation coefficient

  • Maximum likelihood (ML) and restricted maximum likelihood (REML) estimation approaches

  • General process for fitting and comparing multilevel models

Data: Music performance anxiety

Today’s data come from the study by Sadler and Miller (2010) of the emotional state of musicians before performances. The data set contains information collected from 37 undergraduate music majors who completed the Positive Affect Negative Affect Schedule (PANAS), an instrument produces a measure of anxiety (negative affect) and a measure of happiness (positive affect). This analysis will focus on negative affect as a measure of performance anxiety.

The primary variables we’ll use are

  • id: unique musician identification number
  • na: negative affect score on PANAS (the response variable)
  • perform_type: type of performance (Solo, Large Ensemble, Small Ensemble)
  • instrument: type of instrument (Voice, Orchestral, Piano)

Fit model in R

library(lme4)
music_model <- lmer(na ~ orchestra + large_ensemble +
       orchestra:large_ensemble + (large_ensemble|id),
       REML = TRUE, data = music)

  • na ~ orchestra + large_ensemble + orchestra:large_ensemble: Represents the fixed effects

  • (large_ensemble|id): Represents the error terms and associated variance components

    • Specifies two error terms: \(u_i\) corresponding to the intercepts, \(v_i\) corresponding to effect of large ensemble
    • Use (1|id) for models with random intercepts and all other effects fixed.

Model results

effect group term estimate std.error statistic
fixed NA (Intercept) 15.930 0.641 24.833
fixed NA orchestra1 1.693 0.945 1.791
fixed NA large_ensemble1 -0.911 0.845 -1.077
fixed NA orchestra1:large_ensemble1 -1.424 1.099 -1.295
ran_pars id sd__(Intercept) 2.378 NA NA
ran_pars id cor__(Intercept).large_ensemble1 -0.635 NA NA
ran_pars id sd__large_ensemble1 0.672 NA NA
ran_pars Residual sd__Observation 4.670 NA NA

Interpretation

Select the best interpretation for orchestra1:large_ensemble1.

  1. For students who play an orchestral instrument, the mean performance anxiety is expected to be 1.424 points lower for large ensemble performances compared to solo and small ensembles.

  2. The mean decrease in performance anxiety from large ensemble performances versus solos or small ensembles is expected to be 1.424 points greater for students who play orchestral instruments than the expected decrease for soloists and pianists.

  3. The mean performance anxiety for students who play orchestral instruments in large ensembles is expected to be -1.424 points.

Interpretation

Select the best interpretation for sd__(Intercept).

  1. The estimated standard deviation of performance anxiety score for students playing in solos and small ensembles is 2.378 points.

  2. The estimated standard deviation of performance anxiety score for vocalists and pianists is 2.378 points.

  3. The estimated standard deviation of performance anxiety score for students playing in solos and small ensemble is 2.378, after adjusting for instrument.

Inference for fixed effects

  • Notice the R model output has test statistic but no p-values for each coefficient

    • Exact \(t\) distribution under the null hypothesis (no fixed effects) and the associated degrees of freedom are not known

  • We can generally conclude coefficients with test statistic with absolute value greater than 2 are statistically significant

  • Some software will produce p-values by making several assumptions, large sample results , or approximate p-values

  • We will introduce a parametric bootstrap approach in the next chapter.

Unconditional means model

The unconditional means model (also known as random intercepts model) is the multilevel model with no predictors at either level

These models are used to estimate between and within group variability

Level One:

\[Y_{ij} = a_i + \epsilon_{ij} \hspace{10mm} \epsilon_{ij} \sim N(0, \sigma^2)\]

Level Two:

\[ a_i = \alpha_0 + u_i \hspace{10mm} u_i \sim(N, \sigma_u^2) \]

Fitting the unconditional means model

uncond_means_model <- lmer(na ~ 1 + (1 | id), 
                           REML = TRUE, data = music)

tidy(uncond_means_model) |> kable(digits = 3)
effect group term estimate std.error statistic
fixed NA (Intercept) 16.237 0.428 37.943
ran_pars id sd__(Intercept) 2.225 NA NA
ran_pars Residual sd__Observation 4.739 NA NA

Intraclass correlation coefficient

The intraclass correlation coefficient \(\rho\) is

\[ \rho = \frac{\text{Between group variability}}{\text{Total variability}} = \frac{\sigma^2_u}{\sigma^2_u + \sigma^2} \]

In this analysis, \(\hat{\rho} = 0.182\) . This value means…

  • About 18.2% of the variability in performance anxiety can be explained by musician to musician differences

  • The correlation of performance anxiety scores within a musician is 0.182

Note

\(\hat{\rho}\) is calculated based on the variance components from the unconditional means model.

Interpreting \(\rho\)

Which of the following values indicates the individual observations are essentially independent?

  1. \(\rho \approx 1\)
  2. \(\rho \approx 0\)

  • When \(\rho \approx 0\) , the effective sample size (how many pieces of independent information we have) approaches \(n\), the sample size of the data

    • Accounting for multilevel structure of the data is less important when modeling

  • When \(\rho \approx 1\), the effective sample size is close to the number of groups

    • Accounting for multilevel structure of the data is very important when modeling

Estimation

ML and REML

Maximum Likelihood (ML) and Restricted (Residual) Maximum Likelihood (REML) are the two most common methods for estimating the fixed effects and variance components

Maximum Likelihood (ML)

  • Jointly estimate the fixed effects and variance components using all the sample data
  • Can be used to draw conclusions about fixed and random effects

Caution

Issue: Fixed effects are treated as known values when estimating variance components

  • Tends to underestimate variance components (biased estimate), particularly when the sample size is small.

ML and REML

Restricted Maximum Likelihood (REML)

  • Estimate the variance components using the sample residuals not the sample data
  • It is conditional on the fixed effects, so it accounts for uncertainty in fixed effects estimates.
    • This results in unbiased estimates of variance components.

Illustration of ML vs. REML

ML

\[ s^2 = \frac{\sum_{i=1}^n(x_i - \bar{x})^2}{n} \]

Treats \(\bar{x}\) as known value (no loss in degrees of freedom)

REML

\[ s^2 = \frac{\sum_{i=1}^n(x_i - \bar{x})^2}{n-1} \]
Treats \(\bar{x}\) as estimated value (accounts for loss of degree of freedom)


Note: These differences are small when \(n\) is large.

ML or REML?

  • Research has not determined one method absolutely superior to the other

  • REML (REML = TRUE; default in lmer) is preferable when

    • the number of parameters is large, or
    • primary objective is to obtain estimates of the model parameters

  • ML (REML = FALSE) must be used if you want to compare nested fixed effects models using a likelihood ratio test (e.g., a drop-in-deviance test).

    • For REML, the goodness-of-fit and likelihood ratio tests can only be used to draw conclusions about variance components

  • In most cases, there is little difference between the results from ML and REML


Comparing ML and REML results

ML
REML
Term Estimate SE Estimate SE
(Intercept) 15.924 0.623 15.930 0.641
orchestra1 1.696 0.919 1.693 0.945
large_ensemble1 -0.895 0.827 -0.911 0.845
orchestra1:large_ensemble1 -1.438 1.074 -1.424 1.099
sd__(Intercept) 2.286 NA 2.378 NA
cor__(Intercept).large_ensemble1 -1.000 NA -0.635 NA
sd__large_ensemble1 0.385 NA 0.672 NA
sd__Observation 4.665 NA 4.670 NA

Strategy for building multilevel models

  • Conduct exploratory data analysis for Level One and Level Two variables

  • Fit model with no covariates to assess variability at each level

  • Create Level One models. Start with a single term, then add terms as needed.

  • Create Level Two models. Start with a single term, then add terms as needed. Start with equation for intercept term.

  • Begin with the full set of variance components, then remove covariance and variance terms as needed.


Note

Alternate model building strategies in BMLR Section 8.6

Saddler and Miller (2010) strategy

On pg. 283 - 284 of Sadler and Miller (2010), the authors describe their model building process. Try to pick out the steps of the model building process.

Click here to access the paper in Canvas.

06:00

Application exercise

📋 sta310-sp24.netlify.app/ae/lec-13-multilevel-contd

References

🔗 STA 310 - Spring 2024

Roback, Paul, and Julie Legler. 2021. Beyond multiple linear regression: applied generalized linear models and multilevel models in R. CRC Press.
Sadler, Michael E, and Christopher J Miller. 2010. “Performance Anxiety: A Longitudinal Study of the Roles of Personality and Experience in Musicians.” Social Psychological and Personality Science 1 (3): 280–87.
Singer, Judith D, and John B Willett. 2003. Applied Longitudinal Data Analysis: Modeling Change and Event Occurrence. Oxford university press.