Which ICC estimator to use under the fixed effects model?

Question

In the random effects model $$y_{ij} = \mu + \alpha_j + x_{ij} + \epsilon_{ij}$$ the intra-class correlation coefficient is given by $$ICC = \frac{\sigma_{\alpha}}{\sigma_{\alpha}+\sigma_{\epsilon}}$$ and can be estimated by plugging in the variance estimates of the random effects model.

When the number of clusters / groups is small, let's say below 20, we typically do not want to use a random effects model and model $\alpha_j$ as fixed effects instead.

What is an accepted estimator for the ICC under the fixed efects model? I am searching the literature but cannot seem to find it.

Answer 1

When the number of clusters / groups is small, let's say below 20, we typically do not want to use a random effects model and model $\alpha_j$ as fixed effects instead.

I think most people would disagree that 20 is too few clusters to use a random effects model. While there is no hard and fast rule, 6 seems to be a reasonable consensus.

What is an accepted estimator for the ICC under the fixed efects model? I am searching the literature but cannot seem to find it.

On way to do this is to estimate a model without the grouping factor (call it m0) and then another model with the grouping factor as a fixed effect (call it m1). Then compute the difference in the residual variance and divide it by the residual variance of m0. The idea behind this is that the fixed effects for the grouping factor absorb a certain amount of variance in the response. When they are normally distributed (as is the assumption in linear mixed models) this should be the same amount of variance as that estimated by the random intercepts.

Here is how it can be done in R:

We simulate clustered data with an expected ICC of 0.8 (variance of random intercepts of 4, and residual variance of 1):

> set.seed(2)
> dt <- expand.grid(hospID = 1:10, patientID = 1:20)
> dt$Y <- 1
> X <- model.matrix(~ 1, data = dt)
> myFormula <- "Y ~ 1 + (1 | hospID)"
> foo <- lFormula(eval(myFormula), dt)
> Z <- t(as.matrix(foo$reTrms$Zt))  # design matrix for random effects
> betas <- 10    # fixed effects (intercept only in this case)
> b <- rnorm(10, 0, 2)  # random effects (standard deviation of 2, variance of 4)

> dt$Y <- X %*% betas + Z %*% b + rnorm(nrow(dt))

Now we fit the linear mixed model:

> library(lme4)
> (lm0 <- lmer(eval(myFormula), dt)) %>% summary()
Random effects:
 Groups   Name        Variance Std.Dev.
 hospID   (Intercept) 4.011    2.003   
 Residual             1.188    1.090   
Number of obs: 200, groups:  hospID, 10

And we see the estimated variance components are as expected.

Now we fit the models m0 and m1 as described above and compute the ICC from the mixed model and also from the linear models:

> m0 <- lm(Y ~ 1, dt)
> m1 <- lm(Y ~ 1 + as.factor(hospID), dt)

> dt.vc <- as.data.frame(VarCorr(lm0))  # extract the variance components
> (ICC.lmm <- dt.vc[1, 4] / (dt.vc[1, 4] + dt.vc[2, 4]))
[1] 0.7715357
> (ICC.lm <- (var(residuals(m0)) - var(residuals(m1))) / var(residuals(m0)) )
[1] 0.7645219

and these seem to agree well. You can change the seed, change the simulated variances add other fixed effect etc as you see fit.

Answer 2

What is an accepted estimator for the ICC under the fixed efects model? I am searching the literature but cannot seem to find it.

That's probably because there's not much to say about the ICC in a standard fixed effects model -- it is necessarily 0 due to the assumptions that the group/cluster effects are "fixed" and thus not random variables (so their variance is 0) and that the errors are uncorrelated.

Mathematical details

Here is some brief background on where the ICC comes from and what it means, with some text repurposed from my answer HERE.

The fixed effects model, as you've written it, is $$ y_{ij} = \mu + \alpha_j + x_{ij} + \epsilon_{ij}, $$ where the intercepts $\alpha_j$ are defined/assumed to be fixed and thus have variance = 0, and the residuals $\epsilon_{ij}$ have variance $\sigma^2_\epsilon$ (in your question you omitted the square on this term, but here I've added it to be more in line with the usual notation).

Now, the correlation between two random variables $x$ and $y$ is defined as $$ \text{corr} = \frac{\text{cov}(x, y)}{\sqrt{\text{var}(x)\text{var}(y)}}. $$

So to find the intra-class correlation, we use the correlation formula and let our two random variables be two observations (say $i = 1, 2$) drawn from the same $j$ group,

$$ \begin{aligned} ICC &= \frac{\text{cov}(\mu + \alpha_j + x_{1j} + \epsilon_{1j}, \mu + \alpha_j + x_{2j} + \epsilon_{2j})}{\sqrt{\text{var}(\mu + \alpha_j + x_{1j} + \epsilon_{1j}) \text{var}(\mu + \alpha_j + x_{2j} + \epsilon_{2j})}} \\ &= \frac{\text{cov}(\epsilon_{1j}, \epsilon_{2j})}{\sqrt{\text{var}(\epsilon_{1j}) \text{var}(\epsilon_{2j})}} \\ &= \frac{0}{\sigma^2_\epsilon} \\ &= 0, \end{aligned} $$

where the numerator simplifies to 0 due to the assumptions that the group/cluster effects are "fixed" and thus not random variables (used on line 2 above) and that the errors are uncorrelated (used on line 3 above).