Introduction

Sirolimus, also known as rapamycin, is an immunosuppressant administered following solid organ transplantation to prevent allograft rejection. Sirolimus exerts its effect by binding to FK-binding protein-12 (FKBP-12) and inhibiting the activation of the mammalian target of rapamycin (mTOR).1 This inhibition impedes the progression from the G1 to S phase of the cell cycle.2,3 Owing to this distinctive mechanism of action, high efficacy and favorable safety profile, sirolimus is widely used in the clinic.4

Sirolimus is rapidly absorbed following oral administration, with peak plasma concentrations typically observed within 0.5 to 3 h.5,6 The drug undergoes extensive distribution throughout various tissues and organs, with a notable preference for partitioning into erythrocytes in whole blood.7 Sirolimus is a substrate for cytochrome P450 (CYP) 3A4, CYP3A5, and efflux transporter P-glycoprotein and primarily eliminated through bile and feces.3,5,6

The metabolism of sirolimus shows a wide between-subject variability (BSV). Previous studies in renal transplant recipients revealed that the apparent clearance (CL/F) varied from 90 to 416 mL/h/kg,8 resulting in considerable variation in drug concentration after administration of the same dose of sirolimus. Additionally, sirolimus has a narrow therapeutic window of 4–8 ng/mL.9,10 Given these challenges, therapeutic drug monitoring (TDM) of sirolimus whole-blood trough concentrations (C0) has been extensively employed to tailor dosages in clinical settings.

However, dosage adjustments based on the TDM approach have certain limitations, primarily because it can only be implemented after treatment initiation. In contrast, population pharmacokinetic (PK) analysis is a powerful tool that can describe the typical PK parameters of the target population, identify sources of variability, and facilitate the design of individualized drug therapies, both at the onset and throughout the course of treatment.

Currently, there have been several population PK studies of sirolimus in both adult and pediatric patients.11–15 However, to the best of our knowledge, only one population PK study has been conducted for sirolimus in adult liver transplant recipients;10 yet, the limited sample size of that study constrained the breadth and applicability of its findings. Therefore, recognizing this knowledge gap, in our study, we aimed to explore the effects of various covariates on the PKs of sirolimus in adult liver transplant recipients and develop dosing regimen recommendations using a modeling and simulation method.

Methods Patient Characteristics

For this study, we retrospectively collected data from patients who received sirolimus tablets following liver transplantation at the Tongji Hospital (Tongji Medical College, Huazhong University of Science and Technology) between January 2018 and August 2024. All organs were donated voluntarily with written informed consent, and these donations were conducted in accordance with the Declaration of Istanbul. Because of its retrospective nature, the study was approved by Tongji Hospital Ethics Committee (number: TJ-IRB202409087) without the need for written informed consent. Additionally, the study was conducted in accordance with the Declaration of Helsinki (2013).

The following clinical information was collected from the medical records: (1) dosing and sampling information, including sirolimus dosing regimens, date and time of each administration and sampling, drug concentrations measured by TDM, and sirolimus daily dose; (2) demographic information, including postoperative day, age, sex, height, weight, and body mass index; (3) physiological index, including hemoglobin, hematocrit (HCT), alanine aminotransferase, aspartate aminotransferase, total protein, albumin, total bilirubin, alkaline phosphatase, γ-glutamyl transpeptidase, blood urea nitrogen, serum creatinine, and creatinine clearance; (4) concomitant medications, including prednisone acetate, mycophenolic acid (MPA), esomeprazole, ganciclovir, Wuzhi capsules (a traditional Chinese herb, which is usually prescribed to treat drug-induced liver dysfunction), and amlodipine. All authors have confirmed that no identifiable personal data were included in the analysis to ensure the confidentiality and privacy of patient information.

Concentration Measurement

The C0 of sirolimus was measured by the enzyme multiplied immunoassay technique using an Architect i1000 Automatic Chemiluminescence Immunoassay Analyzer (Abbott Diagnostics Inc., IL, USA) and an ARCHITECT Sirolimus Reagent Kit. The intra- and inter-day coefficients of variation were <10%. The assay range was 2–30 ng/mL with a detection limit of 0.3 ng/mL.

Model Development

The data analysis was performed using a nonlinear mixed effects modeling program (NONMEM, v7.5, Icon Inc., PA, USA) with first-order conditional estimation method with η-ε interaction (FOCE-I) estimation method. Perl-speaks-NONMEM (PsN, v5.2.6) were employed to establish and evaluate the population PK model. All statistical analyses and data visualization were conducted using R v4.2.1.

Base Model

A one-compartment model with first-order absorption and elimination was used as the structural model. As only C0 data were collected, the value of absorption rate constant (Ka) was fixed at 0.75 h−1, as reported in previous studies.10 The BSV of the PK parameters was characterized using exponential models (Eq. 1).

(1)

where Pi is the PK parameter estimation of the ith subject, TV(P) is the typical population value for the PK parameter, and represents the variability between subjects, which is assumed to be normally distributed with a mean of 0 and a variance of ω2.

Additive (Eq. 2), proportional (Eq. 3), and combined (proportional and additive) models (Eq. 4) were evaluated to describe the residual unexplained variability (RUV).

(2)

(3)

(4)

where Y and F represent the observed and model-predicted sirolimus concentration, respectively, and represents residual errors and is assumed to follow a normal distribution with a mean of 0 and a variance of σ2.

Base model selection involved the evaluation of multiple factors, including parameter estimation precision, objective function value (OFV), Akaike information criterion,16 Bayesian information criterion,17 and condition number.

Covariate Model

Correlations between demographic information and physiological indices were tested to avoid collinearity effects. Only one covariate model was used for the covariates with a strong correlation. Therefore, the following covariates were investigated: (1) continuous covariates, including age, daily dose, postoperative day, weight, HCT, alanine aminotransferase, and creatinine clearance, and (2) categorical covariates, including sex and concomitant medications.

Continuous covariates were normalized to median population values using a power formula, as shown in Eq. 5.

(5)

where is the continuous covariate value of the ith individual; is the median value of the continuous covariate; is the coefficient of the influence of covariates on the PK parameters.

Categorical covariates were integrated using a proportional model, as shown in Eq. 6.

(6)

where is the categorical covariate with value of 0 or 1.

Potential covariates were evaluated using a stepwise forward inclusion and backward elimination approach. A covariate was considered significant and finally included in the model if its addition decreased the OFV by > 3.84 (P < 0.05; df = 1), and its removal increased the OFV by > 6.63 (P < 0.01; df = 1).

Model Evaluation

The final population PK model was assessed using goodness-of-fit (GOF) plots to evaluate how well the observations aligned with the model predictions. The stability of the final parameter estimates was investigated using a non-parametric bootstrap method involving 1000 resampling iterations. A prediction-corrected visual predictive check (pc-VPC) was performed to assess the predictive performance of the final model. Additionally, normalized prediction distribution error (NPDE) plots were analyzed to determine whether the model followed a normal distribution with a mean of 0 and variance of 1. For a superior visual assessment, four plots of normalized prediction distribution (NPD) were generated instead of NPDE.18

Dosing Regimen Recommendation

Monte Carlo simulations were conducted to estimate the steady-state C0 values for nine commonly used sirolimus dosing regimens under the influence of different covariates based on the final population PK model. The nine dosing regimens evaluated were: (1) 0.5 mg qd (once every day); (2) 1 mg qd; (3) 1.5 mg qd; (4) 2 mg qd; (5) 1 mg qod (once every other day); (6) 2 mg qod; (7) 0.5 mg qd alternating with 1 mg qd; (8) 1 mg qd alternating with 1.5 mg qd; (9) 1 mg qd alternating with 2 mg qd. The recommended sirolimus dosing regimens for patients with varying levels of covariates were determined based on the simulated steady-state C0 range, ensuring that most concentrations were within the therapeutic window of 4–8 ng/mL.9,10

Results Patient Characteristics

A total of 216 whole-blood sirolimus concentrations obtained from 103 adult patients were available for the population PK analysis. These patients were all liver transplant recipients and were prescribed sirolimus to prevent rejection. The median age, weight, and body mass index were 51 years, 68.0 kg and 23.8 kg/m2, respectively. The sirolimus dosing regimens used in the study population included 0.5 mg qd, 1 mg qd, 2 mg qd, 1 mg qod, 2 mg qod, 0.5 mg qd alternating with 1 mg qd, and 1 mg qd alternating with 2 mg qd. The demographic characteristics of the enrolled patients are summarized in Table 1.

Table 1 Baseline Demographic Characteristics of the Population Included in the Population Pharmacokinetic Model

Model Development

A one-compartment model with first-order elimination best described the data. The BSV was estimated for both the CL/F and apparent volume of distribution (V/F) values. A combined error model was used to characterize the RUV. The correlation between CL/F and V/F was calculated.

A statistically significant decrease in OFV (32.9) was observed when HCT was included in the model. No covariates were removed during backward elimination (P > 0.01). Therefore, after the stepwise procedure, the HCT was retained in the final model. The final population PK model was as follows (Eqs.7–9):

(7)

(8)

(9)

Table 2 presents parameter estimates for the final model. The relative standard error of all the parameters was < 30%, indicating that the parameters were estimated with good precision.

Table 2 Population Pharmacokinetic Parameter Estimates and Bootstrap Results of the Final Model

Model Evaluation

The GOF plots showed that the final model provided a good fit to the observed data (Figure 1). The observed concentrations versus the population and individual predictions were evenly distributed on both sides of the identity line (Figure 1A and B). Most of the conditional weighted residuals were densely scattered around 0 and were in the ± 2 range (Figure 1C and D).

Figure 1 Goodness-of-fit plots of the final model. (A) Observed versus individual predicted concentration; (B) observed versus population predicted concentration (PRED); (C) conditional weighted residuals (CWRES) versus PRED; (D) CWRES versus time after dose. The black solid lines are the identity lines, and the red solid lines are the loess smooth lines.

Bootstrap analysis was successful in 97.5% of the 1000 runs. In addition, the parameter estimates of the final models were within the 95% confidence intervals (CI) of the bootstrap results, demonstrating the robustness of the model (Table 2).

The pc-VPC results of the final model indicated a good fit between the predicted and observed concentrations (Figure 2). Notably, the 5th, 50th, and 95th percentiles of the observed concentrations were generally within the corresponding 90% CIs.

Figure 2 Prediction-corrected visual predictive check of the final model. The dots represent observed concentrations; solid lines represent the median (red), 5th, and 95th percentiles (blue) of the observations. The shaded areas represent the 90% confidence intervals for the median (red) and the 5th and 95th percentiles (blue) of the simulated values.

No such trend was observed in the NPD plots (Figure 3A–D). The P values of the Wilcoxon signed-rank, Fisher’s variance, Shapiro–Wilk, and global tests were all > 0.5, suggesting that the final population PK model could describe the observed data well.

Figure 3 Normalized prediction distribution (NPD) plots of the final model. (A) Quantile-quantile plot of the distribution of the NPD versus theoretical normal distribution; (B) histogram of the distribution of the NPD with the density of the standard normal distribution overlaid; (C) NPD versus time after first dose; (D) NPD versus population prediction.

Dosing Regimen Recommendation

Because HCT was included as a covariate in the final population PK model, the steady-state C0 values were simulated for the nine commonly used sirolimus dosing regimens on the background of different HCT levels set to the 10th, 50th, and 90th percentiles of the HCT range observed in the study population, namely 28%, 38%, and 46%, respectively. The simulation results are shown in Figure 4.

Figure 4 Boxplot of the distributions of simulated steady-state sirolimus whole-blood trough concentrations for patients with different hematocrit levels of 28%, 38%, and 46%. The dashed horizontal lines represent the therapeutic range of sirolimus (4–8 ng/mL).

Monte Carlo simulations demonstrated that as the HCT levels increased, the daily dose with the same dosing interval required to achieve the therapeutic window gradually decreased. For patients with a low HCT level of 28%, dosing regimens of 1.5 mg qd as well as 1 mg qd alternating with 1.5 mg qd are recommended. For patients with a normal HCT level of 38% or 46%, the recommended dosing regimens are 1 mg qd, 2 mg qod, or 0.5 mg qd alternating with 1 mg qd.

Discussion

To the best of our knowledge, in this study we developed a population PK model of sirolimus in liver transplant recipients based on the largest sample size to date. A one-compartment model with first-order elimination adequately described the PKs of sirolimus. Furthermore, based on the analysis of our model, several clinical dosing regimen options were recommended for patients with different HCT levels.

The PK parameters estimated in this study using sparse sampling data were consistent with previously reported values. Zhang et al10 reported a CL/F of 6.97 L/h and a V/F of 457.85 L in adult liver recipients, which were comparable to our findings (CL/F: 7.09 L/h, V/F: 496 L).

Additionally, other population PK analyses of sirolimus have been conducted in various transplant populations, including renal transplant14,19–21 and heart transplant recipients.15 Our estimated CL/F was the same as that of heart transplant recipients,15 but slightly lower than the range reported for renal transplant recipients (8.91–14.4 L/h).19,21 This may be owing to the extensive metabolism of sirolimus by CYP3A in the liver. However, because in liver transplant patients, the function of the transplanted liver did not fully recover, which may have resulted in a lower CL/F.

Regarding V/F, our estimate fell within the range observed in the renal transplant recipients (322–727 L),20,21 yet it was significantly lower than the 1350 L reported for the heart transplant recipients.15 This difference may stem from the fact that all samples were C0, which makes the accurate estimation of V/F challenging.

HCT was the only significant covariate identified for CL/F in the final population PK model. An increase in HCT was associated with a decrease in CL/F, which is consistent with the results of previous studies in patients with advanced cancer.22 Given that sirolimus is extensively bound to red blood cells (approximately 95%),23 higher HCT levels may reduce the concentration of unbound free sirolimus and thereby lead to lower CL/F.

The influence of concomitant medication on the PKs of sirolimus was also investigated. In the forward inclusion selection, co-administration of MPA or Wuzhi capsules were included in the model as binary covariates; however, during the backward elimination process, these covariates were removed because they failed to achieve statistical significance and poor estimation of the influence coefficients for the covariates on the PK parameters.

Wuzhi capsules are potent inhibitors of CYP3A4.24 Because sirolimus is primarily metabolized by CYP3A4, the concomitant use of CYP3A4 inhibitors can affect its metabolism. Moreover, there is a potential interaction between MPA and sirolimus. Dösch et al25 reported that the dose-adjusted MPA C0 in the sirolimus/MPA group was lower than that in the cyclosporin A/MPA group. Therefore, further studies with additional information on concomitant medications are required to optimize sirolimus dosing regimens.

In clinical practice, adult liver transplant recipients typically follow a complex sirolimus treatment plan because of its narrow therapeutic window. The simulation results indicated that when the HCT levels were normal, the recommended dosing regimen included an additional dosing interval option, which was every other day. This may be attributed to the long elimination half-life of sirolimus (62 h), which allows C0 to remain within the therapeutic window in these patients. Thus, an optimized dosing regimen can help clinicians avoid frequent dosing adjustments and reduce the occurrence of adverse effects associated with excessive concentrations.

This study had several limitations. First, only trough concentrations were obtained, resulting in a fixed Ka parameter and less reliable estimation of V/F. Therefore, further prospective studies are warranted to investigate the absorption and distribution phase of sirolimus. Second, although this study had the largest sample size of sirolimus concentrations in adult liver transplant recipients, the data remain insufficient to fully evaluate the influence of concomitant medications on the PKs of sirolimus. Third, because all data were collected at a single center, the extrapolation performance of this model requires further investigation.

Conclusion

In conclusion, in this study we developed a population PK model for sirolimus in adult liver transplant recipients based on the largest sample size to date. HCT was identified as a significant factor that influenced the CL/F of sirolimus. We recommend that the dosing regimens be tailored to various HCT levels across the nine frequently prescribed regimens.

Data Sharing Statement

The original data that support the findings of this study can be available from the corresponding authors upon reasonable request.

Acknowledgments

The authors thank Editage (editage.cn) for assistance in English language editing.

Funding

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

Disclosure

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

References

1. Sehgal S. Sirolimus: a new immunosuppressive agent: a historical perspective and immunosuppressive profile. In: Principles of Drug Development in Transplantation and Autoimmunity. Austin, TX: RG Landes; 1996:271–282.

2. Stenton SB, Partovi N, Ensom MH. Sirolimus: the evidence for clinical pharmacokinetic monitoring. Clin Pharmacokinet. 2005;44(8):769–786. doi:10.2165/00003088-200544080-00001

3. MacDonald A, Scarola J, Burke JT, Zimmerman JJ. Clinical pharmacokinetics and therapeutic drug monitoring of sirolimus. Clin Ther. 2000;22(Suppl B):B101–21. doi:10.1016/s0149-2918(00)89027-x

4. Sehgal SN. Sirolimus: its discovery, biological properties, and mechanism of action. Transplant Proc. 2003;35(3 Suppl):7s–14s. doi:10.1016/s0041-1345(03)00211-2

5. Mahalati K, Kahan BD. Clinical pharmacokinetics of sirolimus. Clin Pharmacokinet. 2001;40(8):573–585. doi:10.2165/00003088-200140080-00002

6. Augustine JJ, Bodziak KA, Hricik DE. Use of sirolimus in solid organ transplantation. Drugs. 2007;67(3):369–391. doi:10.2165/00003495-200767030-00004

7. Meier-Kriesche HU, Kaplan B. Toxicity and efficacy of sirolimus: relationship to whole-blood concentrations. Clin Ther. 2000;22(Suppl B):B93–100. doi:10.1016/s0149-2918(00)89026-8

8. Zimmerman JJ, Kahan BD. Pharmacokinetics of sirolimus in stable renal transplant patients after multiple oral dose administration. J Clin Pharmacol. 1997;37(5):405–415. doi:10.1002/j.1552-4604.1997.tb04318.x

9. Masters JC, Shah MM, Feist AA. Drug interaction between sirolimus and ranolazine in a kidney transplant patient. Case Rep Transplant. 2014;2014:548243. doi:10.1155/2014/548243

10. Zhang Y, Zhang X, Zou Y, Sun Y, Li X. Population pharmacokinetics of sirolimus in Chinese adult liver transplant recipients: a retrospective study. Xenobiotica. 2021;51(12):1408–1415. doi:10.1080/00498254.2022.2025628

11. Shen G, Moua KTY, Perkins K, et al. Precision sirolimus dosing in children: the potential for model-informed dosing and novel drug monitoring. Front Pharmacol. 2023;14:1126981. doi:10.3389/fphar.2023.1126981

12. Li S, Zhan M, Wu S, et al. Population pharmacokinetic analysis and dosing optimization of sirolimus in children with tuberous sclerosis complex. J Clin Pharmacol. 2022;62(8):948–959. doi:10.1002/jcph.2033

13. Liu B, Zhang X, Zhao Y, et al. Model-Informed individualized dosage regimen of sirolimus in pediatric patients with intractable lymphatic malformations. Eur J Pharm Sci. 2024;200:106837. doi:10.1016/j.ejps.2024.106837

14. Djebli N, Rousseau A, Hoizey G, et al. Sirolimus population pharmacokinetic/pharmacogenetic analysis and Bayesian modelling in kidney transplant recipients. Clin Pharmacokinet. 2006;45(11):1135–1148. doi:10.2165/00003088-200645110-00007

15. Zahir H, Keogh AM, Akhlaghi F. Apparent clearance of sirolimus in heart transplant recipients: impact of primary diagnosis and serum lipids. Ther Drug Monit. 2006;28(5):614–622. doi:10.1097/01.ftd.0000246765.05248.fa

16. Yamaoka K, Nakagawa T, Uno T. Application of Akaike’s information criterion (AIC) in the evaluation of linear pharmacokinetic equations. J Pharmacokinet Biopharm. 1978;6(2):165–175. doi:10.1007/bf01117450

17. Vrieze SI. Model selection and psychological theory: a discussion of the differences between the Akaike information criterion (AIC) and the Bayesian information criterion (BIC). Psychol Methods. 2012;17(2):228–243. doi:10.1037/a0027127

18. Nguyen TH, Mouksassi MS, Holford N, et al. Model evaluation of continuous data pharmacometric models: metrics and graphics. CPT Pharmacometrics Syst Pharmacol. 2017;6(2):87–109. doi:10.1002/psp4.12161

19. Ferron GM, Mishina EV, Zimmerman JJ, Jusko WJ. Population pharmacokinetics of sirolimus in kidney transplant patients. Clin Pharmacol Ther. 1997;61(4):416–428. doi:10.1016/s0009-9236(97)90192-2

20. Golubović B, Vučićević K, Radivojević D, Kovačević SV, Prostran M, Miljković B. Exploring sirolimus pharmacokinetic variability using data available from the routine clinical care of renal transplant patients - population pharmacokinetic approach. J Med Biochem. 2019;38(3):323–331. doi:10.2478/jomb-2018-0030

21. Shi HQ, Yang J, Zhang LQ, et al. The population pharmacokinetics of sirolimus and CYP3A5*3 polymorphism in Chinese renal transplant patients. Int J Clin Exp Med. 2016;9(3):5854–5866.

22. Wu K, Cohen EE, House LK, et al. Nonlinear population pharmacokinetics of sirolimus in patients with advanced cancer. CPT Pharmacometrics Syst Pharmacol. 2012;1(12):e17. doi:10.1038/psp.2012.18

23. Methaneethorn J, Art-Arsa P, Kosiyaporn R, Leelakanok N. Predictors of sirolimus pharmacokinetic variability identified using a nonlinear mixed effects approach: a systematic review. J Popul Ther Clin Pharmacol. 2022;29(4):e11–e29. doi:10.47750/jptcp.2022.940

24. Li ZR, Li RD, Niu WJ, et al. Population pharmacokinetic modeling combined with machine learning approach improved tacrolimus trough concentration prediction in Chinese adult liver transplant recipients. J Clin Pharmacol. 2023;63(3):314–325. doi:10.1002/jcph.2156

25. Dösch AO, Ehlermann P, Koch A, Remppis A, Katus HA, Dengler TJ. A comparison of measured trough levels and abbreviated AUC estimation by limited sampling strategies for monitoring mycophenolic acid exposure in stable heart transplant patients receiving cyclosporin A-containing and cyclosporin A-free immunosuppressive regimens. Clin Ther. 2006;28(6):893–905. doi:10.1016/j.clinthera.2006.06.015