#### Introduction ####

#Meta Analysis Example
#Crash Course in Meta-Analysis Worshop
#Rob Cribbie
#CPA 2019, Halifax, Nova Scotia


#Difference between Males and Females 
#in Maladaptive Perfectionism

#Effect sizes are reported as Hedge's g
#(which is almost identical to Cohen's d)

dat<-read.csv(file.choose())
dat

#We will use the package "metafor" which 
#has some excellent meta-analysis tools
library(metafor)

#### Fixed Effects Model ####

#Fixed Effects Model (method="FE")

ma1<-rma(data=dat,
         yi=g, sei=se_g, 
         method="FE")
summary(ma1)


#### Random Effects Model ####

#Random Effects Model
#Restricted Maximum Likelihood 
#Estimation (method="REML")

ma2<-rma(data=dat,
         yi=g, sei=se_g, 
         method="REML")
summary(ma2)


#### Diagnostics (Random Effects Model) ####

#Forest Plot for the Standardized Mean Difference
forest(ma2,level=95)

#Identify Influential Cases
ic<-influence(ma2)
plot(ic)
(cd_cut<-4/(length(dat$g)))
ic$inf$cook.d
which(ic$inf$cook.d > cd_cut)

#Influential Cases
#One at a time leave out a study
#Plot the effect sizes from each analsysis
plot(leave1out(ma2)$estimate)

#Publication Bias - Funnel Plot
funnel(ma2)
regtest(ma2) #Test for Funnel Plot Asymmetry

#Four helpful plots with one command
plot(ma2)

#### Moderator Analyses ####

library(ggplot2)

#Country as a Moderator of Effect Size
boxplot(dat$g ~ dat$country)

#Random Effects Model
#Restricted Maximum Likelihood Estimation (method="REML")
#With moderator "Country"

ma4<-rma(data=dat,
         yi=g, sei=se_g, method="REML",
         mods=~country)

summary(ma4)
tapply(dat$g,dat$country,mean)


#Forest Plot with Fitted Values
forest(ma4)
forest(ma4,order="obs")




#### Let's redo the analysis without Study 6 ####
# Recall: Study 6 g=3.97
dat2<-dat[-6,]
dat2

#Random Effects Model
#Restricted Maximum Likelihood Estimation (method="REML")

ma3<-rma(data=dat2,
         yi=g, sei=se_g, method="REML")

summary(ma3)

#Forest Plot for the Standardized Mean Difference
forest(ma3,level=95)

# Funnel Plot to Explore Publication Bias
funnel(ma3)
regtest(ma3) #Test for Funnel Plot Asymmetry

#Explore the Effect of Publication Bias
#By Imputing Missing Studies
trimfill(ma3, estimator="R0")

#### Moderator Analyses ####

#Country as a Moderator of Effect Size
boxplot(dat2$g ~ dat2$country)

#Random Effects Model
#Restricted Maximum Likelihood Estimation (method="REML")
#With moderator "Country"

ma4<-rma(data=dat2,
         yi=g, sei=se_g, method="REML",
         mods=~country)

summary(ma4)
tapply(dat2$g,dat2$country,mean)

 
#Forest Plot with Fitted Values for Each Country
forest(ma4)
forest(ma4,order="obs")

plot(ma4)

# Likelihood Ratio Test for Moderator "Country"
# Must use ML, not REML

ma5a<-rma(data=dat2,
          yi=g, sei=se_g, method="ML",
          mods=NULL)
ma5b<-rma(data=dat2,
          yi=g, sei=se_g, method="ML",
          mods=~country)

anova(ma5a,ma5b)


#### The End! ####