General Notes / Setup

The book Practical Meta-Analysis by Lipsey and Wilson (2001) is an excellent introduction to meta-analysis and covers (in less than 250 pages) the most important aspects of conducting a meta-analysis. This document provides the R code to reproduce the analyses conducted in chapter 7 (which is focused on the computational aspects of a meta-analysis, including the fixed-, random-, and the mixed-effects model) using the metafor package. To read more about the package, see the package website and the package documentation.

The package can be installed with:

install.packages("metafor")

Once the package is installed, we can load it with:

library(metafor)

1) Data

The data used for the analyses are provided in the book on page 130 in Table 7.1. We can create the same dataset with:

dat <- data.frame(
   id = c(100, 308, 1596, 2479, 9021, 9028, 161, 172, 537, 7049),
   yi = c(-0.33, 0.32, 0.39, 0.31, 0.17, 0.64, -0.33, 0.15, -0.02, 0.00),
   vi = c(0.084, 0.035, 0.017, 0.034, 0.072, 0.117, 0.102, 0.093, 0.012, 0.067),
   random = c(0, 0, 0, 0, 0, 0, 1, 1, 1, 1),
   intensity = c(7, 3, 7, 5, 7, 7, 4, 4, 5, 6))
dat
##      id    yi    vi random intensity
## 1   100 -0.33 0.084      0         7
## 2   308  0.32 0.035      0         3
## 3  1596  0.39 0.017      0         7
## 4  2479  0.31 0.034      0         5
## 5  9021  0.17 0.072      0         7
## 6  9028  0.64 0.117      0         7
## 7   161 -0.33 0.102      1         4
## 8   172  0.15 0.093      1         4
## 9   537 -0.02 0.012      1         5
## 10 7049  0.00 0.067      1         6

The yi values are standardized mean differences based on studies examining the effectiveness of 'challenge programs' to treat juvenile delinquency (presumably with positive values indicating that delinquent behaviors were reduced in the group of juveniles participating in such a challenge program compared to those who did not). The vi values are the corresponding sampling variances. The dummy variable random indicates whether juveniles were randomly assigned to the two conditions (1 = yes, 0 = no). Finally, intensity is a variable coded for the purposes of the meta-analysis which indicates (based on assessment of the coder) the intensity of the challenge program (coded from 1 = 'very low' to 7 = 'very high'). These last two variables are potential moderators of the treatment effectiveness.


2) Fixed-Effects Model

We can fit a fixed-effects model to these data with:

res.fe <- rma(yi, vi, data=dat, method="FE")
res.fe
## Fixed-Effects Model (k = 10)
## 
## I^2 (total heterogeneity / total variability):   39.04%
## H^2 (total variability / sampling variability):  1.64
## 
## Test for Heterogeneity:
## Q(df = 9) = 14.7640, p-val = 0.0976
## 
## Model Results:
## 
## estimate      se    zval    pval   ci.lb   ci.ub 
##   0.1549  0.0609  2.5450  0.0109  0.0356  0.2742  * 
## 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

These results match what is reported on page 132-133 (see Exhibit 7.3).


3) Random-Effects Model

A random-effects model (using the DerSimonian-Laird estimator for \(\tau^2\)) can be fitted with:

res.re <- rma(yi, vi, data=dat, method="DL")
res.re
## Random-Effects Model (k = 10; tau^2 estimator: DL)
## 
## tau^2 (estimated amount of total heterogeneity): 0.0260 (SE = 0.0328)
## tau (square root of estimated tau^2 value):      0.1611
## I^2 (total heterogeneity / total variability):   39.04%
## H^2 (total variability / sampling variability):  1.64
## 
## Test for Heterogeneity:
## Q(df = 9) = 14.7640, p-val = 0.0976
## 
## Model Results:
## 
## estimate      se    zval    pval    ci.lb   ci.ub 
##   0.1534  0.0858  1.7893  0.0736  -0.0146  0.3215  . 
## 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Again, these results match what is shown in Exhibit 7.3 and on pages 134-135.

A forest plot showing the results of the individual studies and the results from the random-effects model can be created as follows:

par(mar=c(4,4,2,2))

forest(res.re, xlim=c(-2,3), header=c("Study", "SMD [95% CI]"), top=2,
       xlab="Standardized Mean Difference")