2
$\begingroup$

I am working (in R 3.2.3 using lme4 for doing mixed effect modeling) with vowel data from many different subjects. The question I'm interested is whether certain vowels undergo change over time (as in the English word 'I') as a function of the following consonant type, and in particular whether this is the case for all subjects, or only for some. Suppose, for simplicity, that I am only measuring a single vowel in many different words, and that there are no relevant differences between these words. My model would then looks as follows:

measure ~ (1|subject) + timepoint + (timepoint|subject) + 
          following_segment + (following_segment|subject) + 
          timepoint:following_segment + (timepoint:following_segment|subject)

where timepoint is a continuous variable referring to the timepoint in milliseconds of each measurement, and following_segment is a binary variable referring to either a neutral context or a context where I am expecting a change over time.

If I run ranef on this model and look in the final column, corresponding to the (timepoint:following_segment|subject) term, I indeed get either very large numbers (e.g. 500 Hz) or very small numbers (e.g. 100 Hz, just to name two random examples; there is a large amount of variability, because these estimates are correlated with the subject's physiological properties such as jaw size).

I now want to know whether these differences in random slope estimates between subjects are significant, for individual subjects (i.e. 'is this subject sensitive to the following_segment manipulation or not?'). Since I know that t = B/SE, all I believe I would need is the standard error for the individual participant within this specific term of the model, and then I can use pnorm to get a p value (please correct me if this is too simple). Regarding obtaining the estimated SE, this question gave me a starting point in suggesting summary(model)@REmat, which unfortunately returns only NULL. Another question suggests using ranef(model,condVar=TRUE), which unfortunately gives an error that conditional variances not currently available via ranef when there are multiple terms per factor. Is there another way to test whether the individual random slope estimates are significantly different from 0?

Here is a histogram of the random slopes:

enter image description here

The peak around 800 Hz looks interesting; any way I can test for each participant's slope estimates whether they are located around that 800-Hz peak, versus around the big peak near the 0 point?

|cite|improve this question
$\endgroup$
  • $\begingroup$ Although asked in terms of R, I think this is a statistical question & merits staying open here. $\endgroup$ – gung - Reinstate Monica♦ Feb 20 '16 at 13:46
  • $\begingroup$ You say that you get either very large or very small random effects. Are they bimodal or just high variance? Do the high & low groupings (clusters?) correspond to anything (eg, native speakers vs ESL)? Note that random effects are modeled w/ a normal distribution (which goes to infinity in both directions), so the model assumes there have to be people who are 0. $\endgroup$ – gung - Reinstate Monica♦ Feb 20 '16 at 13:53
  • $\begingroup$ High variance: certainly (because the data is correlated with the subject's physiology). Bimodal & are there meaningful clusters: I suppose that's what my research question is - are there indeed people who are/are not sensitive to the manipulation, and if so, does this sensitivity vs nonsensitivity correlate with (in my case) geographical regions? A plot does appear to show some bimodality, but whether this is statistically meaningful and which subjects may be reasonably certainly grouped into which peak, that is exactly what I would like to figure out. $\endgroup$ – LostLinguist Feb 20 '16 at 14:48
  • $\begingroup$ Do you have data on each person's geographical region? $\endgroup$ – gung - Reinstate Monica♦ Feb 20 '16 at 15:01
  • $\begingroup$ Yes. I have 160 subjects hailing from 8 regions, all properly identified in the variable coding (i.e. my data frame contains a field 'region' with 8 levels). I was afraid to enter region as a predictor, however, because it should (by design) be completely collinear with the variable 'subject' (which I do still need to control for individual variation). This is why I am now trying to find some kind of structure within the individual's random slopes, and I was planning to then test posthoc whether such structure, if found, would happen to correlate with the regions the individuals would be from. $\endgroup$ – LostLinguist Feb 20 '16 at 15:39
2
$\begingroup$

The ranef(model,condVar=TRUE) method is not working because you have the random effects written in separate parenthesis blocks. In this specification, no covariances are estimated among the random slopes (because covariances are only included for terms that share a parenthesis block), but covariances are estimated between the the random intercept and each random slope (because each parenthesis block implicitly also includes the intercept). That's fine, but another totally reasonable specification (although with more parameters) is to allow all the random effects to have nonzero covariance. In other words, you can rewrite the model syntax as

measure ~ timepoint*following_segment + (timepoint*following_segment|subject)

For this rewritten model, the ranef(model,condVar=TRUE) method should work.

If you don't want to rewrite the model in this way and include the extra covariance parameters, then another way to implement the tests you want is via bootstrapping. In each iteration, sample with replacement from each subject's data (do not mix observations between subjects), fit your mixed model, and save the array of random effects. After some large number of iterations (1000? 10000?), you will have distributions for each subject's random effects, which you can compare to 0. Note that this could take a while for your complex model!

|cite|improve this answer
$\endgroup$
  • $\begingroup$ I've accepted your answer because this was indeed the problem with the condVar option, thank you! As a follow-up question: would it now be valid to take the conditional variance for the parameters I'm interested in, to take that value's square root, and to divide the individual subject's random-effect estimate by the resulting value, yielding a t-like statistic that can be checked against the normal distribution for significance? (In addition, am I correct in presuming that if I do this, I will need to correct for multiple comparisons, because I have 160 subjects to test?) $\endgroup$ – LostLinguist May 25 '16 at 11:42
  • $\begingroup$ I think there is nothing technically incorrect about doing significance tests on the individual random effects in this way. There's of course going to be a big multiple comparisons issue, but it sounds like you've already thought of that. In my opinion, testing the individual effects would be a little odd -- if we really do have a priori interest in whether individuals' effects differ from 0 (or some other value), then the individuals should arguably be treated as a fixed rather than random factor -- but, like I said, I don't think it's technically incorrect. $\endgroup$ – Jake Westfall May 25 '16 at 13:51
  • $\begingroup$ The issue is that I don't a priori know which subjects, if any, are going to be the ones exhibiting significantly different behavior - I have a guess that it's correlated with their geographical region, but I don't feel that that guess is solid enough to warrant including region as a fixed effect a priori (even though I could do that). In other words, it seems more objective to me to let the data (of the individual subjects) itself tell me which subjects fall into my target group, and to check only post hoc whether that happens to correlate with e.g. region. Does that make statistical sense? $\endgroup$ – LostLinguist May 25 '16 at 14:35
  • $\begingroup$ @LostLinguist Yeah, it does $\endgroup$ – Jake Westfall May 25 '16 at 14:55

Your Answer

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

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