{'fraction_non_alphanumeric': 0.033338587641866334, 'fraction_numerical': 0.01455443463640185, 'mean_word_length': 4.8008229198415115, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Evidence synthesis with reconstructed survival data

We present a general approach to synthesizing evidence of time-to-event endpoints in meta-analyses of aggregate data (AD). Our work goes beyond most previous meta-analytic research by using reconstructed survival data as a source of information. A Bayesian multilevel regression model, called the"meta-analysis of reconstructed survival data"(MARS), is introduced, by modeling and integrating reconstructed survival information with other types of summary data, to estimate the hazard ratio function and survival probabilities. The method attempts to reduce selection bias, and relaxes the presumption of proportional hazards in individual clinical studies from the conventional approaches restricted to hazard ratio estimates. Theoretically, we establish the asymptotic consistency of MARS, and investigate its relative efficiency with respect to the individual participant data (IPD) meta-analysis. In simulation studies, the MARS demonstrated comparable performance to IPD meta-analysis with minor deviation from the true values, suggesting great robustness and efficiency achievable in AD meta-analysis with finite sample. Finally, we applied MARS in a meta-analysis of acute myeloid leukemia to assess the association of minimal residual disease with survival, to help respond to FDA's emerging concerns on translational use of surrogate biomarker in drug development of hematologic malignancies.


Introduction
The time to occurrence of a disease-related event is the most important outcome measure in clinical trials to assess the efficacy of an intervention. Because the clinical and regulatory evaluation of a new intervention usually involves multiple trials, evidence synthesis using meta-analysis that integrates the survival data of such trials is essential (FDA 2018).
In systematic reviews, aggregate data (AD) meta-analysis is the mainstay approach to combining the results of multiple trials. Even though individual patient data (IPD) metaanalysis is considered the gold standard (Stewart & Clarke 1995, Riley et al. 2010, IPD were hard to obtain (Tierney et al. 2015, Marciniak 2017. A recent systematic survey showed that less than 1% of 78,629 meta-analyses included complete or partial IPD, while the vast majority were based on AD only (Nevitt et al. 2017). Thus, it is of utmost importance to develop sensible methodology for AD meta-analysis.
We are concerned in this work with AD meta-analytic methods to synthesize evidence in time-to-event data. For survival outcomes, the standard approach entails a fixed or random effects model on aggregated study-level hazard ratio (HR) estimates (Moher et al. 2012, Higgins et al. 2020. Owing to data censoring in studies of survival outcomes, HRs are useful for describing differences in treatment effects between an intervention group and control group. Abundant open-source and commercial computational software programs are also conveniently available for conducting meta-analysis. However, this standard approach poses several methodological limitations and challenges in analysis and reporting. First, applications are constricted under the assumption of proportional hazards. Recent systematic reviews of randomized controlled trials showed that the proportional hazards assumption could be significantly violated in one-fifth of treatment comparisons (Rulli et al. 2018), while testing of the proportional hazards as-sumption seldom was conducted (Batson et al. 2016). In fact, significant violations, such as survival curve convergences and crossings, are common in clinical studies, which gives rise to concerns upon general model validity of AD meta-analysis for time-to-event outcomes (Kristiansen 2012).
Second, the exclusive use of HR estimates impairs the effectiveness and generalizability of meta-analysis. The strength of systematic review is that it allows evidence to be synthesized from all available clinical information in the literature and allows meta-analysis to combine the information comprehensively. Using HR estimates as the only source of "survival data" simplifies analyses but undermines the generalizability of results by excluding studies that alternatively report other kinds of survival data, such as Kaplan-Meier (KM) curves, log-rank test statistics, median survival times or survival rates. Existing methods to convert those types of survival data to HR estimates involve additional assumptions about time-to-event outcomes (Tudur et al. 2001, Tierney et al. 2007) and, given the various duration of studies, could assign study weights improperly (Riley et al. 2007).
The last limitation is pertinent to AD meta-analytic reporting for intervention effect.
If the proportional hazards assumption is violated in a single study, the interpretation of the overall meta-analytic HR statistics could be a problem. Moreover, it is preeminent to determine survival curves, rather than a hazard ratio, to gauge the risk of the event over time (Pocock et al. 2002). Although various methods have been developed in the literature for combining published survival curves (Begg et al. 1989, Voest et al. 1989, Shore et al. 1990, Reimold et al. 1992, Hunink & Wong 1994, Dear 1994, Arends et al. 2008, they were all rooted in the proportional hazards assumption except Combescure et al. (2014). That approach applies a distribution-free approach to estimating the summary survival curve for single-arm studies, but is subject to a risk of bias when multiple curves from the same study are modeled. In brief, AD meta-analytic reporting that entails comparative trials without the proportional hazards assumption remains unclear.
The main contribution of this paper is to address these challenges and present a novel AD meta-analysis framework for survival-type endpoints. Particularly, a Bayesian multilevel modeling strategy, coined "meta-analysis of reconstructed survival data" (MARS), is proposed to exploit reconstructed survival curves as a richer source of survival information.
In general, survival data can be reconstructed by 1) acquiring the positions of censoring and event points on images of KM survival curves; 2) performing iterative numerical computation to convert 2-dimensional position data into time-to-censoring and time-to-events.
The use of the reconstructed survival data, which can be regarded an alternative to IPD, enables several advantages similar to using IPD in both modeling and reporting. For example, instead of assuming proportional hazards, MARS can estimate the time-varying hazard rates and HRs. We note that KM curves may not be provided in every published clinical study. To conduct a meta-analysis that thoroughly covers the literature, MARS also simultaneously accommodates 2 more types of AD: HR estimates and survival rates at various time points.
Our goal is to deliver effective estimation for the comparison of efficacy of treatments, and thus expand the traditional framework for AD meta-analysis. The results of simulation study suggest that MARS is competent to provide valid estimation of a treatment effect in both cases of constant HR and time-varying HR, and capable of reporting informative summary metrics, such as survival probabilities, median survival times and restricted mean survival times (RMSTs). Moreover, MARS generates HR estimates that are consistent and have decent asymptotic efficiency with respect to IPD meta-analysis.
The outline of the paper is as follows. In Section 2, we present a data set in adult leukemia with disease-free survival (DFS) as the time-to-event outcome of interest. Our work was originally motivated by a larger AD study consisting of pediatric and adult leukemia patients with DFS endpoints. Here we apply a subset of 64 non-transplantation trials to illustrate the concept. In section 3, we introduce the MARS, which can integrate various types of survival information, and discuss the metrics to quantify treatment effects for informative reporting. The theoretical properties of MARS are investigated in section 4. Simulation studies are illustrated in section 5 for finite sample performance and an application of MARS to leukemia data set is demonstrated in section 6. Lastly, section 7 provides a discussion of our work, with suggestions for further research.

Case study: minimal residual disease in leukemia
Responding to the emerging development of targeted and immunological treatments, the U.S. Food and Drug Administration (FDA 2020) published a guideline clarifying for the use of minimal residual disease (MRD) status in the development of drugs for the treatment of hematologic malignancies. An MRD-positive status means that the disease was still detectable after treatment, whereas an MRD-negative status means that no disease was detected after treatment. Potentially, the MRD status can serve as a prognostic biomarker to assess a patient's risk of future relapse for regulatory and clinical uses.
Evidence supporting the clinical validity of MRD status as a surrogate biomarker varies across hematologic cancer types and patient populations (Rach et al. 2016, Berry et al. 2017). To gain a better understanding of the state-of-the-science of the use of MRD for specific hematologic malignancies, the FDA recommended a systematic review and metaanalysis. In particular, to ensure the effective validation of MRD to support the marketing approval of drugs and biological products, the FDA gave special meta-analytic considera-tions to (i) prespecified inclusion criteria to ensure comprehensive literature search of both positive and negative results, to limit selection bias, and (ii) sensitivity analyses (e.g., alternative statistical methods) to demonstrate the robustness of MRD as a surrogate marker.

Model for survival AD meta-analysis
We present a general framework for survival meta-analysis, called the "meta-analysis of reconstructed survival data" (MARS). MARS attempts to establish a one-stage Bayesian multilevel model with separate AD components, such that the most informative summary survival data available in individual studies can be properly integrated to the parameter estimation and uncertainty quantification.
In this section, we first review the general strategy to reconstruct the survival data from KM curves, and discuss the model specification to relax the proportional hazards assumption. If KM curves in a decent quality are unavailable in some clinical studies, we elicit and incorporate other types of survival data, though less informative, to complement a comprehensive AD meta-analysis and diminish selection bias. We provide the hierarchical prior distributions to complete the model specification, and finally present several alternative summary measures from MARS to better inform AD meta-analysis reporting.

Reconstructed survival data
KM curves are the richest source of information for AD meta-analysis. In the absence of IPD, it can be considered as a poor man's IPD if the survival data are properly reconstructed. The overview of a general strategy for data extraction is as follows.
The first step in reconstructing the survival data is to digitize the KM curves using sufficiently large, clear, and high-resolution images. In electronically published articles, images are typically either 'raster images' or 'vector objects'. A raster image consists of pixels, each represented by a dot or square with its own coordinates in a 2-dimensional grid. A graphical digitizing software, such as GraphClick or DigitizeIt, can be employed to convert and store the corresponding coordinates of the change points on the KM curves by clicking on the individual points using a mouse or automatic read-off (Guyot et al. 2012).
Alternatively, Liu et al. (2014) proposed directly reading in the lines of the PostScript file of a vector image, rather than using a digitizing software program to read in the coordinates from a raster image, and reported an improvement with 200 times better resolution.
After the curves have been digitized, the next step is to implement a numerical algorithm to elicit survival data. Early work in this area involved dividing the KM curves into several time intervals, which gives a good representation of event rates over time, while limiting the number of events (e.g., 20%) within a given time interval Parmar et al. (1998). It was later extended and implemented in a spreadsheet calculator by Tierney et al. (2007). Using supplementary R code, Guyot et al. (2012) Parmar et al. 1998, Williamson et al. 2002, and reported that the Guyot method was consistently superior to others.
Although the modeling of reconstructed survival data is lacking in the literature, the general ideas can be analogues to that of IPD without patient-level covariates. In IPD metaanalysis literature for survival outcomes, it was popular to apply a Poisson generalized linear model, which could be implemented with either fixed or random effects and with baseline hazard stratified by trial (Dear 1994, Arends et al. 2008. Crowther et al. (2012) extended this model to allow non-proportional hazards of the treatment effects.
In MARS, we assume the reconstructed survival data follows a piecewise exponential distribution in the AD framework. In theory, the piecewise exponential survival model can be equivalent to a Poisson log-linear model (Holford 1980, Laird & Olivier 1981. To model the piecewise exponential, we first separate the overall time scale into J prespecified disjointed small segments I j = (t j−1 , t j ] for j = 1, 2, . . . , J, where 0 = t 0 < t 1 < · · · < t J < ∞. For each patient, we create a set of pseudo-observations, one for each interval, indicating whether the individual had an event in that interval. The piecewise exponential model incorporates the data by treating the pseudo-observations as following Poisson distributions with means that are equal to the cumulative hazards in the intervals. Thus, the likelihood function of a piecewise exponential model is equivalent to the product of Poisson likelihood functions, one for each combination of an individual and an interval.
Within each interval, the baseline hazard rate and the log HR for group difference are assumed to be constant, denoted by λ j and β j for j = 1, 2, . . . , J. We denote t J the largest survival or censored time. If some study has an extraordinary long follow-up period compared with other studies, the reconstructed survival time could be censored by a reasonable maximum time to make the inference of the meta-analysis meaningful and robust. The hazard rate for the ith subject in the kth study is modeled as where α k allows the baseline hazard rate to vary from study to study and X k,i indicates the treatment assignment. In the motivating case study (described in Section 2), X k,i is the indicator of MRD status. The study-specific log HR θ k (t) is decomposed into the time-dependent effect β j and an inter-study effect µ k The survival function is uniquely defined given the hazard function. If the study-level covariates (moderators) Z k are available, it is also straightforward to accommodate the study-level covariate effects, γZ k , into µ k in (2) for subgroup analysis.

Supplemental survival data
One objective of MARS is to diminish selection bias. In meta-analysis applications, studies with survival information on intervention effect may not all publish KM curves, but report alternative summary statistics. To include all available information from eligible studies in a comprehensive systematic review, we incorporate into MARS with other common types of survival data, ranked by their richness of survival information. Without ambiguity, we denote the reconstructed survival data as the type I data hereinafter.

Hazard ratios (type II data)
In MARS, HR estimates are synthesized as type II data if the KM curves are not available or are provided in low-resolution images, which is often the case in early publications.
When the reported HRs are estimated under the proportional hazards assumption, they are presumed to be time-independent. Estimating a time-varying HR based solely on reported constant HRs is not feasible. However, given the various lengths of follow-up periods, it is natural to regard the observed constant log HRθ k as an estimate of the weighted average of stepwise time-dependent effects β j in addition to the study-specific effect µ k , with the HR aligned with study duration where θ k (t) is defined in equation (2), |I j | = t j −t j−1 , and T k is the duration of the kth study.
The corresponding likelihood provides no information on the baseline hazard function.
In general, type II survival data can consist of the HRs converted from variants of indirect statistics comparing group differences in clinical studies lacking reported HRs.
A typical example is a log-rank test statistics and its p-value, from which the log HR and its variance can be reasonably approximated (Parmar et al. 1998). This allows the totality of evidence to be assessed and improves the efficiency and reliability of metaanalysis for survival-type endpoints (Tudur et al. 2001). Tierney et al. (2007) reviewed the relevant methods that carefully manipulate published or other summary data to obtain HR summary statistics, classified by the type of information presented in trial reports. The resulting converted HRs can then be conveniently used as type II data in the MARS model.

Survival probabilities (type III data)
Estimates of survival probabilities at particular time points (e.g., 1-year or 3-year DFS rates for individual groups) are often reported in randomized controlled clinical studies, but usually ignored and rarely used in meta-analysis even if HR was not published. Compared with HR estimates, survival rates provide limited information and thus are not considered valid surrogate measures for the meta-analysis of survival outcomes (Michiels et al. 2005). Although it can also be used to estimate HR estimates, such estimation requires assumptions in addition to the proportional hazards assumption. Nevertheless, it put an important piece to complete the information synthesis for systematic review when other types of survival data of a study are not available. For example, Berry et al. (2017) used it in a 2-step approach by first estimating the HR after assuming an exponential distribution of the time-to-event outcomes for the totality of evidence.
In MARS, the reported survival rates are treated as the realizations of the survival function derived from equation (1) S When survival rateŜ k (t, X = x) and its standard errorŝe(Ŝ k (t, X = x)) for xth treatment group are provided at a specific time point t for the kth study, we assume that logit(Ŝ k (t, X = x)) follows a normal distribution as to reduce the bias from asymmetry for interval-ranged S k (t, x) on (0, 1).Ŝ k (t, x) = 1 or 0 has been rarely reported, as the observation would be neglected to avoid the degenerate likelihood.

Priors specification and model implementation
For the log HR estimation, β j represents the time-dependent treatment effect in time interval I j and µ k describes the inter-study heterogeneity. We assume that prior of µ k follows a normal distribution µ k ∼ N (µ 0 , σ 2 µ ). To specify the prior for β = (β 1 , β 2 , . . . , β J ), we rely on the multivariate normal prior: where the mean µ β = µ β 1, and the variance-covariance matrix Σ β takes a first-order autoregressive structure The hyperparameter ρ β induces the association between the time segments and thereby the smoothness of the estimated piecewise HR function. For example, if ρ β = 1, this model implies that the time-dependent treatment effects are identical, that is, degenerate to a proportional hazards model.

Summary measures for meta-analytic reporting
We specify several summary measures, moving beyond the HR estimator, for statistical reporting and result dissemination using MARS. Having precise estimates of the difference in survival between treatments is extremely important, especially in oncology (Michiels et al. 2005). In conventional meta-analysis, HR estimates are routinely used to summarize differences in survival endpoints between groups (Moher et al. 2012, Higgins et al. 2020).
If the underlying proportional hazards assumption is violated, the clinical interpretation of the HR estimate becomes difficult. Alternatively, owing to its flexible Bayesian model structure, MARS can generate more informative summary measurements to quantify the between-group difference with exact inference for uncertainty quantification.

Survival probabilities.
Survival probabilities are essential in survival analysis as they provide crucial summary information about the survival experiences of study populations. MARS readily yields the estimation of the overall survival function The Bayesian model can deliver exact inference on parameter estimation without assuming asymptotic normality. Estimated survival probabilities can function similarly as KM estimates of the survival function so that treatment effects summarized in them are intuitive and straightforward. Using MCMC can generate the 95% credible intervals (CrIs) for uncertainty quantification. Moreover, it is appealing to generate graphical representations of meta-analytic survival curves (see, for example, Figure 3A and 3B).

Median survival time.
The median survival time defines as the amount of time t such that S(t) = 0.5, which is equivalent to that the cumulative hazard function H(t) is In MARS, the cumulative hazard is given by a piecewise linear function where t j−1 < t ≤ t j . Subsequently, the median survival times t M S can be estimated solving of the following equation:Ĥ(t M S ) = log(2).

Restricted mean survival time (RMST).
The RMST is a well-established alternative measure of the "life expectancy" from the time of treatment to a specific time horizon of interest (Royston & Parmar 2013). The RMST is valid in survival models under any distribution, including the proportional hazards model. It is as powerful as the HR when the proportional hazards assumption is valid and even more powerful than the HR when the assumption is invalid. Another advantage of the RMST is that it primarily relies on the chosen time period so that it is less subject to follow-up time (Kim et al. 2017).
The RMST is defined as the expected value of survival time up to a pre-specified time In the proposed piecewise exponential model, the RM ST (τ ) can be approximated using the estimated mean survival probabilitiesŜ: where τ = t j andŜ(t i , X) = exp(− j i=1λ i exp (Xβ)|I i |). Using posterior samples from the MCMC algorithm, estimating the RMST and its standard deviation is succinct.

Asymptotic properties
Meta-analysis is widely used in the current big data era in all scientific areas for thorough statistical inference. As both the scale and number of individual studies increase, it is essential to understand the theoretical properties of MARS.

Consistency
To account for the heterogeneity among studies, MARS employs a hierarchical Bayesian structure in meta-analysis. Using a hierarchical Bayesian model provides consistent and accurate estimation (Rouder & Lu 2005). Assume there are K studies and n k subjects in the k th study. We denote n as average sample size. Let β = (β 1 , β 2 , . . . , β J ) T and η = (η 1 , η 2 , . . . , η J ) T be the parameters of piecewise log HR and piecewise log baseline hazard rates, respectively, and denoteβ andη as the estimates generated by MARS. We denote K 1 , K 2 , K 3 , and K as sets of studies that provide reconstructed survival data (type I data), HRs (type II data), and survival rates (type III data), respectively. We have the following result: Theorem 1. Under the regularity conditions (a)-(f ) specified in Appendix A and for k ∈ K 1 , as the sample size n k → ∞ and the number of patients in each interval n kj becomes larger at a constant ratio n kj /n k = c j , j = 1, . . . , J; Then, as the total number of studies K → ∞ and the average sample size n → ∞,  According to Theorem 1, MARS yields consistent estimators for β and η with sufficient reconstructed survival data. As a result, the estimated summary measures, such as the overall survival rates, median survival time and RMSTs in 3.4, are consistent as well.

Relative efficiency
The relative efficiency of conventional fixed-effect and random-effect AD meta-analyses has been rigorously studied and shown to be at most as efficient as IPD meta-analysis and asymptotically equivalent (Lin & Zeng 2010, Chen et al. 2020  cating that the efficiency of MARS is at most the same as that of IPD meta-analysis. The second statement implies that asymptotically, MARS can achieve the full efficiency when the sample size increases at a higher order rate than the study number. Note that the condition of variance for type II data is met when the hazard ratio is constant. In particular, it holds even if hazard ratio is constant within the period of study duration T k and becomes time-varying. This might not be practical. Thus, we encounter acceptable efficiency loss in reality. Corollary. Letβ denote the estimates generated from MARS by using all 3 types of data, for fixed K and n, V ar(β) ≤ V ar(β * ).
The corollary indicates that for a meta-analysis with finite number of studies, including studies that only reported type III data increases efficiency.

Simulation
We empirically evaluated the finite sample performance of the proposed meta-analysis method in 2 simulation studies. First, we compared the MARS to conventional AD metaanalysis method and IPD meta-analysis in a setting in which the proportional hazards assumption is valid (Case 1). Second, we assessed the validity of MARS and its relative efficiency as compared to IPD meta-analysis when the proportional hazard assumption is invalid but two survival curves cross each other (Case 2). The performances of the methods are compared based on 1,000 replications.

Simulation setting
Case 1.
In each replication, we assumed 13 eligible studies and proportional hazards over time. The random effects for each study were generated as follows: where k = 1, 2, . . . , 13 and X's are independent binary variables generated from a Bernoulli distribution with a rate of 0.55. Given α k and µ k , the survival times in the k th study were independently generated from the hazard function where h 0 (t) is a Weibull distribution with shape 0.1 and scale 0.5 and β = −0.6 (exp(β) ≈ 0.55) is the time independent effect. We assumed that survival times were right-censored according to an exponential distribution with a censoring rate of 0.01. The sample size for each study was randomly sampled from 50 to 150. The follow-up time for each study was randomly sampled from 80 to 120 months, at multiples of 5 months only. We equally partitioned the 120 months into 10 segments, each 12 months long.
To mimic the real-world scenario that HRs may not be reported in every clinical study, we assumed that 6 studies reported KM plots, 6 studies reported HR (including 3 studies with KM plots and 3 without), and 4 studies reported survival rates. For 3 studies that had both KM data and HR estimates, only the KM data were used for MARS, and only the HR estimates were used for AD random-effects model. Therefore, we applied MARS to 6 studies of type I data, 3 studies of type II data, and 4 studies of type III data. In contrast, in the conventional AD meta-analysis, 6 studies were included to estimate the HR. Original survival data from all 13 studies were included in an IPD meta-analysis as a gold standard for comparison.
We computed the log HRs using MARS, IPD meta-analysis and the conventional AD meta-analysis. For MARS and IPD meta-analysis, we also computed survival probabilities, median survival times and RMSTs at 1, 3, 5, and 10 years for each group.
In each replication, we included 20 studies, for which α k , µ k and X were generated the same as in case 1. The survival times in each study were generated from the hazard function h k (t) = 0.01 × exp(α k + (µ k + β(t))X) where the log HR β(t) = 1 − 0.2 (t/30−1) is time-varying and for simplicity, the baseline hazard rate was constant over time. The censoring mechanism and the follow-up time were the same as those in case 1.
We assumed that 9 studies reported KM plots, 12 studies reported HR (including 5 studies with KM plots), and 4 studies reported survival rates. Thus, in MARS, we accounted for 9 studies of type I data, 7 studies of type II data, and 4 studies of type III data. In contrast, in the conventional AD meta-analysis, 12 studies were included to estimate the HR. Data from all 20 studies were included in an IPD meta-analysis.
As we did in case 1, we calculated survival rates, median survival times, and RMSTs at 1, 3, 5 and 10 years. Alternative to the point estimate of log HR, which is meaningless when the HR is not constant, we computed the average of the squared error between the log HR function β(t) and the estimate log HR functionβ(t) over time, whereβ(t) would be a stepwise function in MARS.

Results
Case 1. Figure 1 shows the performance of MARS in which the proportional hazard assumption holds. The piecewise HRs estimated by MARS were approximately constant ( Figure 1A), highlighting that MARS is able to identify a proportional hazards situation when the HR is actually constant. Thus, we compared the average piecewise HR estimates generated by MARS with the HR estimates generated by AD meta-analysis and IPD meta-analysis (Table 1). Conventional AD meta-analysis using the available HR estimates from individual studies yielded the log HR estimates closest to the true value of -0.6. Although MARS generated a point estimate slightly smaller than the true value, the coverage probability of the 95% CrIs was non-inferior to those of the other 2 methods. The estimated survival probabilities of both groupsgroup X=0 and group X=1 are shown in Figure 1B and 1D, respectively. MARS provided similar point estimates with slightly larger variability. Specifically, the true survival rates of the X = 0 group at 3, 5, and 10 years were 0.55, 0.46 and 0.33, respectively (table 2). With MARS, the 95% CrIs were more likely to cover the true survival rates at 3 years, 5 years and 10 years, whereas the IPD meta-analysis provided a more precise estimate of those survival rates.
Regarding the median survival time and RMST ( Figure 1C and Table 2), MARS had similar operating characteristics to those of the IPD analysis, and a higher coverage probability due to wider 95% CrIs length. The mean relative efficiency of MARS to IPD metaanalysis was 1.22, indicating that MARS had reasonable efficiency loss compared with IPD meta-analysis.
Case 2. IPD meta-analysis but much smaller than that of conventional AD meta-analysis, which is inappropriate in this situation (Table 1).
MARS had good performance in estimating survival probabilities ( Figure 2B and 2D).
The survival rates MARS generated for both groups were almost identical to those generated from IPD meta-analysis. In long-term follow-up, both methods slightly over-estimated survival rates for group X = 1. This can be explained by the shrinkage effect in Bayesian estimation of time-varying HR function 2A). Even though the HRs were over-estimated in early years, it effected less to the estimated survival rates due to the low hazard rates in that period. The 95% CrIs generated with MARS were a little wider, but also had a higher   (Table 2).
MARS was also able to give valid estimates of median survival time and RMSTs at 5 years with minor deviations from the true values (Table 2).

Application
In our application of MARS to assess the MRD effects in AML dataset, we used time segments of 0-6 months and every 6 months of follow-up thereafter to a maximum of 11 years (22 intervals total).
The meta-analysis survival curves for DFS are shown in the figure 3A, representing the survival probabilities over years. Generally, the DFS was better for subjects who showed MRD negative than those who had positive MRD status. The estimated DFS at 5 years were 65% (95% CrI, 59%-70%) for MRD negative group and 25% (95% CrI, 19%-32%) for MRD positive group. Figure 3B shows the survival curves for DFS synthesizing only reconstructed survival data (type I data) of 16 studies. Survival probabilities of both groups in figure 3B are lower than those in figure 3A, and the 95% credible interval (CrI) was significantly narrower in figure 3B than CrI in figure 3A. Both discrepancies in point estimates and variability reflect the selection bias and extractor bias resulting from using only type 1 data and consequently suggest to conduct complete data analysis for heterogeneity estimation rather than using partial data, especially when the total number of studies included for meta-analysis is not large. The median survival time for MRD negative group is unavailable as the estimated survival probability at the end of study is about 0.56, while the median survival time for MRD positive group is around 1.08 years. Survival probabilities, with 95% CrIs, estimated using only reconstructed survival data.
(C) Estimated RMSTs and corresponding 95% CrIs for MRD-negative patients and MRDpositive patients at various times.
The RMSTs of both groups at 1, 3, 5, and 10 years are shown in Figure 3C. Patients who were MRD-negative had a longer RMST than patients who were MRD-positive. The estimated RMST at 5 years was 3.73 years (95% CrI, 3.51-3.93) for MRD negative group, while that for MRD positive group was 2.06 years (95% CrI, 1.79-2.35).

Discussion
In this work, we present MARS, a novel AD meta-analytic approach to exploit reconstructed survival data as a valuable source of information. Reconstructed survival data, which can be deemed as an AD alternative to IPD, enables theoretical advantages similar to using IPD in both meta-analytic modeling and reporting. To systematically assess previous research studies, we develop a Bayesian multilevel modeling framework consolidating 3 regular types of survival information: 1) the reconstructed survival data elicited from published Kaplan-Meier curves, 2) log HR estimates, and 3) survival rates at specific time points. For the estimation of the HR function and survival probabilities, the piecewise exponential functions are employed with AR(1) correlations to model time-varying effects, the study effects are intrinsically weighted by the length of follow-up study duration, and MCMC allows exact inference for point estimation and uncertainty quantification. We establish the theoretical large sample properties, including consistency and relative efficiency, and provide conditions under which the MARS may coincide in efficiency with the IPD analysis. In finite sample, the MARS demonstrates comparable performance to IPD meta-analysis through simulation studies, under both the constant and time-varying hazard ratio scenarios. As suggested by the motivating MRD example, this research endeavor can congruously promote clinical research and assist regulatory decisions.
Meta-analysis modeling approaches for continuous and binary outcomes are relatively well developed, but less work has been done with time-to-event outcomes (Freeman & Carpenter 2017). Given the high prevalence of clinical studies in which the proportional hazards assumption is violated, the validity of AD meta-analysis is often in question. The considerable popularity in meta-analysis of survival endpoints suggests the critical need for a modeling approach that can relax proportional hazards assumption (Rulli et al. 2018).
The development of the MARS method contributes to the meta-analysis literature in several ways. In particular, MARS relaxes the proportional hazards assumption in individual clinical studies, which is a fundamental restriction in the meta-analysis using HR estimates (Higgins et al. 2020). To derive valid inference for a given population, MARS comprehensively utilizes the most evidence available in the literature, and thereby mitigates the selection bias by including studies that have partial survival information. Sensible summary measurements can be derived conveniently to enrich AD meta-analytic reporting for time-to-event endpoints. For example, survival probabilities can be generated for visual presentation to inform the risk of study group over time. Alternative statistics, such as survival probabilities and RMSTs, can be derived to better inform meta-analysis reporting especially if the proportional hazard assumption is not met. In summary, by addressing various limitations of current approaches, the MARS method effectively extends the classical framework of AD meta-analysis for survival endpoints.
Only 25% of published IPD meta-analyses have had access to all IPD (Nevitt et al. 2017). A generalization of the MARS method can be applied to combine IPD and AD data for survival endpoints (Riley et al. 2007), or process information together with singlearm trials or multiple treatments to incorporate real-world evidence (Begg & Pilote 1991, Zhang et al. 2019). Additionally, one future methodological exploration is to employ alternative strategies for estimating time-varying effects. Jansen (2011), Jansen & Cope (2012) used fractional polynomials (Royston & Altman 1994) to model the time-varying hazard functions in a two-step IPD analysis of time-to-event data. Fractional polynomials are a parsimonious alternative to regular polynomials and provide flexible parameterization for continuous variables. However, the shape of a fractional polynomial may result in unexpected end effects which may jeopardize the overall accuracy. Freeman & Carpenter (2017) applied the Royston-Parmar model with the baseline log-cumulative hazard modeled by restricted cubic splines. However, their formulation failed to guarantee the monotonicity of the cumulative hazard function. Possibly, one solution is to adopt the extended multiresolution hazard model (Hagar et al. 2017) for reconstructed survival data.