I have always been interested in understanding the following point: During the parameter estimation process, why are Mixed Effects Models able to "implicitly" perform shrinkage/regularization (I have heard this being referred to as Empirical Bayes) in the background ? Apparently this shrinkage process is seen as a good thing, as extreme/large values of parameter estimates are believed to be unlikely and caused by outliers more than actual trends ... thus shrinking large parameter estimates makes the model more realistic by shrinking the effects of outliers on parameter estimates.

I tried to manually write down the Likelihood Functions for such a problem and see if I can pinpoint the aspects where shrinkage/regularization is taking place.

Model and Likelihood: Suppose we take a simple mixed effects model with a single random intercept for group $j$, the model is: $$ y_{ij} = \beta_0 + b_j + \epsilon_{ij} $$ The random intercept $b_j$ is assumed to follow a normal distribution with mean 0 and variance $\sigma_b^2$: $$ b_j \sim N(0, \sigma_b^2) $$ The likelihood function for this model is the product of two parts: the likelihood of the observed data given the fixed and random effects, and the likelihood of the random effects given their distribution: $$ L(\beta, b, \sigma^2 | y) = \prod_{j=1}^{J} \left[ \frac{1}{\sqrt{2\pi\sigma^2}} \exp \left( -\frac{(y_{ij} - \beta_0 - b_j)^2}{2\sigma^2} \right) \right] \times \prod_{j=1}^{J} \left[ \frac{1}{\sqrt{2\pi\sigma_b^2}} \exp \left( -\frac{b_j^2}{2\sigma_b^2} \right) \right] $$

Penalty Term: Looking at the second term:

$$ \prod_{j=1}^{J} \left[ \frac{1}{\sqrt{2\pi\sigma_b^2}} \exp \left( -\frac{b_j^2}{2\sigma_b^2} \right) \right] $$

  • The second part of this likelihood function seems to be related to shrinkage. This term looks like a penalty that discourages the estimates of the random effects from deviating too far from 0. The larger the variance $\sigma_b^2$, the less penalty there is for deviating from 0, and the less shrinkage there is. On the other hand, the smaller the variance, the more penalty there is for deviating from 0, and the more shrinkage there is.

  • The $-\frac{b_j^2}{2\sigma_b^2}$ term in the exponent looks like a penalty that increases as $b_j$ deviates further from 0. When trying to maximize the likelihood function, this penalty discourages large values of $b_j$, effectively "shrinking" the estimates towards 0. The amount of shrinkage seems to depends on the variance $\sigma_b^2$: a larger variance means less shrinkage, while a smaller variance means more shrinkage.

  • It appears that the fixed effects part of the model (i.e. first part of the likelihood) does not contain any random effects, and therefore does not include a penalty term in its likelihood function. As a result, there is no shrinkage towards 0. This is also true in a fixed effects model with interactions: no shrinkage

Case by Case Analysis: I also tried to consider some cases to see what happens when $b_j$ deviates from 0:

  • If $b_j$ is far from 0 (either positively or negatively), then $b_j^2$ is a large positive number. This makes the fraction $\frac{b_j^2}{2\sigma_b^2}$ large, which in turn makes the entire exponent $-\frac{b_j^2}{2\sigma_b^2}$ a large negative number. The exponential of a large negative number is close to 0, so the overall likelihood is decreased. This is the penalty for deviating far from 0.
  • If $b_j$ is close to 0, then $b_j^2$ is close to 0, which makes the fraction $\frac{b_j^2}{2\sigma_b^2}$ close to 0, and the entire exponent $-\frac{b_j^2}{2\sigma_b^2}$ close to 0. The exponential of 0 is 1, so the overall likelihood is not decreased. There is no penalty for being close to 0.

Similar analysis for $\sigma_b^2$:

  • If $\sigma_b^2$ is large, then the fraction $\frac{b_j^2}{2\sigma_b^2}$ is small, even for large values of $b_j$. This makes the exponent $-\frac{b_j^2}{2\sigma_b^2}$ close to 0, and the overall likelihood is not decreased much. There is less penalty for deviating from 0, and therefore less shrinkage.
  • If $\sigma_b^2$ is small, then the fraction $\frac{b_j^2}{2\sigma_b^2}$ is large, even for small values of $b_j$. This makes the exponent $-\frac{b_j^2}{2\sigma_b^2}$ a large negative number, and the overall likelihood is decreased significantly. There is more penalty for deviating from 0, and therefore more shrinkage.

Is my statistical reasoning correct here?

  • Note: I deliberately left out any mention of "partial pooling" here. I am planning on posting another question on the mechanics of partial pooling.