Population pharmacokinetic/pharmacodynamic joint modeling of ixazomib efficacy and safety using data from the pivotal phase III TOURMALINE-MM1 study in multiple myeloma patients.

Ixazomib is an oral proteasome inhibitor approved in combination with lenalidomide and dexamethasone for the treatment of relapsed/refractory multiple myeloma (MM). Approval in the United States, Europe, and additional countries was based on results from the phase III TOURMALINE-MM1 (C16010) study. Here, joint population pharmacokinetic/pharmacodynamic time-to-event (TTE) and discrete time Markov models were developed to describe key safety (rash and diarrhea events, and platelet counts) and efficacy (myeloma protein [M-protein] and progression-free survival [PFS]) outcomes observed in TOURMALINE-MM1. Models reliably described observed safety and efficacy results; prior immunomodulatory drug therapy and race were significant covariates for diarrhea and rash events, respectively, whereas M-protein dynamics were sufficiently characterized using TTE models of relapse and dropout. Moreover, baseline M-protein was identified as a significant covariate for observed PFS. The developed framework represents an integrated approach to describing safety and efficacy with MM therapy, enabling the simulation of prospective trials and potential alternate dosing regimens.

WHAT IS THE CURRENT KNOWLEDGE ON THE TOPIC?

Although proteosome inhibitors, including ixazomib, are clinically well‐studied and approved for patient treatment, integrated safety/efficacy models remain an open challenge, given the underlying complexities in disease dynamics.

WHAT QUESTION DID THIS STUDY ADDRESS?

How can pharmacometric modeling be used to accurately describe complex safety and efficacy signals as a consequence of ixazomib treatment?

WHAT DOES THIS STUDY ADD TO OUR KNOWLEDGE?

The current study demonstrated the use of pharmacometrics modeling to integrate multiple key safety and efficacy end points with ixazomib pharmacokinetics, enabling a knowledge‐based framework of trial outcomes.

HOW MIGHT THIS CHANGE CLINICAL PHARMACOLOGY AND TRANSLATIONAL SCIENCE?

The models presented in this work can enable the simulation and analysis of prospective clinical trial designs with ixazomib; the integrated efficacy/safety framework may be potentially applied to other therapies in multiple myeloma.

INTRODUCTION

Ixazomib is the first oral proteasome inhibitor developed for hematologic malignancies following the approval of parenteral bortezomib and carfilzomib. Ixazomib is administered as a stable citrate ester, ixazomib citrate, which rapidly hydrolyzes under physiological conditions to the biologically active form of the drug, ixazomib. , , Ixazomib is approved in combination with oral lenalidomide and dexamethasone (LenDex) in the United States, the European Union, and other countries for the treatment of multiple myeloma (MM) in patients who previously received at least one line of therapy, based on results from a global phase III, randomized, double‐blind, placebo‐controlled trial, TOURMALINE‐MM1 (C16010). Additional studies are investigating ixazomib as a single‐agent therapy versus placebo as maintenance therapy in patients with newly diagnosed MM following autologous stem cell transplant (ASCT) in the TOURMALINE‐MM3 study (C16019), and in patients newly diagnosed with MM not treated with ASCT in the TOURMALINE‐MM4 study (C16021). Based on these latter two studies, which met their primary progression‐free survival (PFS) end point and are currently pending overall survival results, ixazomib was approved as maintenance therapy in Japan for the treatment of both post‐transplant and transplant‐ineligible patients with MM.

Model‐informed drug development approaches are increasingly prevalent in oncology drug development and regulatory strategy, where they inform a variety of decisions, including trial design, dose regimen selection, benefit–risk assessment, and product labeling. , Ixazomib clinical development (first approval granted ~ 6 years after the initial first‐in‐human study) was furthered by results from extensive clinical pharmacology studies and model‐informed analyses, which directly contributed to dosing regimen selection, QTc risk assessment, and product labeling. , ,

In this current analysis, we present an integrated platform model that describes the pharmacokinetic/pharmacodynamic (PK/PD) relationships among ixazomib dosing, systemic exposure, and clinically relevant safety (diarrhea, rash, and platelet count) and efficacy (myeloma protein [M‐protein] and PFS) end points in the TOURMALINE‐MM1 study. The study demonstrated improved PFS with ixazomib treatment versus placebo, in combination with LenDex (20.6 vs. 14.7 months; hazard ratio 0.74, p = 0.012). Moreover, ixazomib was associated with minimal additional toxicity compared to the LenDex regimen and resulted in no adverse impacts on patients' quality of life. The modeling framework presented here integrates a previously developed population PK model , with time‐independent Markov models describing diarrhea and rash as safety end points; platelet counts were captured by a semi‐physiological PK/PD model. Finally, a two‐population PK/PD model was developed to describe M‐protein dynamics and integrated with time‐to‐event (TTE) models describing the time to relapse and dropout to reliably capture observed PFS events. Collectively, these models robustly captured PK, safety, and efficacy results from the TOURMALINE‐MM1 study. Notably, our results indicate that increased exposure led to deeper M‐protein responses that may be indicative of improved PFS ; however, this potential benefit was counterbalanced by greater safety risk, supporting the clinically approved dosing regimen. The joint modeling scheme presented here represents a holistic yet modular approach to modeling safety and efficacy outcomes in the context of ixazomib treatment; future clinical trials may evaluate the benefit/risk of potential alternate dosing regimens.

METHODS

Studies

The TOURMALINE‐MM1 (C16010) study was a phase III, randomized, double‐blind, placebo‐controlled study that compared the effects of ixazomib versus placebo, in combination with LenDex, in adult patients with relapsed or refractory MM. A total of 722 patients were randomized 1:1 to receive either 4 mg oral ixazomib or matching placebo on days 1, 8, and 15 of 28‐day cycles. All patients received 25 mg oral lenalidomide (10 mg for patients with creatinine clearance ≤60 or ≤50 ml/min based on local prescribing guidelines) on days 1–21 and 40 mg oral dexamethasone on days 1, 8, 15, and 22. Patients received treatment until progressive disease or unacceptable toxicity, whichever occurred first. The trial was conducted according to International Conference on Harmonization Good Clinical Practice guidelines and appropriate ethical and regulatory requirements, including institutional review board approval.

Data collection

Systemic ixazomib concentrations were quantified from patient plasma samples using a validated liquid chromatography/tandem mass spectrometry assay with a range of 0.5 to 500 ng/ml. Serum M‐protein measurements were obtained prior to treatment in each cycle and during follow‐up, via electrophoresis with immunofixation. Platelets were measured through the first six treatment cycles using volume conductivity scanning methods. Finally, adverse events (AEs) were recorded from the first dose of study drug through 30 days after the last dose. Individual patient ixazomib exposures over time were obtained by numerical integration of the weekly area under the plasma concentration–time curve (AUC) predicted by the previously established population PK model. Software used in the analysis is described in Supplement 3 in Appendix S1.

Pharmacokinetic/pharmacodynamic modeling

The integrated modeling framework (Figure 1a) included the impact of ixazomib on efficacy (M‐protein, dropout, relapse, and PFS) and safety (diarrhea, rash, and platelet dynamics) end points. Ixazomib PK was modeled via Bayesian re‐estimation of individual parameters using a previously published population PK model and used to predict the ixazomib concentrations in the PK/PD analyses. Given the impact of patient dropout on efficacy (PFS) and that relapse was not reliably described by a single model (due in part to complex per protocol definitions), TTE models were used to describe these phenomena and link them to ixazomib PK/PD. Covariate evaluation for all models was performed using a stepwise procedure in Pearl‐speaks‐NONMEM , ; statistical significance was assessed at p = 0.01 and p = 0.001 during univariate forward inclusion and backward elimination steps, respectively. Continuous covariates were evaluated using linear, exponential, or power function models; categorical covariates were evaluated as linear predictors. Evaluation of the final models was based on a spectrum of statistical tests (e.g., likelihood ratio test, precision of parameter estimates, model condition number, and η and ε‐shrinkage) and graphical analyses (e.g., goodness‐of‐fit plots). Visual predictive checks evaluated the predictive performance of the models based on 250 independent replicates. , A nonparametric bootstrap procedure (250 replicates with replacement) was performed to evaluate the stability of the final models and to obtain confidence intervals (CIs) of model parameters.

FIGURE 1

Model‐informed drug development framework. E IXA, the drug effect for ixazomib; E LENDEX, drug effect for LenDex; k in, zero‐order production rate; k out, degradation rate of circulating cells; K‐PD, kinetic‐pharmacodynamic; k prp, maturation rate; MIDD, model‐informed drug development; PFS, progression‐free survival; pGx → Gy, probability of transition from state x to state y; PK, pharmacokinetics; PK/PD, pharmacokinetic/pharmacodynamic; R, M‐protein concentration relative to Y BL; R+, increasing M‐protein concentrations from resistant population relative to Y BL; TTE, time‐to‐event; Y, observed M‐protein concentration; Y BL, baseline M‐protein concentration.

M‐protein model development

M‐protein pharmacokinetic‐pharmacodynamic model

A modified two‐population indirect response model was developed to describe M‐protein dynamics (Figure 1b); baseline sensitive (R) and resistant (R+, with exponentially increasing M‐protein) populations were defined. Based on a comparison of various model types, , , including tumor growth inhibition and multi‐clone models, an adapted type 1 indirect response model (ixazomib concentration‐dependent inhibition of the zero‐order production rate of M‐protein) was used as most suitable to describe TOURMALINE‐MM1 data. Moreover, prior to performing the PK/PD analyses, the lowest M‐protein concentration (t nadir) was obtained via curve fitting to the observed M‐protein data. Specifically, various combinations of exponential functions were used to describe qualitatively different patient M‐protein profiles (e.g., observed M‐protein initially decreased and then increased, was consistently zero, immediately decreased to zero in most/all post‐dose assessments, or immediately increased from the baseline value). In ~ 4.71% of patients (P tnadir,0), t nadir could not be reliably estimated, which did not affect subsequent analysis; the subsequent PK/PD model was fit to the remaining patients. The minimum value of the model predicted concentration plus 0.5 g/L was defined as the nadir and the corresponding time was denoted t nadir; the 0.5 g/L increment was included to reflect that M‐protein data were reported as integer values; therefore, the midpoint between two integers was used to minimize potential bias. Moreover, although International Myeloma Working Group criteria generally require 1 g/L for disease to be considered measurable, a 0.5‐g/L threshold is often used in MM studies. The observed M‐protein concentration was defined as the sum of the sensitive and baseline‐resistant populations, expressed relative to the observed baseline M‐protein concentration, Y BL. Y and Y BL denote the observed and baseline M‐protein concentrations, respectively. R and R + denote the M‐protein concentration relative to Y BL and increasing resistant M‐protein concentrations relative to Y BL, respectively, and k is the individual first order growth rate of the drug‐sensitive subpopulation. Analogously, k represents the growth rate of resistant subpopulation (and therefore informs the individual rate of M‐protein increase post‐nadir); R+ is set to zero prior to treatment initiation, under the assumption that drug resistance emerges during the course of treatment. Y ss, represents the individual patient nadir (steady‐state M‐protein) relative to Y BL in the absence of ixazomib concentration. Maximum fractional inhibition (Imax) and the concentration associated with 50% of maximal inhibition (IC50) define the ixazomib concentration (C) dependent reduction in Y ss leading to R ss, the corresponding steady‐state nadir as a function of ixazomib concentration. The individual patient k , Y ss,, and k were assumed to be log‐normally distributed. The residual error was described via a combined additive and proportional model, where the additive component was fixed to a value of 0.5 g/L (reflecting that the raw data were reported in integer values, noted in the derivation of t nadir).

A covariate evaluation was conducted to assess the significance of all variables listed in Table 1 on random model parameters (k , k , and Y ss) was carried out using standard forward addition (p < 0.01) and backward elimination (p < 0.001) of relationships identified as significant in the initial screening step (p < 0.01). Cytogenic risk category and previous treatment with an immunomodulatory drug were identified as covariates on k R and Y ss, respectively; both were included in a linear form: k and Y ss, are typical parameter estimates of θ kR and θ Yss, for patient i, respectively, and θ coef are the coefficients that describe the covariate effects.

TABLE 1

Summary of demographics and baseline characteristics for the efficacy dataset

	Ixazomib (N = 240)	Placebo (N = 227)	Overall (N = 467)
Age, years
Median [min, max]	66.0 [40.0, 91.0]	66.0 [42.0, 89.0]	66.0 [40.0, 91.0]
Creatinine clearance (ml/min)
Median [min, max]	82.1 [22.9, 231]	79.8 [26.5, 203]	80.8 [22.9, 231]
Hematocrit
Median [min, max]	0.350 [0.200, 0.480]	0.340 [0.210, 0.500]	0.350 [0.200, 0.500]
Hemoglobin (g/L)
Median [min, max]	113 [68.0, 148]	111 [71.0, 167]	113 [68.0, 167]
Baseline M‐protein (g/L)
Median [min, max]	23.0 [10.0, 81.0]	23.0 [10.0, 102]	23.0 [10.0, 102]
Baseline platelet count (109/L)
Median [min, max]	201 [75.0, 367]	194 [44.0, 666]	197 [44.0, 666]
Cytogenetic risk
High risk	57 (23.8%)	40 (17.6%)	97 (20.8%)
Standard	125 (52.1%)	132 (58.1%)	257 (55.0%)
Not available	58 (24.2%)	55 (24.2%)	113 (24.2%)
Eastern Cooperative Oncology Group score
0	119 (49.6%)	104 (45.8%)	223 (47.8%)
1	107 (44.6%)	112 (49.3%)	219 (46.9%)
2	14 (5.8%)	11 (4.8%)	25 (5.4%)
International Staging System stage
I or II	205 (85.4%)	200 (88.1%)	405 (86.7%)
III	35 (14.6%)	27 (11.9%)	62 (13.3%)
Prior immunomodulatory drug therapy
Naive	110 (45.8%)	96 (42.3%)	206 (44.1%)
Exposed	130 (54.2%)	131 (57.7%)	261 (55.9%)
Prior proteasome‐inhibitor therapy
Naive	79 (32.9%)	71 (31.3%)	150 (32.1%)
Exposed	161 (67.1%)	156 (68.7%)	317 (67.9%)
Prior lines of therapy
1	144 (60.0%)	135 (59.5%)	279 (59.7%)
2 or 3	96 (40.0%)	92 (40.5%)	188 (40.3%)
Race
White	214 (89.2%)	195 (85.9%)	409 (87.6%)
Black or African native	5 (2.1%)	2 (0.9%)	7 (1.5%)
Asian	15 (6.2%)	19 (8.4%)	34 (7.3%)
Not reported	3 (1.2%)	8 (3.5%)	11 (2.4%)
Other	3 (1.2%)	3 (1.3%)	6 (1.3%)
Sex
Male	136 (56.7%)	127 (55.9%)	263 (56.3%)
Female	104 (43.3%)	100 (44.1%)	204 (43.7%)
Baseline rash grade
No rash (Grade 0)	238 (99.2%)	224 (98.7%)	462 (98.9%)
Grade 1 rash	1 (0.4%)	3 (1.3%)	4 (0.9%)
Grade 2 rash	1 (0.4%)	0 (0%)	1 (0.2%)
Grade 3 rash	0 (0%)	0 (0%)	0 (0%)
Worst rash grade
No rash (Grade 0)	164 (68.3%)	173 (76.2%)	337 (72.2%)
Grade 1 rash	48 (20.0%)	32 (14.1%)	80 (17.1%)
Grade 2 rash	16 (6.7%)	17 (7.5%)	33 (7.1%)
Grade 3 rash	12 (5.0%)	5 (2.2%)	17 (3.6%)
Baseline diarrhea grade
No diarrhea (Grade 0)	237 (98.8%)	224 (98.7%)	461 (98.7%)
Grade 1 diarrhea	3 (1.2%)	3 (1.3%)	6 (1.3%)
Grade 2 diarrhea	0 (0%)	0 (0%)	0 (0%)
Grade 3 diarrhea	0 (0%)	0 (0%)	0 (0%)
Worst diarrhea grade
No diarrhea (Grade 0)	131 (54.6%)	141 (62.1%)	272 (58.2%)
Grade 1 diarrhea	60 (25.0%)	46 (20.3%)	106 (22.7%)
Grade 2 diarrhea	34 (14.2%)	35 (15.4%)	69 (14.8%)
Grade 3 diarrhea	15 (6.2%)	5 (2.2%)	20 (4.3%)

Dropout time‐to‐event model

Overall dropout was described via an empirical hazard model, where the hazard was described as the sum of the baseline hazard and a hazard dependent on steady‐state M‐protein, modulated by an onset function of time. An additional term described the increased M‐protein dependent hazard when M‐protein increased above a threshold relative to baseline (Equations 8 and 9). This model formulation was preferred to conventional dropout models (e.g., constant hazard, Gompertz, or Weibull models), which were unable to robustly describe observed M‐protein dynamics; a data‐driven empirical formulation was selected. λ, λ 0, λ RSS and λ M‐protein represent the hazard, baseline hazard, hazard due to R ss, and hazard due to M‐protein, respectively. The γ denotes the time‐dependent modulation in dropout (increased dropout over time), which starts increasing at time t 0. R ss and Y BL are defined in Equation 5. The α BL is the M‐protein threshold value relative to baseline that defines whether the additional hazard due to λ M‐protein is invoked. ET50 represents the time at which the hazard attains half its maximum value.

Relapse time‐to‐event model

Preliminary analysis indicated that time to relapse could not be adequately described by a normal or log‐normal distribution; therefore, it was not included directly as a random parameter in the M‐protein PK/PD model. Instead, a log‐logistic accelerated failure time (AFT) TTE model described relapse. In contrast to a proportional hazard model where covariate effects are assumed to have a multiplicative effect on the hazard, AFT models assume that the covariate modifies the disease duration by a constant value. Ixazomib concentration and k R were included as predictors of AFT. λ, λ 0, k , t and C are defined previously (Equation 8). α and β are the acceleration factor and shape parameter of the AFT component of the model. α 0 and α IXA are baseline α and the coefficient describing the ixazomib concentration effect on α. α kR and α kR,0 parameterize the covariate effect of k .

Progression‐free survival model development

A TTE model was used to quantify PFS. Exposure to ixazomib (AUCweekly) was a linear predictor of the log‐hazard. A maximal effect at high drug concentrations (Emax) relationship to modulate the onset of the hazard with respect to time was also incorporated. Covariate relationships among k , Y ss, and M‐protein and the hazard were also included: λ, λ 0, k , Y ss, M‐protein, and ET50 are defined previously (Equations 4, 5, and 8). λIXA represents the hazard associated with exposure to ixazomib. The t represents the time since first exposure to ixazomib. The θ kR, θ , and θ M‐protein represent coefficients that describe the covariate effects of k , Y ss, and M‐protein on the baseline hazard, respectively. Each covariate effect was centered on the median value of that covariate in the analysis dataset (i.e., the denominator in each relationship).

Ixazomib safety models

Discrete time Markov models (DTMM) were used to quantify the grade‐to‐grade transition probabilities of diarrhea and rash (see Supplement 2.2 in Appendix S1). As opposed to logistic regression‐based analyses, DTMM formulations of these end points enabled the models to describe longitudinal characteristics (e.g., times of onset/resolution). Transitions between grades of diarrhea and rash are illustrated in Figure 1c. Specifically, longitudinal data were treated as ordered categorical (i.e., grades 0, 1, 2, and 3); the probability of transitioning to any other grade was dependent solely on the current grade. Transition probabilities were described on the logit scale to constrain the estimated probabilities to values between zero and one under a constant proportional odds assumption. P represents the probability of transition from state i to state j or better for i = 0, 1, 2, and 3 and j = 1, 2, and 3. B + η defines the logit‐transformed individual intercept of the probability of transitioning from state i to at least state 1. The function f represents the effect of a predictor, X, such as ixazomib exposure (AUCweekly) or other covariate(s) on the logit of the transition probability, and β denotes the corresponding coefficient.

Transitions between different AE severity grades (0, 1, 2, and 3) were considered in the analysis, totaling 16 different possible transitions. Likelihood ratio tests were used to assess the appropriateness of the proportional odds assumption and for model selection.

Random effects were modeled as uncorrelated by a diagonal Ω matrix. Alternative random effect variance–covariance structures for Ω, including partial and full block structures, were considered as part of model development.

Discrete time Markov model of diarrhea AE

The DTMM of diarrhea AE was based on constant transition probabilities, with the exception of p 1|0, where a time‐dependent change described the slow transition from grade 0 to grade 1, similar to previous modeling approaches. A linear effect of ixazomib exposure affected only the transitions from grade 0. A time‐dependent model component described the slow increase in grade 1 AE prevalence over time. Random effects were included on all grades. In a model‐based covariate evaluation, prior immunomodulatory drug therapy (IMiD; exposed vs. naive) was identified as a covariate on transitions from grade 3. Details of the DTMM for rash are shown in Supplement 2.2.2 in Appendix S1.

Pharmacokinetic–pharmacodynamic model of platelet count

A semimechanistic model was developed to describe the observed platelet count data, based on an adapted Friberg myelosuppression model. The effect of LenDex was described by a hypothetical kinetic‐pharmacodynamic (K‐PD) model of lenalidomide only, representing the combined background effect (ELENDEX) of lenalidomide and dexamethasone, and was driven by individual observed lenalidomide dosing histories.

The model structure (Figure 1d) is summarized in Equations (15), (16), (17), (18): k In represents the zero‐order production rate, E IXA and E LENDEX are the drug effects of ixazomib and LenDex, respectively. T 1 and T 2 represent the amount in the transit compartments, and k prp represents the maturation rate. Propl and Circ represent proplatelet and circulating cells, respectively, and it is assumed that the upstream proliferating stem cell population is at steady‐state. Degradation of circulating cells is described by the rate constant k out. The system of differential equations was in steady‐state without treatment and normalized to baseline.The individual baseline used to predict the platelet count versus time profile was estimated as a population typical value with a linear covariate effect of the individual observed baseline and individual variability with a variance identical to the residual unexplained variability (Equation 19). is the population typical baseline, W is the coefficient of the linear covariate effect of the individual observed baseline, normalized by the median, . RV and RV are the variances of the proportional and additive residual error, respectively, and is an individual random effect with mean 0 and variance 1, similar to the approach reported by Dansirikul et al. The effect of treatment with ixazomib (EIXA) was described by a linear effect of the plasma concentration (CpIXA) and AUC (scaled by kIxa) predicted by the final ixazomib population PK model. Although the original model formulation by Friberg et al. used only concentration, prior model evaluation indicated that a combined effect yielded a statistical improvement in model fit (objective function value difference = −197). Both effects resulted in a first order depletion from the pool of precursor cells (Equation 20). The effect of LenDex (E LENDEX) similarly resulted in a first order depletion from the precursor pool Prol (Equation 21). The effect of LenDex was driven by a (hypothetical) lenalidomide effect‐site concentration (C eff,LEN). The combined treatment effect of ixazomib and LenDex was assumed to result in an additive effect.

RESULTS

Model‐informed drug development framework

Analysis dataset

Data from 722 patients contributed to this analysis; 361 patients received ixazomib and 359 received placebo in combination with LenDex (2 patients were randomized but not dosed). These 720 patients were defined as the safety population and informed the exposure‐safety analyses. Consistent with study inclusion criteria, eligible patients were to have at least one of the following measures: serum M‐protein ≥1 g/L, urine M protein ≥200 mg/24 h, or elevated involved serum free light chain ≥100 mg/L. In addition, patients who did not provide at least 3 M‐protein observations or did not have a baseline M‐protein concentration of at least 10 g/L were excluded from the analysis (16 each in the ixazomib and placebo groups), leaving a total of 467 patients (240 and 227 in the ixazomib and placebo groups, respectively) in the analysis dataset. These patients were defined as the exposure–efficacy population and informed the exposure‐efficacy analyses. A summary of the key demographics and baseline characteristics of patients included in the exposure‐efficacy analysis is provided in Table 1. Corresponding summaries for the safety population (all treated patients) are presented in Supplement 2.1 Table S4 in Appendix S1.

Ixazomib efficacy

The final M‐protein model and its supporting TTE formulations of dropout and relapse generally described the observed data well. Figure 2 presents a visual predictive check (VPC) of M‐protein dynamics (A) and Kaplan–Meier plots of patient dropout (B), relapse (C), and PFS (D). Final model parameter estimates as well as CIs based on bootstrap analyses are summarized in Table 2. A slight tendency to overpredict dropout was observed, corresponding to a similar slight overprediction of PFS; this was likely due to the fact that both dropout and relapse were impacted by additional factors in addition to M‐protein dynamics. Diagnostic plots of the final M‐protein PK/PD model are presented in Supplement 1.1 Figure S1 in Appendix S1. The VPCs of overall dropout, stratified into nine percentile groups of M‐protein and relative steady‐state M‐protein, respectively, are presented in Supplement 1.2 Figures S2 and S3 in Appendix S1. The VPCs of the time to relapse, stratified into nine percentile groups of M‐protein elimination rate, are presented in Supplement 1.3 Figure S4 in Appendix S1. Finally, VPCs based on the model of PFS, stratified into nine percentile groups of M‐protein, relative steady‐state M‐protein, and M‐protein elimination rate are presented in Supplement 1.4 Figures S5, S6, and S7 in Appendix S1. As illustrated, the tendency of these models to slightly overpredict dropout and PFS emerged at high relative steady‐state M‐protein levels and high M‐protein elimination rates, suggesting that deeper, but not necessarily faster, M‐protein responses may have been indicative of improved survival and a lower likelihood of patient dropout.

FIGURE 2

Visual predictive checks of model components describing the efficacy of ixazomib. (a) Solid and dashed black lines represent the median, and 5th and 95th percentiles of the observations. The corresponding shaded areas denote the 95% confidence interval (CI), obtained from the simulations. Simulated M‐protein concentrations were rounded to integer values to reflect the format of the observed data. (b–d) Shaded areas denote the 95% CI, obtained from the simulations. PFS, progression‐free survival; VPC, visual predictive check.

TABLE 2

Parameter estimates for all models in the exposure‐efficacy analyses

Parameter	Final parameter estimate (%RSE)	Untransformed parameter value [bootstrap 95% CI]
M‐protein PK/PD Model
k R	−6.70 (0.65%)	0.206 [0.187; 0.225]/week
Y SS	−1.95 (2.3%)	14.3 [11.3; 18.8] %
K L	−9.78 (2.5%)	0.00951/week
Ixazomib
Imax	0.758 (2.1%)	75.8 [54.8; 116] %
IC50	1.19 (2.4%)	3.29 [1.48; 9.16] ng/ml
Covariates
k R (BCYABCAT)	0.590 (2.7%)	59.0 [31.3; 96.6] % change in high risk vs. other
Y ss (PIMID)	−0.427 (3.3%)	−42.7 [−59.4; −26.5] % change in naive vs. experienced
Interindividual variability
IIV k R	0.655 (15%)	81.0 [73.0; 91.2] %CV
IIV Y ss	2.39 (15%)	155 [141; 170] %CV
IIV K L	1.34 (28%)	116 [98.4; 131] %CV
Residual variability
Prop. Error	0.218 (3.5%)	21.8 [20.5; 23.7] %CV
Add. Error	0.500 (fix)	0.500 [0.500; 0.500] g/L
Dropout model
λ 0	−10.1 (1.6%)	0.00697 [0.00509; 0.0109]/week
λ RSS	0.0000125 (19%)	0.00210 [0.000978; 0.00308] L/(g.week)
ET50	8.14 (3.6%)	20.4 [8.89; 52.9] week
t 0	6.86 (1.4%)	5.68 [4.57; 9.89] week
λ M‐protein	−6.60 (4.5%)	0.228 [0.140; 0.462]/week
α BL	0.148 (0.21%)	116 [74.8; 123] % of M‐protein baseline
Relapse model
λ 0	−1.23 (24%)	49.0 [26.9; 103]/week
α	10.9 (3.9%)	321 [161; 912] week
β	0.938 (15%)	2.56 [1.89; 3.73]
P tnadir,0 a	−1.67 (6.0%)	4.71 [2.96; 6.59] %
α IXA	0.0509 (62%)	0.0509 [0.00425; 0.123] ml/ng
α kR	−2960 (19%)	−17.6 [−26.0; −11.8] week
α kR,0	0.0509 (39%)	0.0509 [0.0183; 0.0870]
PFS model
λ 0	−9.21 (1.8%)	0.0168 [0.0123; 0.0254]/week
λ IXA	−0.00108 (17%)	0.999 [0.998; 0.999] ml/(μg.h)
ET50	7.82 (5.1%)	14.8 [6.95; 34.7] week
Covariates
λ 0 (k R )	0.798 (15%)	7.90 [5.11; 10.8] % change with 10% increase
λ 0 (Y ss)	0.680 (11%)	6.69 [5.47; 8.61] % change with 10% increase
λ 0 (M‐protein)	0.784 (19%)	7.76 [4.84; 11.1] % change with 10% increase

Note: Final parameter estimates denote NONMEM values based on parameterization in model control streams; where applicable, corresponding untransformed values express parameters in actual units, for completeness.

Abbreviations: %CV, percent coefficient of variation; BCYABCAT, cytogenic risk factor; CI, confidence interval; ET50, time of 50% of maximal effect; IC50, half‐maximal inhibitory concentration; IIV, interindividual variability; Imax, maximum inhibition; IXA, ixazomib; K , rate of M protein increase post‐nadir; kR, rate constant describing changes in M protein; k , intercept for effect of k ; PFS, progression‐free survival; PIMID, prior immunomodulatory drug treatment; PK/PD, pharmacokinetic/pharmacodynamic; P tnadir,0, fraction of patients for whom tnadir could not reliably be estimated; RSE, relative standard error; RSS, relative steady‐state M‐protein with ixazomib; t 0, dropout hazard lag time; YSS, relative steady‐state M‐protein without ixazomib; α, acceleration factor, β; shape parameter, λ 0, baseline hazard.

P tnadir,0 was Φ transformed, where Φ represents the cumulative distribution of a normal distribution.

Ixazomib safety

In general, the DTMM and mechanistic PK/PD models used to describe diarrhea incidence and platelet dynamics, respectively, described the observed data well. Table 3 presents the parameter values for both models in addition to CIs from the corresponding bootstrap analyses. The VPCs illustrate the final DTMM of diarrhea incidence for the first 1.5 years of treatment (Figure 3); corresponding plots for placebo patients and for the full study duration are provided in Supplement 2.2.1 in Appendix S1. The model‐predicted transition probability matrix P for diarrhea incidence is presented in Supplement 2.2.1 Equation 13 in Appendix S1. Increasing ixazomib exposure increased the transition from grade 0 diarrhea to all other grades (Supplement 2.2.1, Figure S10 in Appendix S1). Prior IMiD therapy was a significant covariate resulting in a decreased probability of transition from grade 3 diarrhea. After 1 year of treatment, 1.2% of ixazomib‐treated patients were predicted to have grade 3 diarrhea, as opposed to 0.1% of IMiD‐naïve patients. A similar effect was predicted for patients in the placebo arm (0.5% vs. 0.1% of IMiD‐treated and IMiD‐naïve patients predicted to have grade 3 diarrhea, respectively).

TABLE 3

Parameter estimates of safety models

Parameter	Estimate (% RSE)	Untransformed parameter value [bootstrap 95% CI]
Diarrhea model
Explicit time effects
P 1|0T a	2.25 (70%)	9.52 [3.58; 5.21 × 106] (log OR)
K 1|0T a	−10.4 (20%)	0.00532 [7.53 × 10−9; 0.0261]/week
P 0|1I b	0.933 (26%)	0.715 [0.434; 1.01] ml/(h.μg) (log OR)
Ixazomib effect		0.933 [0.391; 1.37] (log OR)
SLP0 c	0.000715 (22%)	0.715 ml/(h.μg) (log OR)
Covariates
B 3x (PIMID)	−2.87 (31%)	−2.87 naive vs experienced (log OR)
Parameters describing the transition probability from state i to state j, p i|j
B 01 d	−5.28 (2.5%)	‐‐
B 02 e	0.211 (34%)	‐‐
B 03	0.657 (16%)	‐‐
B 11	0.677 (57%)	‐‐
B 12	1.99 (3.5%)	‐‐
B 13	0.541 (40%)	‐‐
B 21	1.94 (11%)	‐‐
B 22	−1.92 (16%)	‐‐
B 23	2.21 (3.4%)	‐‐
B 31	3.63 (19%)	‐‐
B32	−2.15 (35%)	‐‐
B 33	−1.78 (44%)	‐‐
Interindividual variability
η 0	1.86 (14%)	136 [119; 158] %CV
η 1	3.72 (20%)	193 [162; 232] %CV
η 2	3.24 (23%)	180 [140; 224] %CV
η 3	2.74 (51%)	165 [74.9; 261] %CV
Platelet PK/PD model
k prp	−3.20 (1.0%)	6.88 [6.34; 7.37]/week
k In	−3.21 (3.3%)	6.78 [5.13; 7.66]/week
BL	203 (0.95%)	203 [198; 206] e9/L
W	0.837 (4.5%)	83.7 [76.9; 91.1] %
Ixazomib effect
slpIXA	0.000859 (25%)	0.144 [0.0707; 0.221]/week/(ng/ml)
k IXA	0.0000818 (33%)	0.0137 [0.00718; 0.0278]/week
LenDex effect
k LEN	0.0483 (23%)	8.11 [5.66; 23.7]/week
slpLEN	0.0134 (5.8%)	2.24 [1.99; 2.62]/week/(len conc.)
Residual variability
Prop. Error	−0.160 (5.8%)	16.1 [13.9; 18.0] %cv
Add. Error	28.1 (6.4%)	28.1 [23.8; 31.5] e9/L

Note: Final Parameter Estimates denote NONMEM values based on parameterization in model control streams; where applicable, corresponding untransformed values express parameters in actual units, for completeness.

Abbreviations: %CV, percent coefficient of variation; BL, typical platelet baseline; CI, confidence interval; kin, platelet production rate; k IXA, ixazomib accumulation rate; k LEN, K‐PD rate, kprp; platelet maturation rate; len conc, lenalidomide concentration predicted by K‐PD model; log OR, log‐odd‐ratio; PIMID, prior immunomodulatory drug treatment; PK/PD, pharmacokinetic/pharmacodynamic; RSE, relative standard error; SLP0, linear effect of ixazomib concentration on diarrhea; slpIXA, ixazomib effect on platelet count; slpLEN, lenalidomide effect on platelet count; W, coefficient of linear covariate effect of the individual observed baseline; η, individual random effect.

P 1|0 T, K 1|0 T: parameters describing slow transition from grade 0 to 1 over time.

Rapid diarrhea onset in week 1.

Ixazomib treatment effect on grade 0 diarrhea.

Bx1 logit probability of transitioning from grade x to grade 1.

Bx2, Bx3: log difference of logit transition probabilities. Probability matrix is included in the Supplement 2.2.

FIGURE 3

Visual predictive check of diarrhea for patients treated with ixazomib. Note: Visual predictive check based on 250 replicated trials. Shaded areas indicate 95% confidence interval (CI).

The observed platelet count dynamics after ixazomib treatment were adequately described by a semimechanistic model, including a K‐PD model component (Figure 1d, Figure S13), in which both ixazomib and LenDex effects resulted in a decrease in the platelet precursor pool, P; over the course of a treatment cycle, the platelet counts declined and recovered approximately back to baseline. The full platelet PK/PD model provided parameter estimates with good precision (relative standard error <33%; Table 3). Moreover, cumulative ixazomib exposure (absolute and cumulative concentrations) resulted in a moderately declining trend in platelet count (on the timescale of years; Figure 4; VPCs for placebo patients are included in Supplement 2.3 in Appendix S1), consistent with the known inhibitory effects of proteosome inhibitors such as ixazomib and bortezomib.

FIGURE 4

Visual predictive check of platelet count in patients treated with ixazomib. Note: Black solid (dashed) lines = observed median (10th and 90th percentiles); Red (blue) shaded areas = 95% confidence interval (CI) of model prediction of median (10th and 90th percentiles). Horizontal dashed reference lines at recovery level 75 × 109/L and thrombocytopenic level 50 × 109/L.

DISCUSSION

These analyses report an integrated model framework to describe the efficacy (M‐protein and PFS) and safety (diarrhea, rash, and platelet count) data observed in patients randomized to ixazomib or placebo, each in combination with LenDex, in the TOURMALINE‐MM1 study (C16010). In all cases, a previously published ixazomib population PK model was utilized to predict individual ixazomib exposures, which were subsequently used as effects on various PD end points.

A modified two‐population model was developed to describe the sum of drug‐sensitive and drug‐resistant M‐protein populations, expressed relative to the baseline M‐protein concentration. The drug‐sensitive M‐protein population was modeled using an indirect response model to quantify the kinetics of M‐protein in response to ixazomib exposure, whereas the increasing resistant M‐protein population was described by an exponentially increasing function following the onset of relapse. The reported PK/PD model in combination with TTE models of relapse and dropout described the observed M‐protein data well. From an efficacy perspective, this analysis suggested that increased ixazomib exposures may lead to deeper, but not necessarily faster, M‐protein response dynamics. These observations are consistent with previously published data that indicated that duration of response was prolonged in patients who responded more slowly and/or deeply. , Furthermore, faster M‐protein reduction (i.e., higher k values) were associated with shorter times to relapse in the TTE model. Taken together, these results suggest that response depth, rather than speed, is more indicative of favorable clinical outcomes. Furthermore, M‐protein dynamics predicted control‐arm PFS well; informing these predictions with ixazomib exposure in the treatment arm further captured PFS trends, which were complicated by composite protocol‐defined relapse criteria beyond serum M‐protein alone. Additionally, the per protocol definition of PFS included death events due to any cause; time to progression (TTP) that right‐censors patients who died before progression may have been better predicted by M‐protein. However, TTP was not among the key primary or secondary end points of the TOURMALINE‐MM1 study and was not considered in the model‐based analysis. Nonetheless, despite the complex definitions of PFS precluding direct prediction via M‐protein, modeling accurately captured the observed data; notably, PFS was slightly overpredicted at high M‐protein baseline levels and elimination rates, consistent with the premise that deeper, if slower, M‐protein dynamics are clinically favorable.

Intuitively, these findings suggest that higher ixazomib doses or an increased dosing frequency may result in deeper M‐protein reductions. Indeed, a phase II trial comparing the current 4‐mg once weekly (q.w.) ixazomib dose to 5.5 mg q.w. (maximum tolerated dose), in combination with dexamethasone, found that while overall response rates trended higher in the 5.5‐mg q.w. arm; the higher dose also resulted in higher proportions of dose modifications and AEs. Similarly, a phase I/II study of twice weekly ixazomib 3 mg dosing demonstrated clinical activity but with greater safety risks compared to q.w. dosing. In totality, these analyses support the approved clinical dose; in addition, the aforementioned correlative M‐protein features are likely relevant for further evaluation.

The analysis of ixazomib safety was described by the relationship between ixazomib exposure (weekly AUC) and time course of AEs (platelet count, diarrhea, and rash) in the TOURMALINE‐MM1 population through the development of PK/PD models. A semimechanistic model was successfully applied to describe platelet count; exposures of ixazomib and LenDex were implemented as linear effects on the platelet precursor elimination rate. Whereas the model adequately described longitudinal platelet dynamics, the effects of ixazomib and LenDex were highly correlated; as such, the developed model may not be applicable to predicting ixazomib monotherapy without additional validation. Moreover, although a linear effect of ixazomib exposure on decreasing platelet counts was sufficient to describe the observed data, a hyperbolic form may be more appropriate with extended longitudinal data. Finally, we note that proteasome inhibition is generally considered to inhibit the distal step of platelet budding from the megakaryocyte. However, mechanistic resolution and attribution of the effects of ixazomib and LenDex on different steps in the platelet biogenesis cascade was not within the scope of this analysis, as all contributing data were for the triplet combination therapy. Separating these potential effects would likely require single agent data across multiple dosing regimens. Nevertheless, model performance was adequate, as verified by visual predictive checks, supporting application for benefit/risk simulations with the models developed for efficacy outcomes.

The DTMM model used to describe transitions between AE grades for rash and diarrhea. These models incorporated linear effects of ixazomib exposure on specific transition probabilities. In the final DTMM models, ixazomib exposure was found to increase only the probability of transition from the “no AE” state to AE presence. Ixazomib patients with prior IMiD exposure were predicted to have higher prevalence of diarrhea, whereas Asian ixazomib‐treated patients were predicted to have higher prevalence of rash (see Supplement 2.2.2 in Appendix S1). The former finding in both arms of the TOURMALINE‐MM1 study is consistent with previous literature evidence highlighting diarrhea as a common and potentially worsening side effect of several IMiDs in oncology, including thalidomide and lenalidomide. ,

These comprehensive pharmacometrics analyses demonstrate how modeling and simulation can be applied to integrate available clinical data representing multiple efficacy and safety variables into a quantitative knowledge‐ and hypothesis‐generating platform in oncology drug development. The analysis confirmed and quantified the characteristics of the M‐protein response to ixazomib treatment, as well as the relationships between drug exposure and various safety end points. By integrating M‐protein biomarker modeling with DTMM models of safety end points, the joint framework described here represents a unified simulation tool to more fully understand ixazomib clinical outcomes.

AUTHOR CONTRIBUTIONS

J.K.S., P.M.D., M.J.H., K.V., R.L., and N.G. wrote the manuscript. J.K.S., M.J.H., K.V., and N.G. designed the research. J.K.S., P.M.D., M.J.H., and K.V. performed the research. P.M.D. analyzed the data.

CONFLICT OF INTEREST

J.K.S., M.J.H., R.L., and N.G. are employees of TDCA. K.V. is a former employee of TDCA. P.M.D. is an employee of Certara USA, Inc., Princeton, New Jersey, USA, who conducted analytical work for TDCA.

Appendix S1

Click here for additional data file.