# Install and load required packages
library(meta)
## Loading required package: metadat
## Loading 'meta' package (version 8.0-1).
## Type 'help(meta)' for a brief overview.
# Input data
study <- c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10)
corr <- c(0.45, 0.3, 0.55, 0.4, 0.5, 0.35, 0.6, 0.25, 0.48, 0.38)
n <- c(150, 120, 200, 180, 160, 130, 140, 170, 110, 190)
# Convert correlations to Fisher's z for meta-analysis (better normality)
meta_data <- data.frame(study, corr, n)
meta_data$z <- 0.5 * log((1 + meta_data$corr) / (1 - meta_data$corr)) # Fisher's z
meta_data$se_z <- 1 / sqrt(meta_data$n - 3) # Standard error of Fisher's z
# Perform random-effects meta-analysis (on Fisher's z scale)
meta_result <- metagen(
TE = z,
seTE = se_z,
studlab = study,
data = meta_data,
sm = "ZCOR", # Fisher's z transformed correlations
comb.fixed = FALSE, # Force random-effects model
comb.random = TRUE,
method.tau = "REML" # Restricted maximum-likelihood estimator
)
## Warning: Use argument 'common' instead of 'comb.fixed' (deprecated).
## Warning: Use argument 'random' instead of 'comb.random' (deprecated).
# Print summary
summary(meta_result)
## COR 95%-CI %W(random)
## 1 0.4500 [0.3123; 0.5692] 10.0
## 2 0.3000 [0.1276; 0.4548] 9.1
## 3 0.5500 [0.4452; 0.6399] 11.0
## 4 0.4000 [0.2695; 0.5161] 10.6
## 5 0.5000 [0.3738; 0.6080] 10.2
## 6 0.3500 [0.1892; 0.4925] 9.4
## 7 0.6000 [0.4821; 0.6966] 9.7
## 8 0.2500 [0.1034; 0.3860] 10.4
## 9 0.4800 [0.3217; 0.6122] 8.8
## 10 0.3800 [0.2512; 0.4955] 10.8
##
## Number of studies: k = 10
##
## COR 95%-CI z p-value
## Random effects model 0.4324 [0.3604; 0.4993] 10.61 < 0.0001
##
## Quantifying heterogeneity (with 95%-CIs):
## tau^2 = 0.0123 [0.0023; 0.0563]; tau = 0.1108 [0.0482; 0.2372]
## I^2 = 65.1% [31.6%; 82.2%]; H = 1.69 [1.21; 2.37]
##
## Test of heterogeneity:
## Q d.f. p-value
## 25.81 9 0.0022
##
## Details of meta-analysis methods:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-Profile method for confidence interval of tau^2 and tau
## - Calculation of I^2 based on Q
## - Fisher's z transformation of correlations
# Back-transform Fisher's z to correlations
summary(meta_result, backtransf = TRUE)
## Warning: Additional arguments provided in '...' are ignored.
## COR 95%-CI %W(random)
## 1 0.4500 [0.3123; 0.5692] 10.0
## 2 0.3000 [0.1276; 0.4548] 9.1
## 3 0.5500 [0.4452; 0.6399] 11.0
## 4 0.4000 [0.2695; 0.5161] 10.6
## 5 0.5000 [0.3738; 0.6080] 10.2
## 6 0.3500 [0.1892; 0.4925] 9.4
## 7 0.6000 [0.4821; 0.6966] 9.7
## 8 0.2500 [0.1034; 0.3860] 10.4
## 9 0.4800 [0.3217; 0.6122] 8.8
## 10 0.3800 [0.2512; 0.4955] 10.8
##
## Number of studies: k = 10
##
## COR 95%-CI z p-value
## Random effects model 0.4324 [0.3604; 0.4993] 10.61 < 0.0001
##
## Quantifying heterogeneity (with 95%-CIs):
## tau^2 = 0.0123 [0.0023; 0.0563]; tau = 0.1108 [0.0482; 0.2372]
## I^2 = 65.1% [31.6%; 82.2%]; H = 1.69 [1.21; 2.37]
##
## Test of heterogeneity:
## Q d.f. p-value
## 25.81 9 0.0022
##
## Details of meta-analysis methods:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-Profile method for confidence interval of tau^2 and tau
## - Calculation of I^2 based on Q
## - Fisher's z transformation of correlations
# Forest plot
forest(meta_result,
leftcols = c("studlab", "corr", "n"),
leftlabs = c("Study", "Correlation", "Sample Size"),
xlab = "Correlation Coefficient (95% CI)",
smlab = "Random-Effects Meta-Analysis")
library(metafor)
## Warning: package 'metafor' was built under R version 4.4.3
## Loading required package: Matrix
## Loading required package: numDeriv
##
## Loading the 'metafor' package (version 4.8-0). For an
## introduction to the package please type: help(metafor)
# Meta-analysis using metafor (Fisher's z scale)
meta_rma <- rma(
yi = z,
sei = se_z,
data = meta_data,
method = "REML",
measure = "ZCOR"
)
# Print results
summary(meta_rma)
##
## Random-Effects Model (k = 10; tau^2 estimator: REML)
##
## logLik deviance AIC BIC AICc
## 5.0683 -10.1366 -6.1366 -5.7421 -4.1366
##
## tau^2 (estimated amount of total heterogeneity): 0.0123 (SE = 0.0090)
## tau (square root of estimated tau^2 value): 0.1108
## I^2 (total heterogeneity / total variability): 65.04%
## H^2 (total variability / sampling variability): 2.86
##
## Test for Heterogeneity:
## Q(df = 9) = 25.8110, p-val = 0.0022
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## 0.4628 0.0436 10.6113 <.0001 0.3773 0.5483 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Funnel plot (for publication bias)
funnel(meta_rma,
xlab = "Fisher's z",
ylab = "Standard Error",
main = "Funnel Plot for Publication Bias")
# Egger's test for funnel plot asymmetry
regtest(meta_rma)
##
## Regression Test for Funnel Plot Asymmetry
##
## Model: mixed-effects meta-regression model
## Predictor: standard error
##
## Test for Funnel Plot Asymmetry: z = -0.1823, p = 0.8553
## Limit Estimate (as sei -> 0): b = 0.5491 (CI: -0.3835, 1.4817)
# Install and load the 'meta' package if not already loaded
library(meta)
# Leave-one-out analysis (using Fisher's z)
leave1out <- metainf(meta_result, pooled = "random")
print(leave1out)
## Influential analysis (random effects model)
##
## COR 95%-CI p-value tau^2 tau I^2
## Omitting 1 0.4304 [0.3496; 0.5048] < 0.0001 0.0146 0.1207 68.9%
## Omitting 2 0.4448 [0.3711; 0.5130] < 0.0001 0.0118 0.1087 64.9%
## Omitting 3 0.4166 [0.3423; 0.4857] < 0.0001 0.0107 0.1036 60.8%
## Omitting 4 0.4361 [0.3556; 0.5102] < 0.0001 0.0145 0.1205 68.6%
## Omitting 5 0.4243 [0.3449; 0.4977] < 0.0001 0.0136 0.1167 67.4%
## Omitting 6 0.4406 [0.3633; 0.5119] < 0.0001 0.0134 0.1156 67.3%
## Omitting 7 0.4124 [0.3448; 0.4758] < 0.0001 0.0080 0.0894 55.3%
## Omitting 8 0.4521 [0.3865; 0.5132] < 0.0001 0.0081 0.0899 54.7%
## Omitting 9 0.4276 [0.3482; 0.5008] < 0.0001 0.0140 0.1183 68.5%
## Omitting 10 0.4384 [0.3588; 0.5118] < 0.0001 0.0142 0.1190 67.9%
##
## Pooled estimate 0.4324 [0.3604; 0.4993] < 0.0001 0.0123 0.1108 65.1%
##
## Details of meta-analysis methods:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Calculation of I^2 based on Q
## - Fisher's z transformation of correlations
# Plot the results
forest(leave1out,
xlab = "Correlation Coefficient (95% CI)",
lab.e = "Overall")
## Warning: Use argument 'label.e' instead of 'lab.e' (deprecated).
# Fit random-effects model
res <- rma(yi = z, sei = se_z, data = meta_data)
# Influence analysis
influence_res <- influence(res)
print(influence_res)
##
## rstudent dffits cook.d cov.r tau2.del QE.del hat weight dfbs inf
## 1 0.1582 0.0540 0.0033 1.2448 0.0146 25.7414 0.0997 9.9658 0.0540
## 2 -1.1272 -0.3568 0.1244 1.0724 0.0118 22.7853 0.0913 9.1314 -0.3572
## 3 1.3116 0.4628 0.1950 1.0315 0.0107 20.4180 0.1096 10.9569 0.4610
## 4 -0.2912 -0.0984 0.0109 1.2499 0.0145 25.4845 0.1061 10.6065 -0.0986
## 5 0.6457 0.2182 0.0510 1.1920 0.0136 24.5375 0.1020 10.1973 0.2183
## 6 -0.7022 -0.2268 0.0541 1.1668 0.0134 24.4644 0.0944 9.4363 -0.2264
## 7 1.9582 0.6310 0.3107 0.8547 0.0080 17.9061 0.0971 9.7132 0.6342
## 8 -1.8496 -0.6415 0.3170 0.8693 0.0081 17.6467 0.1041 10.4101 -0.6383
## 9 0.4128 0.1306 0.0184 1.1971 0.0140 25.4108 0.0879 8.7943 0.1300
## 10 -0.4752 -0.1629 0.0294 1.2309 0.0142 24.9387 0.1079 10.7882 -0.1635
# Plot influence diagnostics
plot(influence_res)
# Cumulative meta-analysis (using 'meta')
cumulative_meta <- metacum(meta_result, pooled = "random")
forest(cumulative_meta)
# Add a hypothetical continuous moderator (e.g., year of publication)
meta_data$year <- c(2010, 2012, 2015, 2011, 2013, 2014, 2016, 2017, 2018, 2019)
# Meta-regression with 'year' as a predictor
meta_reg <- rma(
yi = z,
sei = se_z,
data = meta_data,
mods = ~ year, # Test if 'year' affects the effect size
method = "REML"
)
# Summarize results
summary(meta_reg)
##
## Mixed-Effects Model (k = 10; tau^2 estimator: REML)
##
## logLik deviance AIC BIC AICc
## 4.0343 -8.0686 -2.0686 -1.8303 3.9314
##
## tau^2 (estimated amount of residual heterogeneity): 0.0147 (SE = 0.0107)
## tau (square root of estimated tau^2 value): 0.1211
## I^2 (residual heterogeneity / unaccounted variability): 68.81%
## H^2 (unaccounted variability / sampling variability): 3.21
## R^2 (amount of heterogeneity accounted for): 0.00%
##
## Test for Residual Heterogeneity:
## QE(df = 8) = 25.7960, p-val = 0.0011
##
## Test of Moderators (coefficient 2):
## QM(df = 1) = 0.0039, p-val = 0.9502
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## intrcpt -1.5599 32.3757 -0.0482 0.9616 -65.0152 61.8953
## year 0.0010 0.0161 0.0625 0.9502 -0.0305 0.0325
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Visualize the relationship
with(meta_data, plot(year, z, pch = 19, xlab = "Year", ylab = "Fisher's z"))
abline(meta_reg, col = "red")
