[R-meta] Queries re: rma.mv
Dear Adelina, See below for my responses. Best, Wolfgang
-----Original Message----- From: R-sig-meta-analysis [mailto:r-sig-meta-analysis-bounces at r-project.org] On Behalf Of Adelina Artenie Sent: Wednesday, 20 April, 2022 19:15 To: r-sig-meta-analysis at r-project.org Subject: [R-meta] Queries re: rma.mv Hello I am trying to do a meta-regression to explore whether incidence rates are changing over time. Some of the estimates are correlated so I am using the rma.mv function. The moderator is a measure of time. I have two questions please: 1. I am trying to understand why is it that the correlated estimates are systematically up-weighted. I would have thought that any correlated data would be down-weighted or, at the very least, have a similar weight as the uncorrelated data. I appreciate that the weight calculation in this case is complex (https://www.metafor-project.org/doku.php/tips:weights_in_rma.mv_models ). I wonder if, conceptually, this would be expected? In my mock example, the correlated estimates have weights of 7-10% whereas the uncorrelated ones have weights of 4%. I see a similar trend in my own data, and I find the difference in the weights to be quite considerable.
I assume you are referring to the weights shown in the forest plot. These are just based on the diagonal of the weight matrix. As you note, the weighting in such models is more complex, so looking at these weights doesn't really allow you to draw such a conclusion. In fact, generally, the effect will be as you expect, namely that correlated estimates will be downweighted. Consider the dat.konstantopoulos2011 dataset and let's make all sampling variances equal to each other to see this more clearly. Also, we will compute the row-sum weights, since these reflect more directly the weight each estimate is getting in this model.
dat <- dat.konstantopoulos2011
dat$vi <- median(dat$vi)
res <- rma.mv(yi, vi, random = ~ 1 | district/school, data=dat)
res
wi <- weights(res, type="rowsum")
data.frame(k = c(table(dat$district)),
weight = tapply(wi, dat$district, mean))
So, as you can see, in clusters with more estimates, each individual estimate gets less weight.
Things will be more complicated when sampling variances differ and also when including moderators in the model.
2. Could I confirm with you that my interpretation of the coefficient for the moderator is correct? In this case, the estimate is 0.0846, which represents the mean/median log IRR (incidence rate ratio) for a one-unit increase in the moderator. Therefore, the IRR is 1.088: for each unit increase in time, the incidence rate increases on average by 8.8% (ignoring the confidence intervals for now).
Correct.
Many thanks in advance,
Adelina
Adelina Artenie, MSc, PhD
CIHR Postdoctoral Research Fellow
Population Health Sciences
Bristol Medical School
University of Bristol
Code
***
ID <- c(1:18)
author <- c("AA", "AA", "BB", "CC", "DD", "EE", "EE", "EE", "FF", "GG", "HH",
"II", "JJ",
"KK", "LL", "MM", "NN", "OO")
cases <- c(8:25)
prs_yrs <-
c(200,150,3000,400,500,100,400,300,1200,456,177,296,664,123,123,432,67,3045)
moderator<-c(1:18)
mydata<-data.frame(ID,author,cases,prs_yrs,moderator)
print(mydata)
my_example <- escalc(measure="IRLN", xi=cases, ti=prs_yrs, data=mydata)
my_example
my_meta_reg <- rma.mv(yi = yi,
V = vi,
slab = author,
data = my_example,
random = ~ 1 | author/ID,
test = "z",
method = "REML",
mods = ~ moderator)
summary(my_meta_reg)
regplot(my_meta_reg, label = T, transf=exp,
ylab = "Incidence rate",
xlab = "Moderator",
col="black",
plim = c(1,NA),
labsize = 0.5,
legend = F)
forest(my_meta_reg, showweights=TRUE)