This page provides tips on how to carry out meta-analysis in R.

Fixed effects meta-analysis

Under the fixed effects model, we assume that all studies in the meta-analysis share the same true effect size. The variation among the studies are assumed to result only from sampling error (i.e., chance).

To estimate the mean effect size, we weight each effect size in the sample according to the magnitude of its sampling variance, which takes into account the different sample sizes of the studies. Calculate each weight as the inverse of the squared standard error,
w = 1/SEz2,
where z is the effect size of interest and SE is its standard error.

Use the weighted.mean command to calculate a weighted mean of the effect sizes using the weights (see the help page for the command).

Calculate the standard error of the weighted mean as the square root of the inverse of the sum of sample weights. This standard error measures uncertainty of the estimate of true effect size.
SEz = √(1/∑w).

To calculate an approximate 95% confidence interval for the weighted mean, use the normal approximation. The critical value (1.96) is obtained as

qnorm(1 - 0.05/2)

The limits of the 95% confidence interval are then calculated as the estimated weighted mean – 1.96 SEz and the mean + 1.96 SEz.

Random effects meta-analysis

Under the random effects model we assume that each study estimates a system-specific effect size. There is no “one true” effect size under this model, only a mean effect size. This is more realistic than the fixed effects model for most meta-analyses in ecology and evolution.

To fit the random effects model, estimate the variance among the system-specific effect sizes, τ2 (“tau squared”). One way to estimate it involves calculating the heterogeneity among the observed effect sizes (Q), and then “correcting” by subtracting the within-study sampling variance. The correction is needed because the variance among the observed effect sizes among studies is inflated by within-study sampling error.

First, calculate Q, the weighted heterogeneity among the observed effect sizes z, as
Q = ∑w(zz)2

Then calculate τ2 by subtraction, being careful not to allow a negative value (since τ2 is a variance, which can’t be negative). In the following, N refers to the number of studies in the meta-analysis.
If Q < N – 1, then set τ2 = 0.
Otherwise, calculate
τ2 = (Q – (N – 1))/(∑w – ∑w2/∑w).

τ2 is then used to calculate new weights for the effect sizes of each study under the random effects model, taking account of both within-study sampling variance and between-study variance in effect sizes. SEz2 is the squared standard error of each effect size z (same as we computed in the fixed effects section above).
w’ = 1/(SEz2 + τ2)

To calculate the weighted mean effect size z under the random effects model, use the same procedure as that used before for the fixed effects model except use the new weights w’ calculated here.

The formula for the standard error (SE) of the mean effect size z is the same as that in the fixed-effects model, except use the new weights, w’.

The 95% confidence interval for the mean effect size under the random effects model is calculated in the same way as under the fixed effects model, except that you would use the standard error for the random effects mean.

Use the metafor package

The metafor package include functions to fit fixed- and random-effects models for meta-analysis. The models can be specified using a formula method analogous to that used in linear models (e.g., lm). The formulas make it easy to include covariates in the model (“moderator variables” in meta-analysis jargon).

Begin by loading the metafor package in R. The variables of interest should be available in a dataframe (here named metadata). Each row in the data frame is a distinct study. Variables in the data frame include include the effect size to be used in the meta-analysis (e,g., Cohen’s d) and a variable containing the standard error of each effect size (here labeled se).

The estimated grand mean in a random effects model is obtained using the rma.uni command. The default method for estimating τ2 is REML. Use the DerSimonian-Laird estimator instead (method = "DL") if you want to obtain the same results as those from your “by hand” calculations based on the equations above.

z <- rma.uni(d ~ 1, vi = se^2, data = metadata, method = "DL")

To make a forest plot, use the forest command. Each row of the plot can be labeled using the study identifier, for example if it is included as a column of the data frame named reference, as follows. The command to create a funnel plot is also given.

z <- rma.uni(d ~ 1, vi = se^2, data = metadata, method = "DL")
forest(z, slab = metadata$reference)

To fit a random effects meta-analysis model that includes a covariate A ('moderator variable"), expand the formula as in an ordinary linear model. The coefficients in the summary table can be interpreted in exactly the same way as the coefficients from an ordinary linear model.

z <- rma.uni(d ~ A, vi = se^2, data = metadata, method = "DL")