I am trying to get my head around the difference between a fixed effect versus a random effect. To do this, I am looking at the Student GPA example seen here (from m-clark.github). This example is "assessing the factors that predict college grade point average (GPA). Each of 200 students are assessed for six occasions (each semester for the first three years), so we have observations clustered within students."
Note that the data for this example can be downloaded here.
The example starts by fitting a standard linear regression, with only occasion as a fixed effect.
load('data/gpa.RData')
gpa_lm = lm(gpa ~ occasion, data = gpa)
summary(gpa_lm)
This produces the following output:
Call:
lm(formula = gpa ~ occasion, data = gpa)
Residuals:
Min 1Q Median 3Q Max
-0.90553 -0.22447 -0.01184 0.26921 1.19447
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.599214 0.017846 145.65 <2e-16 ***
occasion 0.106314 0.005894 18.04 <2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.3487 on 1198 degrees of freedom
Multiple R-squared: 0.2136, Adjusted R-squared: 0.2129
F-statistic: 325.3 on 1 and 1198 DF, p-value: < 2.2e-16
Subsequently, the example demonstrates the fitting of a mixed effects model, where occasion is again included as a fixed effect, but this time student is included as a random (intercept) effect:
library(lme4)
gpa_mixed = lmer(gpa ~ occasion + (1 | student), data = gpa)
summary(gpa_mixed)
Producing the following output:
Linear mixed model fit by REML ['lmerMod']
Formula: gpa ~ occasion + (1 | student)
Data: gpa
REML criterion at convergence: 408.9
Scaled residuals:
Min 1Q Median 3Q Max
-3.6169 -0.6373 -0.0004 0.6361 2.8310
Random effects:
Groups Name Variance Std.Dev.
student (Intercept) 0.06372 0.2524
Residual 0.05809 0.2410
Number of obs: 1200, groups: student, 200
Fixed effects:
Estimate Std. Error t value
(Intercept) 2.599214 0.021696 119.8
occasion 0.106314 0.004074 26.1
Correlation of Fixed Effects:
(Intr)
occasion -0.469
Now I can clearly see that in both models, the estimate for the global intercept remains the same (2.599). And I understand that in the mixed effects model each student gets their own intercept estimate; I can see these by running the following code:
coef(gpa_mixed)$student[1:5,]
Resulting in the following output (only the first 5 students shown):
(Intercept) occasion
1 2.528319 0.1063143
2 2.383636 0.1063143
3 2.687471 0.1063143
4 2.412573 0.1063143
5 2.629598 0.1063143
My question, hence, is what is the difference between the above mixed effects model, which includes student as a random effect, versus a fixed effects model that includes student as a fixed effect?
In both cases an intercept effect is estimated for each student, so what is it about the two models that is different?
To exemplify, here is the fixed effects model:
gpa_fixed = lm(gpa ~ occasion + student, data = gpa)
summary(gpa_fixed)
This model produces the following summary – note that again only students 1-5 are shown, with the estimate for student 1 being represented by "(Intercept)":
Call:
lm(formula = gpa ~ occasion + student, data = gpa)
Residuals:
Min 1Q Median 3Q Max
-0.8676 -0.1474 0.0072 0.1438 0.7325
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.518e+00 9.892e-02 25.450 < 2e-16 ***
occasion 1.063e-01 4.074e-03 26.096 < 2e-16 ***
student2 -1.667e-01 1.392e-01 -1.198 0.231301
student3 1.833e-01 1.392e-01 1.318 0.187967
student4 -1.333e-01 1.392e-01 -0.958 0.338196
student5 1.167e-01 1.392e-01 0.838 0.401995
From this I can see that the estimate for the slope effect of occasion remains unchanged (0.10634 in both the fixed effects and the mixed effects model), but the estimates for the intercept effect for each student are slightly different.
For example, in the mixed effects model, student 1 has an intercept estimate of 2.528319
, whereas in the fixed effects model the intercept estimate for student 1 is 2.518
.
Similarly, in the mixed effects model, student 2 has an intercept estimate of 2.383636
, whereas in the fixed effects model the intercept estimate for student 2 is 2.518-0.1667=2.3513
.
Admittedly these values are very similar, but they are not identical.
Are these values different only because of the way that they are estimated? i.e. partial pooling/shrinkage in the case of the mixed effects model versus least-squares in the fixed effects model?
My understanding is that when you include student as a random effect, the estimate for each student's intercept is drawn from a normal distribution centered on the global intercept estimate. Whereas when you include student as a fixed effect, each student's intercept can be estimated with whatever value is best (i.e. best according to least squares) without any constraint as to what distribution that estimate is drawn from. Is my understanding here correct?
Why would it matter, from a statistical point of view, if I modeled student as a fixed effect, rather than as a random effect? Is it for example because with the fixed effect you cannot "borrow strength" from other groups? And because with the fixed effect you are having to estimate a greater number parameters than with the random effect, which only considers variance around the global intercept estimate?
Thank you in advance!
Edit: I also have another question - how exactly do you interpret the value of the random effect?
Random effects:
Groups Name Variance Std.Dev.
student (Intercept) 0.06372 0.2524
Residual 0.05809 0.2410
Number of obs: 1200, groups: student, 200
The value for the variance effect size is 0.06372
. Does this mean that, when variance of all the individual student intercept estimates is 0.06372
? I tried that using var(ranef(gpa_mixed)$student)
but this returns a value of 0.0553
. If anyone could shed some light on how to intuitively interpret the random effect estimates, that would really help!