I'm new to mixed modelling and I'm confused as to whether its appropriate to use a random effect in an analysis I'm doing. Any advice would be appreciated.

My study is testing how well a newly developed index of mammal abundance can predict the value of an established but more labour intensive index. I've been measuring these indices in multiple forest patches, with multiple plots in each forest patch.

Because I'm not directly interested in the effect of forest patches, and because my sample plots are nested within forest patches, I've been using forest patch as a random effect. However, I've got a couple of questions about this:

  1. First, I know that random effects allow you to generalise your results across all possible levels of the random factor, not just the ones you sampled. But it seems to me that to make this kind of inference your levels would have to be randomly sampled? My forest patches were not randomly sampled, so can I still use them as a random effect?
  2. Second, Ive read that you can test whether it is necessary to have a random effect by doing, e.g., a likelihood ratio test to compare models with and without the effect. I've done this, and it suggests that the random effect model does not explain the data as well as a fixed effects only model. My issue with this is that my plots are still nested within forest patches, and so presumably not independent. So, can I use this LRT approach to justify excluding the random effect, or do I still need to include it to account for nestedness? And if I do end up removing the random effect, is there a way to verify that plots within forest patches can be considered independent?
  • $\begingroup$ I think the main conceptual thing with random effects is that they should all be about the same magnitude and be exchangeable - this is what makes the inference about non-sampled random effects possible. Also, you should be careful using LR tests for random effects as there may be considerable uncertainty regarding the variance component even if the ML/REML estimate is zero or close to zero. $\endgroup$ Commented Jan 28, 2013 at 22:38
  • $\begingroup$ Great thanks for that.So would there be some way to work out if its necessary to keep the random effect in? $\endgroup$
    – jay
    Commented Jan 29, 2013 at 23:30
  • $\begingroup$ From what I understand, it is best to not do LR comparisons with models fit by REML. In R's lmer, for example, you should set REML=FALSE when doing LRT. (It defaults to TRUE, which is otherwise better.) $\endgroup$
    – Wayne
    Commented Feb 8, 2013 at 15:56

2 Answers 2


As I understand, you have a simple nested observational design (plots within patches) and your interest is in a correlation/regression between two continuous variables (the two indices). Your sample size is m patches x n plots = N pairs of observations (or the appropriate sumatory if unbalanced). No proper randomization was involved, but maybe you can/should/want to consider that (1) the patches were "randomly" selected from all the patches of this kind or in some area, and then (2) the plots were "randomly" selected within each patch.

If you ignore the random factor Patch, you may be pseudoreplicating by considering that you have randomly selected N plots "freely", without constraining them to be (in number or type) in those (previously) selected patches.

So, your first question: yes, that is what a random factor allows. The validity of such inference depends on the validity of the assumption that haphazard selection is equivalent to random selection of patches (e.g., that your results would not be different if a different set of forest patches was selected). That puts a limit also on your space of inference: the kind of forest or geographical area up to which your results extend depends on the maximal (imaginary) population of patches from where your sample is a credible "random" sample. Maybe your observations are a "reasonable random" sample of the mammals of the forest patches in your region but would be a suspiciously aggregated sample of the mammals of the whole continent.

The second one: the test will depend on "the degree of pseudoreplication", or the evidence in your sample that plots "belong" to patches. This is, how much variation there is among patches and among plots within patches (search for intraclass correlation). In an extreme, only variation among patches is present (plots within a patch are all the same) and you have "pure pseudoreplication": your N should be the number of patches, and sampling one or many plots from each of them does not provide new information. On the other extreme, all variation happens between plots, and there is no extra variation explained by knowing to which forest patch each plot belongs (and then the model without the random factor would appear more parsimonious); you have "independent" plots. NONE of the extremes are very likely to happen... particularly for biological variables observed on the ground, if only because of spatial autocorrelation and geographical distributions of the mammals. I personally prefer to keep factors by design anyway (e.g., even when patches is not a relevant source of variation IN THIS SAMPLE) to sustain the "experimental-observational" analogy explained above; remember: not having evidence in your sample to reject the null hypothesis stating that variation among patches is zero does not mean that variation is really zero in the population.


Random effects induce heteroskedasticity and correlation in "error terms" in the model

A useful way to look at random effects models is to use some mathematical manipulation to mash them back onto the form of a traditional model without a random effect term. This can be accomplished by absorbing the random effects terms into the error terms of a traditional fixed effect model. We can then examine the properties of the new error term in the traditional fixed effect model form, to see what happens with this variation.

To illustrate this technique, consider a Gaussian linear random effects model with observations $Y_{i,j}$ taken over categories $j=1,...,k$. The model can be written as:

$$Y_{i,j} = \beta_0 + \beta_j + u_j + \varepsilon_{i,j} \quad \quad \quad u_j \sim \text{N}(0,\sigma_j^2) \quad \quad \quad \varepsilon_{i,j} \sim \text{N}(0,\sigma^2),$$

where the random effects terms and error terms are mutually independent. Since both the error terms and the random effect terms are random variables in the model, we can combine them into a single alternative error term and write the model in a form that does not show a separate random effect. Specifically, we write the model in the traditional fixed effect form:

$$Y_{i,j} = \beta_0 + \beta_j + \eta_{i,j},$$

where we have defined the quantities:

$$\eta_{i,j} = u_j + \varepsilon_{i,j} \quad \quad \quad \quad \quad \rho_j = \frac{\sigma_j^2}{\sigma_j^2+\sigma^2} \quad \quad \quad \quad \quad \sigma_{j*}^2 = \sigma_j^2+\sigma^2.$$

In this latter form, the error terms $\eta_{i,j}$ are still (jointly) normally distributed with zero mean, but they are no longer homoskedastic and uncorrelated --- they have covariance values:

$$\mathbb{C}(\eta_{i,j},\eta_{i',j'}) = \mathbb{C}(z_j + \varepsilon_{i,j}, z_{j'} + \varepsilon_{i',j'}) = \begin{cases} \sigma_{j*}^2 & & \text{if } i = i' \text{ and } j = j', \\[6pt] \rho_j \sigma_{j*}^2 & & \text{if } i \neq i' \text{ and } j = j', \\[6pt] 0 & & \text{otherwise}. \\[6pt] \end{cases}$$

This form means that there is heteroskedasticity in the model (i.e., with variance $\sigma_{j*}^2$ for observations in category $j$) and the errors within a group are positively correlated (with correlation coefficient $\rho_{j}$).

Should you use a random effects model? As you can see from the above, a Gaussian linear random effects model using categorical predictors is equivalent to a traditional Gaussian linear regression model using categorical predictors, where the latter has heterosedasticity across the categories of observations and positively correlated errors within each category. Consequently, your choice of whether or not to include random effects terms can be framed equivalently as a choice of whether or not to generalise the behaviour of the error terms in the model to allow heteroskedasticity across categories and correlation of error terms within categories.

In the study you are conducting, your forest patches are the categories and your plots within these forest patches are your individual observations. In this case, including a random effect in your model (taken at the level of the forest patches) is equivalent to allowing heteroskedasticity across different forest patches and also having positive correlation of the mammal abundance of different plots within each forest patch. In order to decide whether or not this is appropriate, you merely need to ask yourself if this type of heteroskedasticity/correlation might plausibly occur in this case.

You have also noted that your forest patches were not randomly sampled. This is not necessarily a problem for the random effects model, and it is no more of a problem than to a fixed effects model. In assessing your sampling method you should consider whether your choice of sites was influenced by any of the variables under study, and consider the types of biases this could induce. However, there is nothing inherent in the random effects model (as opposed to the fixed effects model) that presents an analytical difference here.

I would suggest that in the present case this kind of heteroskedasticity/correlation might plausibly occur, owing to the closeness of the plots within a patch and the possible movement behaviour of the mammals under study. Different plots within the same forest patch might plausibly have positively correlated mammal abundance due to the movement of mammals around a forest patch, breeding of mammals from amongst nearby plots within a forest patch, and common conditions/threats to mammals within the same forest patch. You have noted that you tried both the random effects and fixed effects models and conducted a likelihood-ratio test, finding that there was no significant evidence of the presence of random effects. That is perfectly fine, and it is one way to conduct your analysis --- you have an initially plausible model form and then you find that it does not operate better than a simpler model form.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.