Using the dat.raudenbush1985
from week 3. Rerun the random-effects meta-analysis with REML estimation from the earlier exercises.
tester
): build one model with a common \(\tau^2\) estimate and one without. What do you conclude?year
). Center the year variable first. Make a bubble plot illustrating this meta-regression. What do you conclude?tester
) compare it to the meta-regression which contains both tester blindness (tester
) and centered publication year (year
). What do you conclude?There are 19 studies.
library(meta)
library(metafor) # needed later on.
Data<-dat.raudenbush1985
head(Data)
## study author year weeks setting tester n1i n2i yi
## 1 1 Rosenthal et al. 1974 2 group aware 77 339 0.0300
## 2 2 Conn et al. 1968 21 group aware 60 198 0.1200
## 3 3 Jose & Cody 1971 19 group aware 72 72 -0.1400
## 4 4 Pellegrini & Hicks 1972 0 group aware 11 22 1.1800
## 5 5 Pellegrini & Hicks 1972 0 group blind 11 22 0.2600
## 6 6 Evans & Rosenthal 1969 3 group aware 129 348 -0.0600
## vi
## 1 0.0156
## 2 0.0216
## 3 0.0279
## 4 0.1391
## 5 0.1362
## 6 0.0106
I have made the names as in our slides. Note that this is unnecessary duplication (as TE = yi ) but then it maps on nicely to our code.
library(dplyr)
## Warning: package 'dplyr' was built under R version 3.5.2
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
Data <- Data %>% mutate(TE=yi, seTE=sqrt(vi))
As in the solution of exercise 3, let’s combine, author and year.
Let’s bracket the year as per convention. Here I rely on base R and this snippet.
Data$year_b <- paste0("(", format(unlist(Data[,3])),")")
Then we combine as done in here. Here I have opted for a ‘tidyverse’ solution.
library(tidyverse)
Data <-Data %>% unite(author_year, c(author, year_b), sep = " ", remove = FALSE)
Let’s redo our model but now with author_year
.
require(meta)
model_reml2<-metagen(TE,
seTE,
data=Data,
studlab=paste(author_year),
comb.fixed = FALSE,
comb.random = TRUE,
method.tau = "REML",
hakn = FALSE,
prediction=TRUE,
sm="SMD")
model_reml2
## SMD 95%-CI %W(random)
## Rosenthal et al. (1974) 0.0300 [-0.2148; 0.2748] 7.7
## Conn et al. (1968) 0.1200 [-0.1681; 0.4081] 6.6
## Jose & Cody (1971) -0.1400 [-0.4674; 0.1874] 5.7
## Pellegrini & Hicks (1972) 1.1800 [ 0.4490; 1.9110] 1.7
## Pellegrini & Hicks (1972) 0.2600 [-0.4633; 0.9833] 1.7
## Evans & Rosenthal (1969) -0.0600 [-0.2618; 0.1418] 9.1
## Fielder et al. (1971) -0.0200 [-0.2218; 0.1818] 9.1
## Claiborn (1969) -0.3200 [-0.7512; 0.1112] 4.0
## Kester (1969) 0.2700 [-0.0515; 0.5915] 5.8
## Maxwell (1970) 0.8000 [ 0.3081; 1.2919] 3.3
## Carter (1970) 0.5400 [-0.0519; 1.1319] 2.4
## Flowers (1966) 0.1800 [-0.2569; 0.6169] 3.9
## Keshock (1970) -0.0200 [-0.5864; 0.5464] 2.6
## Henrikson (1970) 0.2300 [-0.3384; 0.7984] 2.6
## Fine (1972) -0.1800 [-0.4918; 0.1318] 6.0
## Grieger (1970) -0.0600 [-0.3874; 0.2674] 5.7
## Rosenthal & Jacobson (1968) 0.3000 [ 0.0277; 0.5723] 7.0
## Fleming & Anttonen (1971) 0.0700 [-0.1139; 0.2539] 9.7
## Ginsburg (1970) -0.0700 [-0.4112; 0.2712] 5.4
##
## Number of studies combined: k = 19
##
## SMD 95%-CI z p-value
## Random effects model 0.0837 [-0.0175; 0.1849] 1.62 0.1051
## Prediction interval [-0.2256; 0.3930]
##
## Quantifying heterogeneity:
## tau^2 = 0.0188; H = 1.41 [1.08; 1.84]; I^2 = 49.8% [14.7%; 70.4%]
##
## Test of heterogeneity:
## Q d.f. p-value
## 35.83 18 0.0074
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
Below we perform the models by group.
Interestingly, it seems that the effect is stronger in the blinded group than in the aware group - perhaps they were better studies? However, there is no significant difference between the groups ( Q=.84, p =.367).
tester_subgroup_common<-update.meta(model_reml2,
byvar=tester,
comb.random = TRUE,
comb.fixed = FALSE,
tau.common=TRUE)
tester_subgroup_common
## SMD 95%-CI %W(random) tester
## Rosenthal et al. (1974) 0.0300 [-0.2148; 0.2748] 7.7 aware
## Conn et al. (1968) 0.1200 [-0.1681; 0.4081] 6.6 aware
## Jose & Cody (1971) -0.1400 [-0.4674; 0.1874] 5.7 aware
## Pellegrini & Hicks (1972) 1.1800 [ 0.4490; 1.9110] 1.7 aware
## Pellegrini & Hicks (1972) 0.2600 [-0.4633; 0.9833] 1.7 blind
## Evans & Rosenthal (1969) -0.0600 [-0.2618; 0.1418] 9.1 aware
## Fielder et al. (1971) -0.0200 [-0.2218; 0.1818] 9.1 blind
## Claiborn (1969) -0.3200 [-0.7512; 0.1112] 4.0 aware
## Kester (1969) 0.2700 [-0.0515; 0.5915] 5.8 aware
## Maxwell (1970) 0.8000 [ 0.3081; 1.2919] 3.3 blind
## Carter (1970) 0.5400 [-0.0519; 1.1319] 2.4 blind
## Flowers (1966) 0.1800 [-0.2569; 0.6169] 3.9 blind
## Keshock (1970) -0.0200 [-0.5864; 0.5464] 2.6 blind
## Henrikson (1970) 0.2300 [-0.3384; 0.7984] 2.6 blind
## Fine (1972) -0.1800 [-0.4918; 0.1318] 6.0 aware
## Grieger (1970) -0.0600 [-0.3874; 0.2674] 5.7 blind
## Rosenthal & Jacobson (1968) 0.3000 [ 0.0277; 0.5723] 7.0 aware
## Fleming & Anttonen (1971) 0.0700 [-0.1139; 0.2539] 9.7 blind
## Ginsburg (1970) -0.0700 [-0.4112; 0.2712] 5.4 aware
##
## Number of studies combined: k = 19
##
## SMD 95%-CI z p-value
## Random effects model 0.0837 [-0.0175; 0.1849] 1.62 0.1051
## Prediction interval [-0.2256; 0.3930]
##
## Quantifying heterogeneity:
## tau^2 = 0.0188; H = 1.41 [1.08; 1.84]; I^2 = 49.8% [14.7%; 70.4%]
##
## Quantifying residual heterogeneity:
## tau^2 = 0.0248; H = 1.38; I^2 = 47.5%
##
## Test of heterogeneity:
## Q d.f. p-value
## 35.83 18 0.0074
##
## Results for subgroups (random effects model):
## k SMD 95%-CI Q tau^2 I^2
## tester = aware 10 0.0463 [-0.0945; 0.1872] 22.20 0.0248 59.5%
## tester = blind 9 0.1488 [-0.0198; 0.3174] 12.96 0.0248 38.3%
##
## Test for subgroup differences (random effects model):
## Q d.f. p-value
## Between groups 0.84 1 0.3606
## Within groups 35.16 17 0.0059
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2 (assuming common tau^2 in subgroups)
This model allows an estimate of variance \(\tau^2\) for each group. The conclusion is the same there is no evidence for a significant difference between these groups ( Q=.81, p =.367).
tester_subgroup_sep<-update.meta(model_reml2,
byvar=tester,
comb.random = TRUE,
comb.fixed = FALSE,
tau.common=FALSE)
tester_subgroup_sep
## SMD 95%-CI %W(random) tester
## Rosenthal et al. (1974) 0.0300 [-0.2148; 0.2748] 7.7 aware
## Conn et al. (1968) 0.1200 [-0.1681; 0.4081] 6.6 aware
## Jose & Cody (1971) -0.1400 [-0.4674; 0.1874] 5.7 aware
## Pellegrini & Hicks (1972) 1.1800 [ 0.4490; 1.9110] 1.7 aware
## Pellegrini & Hicks (1972) 0.2600 [-0.4633; 0.9833] 1.7 blind
## Evans & Rosenthal (1969) -0.0600 [-0.2618; 0.1418] 9.1 aware
## Fielder et al. (1971) -0.0200 [-0.2218; 0.1818] 9.1 blind
## Claiborn (1969) -0.3200 [-0.7512; 0.1112] 4.0 aware
## Kester (1969) 0.2700 [-0.0515; 0.5915] 5.8 aware
## Maxwell (1970) 0.8000 [ 0.3081; 1.2919] 3.3 blind
## Carter (1970) 0.5400 [-0.0519; 1.1319] 2.4 blind
## Flowers (1966) 0.1800 [-0.2569; 0.6169] 3.9 blind
## Keshock (1970) -0.0200 [-0.5864; 0.5464] 2.6 blind
## Henrikson (1970) 0.2300 [-0.3384; 0.7984] 2.6 blind
## Fine (1972) -0.1800 [-0.4918; 0.1318] 6.0 aware
## Grieger (1970) -0.0600 [-0.3874; 0.2674] 5.7 blind
## Rosenthal & Jacobson (1968) 0.3000 [ 0.0277; 0.5723] 7.0 aware
## Fleming & Anttonen (1971) 0.0700 [-0.1139; 0.2539] 9.7 blind
## Ginsburg (1970) -0.0700 [-0.4112; 0.2712] 5.4 aware
##
## Number of studies combined: k = 19
##
## SMD 95%-CI z p-value
## Random effects model 0.0837 [-0.0175; 0.1849] 1.62 0.1051
## Prediction interval [-0.2256; 0.3930]
##
## Quantifying heterogeneity:
## tau^2 = 0.0188; H = 1.41 [1.08; 1.84]; I^2 = 49.8% [14.7%; 70.4%]
##
## Quantifying residual heterogeneity:
## H = 1.44 [1.10; 1.88]; I^2 = 51.6% [17.1%; 71.8%]
##
## Test of heterogeneity:
## Q d.f. p-value
## 35.83 18 0.0074
##
## Results for subgroups (random effects model):
## k SMD 95%-CI Q tau^2 I^2
## tester = aware 10 0.0466 [-0.0954; 0.1887] 22.20 0.0256 59.5%
## tester = blind 9 0.1474 [-0.0193; 0.3140] 12.96 0.0236 38.3%
##
## Test for subgroup differences (random effects model):
## Q d.f. p-value
## Between groups 0.81 1 0.3672
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
Data <- Data %>% mutate(year_cent = year-mean(year))
We need to rerun our model as we have added our new moderator variable, note that this overrides our previous model.
require(meta)
model_reml2<-metagen(TE,
seTE,
data=Data,
studlab=paste(author_year),
comb.fixed = FALSE,
comb.random = TRUE,
method.tau = "REML",
hakn = FALSE,
prediction=TRUE,
sm="SMD")
model_reml2
## SMD 95%-CI %W(random)
## Rosenthal et al. (1974) 0.0300 [-0.2148; 0.2748] 7.7
## Conn et al. (1968) 0.1200 [-0.1681; 0.4081] 6.6
## Jose & Cody (1971) -0.1400 [-0.4674; 0.1874] 5.7
## Pellegrini & Hicks (1972) 1.1800 [ 0.4490; 1.9110] 1.7
## Pellegrini & Hicks (1972) 0.2600 [-0.4633; 0.9833] 1.7
## Evans & Rosenthal (1969) -0.0600 [-0.2618; 0.1418] 9.1
## Fielder et al. (1971) -0.0200 [-0.2218; 0.1818] 9.1
## Claiborn (1969) -0.3200 [-0.7512; 0.1112] 4.0
## Kester (1969) 0.2700 [-0.0515; 0.5915] 5.8
## Maxwell (1970) 0.8000 [ 0.3081; 1.2919] 3.3
## Carter (1970) 0.5400 [-0.0519; 1.1319] 2.4
## Flowers (1966) 0.1800 [-0.2569; 0.6169] 3.9
## Keshock (1970) -0.0200 [-0.5864; 0.5464] 2.6
## Henrikson (1970) 0.2300 [-0.3384; 0.7984] 2.6
## Fine (1972) -0.1800 [-0.4918; 0.1318] 6.0
## Grieger (1970) -0.0600 [-0.3874; 0.2674] 5.7
## Rosenthal & Jacobson (1968) 0.3000 [ 0.0277; 0.5723] 7.0
## Fleming & Anttonen (1971) 0.0700 [-0.1139; 0.2539] 9.7
## Ginsburg (1970) -0.0700 [-0.4112; 0.2712] 5.4
##
## Number of studies combined: k = 19
##
## SMD 95%-CI z p-value
## Random effects model 0.0837 [-0.0175; 0.1849] 1.62 0.1051
## Prediction interval [-0.2256; 0.3930]
##
## Quantifying heterogeneity:
## tau^2 = 0.0188; H = 1.41 [1.08; 1.84]; I^2 = 49.8% [14.7%; 70.4%]
##
## Test of heterogeneity:
## Q d.f. p-value
## 35.83 18 0.0074
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
There is no suggestion that publication year is a viable predictor of the effect size (Q = .41, p= .52).
metareg_pub_year<-metareg(model_reml2,year_cent)
metareg_pub_year
##
## Mixed-Effects Model (k = 19; tau^2 estimator: REML)
##
## tau^2 (estimated amount of residual heterogeneity): 0.0201 (SE = 0.0167)
## tau (square root of estimated tau^2 value): 0.1417
## I^2 (residual heterogeneity / unaccounted variability): 42.84%
## H^2 (unaccounted variability / sampling variability): 1.75
## R^2 (amount of heterogeneity accounted for): 0.00%
##
## Test for Residual Heterogeneity:
## QE(df = 17) = 34.7475, p-val = 0.0067
##
## Test of Moderators (coefficient 2):
## QM(df = 1) = 0.4081, p-val = 0.5229
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## intrcpt 0.0851 0.0524 1.6239 0.1044 -0.0176 0.1879
## year_cent -0.0187 0.0293 -0.6388 0.5229 -0.0761 0.0387
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
bubble.metareg(metareg_pub_year,
xlab = "Publication Year (centered)",
ylab = "SMD",
col.line = "hotpink",
studlab = TRUE)
We’ll use set.seed()
to ensure we get the same results every single time.
The permutation test corroborates our findings from the meta-regression. There is no substantial evidence for an effect of publication year.
set.seed(1981)
permutest(metareg_pub_year, iter=5000, progbar=F)
##
## Test of Moderators (coefficient 2):
## QM(df = 1) = 0.4081, p-val* = 0.5482
##
## Model Results:
##
## estimate se zval pval* ci.lb ci.ub
## intrcpt 0.0851 0.0524 1.6239 0.4716 -0.0176 0.1879
## year_cent -0.0187 0.0293 -0.6388 0.5482 -0.0761 0.0387
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Here we make a switch to metafor
.
As can be expected this corroborates what we previously found when examining subgroups. Note that I have opted for ML estimation with the Knapp-Hartung estimation. As model comparisons aren’t well defined with REML methods, we switched to ML estimation in order to allow
metareg_tester<-rma(yi=TE,
sei=seTE,
data=Data,
method = "ML",
mods = ~ tester,
test="knha")
metareg_tester
##
## Mixed-Effects Model (k = 19; tau^2 estimator: ML)
##
## tau^2 (estimated amount of residual heterogeneity): 0.0133 (SE = 0.0126)
## tau (square root of estimated tau^2 value): 0.1155
## I^2 (residual heterogeneity / unaccounted variability): 32.78%
## H^2 (unaccounted variability / sampling variability): 1.49
## R^2 (amount of heterogeneity accounted for): 0.00%
##
## Test for Residual Heterogeneity:
## QE(df = 17) = 35.1572, p-val = 0.0059
##
## Test of Moderators (coefficient 2):
## F(df1 = 1, df2 = 17) = 0.5293, p-val = 0.4768
##
## Model Results:
##
## estimate se tval pval ci.lb ci.ub
## intrcpt 0.0420 0.0786 0.5348 0.5997 -0.1237 0.2077
## testerblind 0.0894 0.1229 0.7275 0.4768 -0.1699 0.3487
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
metareg_pub_tester<-rma(yi=TE,
sei=seTE,
data=Data,
method = "ML",
mods = ~ tester + year_cent,
test="knha")
metareg_pub_tester
##
## Mixed-Effects Model (k = 19; tau^2 estimator: ML)
##
## tau^2 (estimated amount of residual heterogeneity): 0.0071 (SE = 0.0099)
## tau (square root of estimated tau^2 value): 0.0842
## I^2 (residual heterogeneity / unaccounted variability): 20.17%
## H^2 (unaccounted variability / sampling variability): 1.25
## R^2 (amount of heterogeneity accounted for): 43.53%
##
## Test for Residual Heterogeneity:
## QE(df = 16) = 33.8848, p-val = 0.0056
##
## Test of Moderators (coefficients 2:3):
## F(df1 = 2, df2 = 16) = 0.4516, p-val = 0.6445
##
## Model Results:
##
## estimate se tval pval ci.lb ci.ub
## intrcpt 0.0390 0.0762 0.5122 0.6155 -0.1225 0.2006
## testerblind 0.0826 0.1190 0.6943 0.4975 -0.1697 0.3349
## year_cent -0.0225 0.0329 -0.6844 0.5035 -0.0922 0.0472
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(metareg_tester,metareg_pub_tester)
##
## df AIC BIC AICc logLik LRT pval QE tau^2
## Full 4 12.6903 16.4681 15.5475 -2.3452 33.8848 0.0071
## Reduced 3 11.4095 14.2428 13.0095 -2.7048 0.7192 0.3964 35.1572 0.0133
## R^2
## Full
## Reduced 46.7824%
Both the AIC/BIC (as well as the Likelihood ratio test (LRT)) do not suggest superiority of one model over another.
We reverted back to the REML model estimation. This model does not support an interaction effect between publication year and tester blindness. Even if it had, we should be very wary as we only have 19 studies in this meta-analysis.
metareg_interaction<-rma(yi=TE,
sei=seTE,
data=Data,
method = "REML",
mods = ~ tester * year_cent
)
metareg_interaction
##
## Mixed-Effects Model (k = 19; tau^2 estimator: REML)
##
## tau^2 (estimated amount of residual heterogeneity): 0.0338 (SE = 0.0242)
## tau (square root of estimated tau^2 value): 0.1840
## I^2 (residual heterogeneity / unaccounted variability): 54.60%
## H^2 (unaccounted variability / sampling variability): 2.20
## R^2 (amount of heterogeneity accounted for): 0.00%
##
## Test for Residual Heterogeneity:
## QE(df = 15) = 33.7799, p-val = 0.0037
##
## Test of Moderators (coefficients 2:4):
## QM(df = 3) = 1.0372, p-val = 0.7922
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## intrcpt 0.0495 0.0784 0.6310 0.5281 -0.1042 0.2031
## testerblind 0.1084 0.1217 0.8907 0.3731 -0.1302 0.3470
## year_cent -0.0127 0.0399 -0.3186 0.7501 -0.0910 0.0655
## testerblind:year_cent -0.0101 0.0737 -0.1364 0.8915 -0.1546 0.1345
##
## intrcpt
## testerblind
## year_cent
## testerblind:year_cent
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
The example is from here.
Note that throughout I have varied the rounding when I reported, you should make your own decisions on how precise you believe your results to be.
Please see the slides for further reading but a good place to start is Chapter 8 on models and approaches to inference in Koricheva, J., Gurevitch, J., & Mengersen, K. (2013). Handbook of Meta-Analysis in Ecology and Evolution. Princeton, NJ: Princeton University Press.
Cited literature
Raudenbush, S. W. (1984). Magnitude of teacher expectancy effects on pupil IQ as a function of the credibility of expectancy induction: A synthesis of findings from 18 experiments. Journal of Educational Psychology, 76(1), 85–97.
Raudenbush, S. W. (2009). Analyzing effect sizes: Random effects models. In H. Cooper, L. V. Hedges, & J. C. Valentine (Eds.), The handbook of research synthesis and meta-analysis (2nd ed., pp. 295–315). New York: Russell Sage Foundation.
sessionInfo()
## R version 3.5.1 (2018-07-02)
## Platform: x86_64-apple-darwin15.6.0 (64-bit)
## Running under: macOS Sierra 10.12.6
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/3.5/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/3.5/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_GB.UTF-8/en_GB.UTF-8/en_GB.UTF-8/C/en_GB.UTF-8/en_GB.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] forcats_0.4.0 stringr_1.4.0 purrr_0.3.2 readr_1.3.1
## [5] tidyr_0.8.3 tibble_2.1.3 ggplot2_3.2.1 tidyverse_1.2.1
## [9] dplyr_0.8.3 metafor_2.1-0 Matrix_1.2-14 meta_4.9-6
##
## loaded via a namespace (and not attached):
## [1] tidyselect_0.2.5 xfun_0.9 haven_2.1.1 lattice_0.20-35
## [5] colorspace_1.4-1 generics_0.0.2 vctrs_0.2.0 htmltools_0.3.6
## [9] yaml_2.2.0 rlang_0.4.0 pillar_1.4.2 withr_2.1.2
## [13] glue_1.3.1 modelr_0.1.2 readxl_1.3.1 munsell_0.5.0
## [17] gtable_0.3.0 cellranger_1.1.0 rvest_0.3.2 evaluate_0.14
## [21] knitr_1.24 broom_0.5.1 Rcpp_1.0.2 scales_1.0.0
## [25] backports_1.1.4 jsonlite_1.6 hms_0.5.1 digest_0.6.20
## [29] stringi_1.4.3 grid_3.5.1 cli_1.1.0 tools_3.5.1
## [33] magrittr_1.5 lazyeval_0.2.2 crayon_1.3.4 pkgconfig_2.0.2
## [37] zeallot_0.1.0 xml2_1.2.0 lubridate_1.7.4 assertthat_0.2.1
## [41] rmarkdown_1.11 httr_1.4.0 rstudioapi_0.9.0 R6_2.4.0
## [45] nlme_3.1-137 compiler_3.5.1