library(readxl)
d = read_excel("~/Desktop/archive/Rick/Combining Effect Sizes Basic Data.xlsx")
library(metafor)
## Loading required package: Matrix
## Loading 'metafor' package (version 2.1-0). For an overview
## and introduction to the package please type: help(metafor).
Fixed effects
fixedeffects = rma(Ti, vi, data = d, method = "FE")
summary(fixedeffects)
##
## Fixed-Effects Model (k = 19)
##
## logLik deviance AIC BIC AICc
## -3.4608 35.8252 8.9217 9.8661 9.1570
##
## Test for Heterogeneity:
## Q(df = 18) = 35.8252, p-val = 0.0074
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## 0.0603 0.0365 1.6539 0.0982 -0.0112 0.1319 .
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Random effects
randomeffects = rma(Ti, vi, data = d, method = "HE") #HE is the method that matches the results of the first example that we did in class however, default is REML
summary(randomeffects)
##
## Random-Effects Model (k = 19; tau^2 estimator: HE)
##
## logLik deviance AIC BIC AICc
## -4.8483 38.6001 13.6966 15.5855 14.4466
##
## tau^2 (estimated amount of total heterogeneity): 0.0804 (SE = 0.0448)
## tau (square root of estimated tau^2 value): 0.2835
## I^2 (total heterogeneity / total variability): 75.43%
## H^2 (total variability / sampling variability): 4.07
##
## Test for Heterogeneity:
## Q(df = 18) = 35.8252, p-val = 0.0074
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## 0.1143 0.0792 1.4430 0.1490 -0.0409 0.2696
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Q = fixedeffects$QE
c = sum(1/d$vi) - ((sum((1/d$vi)^2))/sum(1/d$vi))
sig = (Q-(nrow(d)-1))/c
d$weights = 1/(d$vi+sig)
randomeffectsQ = rma(Ti, vi, data = d, weights = weights, method = "GENQ") #matches second example in class
summary(randomeffectsQ)
##
## Random-Effects Model (k = 19; tau^2 estimator: GENQ)
##
## logLik deviance AIC BIC AICc
## -3.7409 36.3852 11.4817 13.3706 12.2317
##
## tau^2 (estimated amount of total heterogeneity): 0.0437 (SE = 0.0260)
## tau (square root of estimated tau^2 value): 0.2091
## I^2 (total heterogeneity / total variability): 62.56%
## H^2 (total variability / sampling variability): 2.67
##
## Test for Heterogeneity:
## Q(df = 18) = 35.8252, p-val = 0.0074
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## 0.0893 0.0649 1.3755 0.1690 -0.0380 0.2166
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
randomeffectsREML = rma(Ti, vi, data = d)
summary(randomeffectsREML)
##
## Random-Effects Model (k = 19; tau^2 estimator: REML)
##
## logLik deviance AIC BIC AICc
## -3.7368 7.4735 11.4735 13.2543 12.2735
##
## tau^2 (estimated amount of total heterogeneity): 0.0188 (SE = 0.0155)
## tau (square root of estimated tau^2 value): 0.1373
## I^2 (total heterogeneity / total variability): 41.85%
## H^2 (total variability / sampling variability): 1.72
##
## Test for Heterogeneity:
## Q(df = 18) = 35.8252, p-val = 0.0074
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## 0.0837 0.0517 1.6203 0.1052 -0.0175 0.1850
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Forest plot
forest(randomeffects) #randomeffects is the name of our previous rma model
Funnel plot
funnel(randomeffects)
Calculating effect sizes
Correlations
ns = c(100, 200, 300) #sample sizes
rs = c(.4, .5, .6) #correlations
es = escalc("COR", ni = ns, ri = rs) #provides effect size and variance
es
## yi vi
## 1 0.4000 0.0071
## 2 0.5000 0.0028
## 3 0.6000 0.0014
#now we can use these in the meta
cormeta = rma(es, yi, vi, method = "FE")
summary(cormeta)
##
## Fixed-Effects Model (k = 3)
##
## logLik deviance AIC BIC AICc
## 3.0148 5.8622 -4.0297 -4.9311 -0.0297
##
## Test for Heterogeneity:
## Q(df = 2) = 5.8622, p-val = 0.0533
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## 0.5482 0.0286 19.1787 <.0001 0.4922 0.6042 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Cohen’s d
n1i = c(64, 34, 78) #sample sizes for group 1
n2i = c(58, 38, 72) #sample sizes for group 2
m1i = c(10.5, 14.6, 12.8) #means of group 1
m2i = c(7.8, 13.1, 12.9) #means of group 2
sd1i = c(2.3, 3.2, 4) #standard deviation of group 1
sd2i = c(2.1, 2.6, 3.9) #standard deviation of group 2
es = escalc("MD", n1i=n1i, n2i=n2i, m1i=m1i, m2i=m2i, sd1i=sd1i, sd2i=sd2i)
es
## yi vi
## 1 2.7000 0.1587
## 2 1.5000 0.4791
## 3 -0.1000 0.4164
Hedge’s g
es = escalc("SMD", n1i=n1i, n2i=n2i, m1i=m1i, m2i=m2i, sd1i=sd1i, sd2i=sd2i)
es
## yi vi
## 1 1.2156 0.0389
## 2 0.5120 0.0575
## 3 -0.0252 0.0267
Log odds ratio
ai = c(14, 34, 23) #frequency in first cell (in class was called celln11)
bi = c(63, 93, 98) #celln12
ci = c(43, 78, 44) #celln21
di = c(125, 154, 76) #celln22
es = escalc("OR", ai=ai, bi=bi, ci=ci, di=di)
es
## yi vi
## 1 -0.4370 0.1186
## 2 -0.3260 0.0595
## 3 -0.9029 0.0896