[R-sig-ME] lmer-model - model ok? (control condition, random effects, log transformation)

Phillip Alday Phillip.Alday at unisa.edu.au
Tue Apr 11 07:00:21 CEST 2017

Hi Diana,

I don't think anybody has replied to your message, so here's my go.

The way your conditions are coded seems weird to me. I would have all five conditions (including control) in a single variable/column "cond":

# NB: I use "dat" and not "data", see PS below
dat$cond <- factor(dat$cond, levels=c("control","test1","test2","test3","test4")

If you haven't changed your options, this will be dummy coding with "control" as the baseline/intercept/implicit level.  If you're interested in pairwise contrasts for each test condition compared to the control, this is definitely the way to go.  Based on your description of your experiment, you don't have nested conditions, so don't use the nesting syntax (which expands in a different way than you think anyway and I've never seen used for fixed effects). If you care about other contrasts, you can compute them post hoc with lsmeans or define your own contrast matrix (more advanced, but actually covered in some intro books like Field et al.)

Also, you have a model with the term "rat3" but you don't describe what that represents.

So your basic model would something like:

m <- lmer(logRT ~ cond + (1|subj) + (1|item_ID), dat, REML=F)

Please note that you should not have any random slopes that are not in your fixed effects. There are exceptions to this rule, but yours is not one of them.  In particular, the by-subject slope for condition should be the same condition predictor as in the fixed effects (assuming that all subjects saw all conditions, but even if not, the variance for between-subjects conditions will just drop to zero because there is no within-subject variance for something measured between subjects). This also holds for your length and frequency  nuisance parameters with the extra note that just because something is a nuisance parameter, doesn't mean it's a random effect. Finally, it's pretty weird to have an interaction without a main effect. So taking all that into account, you have a model like this:

m <- lmer(logRT ~ cond * length * frequency + (1+cond|subj) + (1|item_ID), dat, REML=F)

(Because length and frequency are continuous predictors, your coding scheme doesn't have to take interactions in to account, which is the usual motivation for using sum encoding in factorial designs. However, you may want to center or even scale the continuous predictors.)

Since length and frequency differ between items, I would not include them in the by-item RE-structure. The effect of length and especially frequency can differ between subjects (in part because no two people have the same language experience), but adding those to your by-subject RE-structure would make your model much more complex than you need and probably not converge anyway with typical psycholinguistic sample sizes. 

You may run into a problem with the strong correlation between length and frequency (Zipf's law), depending on how carefully you controlled your experiment. Whether you should then drop length or frequency from your model first is a decision you should ideally make based on theoretical grounds.  It's not clear whether this is self-paced reading or eye-tracking, but that could also play a role in your decision.

Finally, in terms of reaction times and the log transform. It's very common to log transform reaction times not just because of the skew but also because of what the log does -- it provides a self-adjusting "zoom" function. For short reaction times, a small difference is a big deal and the log lets that through (e.g. 300 vs 400ms is a big difference). But for long reaction times, a small difference is not a big deal and the log transforms that away, making it so that only big differences come through (e.g. the difference between 2000 and 21000ms is usually not particularly interesting and that difference will be minuscule on the log scale). The assumption of normality is also on the *conditional* distribution (or equivalently, the residuals) and not on the *marginal* (overall) distribution .


PS: Try not to call your data frame "data" as it shadows the builtin "data()" function and this can lead to really weird errors. For example, if you restart your R session and forget to load the data frame before running the lmer() call, then you'll get a weird error message "object of type 'closure' is not subsettable".

> On 7 Apr 2017, at 17:29, Diana Michl <dmichl at uni-potsdam.de> wrote:
> Dear all,
> I'd be very thankful for a short feedback on whether my lmer-model is 
> good the way it is. Thank you very much in advance!
> I'm fitting log-10 reaction times (RT), reading times of sentences (s) 
> of participants. The sentences are constructed in pairs: 1 test s. - 1 
> control s. (as similar as possible to the test s.). The TEST sentences 
> are in 4 conditions, i.e. they belong to one of four KINDS of sentences.
> Sentences are matched in length (number of letters) and frequency 
> (frequency of words averaged across sentence). Matching across ALL 
> material is attempted but not quite possible, so the matching is done as 
> well as possible within each test and control sentence pair.
> First, I treated the material as 5 conditions: 1 control + 4 test 
> conditions -> "cond5" (with contrast.treatment, i.e. control = dummy).
> "Cond2| subj" means: random intercepts for test/control condition per 
> subject. (Random slopes don't improve fit and aren't that important to me.)
> The converging model with the best fit (AIC/BIC) is:
> m1 <- lmer(logRT ~ cond5 + rat3 + (cond2|subject) + 
> (length:frequency|item_ID), data, REML=F)
> But I'm going back and forth on whether "cond5" is actually ok - or 
> whether I need to compare the sentences ONLY pairwise, like: each of the 
> 4 conditions (or "kinds") only in test vs. control.
> So I fit another model with a nested fixed effect. "cond4/cond2" means: 
> 4 conditions ("kinds") within test or control condition (with contr.sum, 
> so sum contrasts).
> The fit (AIC/BIC) is only a bit worse (910 vs. 905 above):
> m2 <- lmer(logRT ~ cond4/cond2 + rat3 + (cond2|subject) + 
> (length:frequency|item_ID), data, REML=F)
> The results of the two models are the same. The fit is almost the same. 
> My question is: Is one model more legitimate than the other, based on 
> what I describe about the matching?
> And is there something objectionable about how I fit the random effects? 
> (My model doesn't necessarily converge when I change them.)
> And I'm not 100% sure about whether I have to use logarithmized reaction 
> times (=response variable). The raw reaction times have huge variance 
> and are right-skewed. The log distribution is more normal, but still not 
> a normal distribution. The model fit seems much better for the logRT, 
> but the results differ slightly for logRT vs. raw RT.
> Again thank you very much and I'd really appreciate your feedback. I 
> hope I provided sufficient and non-confusing information.
> Best
> Diana
> -- 
> Diana Michl, M.A.
> PhD candidate
> International Experimental
> and Clinical Linguistics
> Universität Potsdam
> www.ling.uni-potsdam.de/staff/dmichl
> 	[[alternative HTML version deleted]]
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

More information about the R-sig-mixed-models mailing list