[R-meta] Problem with sigma-square in 3-level meta-analysis

[Viechtbauer, Wolfgang (NP)]

Like I said, the peak may also be at zero, which I assume is the case here and which is perfectly fine.

Best,
Wolfgang

>-----Original Message-----
>From: R-sig-meta-analysis [mailto:r-sig-meta-analysis-bounces using r-project.org] On
>Behalf Of Röhl, Sebastian via R-sig-meta-analysis
>Sent: Friday, 10 March, 2023 15:43
>To: R Special Interest Group for Meta-Analysis
>Cc: Röhl, Sebastian
>Subject: Re: [R-meta] Problem with sigma-square in 3-level meta-analysis
>
>Dear Michael and Wolfgang,
>
>thank you very much for your answers.
>The profile plot for sigma-2 between studies does not show a peak. It is simply
>decreasing...
>Do you have a suggestion which alternative approach or model could be used if we
>just want to report overall effects?
>
>Best,
>Sebastian
>
>-----Ursprüngliche Nachricht-----
>Von: R-sig-meta-analysis <r-sig-meta-analysis-bounces using r-project.org> Im Auftrag
>von Michael Dewey via R-sig-meta-analysis
>Gesendet: Freitag, 10. März 2023 14:27
>An: R Special Interest Group for Meta-Analysis <r-sig-meta-analysis using r-
>project.org>
>Cc: Michael Dewey <lists using dewey.myzen.co.uk>
>Betreff: Re: [R-meta] Problem with sigma-square in 3-level meta-analysis
>
>Dear Sebastian
>
>Just to reinforce Wolfgang's advice whenever I have had something unexpected
>happen with rma.mv() then profile() has usually been very illuminating. Sadly it
>usually tells me I did not have enough data points but life is like that.
>
>Michael
>
>On 10/03/2023 09:19, Viechtbauer, Wolfgang (NP) via R-sig-meta-analysis
>wrote:
>> Dear Sebastian,
>>
>> Just a general note when posting to this list: Please do not reply to a message
>and then change the subject to something new. There is information stored in the
>email header that is used to create threads in email clients and the mailing list
>archives and when doing this, it messes things up. See
>https://stat.ethz.ch/pipermail/r-sig-meta-analysis/2023-March/thread.html and
>note how your message is considered a reply to the 'margins command' thread,
>which is obviously not what you intended. (I am making this reply in an entirely
>new post, so this should hopefully show up in a proper new thread, although of
>course your post won't show up at the beginning of this thread, so let me just
>put the archive URL here: https://stat.ethz.ch/pipermail/r-sig-meta-
>analysis/2023-March/004447.html).
>>
>> Moving on to your actual questions:
>>
>> (a) The observed variances are a combination of the variance in the underlying
>true effects (apparently r-to-z transformed correlations based on 'zr') and the
>sampling variances and hence cannot be directly used to examine whether something
>has gone wrong here.
>>
>> I would examine profile(overall_mod_cc) to see if the likelihood profile plots
>are properly peaked at the estimates (which could very well also happen around
>0).
>>
>> (b) I wouldn't put it this way, but yes, in general (and this is not specific
>to rma.mv()) it is difficult to estimate variance components when the number of
>levels is small. The REML estimates should still be approximately unbiased, but
>will then tend to high variance themselves.
>>
>> (c) Again, I wouldn't automatically conclude here that something is off. It
>might very well be that the variance in the underlying true effects is due to
>difference between rows within samples and not between samples.
>>
>> Speaking of that -- are different rows for the same level of 'Sample_ID' based
>on the same subjects? In this case, those (r-to-z transformed) correlations are
>not independent and should not be treated as such. rcalc() is a function that
>allows to construction of the appropriate var-cov matrix to account for this (but
>requires information about correlations that may not be available).
>Alternatively, vcalc() could be used to construct a rough approximation.
>>
>> And finally -- if you still think something is off here, a possible way forward
>would be to use a Bayesian model where you can use priors on the variance
>components to push their posterior distributions in the desired direction (more
>or less, depending on how much information there is in the data).
>>
>> Best,
>> Wolfgang
>>
>>> -----Original Message-----
>>> From: R-sig-meta-analysis
>>> [mailto:r-sig-meta-analysis-bounces using r-project.org] On Behalf Of Röhl,
>>> Sebastian via R-sig-meta-analysis
>>> Sent: Friday, 10 March, 2023 8:40
>>> To: R Special Interest Group for Meta-Analysis
>>> Cc: Röhl, Sebastian
>>> Subject: [R-meta] Problem with sigma-square in 3-level meta-analysis
>>>
>>> Dear all,
>>>
>>> We are conducting meta-analyses of correlational effects which are
>>> clustered in studies (or independent samples). Therefore we are using
>>> 3-level hierarchical random effects models (Sampling error within
>>> effect sizes within independent samples).
>>> For one of our sub-analyses which 50 effects of 8 independent
>>> samples, the results of the rma.mv pointed to nearly no
>>> (5.710374e-12) variance component between the samples:
>>>
>>>> overall_mod_cc <- rma.mv(zr, V=var, random = ~ 1| Sample_ID /
>>>> nummer, data =
>>> data_cc, method = "REML")
>>>> overall_mod_cc
>>> Multivariate Meta-Analysis Model (k = 50; method: REML) Variance
>>> Components:
>>>             estim    sqrt  nlvls  fixed            factor
>>> sigma^2.1  0.0000  0.0000      8     no         Sample_ID
>>> sigma^2.2  0.0222  0.1489     50     no  Sample_ID/nummer
>>>
>>> However, when I look directly at the manifest variances for the group
>>> means
>>> (0.225) and all effects (0.392), I conclude that there should
>>> actually be a much larger variance component between the independent samples.
>>> (a) Am I making a thinking error here?
>>> (b) Or could this be because the estimation of variance components
>>> using rma.mv() does not work reliably for the small number of studies (8)?
>>> (c) And if (b) is correct - what model or function could I use
>>> instead for this partial analysis?
>>>
>>> Thanks a lot for your help!
>>>
>>> Best,
>>> Sebastian