[R-meta] Random effects structure

The Berkey code doesn't quite apply here. This is more relevant:

http://www.metafor-project.org/doku.php/analyses:gleser2009#multiple-treatment_studies

But your splitting variable is 'exp' within 'Site', so you need to first create a variable that is the combination of the two. This should do it:

dat$Site.exp <- paste0(dat$Site, ".", dat$exp)

Also make sure that the data are sorted accordingly:

dat <- dat[order(dat$Site, dat$exp),]

Then:

calc.v <- function(x) {
   v <- matrix(x$NPSD^2[1] / (x$NPN[1] * dat$NP[1]^2), nrow=nrow(x), ncol=nrow(x))
   diag(v) <- x$var
   v
}
V <- bldiag(lapply(split(dat, dat$Site.exp), calc.v))

You should double-check that the blocks in V are correct (with 'var' along the diagonal and the correct covariances along the off-diagonal) before trying to fit the model.

Best,
Wolfgang

-----Original Message-----
From: C?sar Terrer [mailto:cesar.terrer at me.com] 
Sent: Wednesday, 25 March, 2020 16:21
To: Viechtbauer, Wolfgang (SP)
Cc: r-sig-meta-analysis at r-project.org
Subject: Re: [R-meta] Random effects structure

Hi Wolfgang,

Thank you for your quick response.

I have never constructed a V matrix. Based on the Berkey tutorial, I have tried to apply the following

dat$v2i <- dat$NPSD^2 / (dat$NPN*dat$NP^2)
V <- bldiag(lapply(split(dat[,c("var", "v2i")], dat$exp), as.matrix))
res <-rma.mv(yi=lnR, V=V, random=~1|Site/exp/obs, mods=~Predator, data=dat)

But it gives me an error that? 'V' must be a square matrix. Is this because the split in V has to be done at the "obs" level? I used "exp" because that indicates the level at which all effect sizes use the same control.

Thank you
Cesar

Hi Cesar,

When the same group is used to compute multiple estimates (i.e., a common control), then this also induces dependency on the sampling errors. For the log response ratio, the covariance is then:

SD_C^2 / (n_C*M_C^2)

where M_C and SD_C are the mean and SD of the common control group and n_C the control group size.

So, you ideally should construct a proper V matrix that includes these covariances and that you can then pass to rma.mv().

But yes, random=~1|Site/exp/obs would be sensible.

Best,
Wolfgang

-----Original Message-----
From: R-sig-meta-analysis [mailto:r-sig-meta-analysis-bounces at r-project.org] On Behalf Of C?sar Terrer
Sent: Wednesday, 25 March, 2020 15:34
To: r-sig-meta-analysis at r-project.org
Subject: [R-meta] Random effects structure

Dear community,

I am conducting a meta-analysis to study the growth rate of bacterial predators as compared to their prey, using the log response ratio. Furthermore, I want to study if this effect varies across different predators. The dataset has the following structure, here showing a subset:

Site CommonControl exp obs Predator lnR var
A Alaska 155 1 1 Bdello -0.6713152 0.03785708
A Alaska 155 1 2 Cytoph -0.0702467 0.05763364
A Alaska 155 1 3 Myxo -0.148982 0.00748768
A Alaska 1510 2 4 Bdello -0.4926361 0.01691187
A Alaska 1510 2 5 Cytoph -0.213787 0.01045785
B Andesite1controlWeek1 9 6 Bdello 0.27873598 0.14129722
B Andesite1controlWeek1 9 7 Cytoph -0.3243682 0.01466085
B Andesite1controlWeek1 9 8 Lyso 1.18302506 0.11663149
B Andesite1controlWeek6 11 9 Bdello -0.8465128 0.03701618
B Andesite1controlWeek6 11 10 Cytoph -0.1559056 0.0283173
B Andesite1controlWeek6 11 11 Lyso -0.8039415 0.04926915

1. There are different sites, thus a potential source of non-independency
2. Within each site,?we use the value for preys in?the denominator multiple times. I guess rows of data using the same denominator (CommonControl) are also potentially correlated and should be also added as a random-effect.

Based on 1., 2., and what I have understood from the?Konstantopoulos (2011) tutorial, I think I should use the following model:

res <-rma.mv(yi=lnR, V=var, random=~1|Site/exp/obs, mods=~Predator, data=data)

Could you please let me know if the structure of random effects seems appropriate, and help me understand why I need to include "obs"?
Thank you.
Cesar