5
$\begingroup$

I have two observations from a person each, where every observation corresponds to a different treatment. The treatments are fixed effects, the persons are random effects. I use the command:

model <- lmer(response ~ treatment + (1|Person))

When I run summary() on this model, I get variances for the residual and the random person effect. I also get the mean values corresponding to what treatment was given. I understand this.

But what's the output from ranef? fixef? coef? I don't understand what these things do. I also know nothing about mixed models other than the very basic definition, so I would like a non-technical explanation as opposed to one that says "oh ranef just gives you the insert technical term here ").

EDIT:

Playing around with things, it seems that the output of 'coef' are sums of the fixef and ranef outputs. The fixed outputs seem to have an obvious explanation, so I guess my question reduces to an inquiry on 'ranef' and 'fitted'. What do these do?

I believe it is especially 'ranef' that I don't undertand. How are these random effects (that's what ranef stands for, no?) when the only estimates we have to work with are parameter estimates for the mean values and a variance estimate for the random effect and a variance estimate for the residual "noise"?

$\endgroup$
11
$\begingroup$
  • fixef() is relatively easy: it is a convenience wrapper that gives you the fixed-effect parameters, i.e. the same values that show up in summary(). Unless you are specifying your model in a very particular way, these are not the "mean values corresponding to what treatment was given" as suggested in your question; rather they are contrasts among treatments. Using R's default setup, the first parameter ("Intercept") is the mean response for the first treatment level, while the remaining parameters are the difference between the mean responses for levels 2 and higher and the mean response for level 1. (From Jake Westfall in comments: "Another way of explaining fixef() is that it returns essentially the same thing as when you call coef() on an lm regression object -- that is, it returns the (mean) regression coefficients.")
  • ranef() gives the conditional modes, that is the difference between the (population-level) average predicted response for a given set of fixed-effect values (treatment) and the response predicted for a particular individual. You can think of these as individual-level effects, i.e. how much does any individual differ from the population? They are also, roughly, equivalent to the modes of the Bayesian posterior densities for the deviations of the individual group effects from the population means (but note that in most other ways lme4 is not giving Bayesian estimates).

    It's not that easy to give a non-technical summary of where the conditional modes come from; technically, they are the solutions to a penalized weighted least-squares estimation procedure. Another way of thinking of them is as shrinkage estimates; they are a compromise between the observed value for a particular group (which is what we would estimate if the among-group variance were infinite, i.e. we treated groups as fixed effects) and the population-level average (which is what we would estimate if the among-group variance were 0, i.e. we pooled all groups), weighted by the relative proportions of variance that are within vs among individuals. For further information, you can search for a non-technical explanation of best linear unbiased predictions (or "BLUPs"), which are equivalent to the conditional modes in this (simple linear mixed model) case ...

  • coef() gives the predicted effects for each individual; in the simple example you give, coef() is basically just the value of fixef() applicable to each individual plus the value of ranef().

I agree with comments that it would be wise to look for some more background material on mixed models:

  • Gelman and Hill's Applied Regression Modeling
  • Pinheiro and Bates Mixed-Effects Models in S and S-PLUS
  • various books by Zuur et al.
  • McElreath's Rethinking Statistics
  • (shameless plug) chapter 13 in Fox et al. Ecological Statistics
$\endgroup$
  • $\begingroup$ Another way of explaining fixef() is that it returns essentially the same thing as when you call coef() on an lm regression object -- that is, it returns the (mean) regression coefficients. $\endgroup$ – Jake Westfall May 23 '16 at 16:06
  • $\begingroup$ So how is the calculation of 'ranef' carried out in practice? Let me give you my estimates: The mean value estimate for treatment A is $95$ and, as you said, the contrast estimate for treatment B is $9.6$ (so B-treatments have mean $95 + 9$). These are the fixed effects of these two treatments. Then, I have two random effects: the residual has variance $18$ and the variance corresponding to the persons (8 different people, 2 observations for each corresponding to each treatment) is $205$. So, how would you obtain 'ranef'? The response from, say, person 1 was 94 (A) and 99 (B). So ranef[1] = ? $\endgroup$ – user116612 May 23 '16 at 16:42
  • $\begingroup$ It's potentially complicated. Can you edit your question to include a reproducible example/some of the output from your model? $\endgroup$ – Ben Bolker May 23 '16 at 16:44
  • $\begingroup$ In the linear case, you point out that these are equivalent to BLUPs of the random effects. In the non-linear (and linear) case, these can also be seen as the mode of the posterior distribution of each individual's (or "cluster"'s) random effect, right? Its not totally clear to me what exactly is output by ranef() and I can't seem to find any documentation of it but the posterior mode is always what I assumed. $\endgroup$ – use_norm_approx Apr 21 '18 at 2:40
2
$\begingroup$

The model fitted with lmer(response ~ treatment + (1|Person)) may be expressed in a matrix form as $$ y = X\beta + Zu + e $$ where $\beta$ is the fixed effect vector, u the random effect vector and e the vector of error terms. The function getME{lme4} may be used to extract both $X$ and $Z$ for a model. For the model in discussion, we apply the typical assumptions of normality, zero means for the random elements, and
$$Var(u) = \sigma_u^2 I, \ Var(e) = \sigma_e^2I\ and \ Cov(u, e') = 0.$$ Therefore, $$Var(y) = \sigma_u^2 ZZ'\ + \sigma_e^2I$$

Since the question was specific about how the "estimated" random effects of $u$ are calculated, I assume there is an understanding of how $\hat\beta, \hat\sigma_u^2,$ and $\hat\sigma_e^2$ are obtained by the method of restricted maximum likelihood (REML), which is the default for lmer().

As mentioned by others, $u$, being random, can be predicted by the BLUP method. Conceptually, the key concept is to find an estimator, not for $u,$ but for the conditional mean of $u$ given $y$, denoted as $E(u|y)$. To do this, observe that $$ \begin{pmatrix} y\\ u \end{pmatrix} = \begin{pmatrix} X\beta\\ 0 \end{pmatrix} + \begin{pmatrix} Z & I\\ I & 0 \end{pmatrix} \begin{pmatrix} u\\ e \end{pmatrix} \sim N\begin{pmatrix} \begin{pmatrix} X\beta \\ 0\\ \end{pmatrix}, & \begin{pmatrix} \sigma_u^2 ZZ' + \sigma_e^2 I & \sigma_u^2 Z\\ \sigma_u^2 Z' & \sigma_u^2I \end{pmatrix} \\ \end{pmatrix} $$

Therefore, by the property of multivariate normal distribution, $$E(u|y) = Z'( ZZ' + \frac{\sigma_e^2}{\sigma_u^2} I)^{-1}(y - X\beta)$$

It then seems natural to plug in the REML estimators $\hat\beta, \hat\sigma_u^2,$ and $\hat\sigma_e^2$ into the above formula for $E(u|y)$ as a predictor for the unobserved $u$. Theoretical works have established that this estimator indeed is "the" BLUP of $u$. I do not know the scripts of lmer() well to claim that it uses this algorithm; however, upon testing, the algorithm did produce numbers consistent with ranef().

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy