Meta-analysis with the R package confMeta

Simulating individual studies

Consider the hypothetical scenario where we have n = 5 individual studies that should be combined into a single confidence interval using a confidence level of 1 - alpha = 0.95. We can simulate these by running the following code

library(confMeta)
library(ggplot2)

set.seed(42)

n <- 5
conf_level <- 0.95
estimates <- rnorm(n, mean = 0.5, sd = 0.2)
SEs <- rgamma(n, shape = 5, rate = 20)
study_names <- c("Study A", "Study B", "Study C", "Study D", "Study E")

Our goal is to combine them. However, this requires the specification of a \(p\)-value function, i.e. a method, that takes the individual studies (argument estimates), their standard errors (argument SEs), and the mean under the null-hypothesis (argument mu) as input and returns the corresponding \(p\)-value at the specified mean value.

The core function of the package is confMeta(). You provide your data and specify which \(p\)-value combination function you want to use. The possible choices provided by the package are:

• Harmonic mean (Function: p_hmean)
• Wilkinson (Function: p_wilkinson)
• Pearson (Function: p_pearson)
• Edgington (Function: p_edgington)
• Weighted Edgington (Function: p_edgington_w)
• Fisher (Function: p_fisher)
• Tippett (Function: p_tippett)
• Stouffer (Function: p_stouffer)

The documentation for this function can be inspected by running

?confMeta
?p_value_functions 

For this example, let’s use Edgington’s method, which is based on the sum of individual \(p\)-values. Thus, we can create the confMeta object as follows

cm_edgington <- confMeta(
  estimates = estimates,
  SEs = SEs,
  study_names = study_names,
  conf_level = conf_level,
  fun = p_edgington,
  fun_name = "Edgington"
)

As the variable cm_edgington now contains the confMeta object, we can inspect it by running the following code

# Use the specific print argument
print(cm_edgington)
#>   Meta-Analysis with p-value combination
#> -----------------------------------------
#> Number of studies:  5 
#> Combination method: Edgington 
#> 
#> Estimate: 0.582,
#> 95% Confidence Interval: [0.415, 0.758]

# See what elements it has
names(cm_edgington)
#>  [1] "estimates"      "SEs"            "w"              "study_names"   
#>  [5] "conf_level"     "p_fun"          "fun_name"       "individual_cis"
#>  [9] "joint_cis"      "gamma"          "p_max"          "p_0"           
#> [13] "aucc"           "aucc_ratio"     "comparison_cis" "comparison_p_0"
#> [17] "heterogeneity"  "table_2x2"


# Check out the combined confidence interval(s)
cm_edgington$joint_cis
#>          lower     upper
#> [1,] 0.4148163 0.7576969

When created, the confMetaobject automatically stores several information:

• The combined confidence interval inverting your chosen \(p\)-value function (identifying the range of \(\mu\) values where the \(p\)-value exceeds \(\alpha\)).
• Individual CIs computing the Wald-type confidence intervals for each individual study.
• Reference Comparisons fitting multiple standard meta-analytic models to your data for comparison. These include fixed-effect, random-effects, Hartung-Knapp and Henmi-Copas.
• Heterogeneity Statistics providing estimates of Cochran’s \(Q\), \(I^2\), and \(\tau^2\).

Weighted Methods

Certain combination rules, such as p_edgington_w or p_stouffer, accept study weights. A common approach is to weight studies by the inverse of their standard errors to give more precise studies greater influence on the combined \(p\)-value.

weights <- 1 / SEs

cm_weighted <- confMeta(
  estimates = estimates, SEs = SEs,
  fun = p_edgington_w, fun_name = "Edgington (Weighted)",
  w = weights, input_p = "two.sided"
)

Heterogeneity adjustments

By default, the \(p\)-value functions are evaluated under a fixed-effect model (heterogeneity = "none"). You can change this by applying either an additive or multiplicative heterogeneity correction.

These corrections can be achieved using estimate_tau2() and estimate_phi(). It is necessary to compute these between-study variance parameters from your data before passing them into the main combination functions.

1. Additive Heterogeneity (Random Effects) The additive model assumes a between-study variance of \(\tau^2\). You can estimate this using the estimate_tau2() function, which is a wrapper of meta::metagen(). Here, we use the Paule-Mandel (“PM”) estimator:

# Estimate tau^2
tau2_est <- estimate_tau2(estimates, SEs, method.tau = "PM")

# Create a random-effects confMeta object using additive heterogeneity
cm_additive <- confMeta(
  estimates = estimates, 
  SEs = SEs,
  study_names = study_names,
  fun = p_edgington, 
  fun_name = "Edgington (Additive RE)",
  heterogeneity = "additive", 
  tau2 = tau2_est,
  input_p = "two.sided"
)

2. Multiplicative Heterogeneity Alternatively, the multiplicative model assumes the variances are scaled by \(\phi\). You can estimate this using the estimate_phi() function:

# Estimate phi
phi_est <- estimate_phi(estimates, SEs)

# Create a confMeta object using multiplicative heterogeneity
cm_multiplicative <- confMeta(
  estimates = estimates, 
  SEs = SEs,
  study_names = study_names,
  fun = p_edgington, 
  fun_name = "Edgington (Multiplicative)",
  heterogeneity = "multiplicative", 
  phi = phi_est,
  input_p = "two.sided"
)

By supplying these adjusted objects to autoplot(), you can visually compare how the choice of heterogeneity model impacts the width of your combined confidence intervals.

autoplot(cm_edgington, cm_additive, cm_multiplicative)