##############################################################
# A simple example:
# Checking the example in Harris Cooper et al.'s book page 388
##############################################################
r <- matrix(c(-0.074, -0.127, 0.324, 0.523, -0.416, -0.414), 1)
n <- 142
computvcov <- r.vcov(n = n, corflat = r,
name = paste("C", c("st", "su", "sv", "tu", "tv", "uv"), sep = ""),
method = "each")
round(computvcov$list.rvcov[[1]], 4)
round(computvcov$ef, 4)
round(computvcov$list.vcov[[1]], 4)
round(computvcov$matrix.vcov, 4)
######################################################
# Example: Craft2003 data
# Preparing covariances for multivariate meta-analysis
######################################################
data(Craft2003)
# extract correlation from the dataset (craft)
corflat <- subset(Craft2003, select=C1:C6)
# transform correlations to z and compute variance-covariance matrix.
computvcov <- r.vcov(n = Craft2003$N, corflat = corflat, method = "average")
# name transformed z scores as y
y <- computvcov$ef
# name variance-covariance matrix of trnasformed z scores as S
S <- computvcov$matrix.vcov
S[1, ]
## fixed-effect model
MMA_FE <- summary(metafixed(y = y, Slist = computvcov$list.vcov))
#######################################################################
# Running random-effects model using package "mixmeta" or "metaSEM"
#######################################################################
# Restricted maximum likelihood (REML) estimator from the mixmeta package
#library(mixmeta)
#mvmeta_RE <- summary(mixmeta(cbind(C1, C2, C3, C4, C5, C6) ~ 1,
# S = S, data = y, method = "reml"))
#mvmeta_RE
# maximum likelihood estimators from the metaSEM package
# library(metaSEM)
# metaSEM_RE <- summary(meta(y = y, v = S))
# metaSEM_RE
##############################################################
# Plotting the result:
##############################################################
obj <- MMA_FE
# obj <- mvmeta_RE
# obj <- metaSEM_RE
# pdf("CI.pdf", width = 4, height = 7)
plotCI(y = computvcov$ef, v = computvcov$list.vcov,
name.y = NULL, name.study = Craft2003$ID,
y.all = obj$coefficients[,1],
y.all.se = obj$coefficients[,2])
# dev.off()
Run the code above in your browser using DataLab