[R-meta] Clarification on ranef.rma.mv()

[Viechtbauer, Wolfgang (SP)]

Hi Luke,

Just to conclude (as far as I am concerned) this thread:

I don't have a good understanding of what taking the SVD of a Cholesky decomposition is doing in the first place. I helped to show how you can obtain this for rma.mv models, but I can't help you with making sense of this.

Best,
Wolfgang

>-----Original Message-----
>From: Luke Martinez [mailto:martinezlukerm using gmail.com]
>Sent: Friday, 17 September, 2021 0:04
>To: Viechtbauer, Wolfgang (SP)
>Cc: R meta
>Subject: Re: [R-meta] Clarification on ranef.rma.mv()
>
>Hi Wolfgang,
>
>A quick follow-up on estimating the proportion of between variance in
>metafor as in lme4:::rePCA.merMod(). We discussed using it in
>correlated random-effects, but can we use that for non-correlated
>(varying intercept) models as well?
>
>Assuming yes, it may be that I've a bug in the code below, but even
>though the "sigma2" for "paper/study" is 0.000 and sigma2 for
>"paper/study/obs" is 0.0037, at the end the POV_S for "paper/study" is
>larger than that for "paper/study/obs", not sure why?
>
>dat <- dat.bornmann2007
>dat <- escalc(measure="OR", ai=waward, n1i=wtotal, ci=maward,
>n2i=mtotal, data=dat)
>dat$paper <- as.numeric(factor(dat$study))
>dat$paper[dat$paper <= 2] <- 1
>fit <- rma.mv(yi, vi, random = ~ 1 | paper/study/obs, data=dat)
>
># Apply the pca approach:
>
>    round(S <- fit$sigma2, 4)
>   #[1]  0.0157  0.0000  0.0037
>    S <- diag(S)
>    colnames(S) <- rownames(S) <- fit$s.names
>    sds <- setNames(svd(chol(S))$d, colnames(S))
>    (pov_S <- round(sds^2 / sum(sds^2), digits = 4))
>
>paper     paper/study    paper/study/obs  ## Does these results make sense?
>0.8077          0.1923          0.0000
>
>On Wed, Sep 15, 2021 at 5:04 PM Luke Martinez <martinezlukerm using gmail.com> wrote:
>>
>> Dear Wolfgang,
>>
>> Thank you so very much! Especially thank you for reminding me
>> regarding the effect of data dependent modifications!
>>
>> All the best,
>> Luke
>>
>> On Wed, Sep 15, 2021 at 1:21 PM Viechtbauer, Wolfgang (SP)
>> <wolfgang.viechtbauer using maastrichtuniversity.nl> wrote:
>> >
>> > Please see below for my comments.
>> >
>> > >-----Original Message-----
>> > >From: Luke Martinez [mailto:martinezlukerm using gmail.com]
>> > >Sent: Tuesday, 14 September, 2021 19:25
>> > >To: Viechtbauer, Wolfgang (SP)
>> > >Cc: R meta
>> > >Subject: Re: [R-meta] Clarification on ranef.rma.mv()
>> > >
>> > >Hi Wolfgang,
>> > >
>> > >Thank you. And since in rma.mv() we can have up to two ~ inner | outer
>> > >random terms, then, I'm assuming to get the proportion of variation
>> > >for the second ~ inner | outer random term, I can do:
>> > >
>> > >sds <- svd(chol(rma.mv_model4$H))$d
>> > >sds^2 / sum(sds^2)
>> >
>> > Correct.
>> >
>> > >I guess one potential problem I'm running into is that what should we
>> > >do if we see that the proportion of explained between-studies variance
>> > >by only one or two levels of a categorical variable is almost zero
>> > >while rest of the levels of that categorical variable make significant
>> > >contributions?
>> > >
>> > >The reason I ask this is that with continuous variables (using struct
>> > >= "GEN"), if a variable's contribution is almost zero, then, you can
>> > >decide not to use that continuous variable in the random part at all
>> > >(that variable altogether is overfitted).
>> > >
>> > >But with categorical variables, when several levels make good
>> > >contributions to the between-studies variance except just one or two
>> > >levels, then, you can't easily decide not to use that whole
>> > >categorical variable in the random part at all.
>> > >
>> > >Do you have any opinion on this dilemma?
>> >
>> > I would choose a random effects structure that is motivated by the structure
>of the data and the possible sources of heterogeneity/variability/dependencies
>that I think may exist in the data. For example, if a slope may vary across
>units, then I would add a random effect for that slope to the model. If it turns
>out that the estimate of the slope variability is very low, dropping that random
>effect (which is the same as assuming that the variance is 0) or not will pretty
>much do the same thing. My preference would be not to change the model, because I
>generally try to avoid making changes to a model (since the consequences of such
>data dependent modifications are hard to predict).
>> >
>> > The same applies to structures like "UN". If a particular tau^2 value is low,
>then no, I would not drop that random effect. You *can* however set particular
>tau^2 values to 0 (the 'tau2' argument of rma.mv() allows you to do that), but
>again, I would personally avoid doing that.
>> >
>> > Best,
>> > Wolfgang