StandardMetaAnalysis {RobustBayesianCopas}

Standard meta-analysis

Description

This function performs maximum likelihood estimation (MLE) of (\theta, \tau) for the standard random effects meta-analysis model,

y_i = \theta + \tau u_i + s_i \epsilon_i,

where y_i is the reported treatment effect for the ith study, s_i is the reported standard error for the ith study, \theta is the population treatment effect of interest, \tau > 0 is a heterogeneity parameter, and u_i and \epsilon_i are independent and distributed as N(0,1).

Usage

StandardMetaAnalysis(y, s, init = NULL, tol=1e-10, maxit=1000)

Arguments

y

An n \times 1 vector of reported treatment effects.

s

An n \times 1 vector of reported within-study standard errors.

init

Optional initialization values for (\theta, \tau). If specified, they must be provided in this order. If they are not provided, the program estimates initial values from the data.

tol

Convergence criterion for the optimization algorithm for finding the MLE. Default is tol=1e-10.

maxit

Maximum number of iterations for the optimization algorithm for finding the MLE. Default is maxit=1000.

Value

The function returns a list containing the following components:

theta.hat

MLE of \theta.

tau.hat

MLE of \tau.

H

2 \times 2 Hessian matrix for the estimates of (\theta, \tau). The square root of the diagonal entries of H can be used to estimate the standard errors for (\theta, \tau).

conv

"1" if the optimization algorithm converged, "0" if algorithm did not converge. If conv=0, then using H to estimate the standard errors may not be reliable.

References

Bai, R., Lin, L., Boland, M. R., and Chen, Y. (2020). "A robust Bayesian Copas selection model for quantifying and correcting publication bias." arXiv preprint arXiv:2005.02930.

Ning, J., Chen, Y., and Piao, J. (2017). "Maximum likelihood estimation and EM algorithm of Copas-like selection model for publication bias correction." Biostatistics, 18(3):495-504.

Examples

############################################
# Example on the antidepressants data set. #
# This is from Section 6.2 of the paper by #
# Bai et al. (2020).                       #
############################################
# Load the full data
data(antidepressants)
attach(antidepressants)

# Extract the 50 published studies
published.data = antidepressants[which(antidepressants$Published==1),]
# Observed treatment effect
y.obs = published.data$Standardized_effect_size
# Observed standard error
s.obs = published.data$Standardized_SE

#################################
# Fit a standard meta-analysis  #
# that ignores publication bias #
#################################

# Set seed
set.seed(123)

SMA.mod = StandardMetaAnalysis(y=y.obs, s=s.obs)

# Point estimate for theta
SMA.theta.hat = SMA.mod$theta.hat
# Use Hessian to estimate standard error for theta
SMA.Hessian = SMA.mod$H
# Standard error estimate for theta
SMA.theta.se = sqrt(SMA.Hessian[1,1]) 
# 95 percent confidence interval 
SMA.interval = c(SMA.theta.hat-1.96*SMA.theta.se, SMA.theta.hat+1.96*SMA.theta.se)

# Display results
SMA.theta.hat
SMA.theta.se
SMA.interval