[R-meta] Plotting the correlation among true/random effects across categories

Hi Yuhang,

1) In the model below, 'outcome' turns into a dummy variable that is 1 for
PD and 0 for AL. So the intercept random effect represents variation in the
treatment effects for AL, while the outcomePD random effect represents
variation in how much the treatment effects differ for PD compared to AL.
The correlation is then the correlation between these two random effects.

2) That depends on what you are interested in. But if you want to know how
much variance there is in a random effect, I would take the variance
estimate (and not compute the variance of the BLUPs) and similarly if you
are interested in the correlation between two random effects.

Best,
Wolfgang

-----Original Message-----
From: Yuhang Hu [mailto:yh342 at nau.edu]
Sent: Friday, 16 June, 2023 17:06
To: Viechtbauer, Wolfgang (NP)
Cc: R Special Interest Group for Meta-Analysis
Subject: Re: [R-meta] Plotting the correlation among true/random effects

across

categories

Hi Wolfgang,

Thank you very much for your response. I wanted to ask two

questions regarding

your responses.

1-- If instead of struct="UN", I use struct="GEN" (below), then what does

the

correlation reported between AL and PD represent?

My understanding is that it represents the correlation between the

'differences'

(not between the means of AL and PD) between AL and PD across the trials,

right?

model <- rma.mv(yi~ outcome, vi, data =  dat.berkey1998,
               random = ~ outcome | trial, struct = "GEN")

2-- You mentioned that the estimated variation by the model (say tau2)

accounts

for an additional piece (i.e., E(var(u_i|y_i))) in random effects (ui)

that is

not present in the BLUPs.

My question is then, is the variation, or correlation obtained by the

model

preferred over (and more reliable than) that obtained by hand-calculating

them

from the BLUPs?

Many thanks,
Yuhang

On Thu, Jun 1, 2023 at 8:18?PM Viechtbauer, Wolfgang (NP)
<wolfgang.viechtbauer at maastrichtuniversity.nl> wrote:
You get the estimated variation in the underlying true means of each

category

(around the fixed effects, which correspond to the expected value of the
underlying true means of each category) and the correlation of that

variation

between the two categories. I have now added the model to the write-up of

this

example here:

https://www.metafor-project.org/doku.php/analyses:berkey1998

so you can see exactly what the code (given on that page) corresponds to

model-

wise.

Best,
Wolfgang

-----Original Message-----
From: Yuhang Hu [mailto:yh342 at nau.edu]
Sent: Thursday, 04 May, 2023 23:39
To: Viechtbauer, Wolfgang (NP)
Cc: R Special Interest Group for Meta-Analysis
Subject: Re: [R-meta] Plotting the correlation among true/random effects

across

categories

Thank you very much, Wolfgang. Regarding my second (unclear) follow-up,

I might

have a conceptual misunderstanding that I hope to clear up.

When we use a "categorical" moderator (with struct="UN") to the left of

| (e.g.,

~outcome | trial_id), conceptually it means that in each unique

trial_id, we fit

"effect_size ~ outcome+0", and obtain the Mean for each outcome category

(say

two

categories). Then in the rma.mv() output, we get the variation in means

of each

category as well as the correlation between the means of the two

categories

across all unique trial_ids.

Is my conceptual understanding correct?

Thanks,
Yuhang

On Wed, May 3, 2023 at 1:59?PM Viechtbauer, Wolfgang (NP)
<wolfgang.viechtbauer at maastrichtuniversity.nl> wrote:
Please see below for responses.

Best,
Wolfgang

-----Original Message-----
From: Yuhang Hu [mailto:yh342 at nau.edu]
Sent: Wednesday, 03 May, 2023 22:48
To: Viechtbauer, Wolfgang (NP)
Cc: R Special Interest Group for Meta-Analysis
Subject: Re: [R-meta] Plotting the correlation among true/random

effects across

categories

Thank you very much, Wolfgang. Two quick follow-ups:

1) To convert these estimated random deviations to true effects, I

should add

the

fixed effect estimates (assuming I use 'outcome + 0' in model formula)

for each

category of outcome to its relevant column, right?

AL = paired[,1] + model$b[1,1]
PD = paired[,2] + model$b[2,1]
plot(PD~AL, pch=21, bg="gray", cex=1.5, lwd=1.2)

Correct.

2) When using categorical variables (with "UN") to the left of |, I

think we

drop

the intercept in the random-effects design matrix, so what is actually

allowed

to

vary across the trials given that eac trial has only one instance of AL

and PD

in

it?

Not entirely sure what you mean by this question. The help file explains

what

the

different structures do:

https://wviechtb.github.io/metafor/reference/rma.mv.html#specifying-random-

effects

Thank you for your time.
Yuhang

On Wed, May 3, 2023 at 12:51?AM Viechtbauer, Wolfgang (NP)
<wolfgang.viechtbauer at maastrichtuniversity.nl> wrote:
You can extract the BLUPs of the random effects and create a

scatterplot based

on

them:

re <- ranef(model)
re

paired <- do.call(rbind, split(re[[1]]$intrcpt, dat.berkey1998$trial))
paired

plot(paired, pch=21, bg="gray", cex=1.5, lwd=1.2)

And before somebody asks why cor(paired) does not yield 0.7752 (or why

the

values

in var(paired) do not match up with the variances as estimated from the

model),

see for example:

https://stats.stackexchange.com/q/69882/1934

Or let me give a more technical explanation based on the standard RE

model:

dat <- escalc(measure="RR", ai=tpos, bi=tneg, ci=cpos, di=cneg,

data=dat.bcg)

res <- rma(yi, vi, data=dat)
res$tau2
var(ranef(res)$pred)

You will notice that the latter is smaller than tau^2. By the law of

total

variance:

tau^2 = var(u_i) = E(var(u_i|y_i)) + var(E(u_i|y_i)).

The conditional means of the random effects (which is what ranef()

provides

estimates of) are E(u_i|y_i) and hence their variance is only part of

the total

variance. Therefore, the estimate of tau^2 and the estimated variance

of the

BLUPs of the random effects will not match up.

In more complex models, this then also affects things like the

correlation

between the BLUPs of the random effects.

Best,
Wolfgang

-----Original Message-----
From: R-sig-meta-analysis [mailto:

r-sig-meta-analysis-bounces at r-project.org]

On

Behalf Of Yuhang Hu via R-sig-meta-analysis
Sent: Wednesday, 03 May, 2023 1:21
To: R meta
Cc: Yuhang Hu
Subject: [R-meta] Plotting the correlation among true/random effects

across

categories

Hello Colleagues,

I was wondering if there is a way to scatterplot the correlation

between

the categories of variable "outcome" (AL and PD) which is estimated to

be

rho = .7752 in my model below?

model <- rma.mv(yi~ outcome, vi, data =  dat.berkey1998,
               random = ~ outcome | trial, struct = "UN")

   rho.AL  rho.PD    AL  PD
AL       1             -   5
PD  0.7752       1    no   -

Thanks,
Yuhang