The use of one-stage meta-analytic method based on individual participant data for binary adverse events under the rule of three: a simulation study

One-stage IPD meta-analysis for rare AE

The general concept of one-stage IPD meta-analysis was described in the introduction. For a binary outcome, it contains two exhaustive and mutually exclusive possible states: occurring or not occurring (Doi & Williams, 2013). The treatment-covariate interaction term is unnecessary under the assumption that all covariates were balanced in groups in RCTs as expected under randomization and assuming missing data does not depend on the covariates differently for some arms of the study. Then the multilevel variance component logistic regression (MVCL) model for one-stage IPD meta-analysis is: $y_{ij}\sim \text{ln}\left(\frac{p_{ij}}{1-p_{ij}}\right)={\text{α}}_j+\text{θ}{x}_{ij}+e_{ij}$ $e_{ij}\sim N\left(0,{\sigma }_e^2\right)$ where $y_{ij}\sim \text{Bin}\left(n_{ij},p_{ij}\right)$ ${\text{α}}_j={\text{α}}_0+u_j$ $u_j\sim N\left(0,{\text{σ}}_u^2\right),\text{var}\left(p_{ij}\right)=p_{ij}\left(1-p_{ij}\right)/n_{ij}$

Here i denotes individuals, j denotes studies included, and y is the binary outcome assigned as 1 (occurring) or 0 (not occurring), which obeys a binomial distribution of n (sample size) and p (probability). p denotes the estimated probability of a rare event occurring (y = 1). α_j is the intercept and consists of the fixed effect term α₀ and the random term u_j, with the variance σ_u² that refers to the between-study variance. θ is the fixed coefficient terms that refers to the occurring of AEs, and x is the matrix of covariates for the group of individuals (1 for the treatment group and 0 for the control group). e_ij is the random error of the individual level, with the variance σ_e² that refers to the within-study variance.

Under this model, the pooled odds ratio (OR) for occurring AEs can be estimated directly in the multilevel regression model as: $\widehat{\text{OR}}=\text{exp}\left(\text{θ}\right)$

The estimation of θ in this model can be achieved by maximum likelihood or restricted maximum likelihood algorithms. The total variance of the model is the sum of between-study variance and within-study variance (σ_e² + σ_u²). The magnitude of heterogeneity in the MVCL one-stage IPD meta-analysis can be expressed as: $I^2=\frac{{\text{σ}}_u^2}{{\text{σ}}_e^2+{\text{σ}}_u^2}$

In the multilevel regression model, this term is also called the variance partition coefficient.

Two-stage IPD meta-analysis for rare AEs

For standard two-stage IPD meta-analysis, the study-specific effect θ_j and the variance were obtained from each study using the regression model or 2 × 2 table.

Then the five classical methods listed above were used to combine θ_j and the variance for rare binary AEs in the second stage. An assumption was made for the θ_j as: ${\text{θ}}_j\sim N\left(\widehat{\text{θ}},{\text{σ}}_j^2+{\tau }^2\right)$

Here τ² is the heterogeneity between studies and σ_j² is the within-study variance for θ_j. When τ² = 0, the model denotes a fixed effect model; otherwise, it is a random effect model. W_j denotes the weight of each study, that: $W_j=1/\left({\text{σ}}_j^2+{\tau }^2\right)$

For the sake of illustration, summary data of the 2 × 2 table of each study is introduced. We denote a_j, b_j, c_j, and d_j as the number of events and the number of non-events in the intervention and control groups. Under the IV method: ${\text{σ}}_{\text{IV}-j}^2=1/a_j+1/b_j+1/c_j+1/d_j$

For the MH method, the variance is: ${\text{σ}}_{\text{MH}-j}^2=\left(a_j+b_j+c_j+d_j\right)/b_j{c}_j$

The pooled estimate of the IV and MH methods is calculated as: $\widehat{\text{θ}}=\frac{\left({{\displaystyle \sum }}^{\text{}}{W}_j{\text{θ}}_j\right)}{\left({{\displaystyle \sum }}^{\text{}}{W}_j\right)}$

Unlike the IV and MH methods, the Peto method only contains a fixed effect model, with the pooled estimate as: $\widehat{{\text{θ}}_{\text{peto}}}=\text{exp}\left\{\left(a_j-E\left[a_j\right]\right)/v_j\right\}$ where E[a_j] is the expected number of events in the intervention group and v_j is the hypergeometric variance of a_j. The calculation of E[a_j] and v_j have been described elsewhere (Harris, 2008).

The magnitude of the heterogeneity of between-study variance on the total variance is then: $I^2=\frac{{\tau }^2}{{\text{σ}}^2+{\tau }^2}$

Note that the within-study variance (σ²) is assumed to be equal across studies.

Simulation parameters setting

The simulation was aimed at comparing the performance of the MVCL one-stage IPD meta-analysis method with five classical methods (based on the two-stage IPD meta-analysis method) for rare binary AEs under different simulation scenarios. In order to generate individual participants’ binary AE data, six key parameters were set, with different values to represent different simulation scenarios. The parameters were: incidence rate of rare AEs in the control group (pc); numbers of patients in control group in each individual study (n); adverse effect size (odds ratio, OR); between-study heterogeneity (tau, tau² = τ²); number of studies in each meta-analysis (m); and the sample size ratio between the treatment and control group (r).

We first assigned pc with four values of 0.05, 0.01, 0.001, and 0.0001, which represented common, uncommon, rare, and very rare AE probabilities according to previous studies (Bhaumik et al., 2012; Rücker et al., 2009; Yusuf et al., 1985; Bradburn et al., 2010). The n can be calculated as $\left(\frac{1}{pc}\right)*3$, based on the rule of three. Then, the sample sizes (n) of the control groups under the rule of three were 60 (pc = 0.05), 300 (pc = 0.01), 3,000 (pc = 0.001), and 30,000 (pc = 0.0001). In order to compare with the standard rule of three, we also considered two other situations to determining the individual study sample size—the rule of two, $\left(\frac{1}{pc}\right)*2$, and the rule of one, $\left(\frac{1}{pc}\right)*1$. The sample sizes (n) of the control group were then 40 (pc = 0.05), 200 (pc = 0.01), 2,000 (pc = 0.001), and 20,000 (pc = 0.0001) under the rule of two; and 20 (pc = 0.05), 100 (pc = 0.01), 1,000 (pc = 0.001), and 10,000 (pc = 0.0001) under the rule of one. The expected events in the control group were 3, 2, and 1 under the rule of three, rule of two, and rule of one, respectively. The group ratio (r) was set as 1:1, and then the incidence rate in the treatment group could be calculated by pc and OR.

We considered the effect sizes as four situations: no effect (OR = 1.0), mild harmful effect (OR = 1.25), moderate harmful effect (OR = 2), and large harmful effect (OR = 5). Theoretically, the beneficial effects would be consistent with the harmful effects on the simulation performance; therefore, we did not set scenarios on beneficial effects.

For τ, we considered three situations, as previous literature demonstrated that 0.1, 0.5, and 1.0 indicated mild, moderate and substantial between-study heterogeneity, respectively (Cheng et al., 2016).

For the number of included studies, Cheng et al. (2016) reviewed the previous literature and set m as 5, because that is the median number of published meta-analyses of RCTs. We agree with this assumption; however, we empirically added two scenarios (10 and 15) to cover more situations that are also commonly used in similar simulation studies (Doi et al., 2015a, 2015b). The key R code for simulations is presented in the Supplementary File.

Simulation scenarios and statistics analysis

We compared the performance of the methods under different scenarios of sample size setting rules, effect sizes, between-study heterogeneity, and numbers of studies in each meta-analysis. This was achieved by setting two types of parameters: fixed and varied. For fixed parameters, a fixed value was assigned in each set of scenario; for varied parameters, different values were assigned according to the aim of comparison. For example, when comparing the performance of the methods under different sample size setting rules, we treated effect size, between-study heterogeneity, and numbers of studies as fixed parameters, and treated the sample sizes as varied parameters. Detailed simulation scenarios can be seen in Table 1.

We compared the statistical performance of the MVCL one-stage method with the five classical meta-analytic methods in different scenarios. The following statistics were used to measure the performance: (1) the percentage bias, which is the percentage difference between the pooled OR and the true value; (2) the mean square error (MSE), which refers to the sample standard deviation of true and observed values; (3) the coverage probability, which indicates the proportion of times that the interval contains the true value; and (4) the average width of the 95% confidence interval, which is a reflection of the precision of effect sizes. Generally, smaller percentage bias, smaller MSE, larger coverage, and smaller average width of the 95% CI indicates better effect estimation.

A total of 52 simulation scenarios were defined according to the above conditions (Table 1). We simulated 1,000 data sets for each scenario, to ensure the accuracy of our simulation results. All of the simulations and statistics were conducted on R software (R Core Team, 2017); the MVCL one-stage IPD meta-analysis was based on the lme4 package and the five classical meta-analytic methods were based on the meta package.

Performance under different sample size setting rules

Table 2 presents the performance of the six methods under three sample size setting rules with fixed values of effect size (OR = 1.25, mild effect), between-study heterogeneity (tau = 0.1, mild heterogeneity), and number of studies (m = 10). Our simulation study suggested that, under the rule of three, the MVCL one-stage IPD method had the lowest percentage of bias (<4%) of OR, regardless of the incidence rate. We then compared the performance of the six methods under the rule of two, in which the expected number of events in the control group decreased to two in a single trial. The bias increased for all of the six methods. However, the one-stage IPD method still generally had the lowest amount of bias, with the percentage of bias ranging from 4.9% to 7.2%. This was not the case under the rule of one when the incidence rate decreased to 0.001—the fixed effect IV, random effect IV, and random effect MH had the lowest biases; while when the incidence rate decreased to 0.0001, except for Peto method, the remaining methods tended to have similar biases (Table 2).

For other performance indicators, generally the MVCL one-stage IPD method had the lowest MSE and the narrowest average width of 95% confidence intervals. However, the MVCL one-stage IPD method did not show better coverage probability compared to the other five methods (Fig. 1).

Performance with different OR

Table 3 presents the performance of the six methods with different magnitudes of OR under the rule of three, with fixed between-study heterogeneity (tau = 0.1, mild heterogeneity) and numbers of studies (m = 10). The results of the simulation showed that the MVCL one-stage IPD meta-analysis generally gave the lowest percentage bias of the estimated OR and the bias remained very low (2.8–4.1%), regardless of the magnitude of the OR. The Peto method ranked the second lowest in bias in nine of 16 scenarios. However, when the incidence rate was 0.001 and 0.0001, the fixed and random effect IV methods and the random effect MH method performed better than the MVCL one-stage IPD when there were moderate and large effects (OR = 2 and OR = 5). We observed unstable results for the Peto method when there was a large effect (OR = 5) in different incidence rates.

Figure 2 presents the average MSE, coverage probability, and average width of the 95% confidence intervals. The MVCL one-stage IPD method had the lowest MSE and the narrowest average width of 95% confidence intervals. All of the six methods had good coverage probabilities (>95%), while the IV method and random effect MH method tended to perform better than the MVCL one-stage IPD. Overall, with increasing OR there was an increasing trend for the MSE and the average width of the 95% confidence intervals, but the coverage remained stable.

Performance with different heterogeneity

Table 4 presents the performance of the six methods with different between-study heterogeneity under the rule of three, with fixed effect size (OR = 1.25, mild effect) and numbers of studies (m = 10). The results showed that the MVCL one-stage IPD method performed well, with lower percentage bias than the other five methods when there was mild or moderate heterogeneity between studies. The Peto method performed better than the IV and MH methods in five of the 12 scenarios. However, when there was substantial between-study heterogeneity, all six methods presented an obviously large amount of bias such the most of the ORs were overestimated by more than 50%. For the IV, MH, and Peto methods, the percentage bias decreased as the heterogeneity decreased.

We observed similar results for other performance indicators—the MVCL one-stage IPD method generated the lowest MSE and the narrowest average width of the 95% confidence intervals, and the IV and random effect MH methods tended to cover more true values. The MSE and the average width of the confidence intervals increased when higher between-study heterogeneity was present (Fig. 3).

Performance with different numbers of studies

Table 5 presents the comparison of the MVCL one-stage IPD method with the other five methods with different numbers of studies included under the rule of three, with fixed effect size (OR = 1.25, mild effect) and between-study heterogeneity (tau = 0.1, mild heterogeneity). Our simulation data showed that the MVCL one-stage IPD method consistently presented lower bias than the other five methods. Furthermore, when the number of included studies increased, the bias decreased with all of the methods. When the number of included studies was 15, the bias remained very low level with all of the methods.

The MSE and average width of confidence intervals of the MVCL one-stage IPD method remained the best of the methods, and the number of studies was inversely associated with the magnitude of MSE and the average width of the 95% confidence intervals. The coverage of the MVCL one-stage method generally performed less well than the IV method and the random effect MH method (Fig. 4).