Dear all, Dr Viechtbauer, Dr Del Ponte,I read what you suggested, and went through the archives of the mailing list, and even though some threads were dealing with similar questions, I didn't really find the answers I was looking for there.I have multiple questions : 1) How exactly is calculated the weight of each estimate in the function rma.mv() ?2) (this question might be answered too when answering the first one but...) Why does that function (rma.mv) attribute more weight to estimates coming from studies that have multiple estimates relatively to ones that come from studies that only have one estimate ?3) I tried to produce a reproducible example of my code, would it be possible to tell me if I got things right ?my dataset is this table (I didn't know if I could attach anything) :
dat
author year estimate std_error p_value r number_of_obs
1 Ballinger et al 1988 -19.3450 5.8220 4.638e-03 -0.6511331 17
2 Matthieu & Norman 2020 -2.7398 0.8103 7.424e-04 -0.0920000 1352
3 Khangura 2011 -1.9610 1.3470 1.536e-01 -0.2298655 40
4 Kutcher 1990 -21.1350 3.9330 1.420e-05 -0.7321321 27
5 Steed 2007 -13.9700 3.2500 7.355e-04 -0.7542836 16
6 Steed 2007 -31.9440 4.7650 2.765e-04 -0.9301781 9
7 Steed 2007 -4.2780 3.0940 1.883e-01 -0.3466844 16
8 Sprague et al 2010 -21.9880 5.3010 7.562e-04 -0.7198487 18
9 Upadhaya et al 2019 -43.6170 3.3440 1.133e-06 -0.9772861 10
10 Khangura et al 2005 -5.3500 NA 5.000e-02 -0.4200000 29
11 Khangura et al 2005 -10.5600 NA 1.000e-03 -0.5700000 29
12 Khangura et al 2005 -9.9700 NA 1.300e-02 -0.4600000 29
13 Khangura et al 2005 -5.4500 NA 1.100e-02 -0.4700000 29
14 Khangura et al 2005 -22.7500 NA 2.800e-02 -0.4200000 29
15 Khangura et al 2005 -16.8300 NA 2.200e-02 -0.4300000 29
16 Khangura et al 2005 -9.2100 NA 3.900e-02 -0.3900000 29
V1 <- c("Ballinger et al", "Matthieu & Norman", "Khangura", "Kutcher", "Steed", "Steed", "Steed", "Sprague et al",
"Upadhaya et al", "Khangura et al", "Khangura et al", "Khangura et al", "Khangura et al", "Khangura et al",
"Khangura et al", "Khangura et al" )
V2 <- c(1988, 2020, 2011, 1990, 2007, 2007, 2007, 2010, 2019, 2005, 2005, 2005, 2005, 2005, 2005, 2005)
V3 <- c(-19.3450, -2.7398, -1.9610, -21.1350, -13.9700, -31.9440, -4.2780, -21.9880, -43.6170, -5.3500,
-10.5600, -9.9700, -5.4500, -22.7500, -16.8300, -9.2100)
V4 <- c(5.8220, 0.8103, 1.3470, 3.9330, 3.2500, 4.7650, 3.0940, 5.3010, 3.3440, NA, NA, NA, NA, NA, NA, NA)
V5 <- c(4.638e-03, 7.424e-04, 1.536e-01, 1.420e-05, 7.355e-04, 2.765e-04, 1.883e-01, 7.562e-04, 1.133e-06,
5.000e-02, 1.000e-03, 1.300e-02, 1.100e-02, 2.800e-02, 2.200e-02, 3.900e-02)
V6 <- c(-0.6511331, -0.0920000, -0.2298655, -0.7321321, -0.7542836, -0.9301781, -0.3466844, -0.7198487,
-0.9772861, -0.4200000, -0.5700000, -0.4600000, -0.4700000, -0.4200000, -0.4300000, -0.3900000)
V7 <- c(17, 1352, 40, 27, 16, 9, 16, 18, 10, 29, 29, 29, 29, 29, 29, 29)dat <- cbind(V1, V2, V3, V4, V5, V6, V7)dat <- as.data.frame(dat)dat$V1 <- as.character(dat$V1)
dat$V2 <- as.integer(as.character(dat$V2))
dat$V3 <- as.numeric(as.character(dat$V3))
dat$V4 <- as.numeric(as.character(dat$V4))
dat$V5 <- as.numeric(as.character(dat$V5))
dat$V6 <- as.numeric(as.character(dat$V6))
dat$V7 <- as.numeric(as.character(dat$V7))str(dat)dat <- dat %>% rename(author = "V1",
year = "V2",
estimate = "V3",
std_error = "V4",
p_value = "V5",
r = "V6",
number_of_obs = "V7")
for (i in 1:nrow(dat)){
if (is.na(dat$std_error[i]) == TRUE ){
dat$std_error[i] <- abs(dat$estimate[i]) / qt(dat$p_value[i]/2,
df=dat$number_of_obs-2, lower.tail=FALSE)
}
}
res <- rma.mv(yi = dat$estimate, V = (dat$std_error)**2, random = ~ 1 | author, data=dat)coef_test(res, vcov="CR2")forest(res, addcred = TRUE, showweights = TRUE, header = TRUE,
order = "obs", col = "blue", slab = dat$author)
funnel(res)
ranktest(res)The for loop that I use (and in which I use a formula that Dr Viechtbauer gave me in a previous answer) to calculate the standard errors for estimates that I didn't calculate myself gives me an error message, but still gives me values.
The error message is :
Warning messages:
1: In data$std_error[i] <- abs(data$estimate[i])/qt(data$p_value[i]/2, :
number of items to replace is not a multiple of replacement length
and it is repeated on 7 lines.I am not entirely sure, but I think this comes from the fact that values are set to NAs before the loop is run and are replaced by values, so replacing an item of length 0 with an item of length 1 I believe.
I hope this example can be run simply with a copy / paste, I think it should.
Did I do things correctly ? If not what should I modify ?Thank you !Norman
----- Mail d'origine -----
De: Viechtbauer, Wolfgang (SP) <wolfgang.viechtbauer at maastrichtuniversity.nl>
?: Norman DAURELLE <norman.daurelle at agroparistech.fr>
Cc: r-sig-meta-analysis <r-sig-meta-analysis at r-project.org>
Envoy?: Thu, 11 Jun 2020 15:33:02 +0200 (CEST)
Objet: RE: [R-meta] weight in rmv metafor
Dear Norman,
To give a simple example: When (some of the) studies contribute multiple estimates, the dataset has a multilevel structure (with estimates nested within studies). A common way to deal with this is to fit a multilevel model with random effects for studies and estimates within studies. Like this:
http://www.metafor-project.org/doku.php/analyses:konstantopoulos2011
However, multiple estimates from the same study are actually often computed based on the same sample of subjects. In that case, the sampling errors are also correlated. The multilevel model does not capture this. For this, one would ideally want to fit a model that also allows for correlated sampling errors. Like this:
http://www.metafor-project.org/doku.php/analyses:berkey1998
However, computing the covariances between the sampling errors within a study is difficult and requires information that is often not available.
We can ignore those correlations and use the multilevel model as a working model that is an approximation to the model that also accounts for correlated sampling errors. After fitting the multilevel model with rma.mv(), one can then use cluster robust inference methods to 'fix things up'.
Quite a bit of this has been discussed at length in previous posts on this mailing list. You might want to search the archives for some of these posts.
Best,
Wolfgang
-----Original Message-----
From: Norman DAURELLE [mailto:norman.daurelle at agroparistech.fr]
Sent: Thursday, 11 June, 2020 15:05
To: Viechtbauer, Wolfgang (SP)
Cc: r-sig-meta-analysis
Subject: Re: [R-meta] weight in rmv metafor
Thank you.
I am not sure I understand exactly what you mean by " if the working model
is only an approximation and doesn't cover all dependencies ".
Could you please explain it ?
For now I used the rma() function to synthesize the available literature
existing on the blackleg - oil seed rape disease-yield relationship, using
slopes as effect-sizes.
the models that gave me the slopes I used in the meta-analysis are all Y = a
+ bX, simple linear regressions with Y being the yield and X being the
diseqse severity.
So my slopes, b, are all negative, and I have obtained a "summary" effect
size through the rma() function.
But I indeed have two studies that for now contribute to most of the effect-
sizes that are included in my meta-analysis.
So why exactly is it necessary to use the rma.mv() function ?
What exactly does the "multivariate" qualificative refer to ?
Thank you,
Norman.
________________________________________
De: "Wolfgang Viechtbauer" <wolfgang.viechtbauer at maastrichtuniversity.nl>
?: "Norman DAURELLE" <norman.daurelle at agroparistech.fr>, "r-sig-meta-
analysis" <r-sig-meta-analysis at r-project.org>
Envoy?: Jeudi 11 Juin 2020 22:34:55
Objet: RE: [R-meta] weight in rmv metafor
Dear Norman,
If you only used rma(), then this is not correct. rma.mv() with an
appropriately specified model (plus clubSandwich::coef_test() if the working
model is only an approximation and doesn't cover all dependencies) would be
more appropriate.
Best,
Wolfgang
-----Original Message-----
From: Norman DAURELLE [mailto:norman.daurelle at agroparistech.fr]
Sent: Thursday, 11 June, 2020 14:13
To: r-sig-meta-analysis
Cc: Viechtbauer, Wolfgang (SP)
Subject: Re: [R-meta] weight in rmv metafor
Hi all,
I read this discussion and one question came to my mind : I also had some
studies that contributed multiple effect sizes in the meta-analysis that I
recently ran thanks to Dr Viechtbauer's advice.
For now I only used the rma function, but should I have used rma.mv because
of these stuides that had multiple effect sizes ?
Thank you !
Norman
________________________________________
De: "James Pustejovsky" <jepusto at gmail.com>
?: "Wolfgang Viechtbauer" <wolfgang.viechtbauer at maastrichtuniversity.nl>
Cc: "r-sig-meta-analysis" <r-sig-meta-analysis at r-project.org>, "Huang Wu"
<huang.wu at wmich.edu>
Envoy?: Mercredi 10 Juin 2020 05:08:09
Objet: Re: [R-meta] weight in rmv metafor
Hi Huang,
I've written up some notes that add a bit of further intuition to the
discussion that Wolfgang provided. The main case that I focus on is a model
that is just a meta-analysis (i.e., no predictors) and that includes random
effects to capture both between-study and within-study heterogeneity. I
also say a little bit about meta-regression models with only study-level
predictors.
https://www.jepusto.com/weighting-in-multivariate-meta-analysis/
Best,
James
On Sun, Jun 7, 2020 at 4:11 PM Viechtbauer, Wolfgang (SP) <
wolfgang.viechtbauer at maastrichtuniversity.nl> wrote:
Of course the weights "impact the estimated fixed effects". But whether
studies with multiple effect sizes tend to receive more weight depends on
various factors, including the variances of the random effects and the
sampling error (co)variances.
A more detailed discussion around the way weighting works in rma.mv
models can be found here:
http://www.metafor-project.org/doku.php/tips:weights_in_rma.mv_models
Note that weights(res, type="rowsum") currently only works in the 'devel'
version of metafor, so follow
https://wviechtb.github.io/metafor/#installation if you want to reproduce
this part as well.
I hope this clarifies things.
Best,
Wolfgang