[R-meta] Weighing factors in rma.mv corresponding to those used in lmer
Viechtbauer Wolfgang (SP)
wolfgang.viechtbauer at maastrichtuniversity.nl
Thu Mar 29 12:20:10 CEST 2018
Let's start with:
res0 <- lmer(y ~ x1 + x2 + (1|PubID), data=Dataset, na.action = "na.omit")
so a model without 'weights'. This model has 'PubID' level random effects and errors, each of which is associated with a variance component (PubID variance and error variance).
Since rma.mv() is geared towards meta-analytic applications, it expects you to specify the var-cov matrix of the sampling errors. So, you would get the same results as above with:
res0m <- rma.mv(y ~ x1 + x2, V = sigma(res0)^2 * diag(nrow(Dataset)), random = ~ 1 | PubID, data=Dataset)
So, sigma(res0)^2 * diag(nrow(Dataset)) is the (diagonal) var-cov matrix of the errors.
Actually, this would be sufficient:
res0m <- rma.mv(y ~ x1 + x2, V = sigma(res0)^2, random = ~ 1 | PubID, data=Dataset)
since if you specify a constant for V, rma.mv() constructs the diagonal var-cov matrix for you with that constant along the diagonal.
But how can we actually get rma.mv() to estimate the error variance (instead of taking that value from 'res0' as estimated by lmer())? With this:
Dataset$ID <- 1:nrow(Dataset)
res0m <- rma.mv(y ~ x1 + x2, V=0, random = list(~ 1 | PubID, ~ 1 | ID), data=Dataset)
So, we set V=0, so the user-supplied error var-cov matrix is just all 0s and instead we add random effects for every row of the dataset - which will then be the error variance.
Now let's go to:
res1 <- lmer(y ~ x1 + x2 + (1|PubID), data=Dataset, weights = n, na.action = "na.omit")
What this means is that the error variances are equal to sigma^2 / n. So, we would get the same results with:
res1m <- rma.mv(y ~ x1 + x2, V=sigma(res1)^2 * diag(1/Dataset$n), random = list(~ 1 | PubID), data=Dataset)
We can also get rma.mv() to estimate the variance component for the errors, but this is a bit more tricky. Essentially, we want to add a var-cov matrix for which we know diag(1/Dataset$n) but not the multiplicative variance component. We can do this via the 'R' argument, which allows you to specify a known matrix with unknown multiplicative variance component (which is then estimated). So, this would be:
Dataset$ID <- 1:nrow(Dataset)
Rmat <- diag(1/Dataset$n)
rownames(Rmat) <- colnames(Rmat) <- Dataset$ID
res1m <- rma.mv(y ~ x1 + x2, V=0, random = list(~ 1 | PubID, ~ 1 | ID), R=list(ID = Rmat), Rscale=FALSE, data=Dataset)
The Rscale=FALSE part is needed since, by default, 'Rmat' would be rescalled to a correlation matrix, but that doesn't make sense here.
res2 and res3 are the same as res1, except that you want to use:
Rmat <- diag(Dataset$SE)
Rmat <- diag(Dataset$SE^2)
Wolfgang Viechtbauer, Ph.D., Statistician | Department of Psychiatry and
Neuropsychology | Maastricht University | P.O. Box 616 (VIJV1) | 6200 MD
Maastricht, The Netherlands | +31 (43) 388-4170 | http://www.wvbauer.com
From: R-sig-meta-analysis [mailto:r-sig-meta-analysis-bounces at r-project.org] On Behalf Of Martineau, Roger
Sent: Tuesday, 27 March, 2018 22:56
To: r-sig-meta-analysis at r-project.org
Subject: [R-meta] Weighing factors in rma.mv corresponding to those used in lmer
Dear list members,
I want to use the rma.mv function and compare the output with that from the lmer function.
I need to compare 3 weighting schemes: (1) the number of replica (n), (2) the inverse of SE and (3) the inverse of variance or 1/SE^2 in each PubID.
(res1 <- lmer(y ~ x1 + x2 + (1|PubID), data=Dataset, weights = n, na.action = "na.omit"))
(res2 <- lmer(y ~ x1 + x2 + (1|PubID), data=Dataset, weights = 1/SE, na.action = "na.omit"))
(res3 <- lmer(y ~ x1 + x2 + (1|PubID), data=Dataset, weights = 1/SE^2, na.action = "na.omit"))
What are the rma.mv scripts corresponding to res1, res2 and res3 ?
Thanks in advance
Roger Martineau, mv Ph.D.
Centre de recherche et de développement
sur le bovin laitier et le porc
Agriculture et agroalimentaire Canada/Agriculture and Agri-Food Canada
2000, Rue Collège / 2000, College Street
Sherbrooke (Québec) J1M 0C8
roger.martineau at agr.gc.ca<mailto:roger.martineau at agr.gc.ca>
More information about the R-sig-meta-analysis