Abstract
Aim : The aim of this study
was to develop a population pharmacokinetic (PPK) model in Chinese
children for intravenous busulfan, and to develop a novel busulfan
dosing regimen to support better area under the concentration-time curve
(AUC) targeting.
Methods : We collected busulfan concentration-time samples from
69 children who received
intravenous
busulfan prior to allogeneic
hematopoietic stem cell transplantation (allo-HSCT). A population
pharmacokinetic model for busulfan was developed by nonlinear mixed
effect modelling and was validated by an external dataset (n=14). A
novel busulfan dosing regimen was developed through simulation on 1000
patients. Limited Sampling Strategy (LSS) was established by the
Bayesian forecasting. Absolute Prediction Error (APE), Mean Absolute
Prediction Error (MAPE) and relative Root Mean Squared Error (rRMSE)
were calculated to evaluate predictive accuracy.
Results : A one-compartment model with first-order elimination
best described the data. GSTA1 genotypes, BSA and AST were found
to be significant covariates of Bu clearance, and BSA had remarkable
impact of the volume. Moreover, recommended dose regimens for children
with different GSTA1 genotypes and BSA were developed with a
perfect AUC targeting. A two-point LSS, two hours and four hours after
dosing, behaved well with acceptable prediction precision.
Conclusion : This study developed a PPK model for busulfan that
firstly incorporated GSTA1 genotypes in Chinese pediatric
population. We recommend a BSA-based dosing for personalizing busulfan
therapy in pediatric population. Additionally, an optimal LSS
(C2h and C4h) provides convenience for
therapeutic drug monitoring (TDM) in the future.
KEYWORDS : busulfan, HSCT, individualized therapy, population
pharmacokinetics, pediatric patients
Introduction
Busulfan (Bu) is a common alkylating agent, which can bind to the
guanine of intracellular DNA to damage its structure and function
(1). Bu-based conditioning schemes
are regarded as the cornerstone of allogeneic hematopoietic stem cell
transplantation (allo-HSCT) because of their myeloablative activity (1,
2). However, IV Bu usually leads
to large variability in pharmacokinetics and the treatment window is
narrow. Higher exposure (expressed as area-under-the-curve of 0 - 6 h,
AUC0-6h)
(>1350 μM·min) or low
AUC0-6h (< 900 μM·min) of Bu may lead to a
higher probability of sinusoidal
obstructive syndrome (SOS) or
acute graft-versus-host disease
(aGVHD) (3-5). The exposure and
clinical outcomes of Bu are largely affected by different conditioning
regimen agents, which has been confirmed in children receiving
Bu/cyclophosphamide (Bu/CTX), Bu/fludarabine (Bu/FLU) and Bu/
cyclophosphamide/etoposide (Bu/CTX/VP16) (6-8). Thus, therapeutic drug
monitoring (TDM) for busulfan is the integral component of HSCT.
Personalized Bu dosing via TDM based on the first-dose pharmacokinetics
(PK) could make contribution to lower the occurrence of toxicity and to
improve rate of engraftment (9, 10).
Bu is metabolized by the formation
of a glutathione conjugate in the liver (11, 12). This reaction is
primarily catalyzed by glutathione S-transferase
(GST)
enzymes, such as GSTA1 ,GSTM1 , GSTT1 and GSTP1 (13).GSTA1 is the predominant
GST enzyme involved in Bu metabolism,
and the activity of GSTM1 is close to the half of GSTA1 .
However, the activity of GSTT1 and GSTP1 is relatively
weak (2, 14). Hence, polymorphisms
in the GSTA1 or GSTM1 isoenzymes would more likely affect
Bu metabolism.
Studies about the relationships between GSTA1 polymorphisms and
Bu PK have yielded inconsistent
results.
Although it has been reported thatGSTA1 haplotypes were not significant influence factors of
intravenous (IV) Bu clearance (15), Yin et al. (16) stated adult
patients with the GSTA1 *A/*B genotype had a significantly
higher AUC, higher peak concentration (Cmax) and lower
clearance (CL). Pediatric patients have an increased Bu clearance when
compared with adults (17), which is partially caused by higher
expression or activity of GST enzymes in children (18, 19).
Research on the relationship between GSTA1 polymorphism and PK in
children is also controversial. Ansari et al. (20) reported thatGSTA1 *A/*A was associated with lower drug concentration and more
extensive metabolism. Conversely, Zwaveling et al. (21) asserted
that GSTA1 polymorphisms had no influence on population PK
parameters of IV Bu in children undergoing HSCT. In addition to age, the
activity of GST enzymes is also affected by race. Minor allele
frequencies (MAF) taken from HapMap were reported to be differing
between Caucasian and Asian populations.
Until
now, the correlation between GST polymorphisms and PK in Chinese
children has not been reported.
It has been suggested that, apart from GSTs polymorphism, the
effect of age, body-weight (WT), primary disease, hepatic function,
renal function and drug interactions (fludarabine and phenytoin) may
partly explain inter-individual variability on Bu PK (22-25). Little is
known about correlates of Bu PK in Chinese children, which may be
potential factors to optimize dosage of Bu.
The aim of this study was thus to build a population pharmacokinetic
(PPK) model with data from pediatric patients, so as to present the PK
feature of Bu, to reveal the variability of PK parameters, and to
identify the potential contribution of covariates on the disposition of
Bu. Furthermore, Maximum A Posteriori (MAP) Bayesian forecasting made
use of a PPK model and limited number of samples to forecast
AUC0-6h and to formulate
an
optimal LSS, which is an alternative monitoring strategy.
Methods
2.1 Patients and treatment regimens
From March 2019 to April 2020, we collected 76 patients received
allo-HSCT. This study was approved by the Beijing’s Children Hospital
and all patients/parents provided informed consent.
In the pre-treatment regimen, intravenous infusion of busulfan
(BUSULFEX, Otsuka Pharmaceutical Co. LTD, Zhejiang, China) typically
started eight days prior to transplantation. The infusion frequency was
once every six hours within 3-4 days and each infusion process last for
two hours. The dosage of Bu was based on weight (1 mg/kg for children
less than 9 kg, 1.2 mg/kg for children within a 9-16 kg range, 1.1 mg/kg
for children between 16-23 kg, 0.95 mg/kg for children between 23-34 kg,
and 0.8 mg/kg for children more than 34 kg) in accordance with the
European
Medicines Agency (EMA) (26). Different conditioning regimens were used
dependent on the primary diseases and the types of donor. Briefly, Bu
combined with cyclophosphamide
(CTX) were the basic components of conditioning chemotherapy. For
patients with malignant diseases, regimens containing
cytarabine
(Ara-C) were commonly applied, while
fludarabine (FLU) was
administrated for myeloablative treatment of patients with non-malignant
tumors. Specific dosage and duration of different conditioning regimens
have been exhibited in Supporting Information Table 1 .
Phenytoin (PHT) started 30 minutes before the initiation of Bu therapy
in order to prevent central nervous system toxicity caused by Bu.
Moreover, cyclosporine (CsA) with or without methotrexate (MTX) /
mycophenolate mofetil (MMF) was given as GVHD prophylaxis.
2.2 Bu determination and genotyping
Blood samples for PK analysis and genotyping were withdrawn from central
venous lines, in heparinized glass tubes, pre-infusion, 0.5, 1, 2, 2.5,
4 and 6 hours after the first infusion. Plasma concentrations of Bu were
determined using a high performance liquid chromatography-tandem mass
spectrometry (HPLC-MS/MS) (1). The lower limit of quantitation was 10
ng/mL and the range of quantitation was from 10 to 10000 ng/ml.
Pre-transplant genomic DNA was isolated and extracted from whole blood
or hemocyte prior to the first Bu infusion. GST genotypes of patients,GSTA1 (rs3957356 and rs3957357, which defines haplotype *Aand *B ) and GSTM1 (rs3754446), were detected with the ABI
3730XL DNA Analyzer (ABICo.; BioSune Biotechnology Co., Ltd, Shanghai,
China). Supporting Information Table 2 displays the primer sets
and Tm used for the genotyping assays.
2.3 PPK analysis
BU plasma concentration-time data was analyzed using a Non-Linear
Mixed-Effects (NLME) model implemented in Phoenix 8.0 (Certara USA Inc.,
Princeton, NJ, USA). Typical PPK parameters and their random
inter-individual
variability (IIV) were estimated using a
First-order
Conditional Estimation method with Extended Least Squares method (FOCE
ELS). One- or two-compartment distribution with first-order elimination
were tested as structural model. Exponential errors follows a log-normal
distribution, assumed to describe the IIV in PK parameters by the
equation Pi = P × \(e^{\eta i}\), wherePi is the individual PK parameter of thei th individual, P is the geometric average population
value, and ηi is the subject-specific random
effect value, normally distributed random variable with a mean of 0 and
a variance of ω2 (27). Additive, proportional,
combined additive and proportional models were evaluated to account for
the intra-individual (residual) variability. The selection of base model
was based on changes in -2 log likelihood (-2LL) and on graphic analyses
of the Goodness-of-fit (GOF).
Based on scatter-plots of individual parameters versus covariates and
clinical plausibility, we selected
candidate
covariates for each PK parameters (28). The influence of potential
covariates [sex, age, body weight, body surface area (BSA), alanine
transaminase (ALT), aspartate
aminotransferase (AST), alkaline phosphatase (ALP), total bilirubin
(TBIL), creatinine (CRE), creatinine clearance (CLcr), GSTA1genotypes (*A/*A, *A/*B, *B/*B ), GSTM5 genotypes, primary
disease (malignant disease, non-malignant disease) and drug interactions
(FLU)] on clearance (CL) and volume (V) were further investigated
using a stepwise procedure. Categorical covariates were coded as number,
and continuous covariates were centered on their median value (27).
During
forward
selection, covariates were defined as significance if the -2LL decreased
by at least 3.84 (p ≤ 0.05) following their inclusion in the
model. During backward elimination, one covariate could remain in the
final model if the -2LL increased by at least 6.63 (p ≤ 0.01)
when removed at a time from the full model.
2.4 Model validation
Shrinkage in individual random effects were evaluated in order to assess
whether the final model could be capable to estimate individual PK
parameters by taking advantage of population typical values and sparse
PK data. Shrinkage values of less than 20% indicate that the individual
data are rich enough to compute the PK parameters, whereas larger
shrinkage values generally mean that individual Bayesian estimates are
biased towards the population mean values(29).
Graphical observation of the final model adopted GOF, including
conditional weighted residuals (CWRES) over population predicted
concentrations (PRED) or time after dose (TAD) and the relationship
between observed (OBS) and PRED or individual predicted value (IPRED)
(30). In addition, the final model was evaluated using visual predictive
check (VPC) and bootstrap analysis. VPC was based on the final
pharmacokinetic estimates and then calculated the 95% confidence
interval (CI) for concentrations by simulating 1,000 individuals.
Simulated concentrations were
compared with the 5th, 50th, 95th percentiles of the observed
concentrations. In the bootstrap analysis (31), the 95% CI of the
parameter estimates were derived from 1000 datasets of 69 subjects
generated by random sampling using the Phoenix NLME.
External validation of the model was performed using an external dataset
to evaluate the predictive performance of the final PPK model. The
external dataset consisted of 81 busulfan concentrations from 14
children undergoing allo-HSCT.
Busulfan
plasma concentrations were predicted by fixing the population PK
parameters to the final estimates of the previously established model
and setting maximum evaluations to
0. From this study, measured concentrations of individual patients
assigned to busulfan were compared with calculated concentrations of
these individual patients at the same time with our PPK model using
their BSA, ALT and GSTA1 genotypes of patients. Differences of
< 20% between calculated and measured concentrations were
allowed (32).
2.5 Dosing Simulations
A new dosage scheme using a simulated dataset of 1000 patients was
designed to achieve a targeted AUC of 1125 μM·min. The BSA
of simulated patients were
distributed from 0.2 m2 to 1.6 m2.
Meanwhile, the genotypes of GSTA1 included *A*A and *A*Band the range of AST was coincident with the clinical dataset. Then,
1000 AUCsimulated values and their variabilities were
generated in Phoenix-NLME. Model-based doses were recommended according
to the following Equation 1 :