[R-sig-ME] var(ranef(Random Effect)) not the same as the variance component
Harold Doran
h@ro|d@dor@n @end|ng |rom c@mb|um@@@e@@ment@com
Wed Sep 9 22:45:21 CEST 2020
Simon
Your original question was, “A very basic question. Generally, `var(ranef(Random Effect))` may not necessarily be the same as the variance component reported for that Random Effect in the model output, correct?”
First, let me define terms since there are a lot of variances in this model. There is the “marginal variance” which is the variance of the random effects reported by lmer. Then there is the “conditional variance of the random effects”, which is the variance of each individual BLUP, then you can take the “observed variance” over all BLUPs.
Now, take the simple model used in my example,
fm1 <- lmer(Reaction ~ 1 + (1|Subject), sleepstudy)
So, this term, var(ranef(fm1)$Subject), is the observed variance of the BLUPs.
The marginal variance of the random effects in this model reported by lmer is 35.75^2 and as Dimitris pointed out, you would not expect the observed variance of the conditional means (i.e., the BLUPs) to be equal to the marginal variance of the random effects.
The random effects are conditional means, E(Y|X) and they have a conditional variance. So:
E(var(Y|X)) = mean(sapply(attr(ranef(fm1)$Subject, "postVar"),function(x) x))
var(E(Y|X) = observed variance of the BLUPs
Combined, their sum is the marginal variance of the random effects. So, concretely, we observe that
> sqrt(var(ranef(fm1)$Subject) + mean(sapply(attr(ranef(fm1, condVar=TRUE)[[1]], "postVar"),function(x) x)))
(Intercept)
(Intercept) 35.75385
Which is the same as the marginal SD of the random effects from lmer.
Hope this helps.
BTW, note to Ben Bolker (if he even read this far 😊 ), I was originally the one who suggested to Doug that there be an extractor called “postVar”. He regretted that term and expressed that to me many times (now I understand why). But, I had to remember how to get the conditional variances of the random effects and the code above is what I came up with since they are an attribute of the ranef extractor. That code is a bit “ugly”, so either I’m doing it wrong (big possibility) or might I suggest a new extractor function can be written that is easier to get those?
From: Simon Harmel <sim.harmel using gmail.com>
Sent: Wednesday, September 9, 2020 4:19 PM
To: Harold Doran <harold.doran using cambiumassessment.com>
Cc: D. Rizopoulos <d.rizopoulos using erasmusmc.nl>; r-sig-mixed-models <r-sig-mixed-models using r-project.org>; Ben Bolker <bbolker using gmail.com>
Subject: Re: [R-sig-ME] var(ranef(Random Effect)) not the same as the variance component
Thank you, Harold.
1) `var(ranef(fm1)$Subject)` is the posterior variance of random effects or their prior variance? (this was the point discussed so far in this thread)
2) Also, what about `mean(sapply(attr(ranef(fm1)$Subject, "postVar"),function(x) x))`, what this expected variance represents in words?
3) What does their sum represent? The total observed variance in random deviations in intercepts of subjects?
Thank you very much
On Wed, Sep 9, 2020 at 1:48 PM Harold Doran <harold.doran using cambiumassessment.com<mailto:harold.doran using cambiumassessment.com>> wrote:
Simon
Here is an example to show what my notation implies with respect to your question:
fm1 <- lmer(Reaction ~ 1 + (1|Subject), sleepstudy)
sqrt(var(ranef(fm1)$Subject) + mean(sapply(attr(ranef(fm1, condVar=TRUE)[[1]], "postVar"),function(x) x)))
From: Simon Harmel <sim.harmel using gmail.com<mailto:sim.harmel using gmail.com>>
Sent: Wednesday, September 9, 2020 11:53 AM
To: D. Rizopoulos <d.rizopoulos using erasmusmc.nl<mailto:d.rizopoulos using erasmusmc.nl>>
Cc: r-sig-mixed-models <r-sig-mixed-models using r-project.org<mailto:r-sig-mixed-models using r-project.org>>; Harold Doran <harold.doran using cambiumassessment.com<mailto:harold.doran using cambiumassessment.com>>; Ben Bolker <bbolker using gmail.com<mailto:bbolker using gmail.com>>
Subject: Re: [R-sig-ME] var(ranef(Random Effect)) not the same as the variance component
First of all, thank you all for your valuable input.
Dimitris,
Thank you I upvoted your answer on CV as well. But please help me understand a few things.
1- By D matrix, you mean the G matrix shown in: https://bookdown.org/marklhc/notes/simulating-multilevel-data.html#linear-growth-model
2- When you say variance components in the output are prior values, can you tell me how these prior values are obtained? I guess from the data itself, but how exactly (do we run individual models first to see how much intercepts and slopes vary & co-vary and take those as prior)?
3- Harold above noted that: "The conditional means of the random effects are E(Y|X) and hence their variance is only one portion of the total variance [i.e., var(y)]." I'm not sure how this directly relates to my question in this thread?
Thank you,
Simon
On Tue, Sep 8, 2020 at 3:47 PM Harold Doran <harold.doran using cambiumassessment.com<mailto:harold.doran using cambiumassessment.com>> wrote:
To add a little notation to this, we can use law of total variance, var(y) = E(var(Y|X)) + var(E(Y|X)). The conditional means of the random effects are E(Y|X) and hence their variance is only one portion of the total variance.
-----Original Message-----
From: R-sig-mixed-models <r-sig-mixed-models-bounces using r-project.org<mailto:r-sig-mixed-models-bounces using r-project.org>> On Behalf Of D. Rizopoulos
Sent: Monday, September 7, 2020 11:02 PM
To: Simon Harmel <sim.harmel using gmail.com<mailto:sim.harmel using gmail.com>>; r-sig-mixed-models <r-sig-mixed-models using r-project.org<mailto:r-sig-mixed-models using r-project.org>>
Subject: Re: [R-sig-ME] var(ranef(Random Effect)) not the same as the variance component
External email alert: Be wary of links & attachments.
Yes, you do not expect these two be the same. The variance components are the prior variances of the random effects, whereas var(ranef(model)) is the variance of the posterior means/modes of the random effects.
Best,
Dimitris
Dimitris Rizopoulos
Professor of Biostatistics
Erasmus University Medical Center
The Netherlands
________________________________
From: R-sig-mixed-models <r-sig-mixed-models-bounces using r-project.org<mailto:r-sig-mixed-models-bounces using r-project.org>> on behalf of Simon Harmel <sim.harmel using gmail.com<mailto:sim.harmel using gmail.com>>
Sent: Tuesday, September 8, 2020 1:22:15 AM
To: r-sig-mixed-models <r-sig-mixed-models using r-project.org<mailto:r-sig-mixed-models using r-project.org>>
Subject: [R-sig-ME] var(ranef(Random Effect)) not the same as the variance component
Hello All,
A very basic question. Generally, `var(ranef(Random Effect))` may not necessarily be the same as the variance component reported for that Random Effect in the model output, correct?
Thank you all,
Simon
[[alternative HTML version deleted]]
_______________________________________________
R-sig-mixed-models using r-project.org<mailto:R-sig-mixed-models using r-project.org> mailing list
https://eur01.safelinks.protection.outlook.com/?url=https%3A%2F%2Fstat.ethz.ch%2Fmailman%2Flistinfo%2Fr-sig-mixed-models&data=02%7C01%7Cd.rizopoulos%40erasmusmc.nl%7Cae4132330fbd412dea9508d85384efba%7C526638ba6af34b0fa532a1a511f4ac80%7C0%7C0%7C637351177705805772&sdata=agUgmCzM5ecsaoGLm8aPX0%2FuHZF1mK%2BXbP%2Fi6KX5UvI%3D&reserved=0
[[alternative HTML version deleted]]
[[alternative HTML version deleted]]
More information about the R-sig-mixed-models
mailing list