Searching strategies and inclusion criteria
The Cochrane Preferred Reporting Items for Systematic Reviews and Meta-Analyses guidelines for meta-analyses were used to collect data and report results accordingly. A comprehensive literature search was performed using the public database of PubMed, Embase, and the Cochrane Library for clinical trials of ER+/HER2- metastatic breast cancer by two independent authors between January 1, 2010, and August 23, 2023. The searching key words included “ER+/HER2- breast cancer”, “CDK4/6 inhibitors” and “PI3K/AKT/mTOR inhibitors”. The details of searching strategies are presented in Table S1.
The inclusion criteria were as follows: (1) single-arm, placebo-controlled, or positive drug-controlled clinical trials; (2) involving adults with ER+/HER2- metastatic breast cancer; (3) participants were treated with regimens containing CDK4/6 or PI3K/AKT/mTOR inhibitors; (4) at least one of the efficacy indicators (OS, PFS, or ORR) was used to evaluate efficacy.
Data extraction and risk-of-bias assessment
Microsoft Excel (version 2021) was used for data collection. The collected information included literature characteristics (authors, publication year, and clinical trial registration number), trial design (groups, administration timing, and sample size), participant characteristics (including age, race, ECOG score, proportion of patients with visceral metastasis, proportion of patients with bone metastases only, and proportion of patients who had previously received treatment), efficacy outcomes (OS, PFS, ORR), and safety outcomes (incidence of grade 3–5 adverse events).
If outcome measures are presented graphically in the literature, data extraction will be conducted using Engauge Digitizer (version 12.1). Two independent researchers will extract data from the studies, and any discrepancies will be resolved by a third researcher. During digitization, the data extraction error should not exceed 2%. If the error surpasses this threshold, digitization will be repeated, and the final value is the average of the two extractions.
Model development
Parametric survival models were used to analyze OS and PFS data. The survival model is associated with a hazard function h(t), which can be interpreted as an instantaneous hazard function determined as follows[12, 13]:
$$S(t)=\exp ( - \int\limits_{0}^{t} {h(t)dt} )$$
1
We evaluated four different hazard functions, namely, the constant, Weibull, Gompertz, and log-normal models (Eqs. 2, 3, 4 and 5). In Eqs. (2–4), parameter λ and β represent the hazard rate at time 0 and the regression coefficient of the hazard rate changing over time, respectively. While in Eq. (5), µ and σ are the mean and standard deviation of the lognormal distribution, respectively, which will determine the shape of the function curve.
Constant \(h(t)=\lambda\) (2)
Gompertz \(h(t)=\lambda \cdot \exp (\beta \cdot \ln (t))\) (3)
Weibull \(h(t)=\lambda \cdot \exp (\beta \cdot t)\) (4)
Lognormal \(h(t)={\frac{{(\sigma t\sqrt {2\pi } )}}{{1 - \phi (z)}}^{ - 1}}{e^{( - \frac{1}{2}{z^2})}}\) (5)
The selection of the four hazard functions was based on the minimum value of the objective function (OFV), parameter plausibility, and diagnosis plots.
The inter-trial variation (Eq. 6) and residuals (Eq. 7) were added to the model parameters in the following form:
$${P_i}={P_{pop}} \times {e^{{\eta _i}}}$$
6
$$Ob{s_{j,i}}=Pre{d_{j,i}}+S{E_{j,i}} \times {\varepsilon _{j,i}}$$
7
$$S{E_{j,i}}=\sqrt {\frac{{Ob{s_{j,i}} \times (1 - Ob{s_{j,i}})}}{{{N_{j,i}}}}}$$
8
In Eq. 6, Pi is the model parameter of the ith study, Ppop is the typical value of the model parameter, ηi is the interstudy variation and ηi follows a normal distribution with a mean of 0 and variance of ω2. In Eq. (7), Obsj,i is the observed survival data at jth time in ith study; Predj,i is the predicted survival data at jth time of ith study; εj,i is the residual error of jth time of ith study and εj,i conforms to a normal distribution with a mean of 0 and variance of σ2. εj,i is weighted by the standard error of the survival data at jth time of ith study (SEj,i), which can be acquired using Eq. (8), where Nj,i is the sample size of jth time in ith study.
Covariate models were used to examine potential factors affecting OS and PFS, including age, race, ECOG score, proportion of patients with visceral metastases, proportion of patients with bone metastases only, proportion of patients who had previously received treatment and the types of endocrine therapy. Missing covariate information was imputed using the median value of the entire study population. Missing covariate proportions > 30% were not assessed.
Continuous variables were introduced using Eq. (9) and Eq. (10), and categorical variables were introduced using Eq. (11):
$${P_{pop}}={P_{typical}}+(COV - CO{V_{_{{median}}}}) \times {\theta _{\operatorname{cov} }}$$
9
$${P_{pop}}={P_{typical}} \times {(\frac{{COV}}{{CO{V_{median}}}})^{{\theta _{\operatorname{cov} }}}}$$
10
$${P_{pop}}={P_{typical}}+COV \times {\theta _{\operatorname{cov} }}$$
11
In Eqs. (9–11), Ppop is the model parameter at different covariable levels, and Ptypical is the typical population value of the model parameter when the categorical covariate is equal to 0 or the continuous covariate is equal to COVmedian. COV is the covariate value, COVmedian is the median value of this covariable in the analysis dataset, and θCOV is the correction coefficient of the covariable for the model parameters.
The stepwise method is commonly used to develop a covariate model. Forward inclusion and backward elimination methods were used for covariate screening. The bound of OFV decrease in the forward method was set at 3.84 (P < 0.05), whereas in the backward method, the bound was set at 6.63 (P < 0.01).
Model evaluation
The goodness-of-fit plots were used to describe the fitness of models. A visual predictive check was performed to evaluate the predictive performance of the model by comparing the 95% confidence interval (CI) of the model predictions using the observed values. The stability of the model was assessed using a bootstrap with 1,000 new datasets acquired from the original dataset to obtain the median of the model parameter distribution and its 95% CI, which were compared with the estimated values of the model parameters obtained from the original dataset.
Model prediction
Having established the final model, Bayesian feedback was used to obtain the individual model parameters and their standard errors. A random-effects model in a single-arm meta-analysis was used to summarize the model parameters for different targeted drugs. Finally, Monte Carlo simulations (1,0000 times) were used to obtain the typical OS and PFS time curves and their 95% CIs for two types of targeted drugs[14].
ORR and grade 3–5 adverse events
To determine the characteristics of ORR and safety, ORR and incidence of common grade 3–5 adverse events were pooled using a random-effect model in a single-arm meta-analysis.
Correlations among ORR, PFS and OS
Investigating the correlation between OS and potential short-term efficacy endpoints (ORR and PFS), can accelerate drug development. Weighted linear regression analysis was used to evaluate correlations among the following efficacy indicators: (1) 0.5-year and 1-year PFS, (2) 2-year OS, (3) ORR. The Pearson correlation coefficient (R2) was used to measure the strength of correlation among them. We considered a prespecified R2 values of ≥ 0.8 to indicate a high correlation, 0.6–0.8 to indicate a moderate correlation, and < 0.6 to indicate a low correlation[12, 15].
Software
Model development and simulation were performed using NONMEM 7.5 (ICON Development Solutions, USA). Statistical analysis、meta-analysis and plotting were performed using R4.3.2 (The R Foundation of Statistical Computing). Quality assessment of literatures was performed using Review Manager (Version 5.4. The Cochrane Collaboration, 2020).