categorical random effects correlation in lme4

Thu, May 23, 2013 2:21 PM

Hi Andrew,

in Baayen, Davidson, & Bates (2008) they discuss a situation very similar to yours and in which they perform similar steps as you have done (full model has a correlation of slopes of -1, then remove correlation between random slopes, and then remove random slopes altogether) and then decide for the smallest model (only random intercepts) as it provides the best account in terms of a LRT .

In your case it now seems to be the case that the full model and intercept only model are identical (as predicted by slide 38 of the presentation you linked).

So you basically only have the options between:
- a model with (0 + factor(white) | primary_ther) [or equivalently (factor(white) | primary_ther)],
- (1|primary_ther) + (0+white|primary_ther) [as I suggested and which should have less random parameters as the model before, which is done by Baayen et al. (2008), and also recommended by Barr et al. (2013) in case of nonconvergence],
- and the random intercept only model (1 | primary_ther).

What does the LRT [i.e., anova()] now say to these three options?

I hope this helps and I really liked the Barr article (also I would not trust their recommendation on how to obtain p-values but resort on either afex or car::Anova(..., test = "F") which are equivalent).

I hope that helps and sorry if you did also this already,
Henrik

Note: Baayen et al use slightly different syntax.

References:
Baayen, R. H., Davidson, D. J., & Bates, D. M. (2008). Mixed-effects modeling with crossed random effects for subjects and items. Journal of Memory and Language, 59(4), 390?412. doi:10.1016/j.jml.2007.12.005

Barr, D. J., Levy, R., Scheepers, C., & Tily, H. J. (2013). Random effects structure for confirmatory hypothesis testing: Keep it maximal. Journal of Memory and Language, 68(3), 255?278. doi:10.1016/j.jml.2012.11.001



Am 22/05/2013 20:22, schrieb Andrew McAleavey:

Hi Malcolm & Henrik & all,

Thanks for the help! I've tried some of these before, but tried even more
now. I think I'm unfortunately still stuck on this, and I think a part of
the reason is that dropping the random intercepts seems weird but I think
it is necessary and helpful with a categorical fixed/random interaction.

I know that removing the random intercept seems weird, but my understanding
is that eliminating the intercept when modeling a categorical random slope
helps interpretation, but doesn't change the model substantively. The
variance of the random slope parameter for "white" changes and the
correlation flips to -1.000 when you add in the intercept, but the model is
actually equivalent. That is, without the intercept in the model, the
random white=0 effect actually is the random intercept for therapists. See
model fit below.

As I understand it, by removing the intercept when using a categorical
random/fixed interaction effect, you get a more interpretable value, namely
variance in the differences between the reference group and the target
group attributable to therapists. I believe that the variance shown when
the intercept is included is too large, since it actually describes all
variance associated with the categorical variable as deviations from an
"average" client who doesn't exist.

When estimating the model suggested by Henrik [(1|primary_ther) +
(0+white|primary_ther)], two things seem odd to me. First, the intercept is
forced to be uncorrelated with the slopes - this is not theoretically
necessary. The other is that there are three random effects being estimated
now, which is too many (with a two level categorical effect, only one
parameter is necessary). Unless I'm wrong about something else, I don't
think I can treat this the way I would a continuous variable. I think this
is redundant and See below again.

So I suppose I'm still stuck with the original question, since I don't
think these solutions hold with a categorical effect. Am I just missing
something else, or maybe the suggestions ought to hold for categorical
effects and it's something else causing difficulties?

Thanks a lot,
Andrew

Here's the model adding back in the random intercepts but including the
white random effect:

print(fm1_ml_int <- lmer(DI ~ first_di + factor(white) +

(factor(white)|primary_ther), rem3post, REML=F), corr=F)
Linear mixed model fit by maximum likelihood
Formula: DI ~ first_di + factor(white) + (factor(white) | primary_ther)
    Data: rem3post
   AIC  BIC logLik deviance REMLdev
  4980 5020  -2483     4966    4983
Random effects:
  Groups       Name           Variance Std.Dev. Corr
  primary_ther (Intercept)    0.041708 0.20423
               factor(white)1     0.018995 0.13782  -1.000
  Residual                         0.509959 0.71411
Number of obs: 2263, groups: primary_ther, 192

Fixed effects:
                Estimate Std. Error t value
(Intercept)       0.45138    0.06007   7.514
first_di            0.48009    0.02360  20.338
factor(white)1  0.01140    0.03298   0.346

Here's the difference test of interest, they appear to be equivalent models:

anova(fm1_ml, fm1_ml_int)

Data: rem3post
Models:
fm1_ml: DI ~ first_di + factor(white) + (0 + factor(white) | primary_ther)
fm1_ml_int: DI ~ first_di + factor(white) + (factor(white) | primary_ther)
                  Df    AIC    BIC  logLik Chisq Chi Df Pr(>Chisq)
fm1_ml       7 4979.9 5019.9 -2482.9
fm1_ml_int  7 4979.9 5019.9 -2482.9     0      0          1

And here's a version of the model with uncorrelated intercept and slopes
(the bigger model has a singular convergence, as anticipated in these
slides:
http://lme4.r-forge.r-project.org/slides/2011-03-16-Amsterdam/2Longitudinal.pdf
):
Formula: DI ~ 1 + (1 | primary_ther) + (0 + factor(white) | primary_ther)
    Data: rem3post_2
   AIC  BIC logLik deviance REMLdev
  5027 5062  -2508     5015    5022
Random effects:
  Groups       Name           Variance          Std.Dev.     Corr
  primary_ther (Intercept)       1.7240e-06  0.001313
  primary_ther factor(white)0  4.6543e-02  0.215738
                     factor(white)1  6.9851e-03  0.083577   0.752
  Residual                            5.1870e-01  0.720212
Number of obs: 2263, groups: primary_ther, 192

Fixed effects:
             Estimate Std. Error t value
(Intercept)  1.48466    0.01791   82.91


On Tue, May 21, 2013 at 1:30 PM, Malcolm Fairbrother <
M.Fairbrother at bristol.ac.uk> wrote:

Dear Andrew,

What if you drop the "0 +" bit? So:

  lmer(DI ~ first_di + factor(white) + (factor(white) | primary_ther),
rem3post, REML=F)

or

  lmer(DI ~ first_di + white + (white | primary_ther), rem3post, REML=F)

Including "0 +" means you're not estimating a random intercept for
"primary_ther", which it sounds like you need/want. Instead, you're
getting two random slopes, which I guess are perfectly correlated
because they're two sides of the same coin (the coin being your binary
dummy variable "white").

If that doesn't solve the problem, it might help for you to post the
results of "str(rem3post)" (i.e., your dataset) as well.

Cheers,
Malcolm

Date: Tue, 21 May 2013 11:41:04 -0400
From: Andrew McAleavey <andrew.mcaleavey at gmail.com>
To: r-sig-mixed-models at r-project.org
Subject: [R-sig-ME] categorical random effects correlation in lme4

Hi,

I'm currently investigating a question of relative effectiveness of
therapists, and the particular question is whether some therapists are
differentially effective with white versus racial/ethnic minority clients
(this is coded as a binary variable called "white" in this data). We have
conceptualized this as a cross-level random effect, so the model has one
random effect for therapist intercept and one effect for the difference

in

effectiveness between their white and nonwhite clients.

I am relatively new to lme4, but I think I have specified the model
correctly (the fixed effects represent client pretreatment severity and

the

nonsignificant fixed effect of binary race; they don't seem to impact the
estimation problem). Here's the model of interest:

print(fm1_ml <- lmer(DI ~ first_di + white + (0 +

factor(white)|primary_ther), rem3post, REML=F), corr=F)

The problem is that the two random effects are appearing to correlate at

= 1.000. I think this is an estimation problem, and probably indicates

that

the random variables aren't accounting for all that much variance. I'm
dubious of interpreting this model, therefore. However, when comparing it
to the random intercepts only model using the LRT, there is a significant
difference, suggesting that even though the explained variance is (very)
small, it may be worth including:

anova(fm1_a_ml, fm1_ml)

Data: rem3post
Models:
fm1_a_ml: DI ~ first_di + factor(white) + (1 | primary_ther)
fm1_ml: DI ~ first_di + factor(white) + (0 + factor(white) |

primary_ther)

                 Df    AIC      BIC      logLik     Chisq    Chi Df
Pr(>Chisq)
fm1_a_ml   5    4982.7  5011.3  -2486.3
fm1_ml      7    4979.9   5019.9  -2482.9  6.7871      2        0.03359 *

My question is basically this: How should I interpret these results?

There

are significant differences between therapists in terms of their relative
effectiveness with white vs. nonwhite clients, but they're just small? Or
is even this not justified? Would it be safer to say that there are

likely

no estimable differences? Am I missing something else?

Thanks a lot,
Andrew McAleavey

Here's the model of interest output:

print(fm1_ml <- lmer(DI ~ first_di + factor(white) + (0 +

factor(white)|primary_ther), rem3post, REML=F), corr=F)
Linear mixed model fit by maximum likelihood
Formula: DI ~ first_di + factor(white) + (0 + factor(white) |

primary_ther)

    Data: rem3post
   AIC  BIC logLik deviance REMLdev
  4980 5020  -2483     4966    4983
Random effects:
  Groups       Name              Variance        Std.Dev.   Corr
  primary_ther factor(white)    0 0.0417050  0.204218
                  factor(white)1     0.0044086     0.066397   1.000
  Residual                            0.5099596     0.714115
Number of obs: 2263, groups: primary_ther, 192

Fixed effects:
                      Estimate Std. Error t value
(Intercept)        0.45138    0.06007   7.514
first_di             0.48009    0.02360   20.338
factor(white)1   0.01140    0.03298   0.346

--
Andrew McAleavey, M.S.
Department of Psychology
The Pennsylvania State University
346 Moore Building
University Park, PA 16802
aam239 at psu.edu

Dipl. Psych. Henrik Singmann
PhD Student
Albert-Ludwigs-Universit?t Freiburg, Germany
http://www.psychologie.uni-freiburg.de/Members/singmann

categorical random effects correlation in lme4

Thread (6 messages)