Plasma metabolomic profiling reveals factors associated with dose-adjusted trough concentration of tacrolimus in liver transplant recipients

Inter- and intrapatient variability of tacrolimus exposure is a vital prognostic risk factor for the clinical outcome of liver transplantation. New factors or biomarkers characterizing tacrolimus disposition is essential for optimal dose prediction in recipients of liver transplant. The aim of the study was to identify potential plasma metabolites associated with the dose-adjusted trough concentration of tacrolimus in liver transplant recipients by using a global metabolomic approach. A total of 693 plasma samples were collected from 137 liver transplant recipients receiving tacrolimus and regular therapeutic drug monitoring. Untargeted metabolomic analysis was performed by ultraperformance liquid chromatography-quadrupole time-of-flight mass spectrometry. Univariate and multivariate analyses with a mixed linear model were conducted, and the results showed that the dose-adjusted tacrolimus trough concentration was associated with 31 endogenous metabolites, including medium- and long-chain acylcarnitines such as stearoylcarnitine (β = 0.222, p = 0.001), microbiota-derived uremic retention solutes such as indolelactic acid (β = 0.194, p = 0.007), bile acids such as taurohyodeoxycholic acid (β = −0.056, p = 0.002), and steroid hormones such as testosterone (β = 0.099, p = 0.001). A multiple linear mixed model including 11 metabolites and clinical information was established with a suitable predictive performance (correlation coefficient based on fixed effects = 0.64 and correlation coefficient based on fixed and random effects = 0.78). These data demonstrated that microbiota-derived uremic retention solutes, bile acids, steroid hormones, and medium- and long-chain acylcarnitines were the main metabolites associated with the dose-adjusted trough concentration of tacrolimus in liver transplant recipients.

1 Introduction

Tacrolimus is considered the primary immunosuppressant in solid organ transplantation and has been used in liver transplant recipients for approximately two decades (Muduma et al., 2016). Pharmacokinetics-guided dosing is recommended to individualize tacrolimus treatment due to the narrow therapeutic window and wide interindividual pharmacokinetic variability (Brunet et al., 2019). Although the dose-interval area under the curve (AUC) of tacrolimus would be best associated with clinical effects, a multiple sampling strategy is not feasible in clinical practice. A limited sampling strategy with trough concentration (C0) is regularly adopted in most center for therapeutic drug monitoring (TDM) (Brunet et al., 2019; Birdwell et al., 2015). Pharmacogenomics is also essential in pharmacokinetics of tacrolimus, with CYP3A5 genotypes used to guide initial tacrolimus dosing (Brunet et al., 2019). Nevertheless, limitations and challenges for tacrolimus use exist. In liver transplantation, guidelines for CYP3A5 genotype and tacrolimus dosing are recommended only when the donor and recipient genotypes are identical (Birdwell et al., 2015). Furthermore, genotypes provide less information and low predictive value for intrapatient variability (IPV) (Shuker et al., 2015), whereas a recent study identified IPV of tacrolimus exposure as a crucial prognostic risk factor for the clinical outcome in solid organ transplantation, including liver transplantation (Del Bello et al., 2018). Although the most important clinical cause influencing such IPV is medication nonadherence, which is modifiable by clinical interventions, other factors, including drug-drug interactions, food intake, herbal or nutritional constituents and gastrointestinal disorders, also have shown significant roles (Gonzales et al., 2020). Pathophysiological changes such as graft function recovery, inflammation, and altered plasma protein concentrations may affect the tacrolimus disposition of recipients after liver transplantation (Ganesh et al., 2017). Thus, the identification of novel factors or biomarkers characterizing tacrolimus disposition is essential for predicting the optimal dose in liver transplant recipients.

Cytochrome P450 3A isoenzymes (CYP3A), mainly CYP3A4 and CYP3A5, are involved in tacrolimus metabolism. Several human studies have been conducted to identify endogenous biomarkers of CYP3A activity. With the advent of serendipity and hypothesis-driven processes, urinary and plasma metabolites were evaluated as endogenous metrics through the main phenotyping validation criteria (Magliocco et al., 2019). The ability to provide insights into phenotypic metabolite changes with wide coverage and high throughput allows metabolomics to be a promising method for biomarker discovery (Tolstikov et al., 2020). Mass spectrum (MS)-based metabolomic methods involving targeted or untargeted approaches have also been used to detect novel CYP3A biomarkers. Urinary cortisol/6β-hydroxycortisol, cortisone/6β-hydroxycortisone, dehydroepiandrosterone (DHEA)/(7β-OH-DHEA plus 16α-OH-DHEA), and plasma cholesterol/4β-hydroxycholesterol are reliable predictive markers of hepatic CYP3A activity (Shin et al., 2013). Phapale et al. conducted a pharmacometabolomic study on the urine of healthy volunteers to identify a predictive metabolic phenotype of individualized tacrolimus pharmacokinetics and identified four metabolites, namely, cortisol, methylguanosine, acetyl-arginine, and phosphoethanolamine, representing steroid-related, nucleotide/purine, amino acid-related, and glycerophospholipid metabolism, respectively, to predict the pharmacokinetic parameters of tacrolimus (Phapale et al., 2010).

With respect to liver transplant recipients, many clinical factors and drugs could alter the metabolomic profile (Oweira et al., 2018). Thus, the biomarkers for tacrolimus disposition identified in healthy controls should be applied in clinical practice with caution, and a study performed in liver transplant recipients may provide more insight views for metabolomics-guided dosing. Additionally, dose-adjusted tacrolimus trough concentration (C0/D) was commonly used as a surrogate marker for dose response or an index of clearance in most of the studies conducted in patient (Riva et al., 2019). Thus, the present study attempted to identify potential plasma endogenous metabolites associated with tacrolimus C0/D in liver transplant recipients by using a global metabolomic approach.

2 Materials and method

2.1 Subjects

The present retrospective study was conducted at the Hepatobiliary Center, the Affiliated Drum Tower Hospital of Nanjing University Medical School. Adult patients receiving liver transplantation between January 2018 and July 2019 were included. All the recipients received tacrolimus-based regimens and TDM. Subjects with incomplete electronic medical records, and patients undergoing multivisceral transplantation or other additional organ transplantation were excluded. Plasma samples for metabolomics study were available residual biosamples for TDM assay from the biobank of our laboratory. Samples with a C0 of tacrolimus lower than the quantitation limit or with a volume not enough for such analysis were excluded. The enrolled subjects were randomly divided into the model development group and a validation group. The model development group data were used to identify potential metabolites associated with C0/D and to construct the multivariate model, whereas validation group data were used to evaluate the results (Figure 1). The study was carried out in accordance with the principles of the Helsinki Declaration and approved by the Institutional Review Board at Nanjing Drum Tower Hospital (No. 2020-053-01). Informed consent was waived due to the deidentified data provided to researchers and residual biosamples used.

FIGURE 1

FIGURE 1. Flow diagram of the study. The enrolled recipients were divided into the model development group (112 recipients and 570 samples) and a validation group (25 recipients and 123 samples). The model development group data were used to identify potential metabolites detected by UHPLC-QTOF-MS and clinical factors associated with C0/D and to construct the multivariate model using LMM with LASSO, whereas validation group data were used to evaluate the model. UHPLC-QTOF-MS: ultrahigh-performance liquid chromatography-quadrupole time-of-flight mass spectrometry; C0/D: dose-adjusted tacrolimus trough concentration; LMM: linear mixed model; LASSO: least absolute shrinkage and selection operator.

2.2 Treatment and TDM

Recipients received a tacrolimus-based immunosuppressive regimen including mycophenolate mofetil and corticosteroids. The initial dose of tacrolimus (Prograf^®, Astellas, Killorglin, Ireland; or Saifukai^®, Huadong Medicine, Hangzhou, China) was generally 2-3 mg orally twice daily, whereas the posttransplantation C0 target levels were 8–12 ng/ml in the first month, 7–10 ng/ml during the next months, and 5–7 ng/ml after the first year. TDM was performed for patients receiving tacrolimus at a fixed dose for about 3 days, and whole blood samples for C0 detection were taken within 1 hour before the morning dose. The concentration of tacrolimus was determined by the enzyme multiplied immunoassay technique (EMIT) using a drug testing system (Viva-E, Siemens Healthcare Diagnostics Inc, Erlangen, Germany) and EMIT^® 2000 tacrolimus assay kit (Siemens Healthcare Diagnostics Inc, Erlangen, Germany). The procedures were conducted according to the manufacturer’s instructions. The calibration range of the assay was 2.0–30 ng/ml, and samples with a C0 level lower than 2.0 ng/ml were excluded from the study. The performance of the tacrolimus assay in our laboratory was verified in an external quality assessment organized by the National Center for Clinical Laboratories. The internal quality assessment during the study was regularly performed, while the average bias was 6.47%, −1.67% and 1.35%, and the average coefficient of variation was 13.44%, 7.69% and 6.30% for low-, medium- and high-quality controls, respectively. An aliquot of 500 μL of whole blood was also centrifuged at 2,500 g for 10 min at 4°C once the sample was received in the clinical laboratory. The plasma and blood cells were stored at −80°C for further analysis.

2.3 Data collection

Patient data around C0 sampling were collected retrospectively from the hospital information system (HIS) and included daily doses of tacrolimus, weight, concomitant drugs (voriconazole, fluconazole, Wuzhi capsule, omeprazole or esomeprazole, caspofungin, and dihydropyridine calcium channel blockers [DHP CCB]), and blood biochemistry parameters. Wuzhi capsule is a Chinese patent medicine, and its abundant active ingredients include deoxyschizandrin, schisantherin A, schisandrol B and schisandrin (Qin et al., 2014). This is a common drug used to increase the exposure to tacrolimus in most transplantation centers in China.

2.4 Untargeted metabolomics

Plasma samples were thawed on ice, and a 40-μL aliquot was mixed with 160 μL of cold acetonitrile by vortexing for 5 min. After centrifugation at 20,000 g for 15 min at 4°C, the supernatant was collected and centrifuged again. The new supernatant was prepared for injection and separation using ultrahigh-performance liquid chromatography (UHPLC) (ExionLC AD UHPLC, ABSciex, Framingham, MA, United States) with a Kinetex C18 (2.6 µm, 100 mm × 2.1 mm, Phenomenex, CA, United States) column at 40°C. All samples were analyzed using a TripleTOF 5,600 + system (ABSciex, Framingham, MA, United States) in positive and negative electrospray ionization modes. Untargeted metabolomics methods and metabolomics data analysis are shown in the Supplementary Material.

2.5 Statistical analyses

The C0/D of tacrolimus was calculated by dividing the trough by the daily dose (D) normalized to actual body weight, whereas the log-transformed C0/D (log₂C0/D) was investigated as the dependent variable in the subsequent univariate and multivariate analyses. The enrolled subjects were randomly divided into the model development group and a validation group using the function of “createDataPartition” in “caret” package in R software, and unsupervised principal component analysis (PCA) was adopted to evaluate the assay quality and the homogeneity of the metabolic data by using the SIMCA 14.1 (Umetrics, Uppsala, Sweden) software package. Considering the log₂C0/D variation at two levels (between subjects and within subjects due to the repeated measures), metabolites associated with log₂C0/D were identified using a linear mixed model, with each feature as the fixed effect and subject as the random effect on log₂C0/D. p values of the fixed effect were adjusted by the false discovery rate (FDR) (“qvalue” package in R software) for multiple comparisons, and power was estimated by a simulation-based method (“simr” package in R software). Mummichog algorithm prediction (Li et al., 2013) was performed to identify enriched pathways.

To build a multivariate model in the model development group, all the metabolites and other factors including demographic information, biochemistry parameters and concomitant drug use with statistically significance in univariate analysis were included, and variable selection and penalization were performed by minimizing Akaike information criteria (AIC) using the method of least absolute shrinkage and selection operator (LASSO) augmented with 10-fold cross-validation. The predictive accuracy of the final model was assessed numerically in the model development and validation group through calculation of the correlation coefficient (R) by the fixed effects and full model, the mean error (ME), the mean absolute error (MAE), the mean relative error (MRE), and the relative root mean squared error (RMSE).

All statistical analyses were performed using STATA MP 16.0 (StataCorp, TX, United States), RStudio 1.3.959 (R Foundation for Statistical Computing, Vienna, Austria), and R-3.6.3 (R Foundation for Statistical Computing, Vienna, Austria). p values were two-sided, and values <0.05 were considered statistically significant.

3 Results

3.1 Subjects

A total of 137 liver transplant recipients and 693 samples were finally include in the study (Supplementary Figure S1). The median age of the patients was 49 years [interquartile range (IQR): 43–57 years], with a body mass index (BMI) of 23.31 (20.64–25.25) kg/m². Most of the recipients were primarily diagnosed with end-stage liver cirrhosis (56.9%) and malignancies (19.7%) and underwent transplantation from a deceased donor (95.6%). Data regarding the concomitant medications were retrieved from HIS. Voriconazole was observed in 22 samples from eight patients, fluconazole in 65 samples from 17 patients, omeprazole or esomeprazole in 68 samples from 43 patients, and DHP CCB in 50 samples from 24 patients. Additionally, Wuzhi capsule was concomitantly administered in 198 samples from 85 subjects. The differences in the demographic information, clinical biochemistry parameters, and tacrolimus dose and concentration between the two groups were statistically nonsignificant (Table 1).

TABLE 1

TABLE 1. Demographic and medical information of liver allograft recipients^a.

3.2 Associations of demographic information and biochemistry parameters with log₂C0/D in model development group

The linear mixed model exhibited that log₂C0/D increased with postoperative days [β = 0.002, standard error (SE) = 0.0003, p < 0.001] and height (β = 0.026, SE = 0.009, p = 0.005). Male subjects exhibited a higher log₂C0/D than females (β = 0.432, SE = 0.138, p = 0.002) (Figure 2. A). Cholinesterase, adenosine deaminase, total protein, globulin, and hemoglobin were positively associated with log₂C0/D. Tacrolimus pharmacokinetics were negatively associated with renal function and positively associated with parameters related to lipids (Figures 2B–F).

FIGURE 2

FIGURE 2. (A) Different distributions of log₂C0/D of tacrolimus between female and male subjects. (B–F) Median spline plots of Z score of biochemical parameters versus log₂C0/D of tacrolimus. All the plotted parameters were identified to be significantly associated with the log₂C0/D of tacrolimus by a linear mixed model.

3.3 Metabolites identification by untargeted metabolomics analysis

The features were subjected to PCA after data normalization, and the results are illustrated in Supplementary Figures S2A,B. Quality control samples were clustered well in both ion modes, indicating platform stability and metabolomics data reliability. The illustration of grouping homogeneity suggested no evident separation trends between the model development and validation groups in either mode. Pearson’s correlation analysis exhibited good agreement of creatinine, urea, and uric acid levels between the levels obtained in the clinical laboratory and intensity in the metabolomics data (creatinine: R = 0.95, p < 0.001; urea: R = 0.94, p < 0.001; uric acid: R = 0.90, p < 0.001; Supplementary Figure S2C). Obvious differences were also observed for voriconazole, fluconazole, and omeprazole at the metabolomics level between samples with and without drug use, indicating high reliability of the metabolomics data (Supplementary Figure S2D).

Considering the inherent limitation of chemical identification in metabolomics, only features with MS2 information were included in the subsequent analysis. The linear mixed model with the subjects as the random effect was performed in the model development group to evaluate the fixed effect of each feature. The Mummichog algorithm using the m/z and p value indicated that the metabolic pathways were significantly enriched in Vitamin D3 (cholecalciferol) metabolism, squalene and cholesterol biosynthesis, bile acid biosynthesis, C21-steroid hormone biosynthesis and metabolism, aminosugars metabolism and so on (Figure 3). The results from the linear mixed model are illustrated as a volcano plot in Figure 4, exhibiting 633 features with significance, of which 112 were identified by accurate m/z and MS2 spectra matching. The power for each identified significant metabolite was evaluated, and the metabolite with a value of power ≥0.80 was retained. Finally, 31 endogenous metabolites were revealed as potential markers of tacrolimus pharmacokinetics (Table 2).

FIGURE 3

FIGURE 3. Enriched metabolic pathways predicted by the mummichog algorithm with a p value ≤0.05. Colors with gradients of the bars represent relatively high to low values of the enrichment factor.

FIGURE 4

FIGURE 4. Volcano plot of the association of 2,110 features with MS2 information and log₂C0/D evaluated by linear mixed model with the fixed effect of each feature and the subjects as the random effect. Horizontal dotted line represents -log10 (p value) of 1.823 controlling FDR at 0.05. Black and white dots represent unidentified metabolites with q value >0.05 and ≤0.05, respectively. Blue dots below the horizontal dotted line represent identified metabolites with q value >0.05, whereas blue dots above the dotted line represent identified metabolites with q value ≤0.05 but power <0.8. Red and yellow dots represent endogenous metabolites and drugs or drug metabolites with q value ≤0.05 and power ≥0.8, respectively. Labels on dots represent the number of metabolites in Table 2.

TABLE 2

TABLE 2. Identified endogenous metabolites with significance using generalized linear mixed model^a.

3.4 Association between concomitant drug use and log₂C0/D and metabolites

The mixed model analysis exhibited higher log₂C0/D in the samples taken after voriconazole or Wuzhi capsule use than the samples not related to drug use (voriconazole: β = 1.878, SE = 0.203, p < 0.001; Wuzhi capsule: β = 0.485, SE = 0.085, p < 0.001) (Figure 5A), confirming significant pharmacokinetic interactions of voriconazole and Wuzhi capsule with tacrolimus. Conversely, DHP CCB and caspofungin exhibited an inhibitory effect on log₂C0/D (DHP CCB: β = -0.611, SE = 0.130, p < 0.001; caspofungin: β = -0.309, SE = 0.130, p = 0.017). The 31 endogenous metabolites with significance were compared between the samples with and without the four drug uses. The Venn diagram demonstrated two overlapping metabolites, namely, carnitine and testosterone, among all the comparisons and five metabolites, phenylalanine, uric acid, 17alpha-hydroxyprogesterone, monoolein, and taurohyodeoxycholic acid, among the comparisons of Wuzhi capsule, DHP CCB, and caspofungin (Figures 5B,C).

FIGURE 5

FIGURE 5. Associations of drug use with log₂C0/D of tacrolimus and the significant metabolites. (A) Wuzhi capsule and voriconazole use exhibited an increased log₂C0/D, whereas dihydropyridine calcium channel blockers and caspofungin indicated a decreased effect by mixed model analysis. n: number of subjects; N: number of samples; *: p < 0.05; ***: p < 0.001. (B) Comparisons of 31 endogenous metabolites between the samples with and without use of Wuzhi capsule (black), dihydropyridine calcium channel blockers (blue), caspofungin (red), and voriconazole (ochre) by linear mixed model. Each dot represents the estimated coefficient of the metabolites, and the vertical lines with caps represent the 95% confidence interval of estimates. (C) Venn diagram of significant metabolites from the comparisons. Wuzhi: Wuzhi capsule; DHP CCB: dihydropyridine calcium channel blockers; VRZ: voriconazole; CAS: caspofungin.

3.5 Multiple regression model analysis

The multivariate model included 11 endogenous metabolites and 12 clinical characteristics as fixed effects (Table 3). A visual inspection of the residual distribution exhibited no obvious deviation from homoscedasticity or normality (Supplementary Figures S3A,B). Figure 6 illustrates the correlations between log₂C0/D and predictions based on the estimated fixed and fixed plus random effects in the model-developed group (R_fixed = 0.64 and R_total = 0.78) and predictions based on the estimated fixed effects in the validation group (R_fixed = 0.57). The predictive performance of the multiple regression model in the model development and validation groups is presented in Table 4. No significant difference of prediction error (%) between the two groups was found (β = 0.871, SE = 1.894, p = 0.645) (Supplementary Figure S4).

TABLE 3

TABLE 3. Parameter estimates of the multiple regression model using LASSO^a.

FIGURE 6

FIGURE 6. Correlations of observed log₂C0/D with predictions based on the estimated fixed (blue dots and fitted line) and fixed + random effects (red dots and line) in the model development group. (A) Predictions based on the estimated fixed effects (blue dots and fitted line) in the validation group. (B) Black lines are the reference line (y = x).

TABLE 4

TABLE 4. Predictive performance of the multiple regression model in model development group and validation group.

4 Discussion

The present study applied the method of plasma metabolomics in biomarker finding of tacrolimus pharmacokinetics in liver transplantation recipients. Results of metabolomics identified steady state C0/D of tacrolimus was associated with several endogenous metabolites mainly including medium long chain acylcarnitines, microbiota-derived uremic retention solutes, bile acids, and steroid hormones. Furthermore, levels of several metabolites were obviously altered according to inhibitor or inducer use. A multiple linear mixed model was established and factors including 11 metabolites and clinical information presented a suitable predictive performance for log₂C0/D.

Our study suggested that biochemistry parameters and metabolites related to renal function, such as creatinine, uric acid, pseudouridine, and eGFR, were associated with tacrolimus log₂C0/D. Studies have reported creatinine as a factor influencing tacrolimus pharmacokinetics. Fukudo et al. treated serum creatinine as a categorical covariate in adult liver transplant recipients, and tacrolimus clearance was decreased by 19% in subjects with creatinine greater than 1 mg/dl (Fukudo et al., 2003). A linear relation described the serum creatinine effect on clearance, with a coefficient of 0.0801 in adult and pediatric liver transplant recipients (Sam et al., 2006). Pseudouridine and uric acid were also in the list of biomarkers for chronic kidney disease progression using a metabonomic approach (Ye and Mao, 2016; Niewczas et al., 2017).

The observed association between tacrolimus disposition and renal function remains to be explained because tacrolimus undergoes almost no renal elimination (Yu et al., 2018). Studies have demonstrated a twofold higher peak tacrolimus concentration after intraintestinal administration in rats with experimental renal dysfunction than in normal controls (Okabe et al., 2000). The absorption rate of tacrolimus in the intestine might be elevated due to the enhanced nonspecific permeability during renal failure. P-glycoprotein (Pgp) function was suppressed in the liver as opposed to that in the kidney in rats with glycerol-induced acute renal failure, resulting in an increase in tacrolimus bioavailability (Okabe et al., 2002; Huang et al., 2000). Interactions of microbiota-derived uremic retention solutes with drug metabolizing enzymes and transporters were considered to account for that (Prokopienko and Nolin, 2018). Endogenous metabolites associated with tacrolimus log₂C0/D in the present study also comprised microbiota-derived uremic retention solutes such as hippuric acid, 3-indoxylsulfate, indolelactic acid, and substances involved in the metabolism of such solutes, including indole-3-carboxaldehyde, tryptophan, and phenylalanine. These results confirmed the significance of tryptophan metabolism by intestinal microorganisms into indole and its derivatives, some of which are ligands of the aryl hydrocarbon receptor (AhR) regulating the expression of genes encoding drug metabolizing enzymes and transporters involved in tacrolimus pharmacokinetics (Agus et al., 2018). Although indoxylsulfate is also considered a ligand to AhR, an inhibitory effect on CYP3A was reported in studies by Hubbard et al. and Prokopienko et al. (Hubbard et al., 2015; Prokopienko and Nolin, 2018). These findings are partially concurrent with the findings of the present study. Additionally, indole-3-carbinol, a dietary indole from cruciferous vegetables, has been reported to reverse multidrug resistance to vinblastine and doxorubicin mediated by Pgp (Arora and Shukla, 2003). The inhibitory effect of indole-3-carbinol on Pgp may account for the positive association between indole-3-carbinol and tacrolimus C0/D in the present study.

Glycocholic acid and taurohyodeoxycholic acid are the primary and secondary conjugated bile acids, respectively. Zhang et al. proposed taurohyodeoxycholic acid as a tertiary bile acid because it is derived from redox modification by the host on the steroid skeleton of secondary bile acid, which is produced by microbial enzymes (Zhang et al., 2019). CYP3A4 is the predominant enzyme involved in the transformation (Chen et al., 2014). However, the apparent Km for hydroxylation by human liver microsomes was high, resulting in limited levels of taurohyodeoxycholic acid. These levels increased when bile acids accumulated in cholestatic liver disease. Thus, taurohyodeoxycholic acid may be an indicator of CYP3A4 activity and biomarker of tacrolimus pharmacokinetics in liver transplant recipients, as the incidence of biliary complications after liver transplantation was up to 50% (Watt and McCashland, 2008). Collectively, the complex relationships among microbiota-derived uremic retention solutes, metabolites involved in tryptophan metabolism, bile acid metabolism influencing the effects of gut bacteria, and tacrolimus log₂C0/D indicate the significance of the microbiome on the liver and interorgan communication through the gut-liver-kidney axis proposed in other studies (Patel and Shawcross, 2018; Nigam and Bush, 2019).

Univariate analysis also revealed a positive association between several medium- and long chain acylcarnitines and tacrolimus C0/D. Kim et al. identified five ω- or (ω-1)-hydroxylated medium-chain acylcarnitines as novel urinary biomarkers for CYP3A activity in healthy subjects (Kim et al., 2018). The previous study demonstrated that hydroxylated medium-chain acylcarnitine levels in urine correlated positively with tacrolimus clearance. However, a negative correlation was observed in our study. Recently, negative effects of a high-fat diet on intestinal permeability were summarized, and such effects were related to modulation of the expression and distribution of tight junctions, stimulations of a shift to barrier-disrupting hydrophobic bile acids, and induction of oxidative stress and apoptosis in intestinal epithelial cells (Rohr et al., 2020). Tomita et al. evaluated the effects of acylcarnitines on the transporting system in Caco-2 cell monolayers and observed that lauroylcarnitine (acylcarnitine 12:0, Car 12:0) and palmitoylcarnitine (acylcarnitine 16:0, Car 16:0) increased influx and decreased efflux in a manner dependent on their concentration and acyl chain lengths (Tomita et al., 2010). The association of the acylcarnitine profile with tacrolimus C0/D may be attributed to the hydroxylation of medium-chain acylcarnitines by CYP3A4, detergent effects on membranes disrupting membrane barriers, and increasing membrane solubility of long chain acylcarnitines. The exact mechanism of these actions must be further investigated.

Significantly, our results indicated that tacrolimus C0/D increased with days after transplantation in univariate and multivariate regression analysis. Tacrolimus is largely bound to erythrocytes (Woillard et al., 2018), and the increased tacrolimus C0/D with postoperative days in our study may be related to the increase of hematocrit over time after liver transplantation. Additionally, increased bioavailability and the variation of drug-drug interaction could also attribute to that (Riva et al., 2019). In fact, the frequency of comedication with Wuzhi capsule was increased when the graft function was stable in our study.

Multiple regression model established in the study included 11 factors of metabolites, nine factors of clinical biochemistry, and three comedications. There is no significant difference of prediction error (%) between the model development and validation group, indicating a successful application of the model to validation dataset. Inter- and intra-individual variability of tacrolimus pharmacokinetics are relatively large, and most published models performed inadequately in adult liver transplant recipients (Cai et al., 2020). The prediction performance for C0/D (or trough) of the multiple regression model in validation group are relative superior among the published models (Cai et al., 2020). Our findings will also provide an insight into candidate covariates for further model development.

The present study has certain limitations. Firstly, some metabolites exhibiting significance in metabolomics data were not confirmed using authentic standards, which limited the confidence levels of the identification. The identification confidence level is at least two according to the method proposed by Schymanski et al. (Schymanski et al., 2014). The high correlations of typical metabolites between levels obtained in the clinical laboratory and intensity in metabolomics data have confirmed the acceptable identification confidence of these approaches. Secondly, we did not incorporate the genotypic information in recipients and donors, which was suggested to play important roles in tacrolimus pharmacokinetics. The primary objective of the present study is to identify potential plasma metabolomic factors, and metabolomics reflects the comprehensive effects of genetic, environmental, and physiological impacts (Pang et al., 2019). Further studies taking genotypes into account may provide more insights about the interactions between gene and metabolites. Thirdly, the present study is repeated cross-sectional designed, and the limitations inherent to this type of descriptive study do not allow it to establish a causal relationship between metabolic factors and tacrolimus C0/D. Therefore, the ability of the metabolites level prior to the drug administration for predicting of tacrolimus C0/D and guiding tacrolimus dosing should be verified in a cohort study with longitudinal design.

5 Conclusion

Thus, the present study focused on metabolomics in liver transplant recipients to identify potential factors associated with tacrolimus dose-adjusted concentration. Microbiota-derived uremic retention solutes, bile acids, steroid hormones, and medium- and long-chain acylcarnitines were the main metabolites identified with significance. Although a multiple regression model incorporating such metabolites was established, further investigation in a larger, more balanced panel is required with additional information on both donor and recipient genotypes and the integration of innovative strategies such as machine learning. Our findings provide valuable insight into metabolomics on individual dosing and guide research on transporters and metabolizing enzymes.

Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding authors.

Ethics statement

The studies involving human participants were reviewed and approved by Institutional Review Board at Nanjing Drum Tower Hospital. Written informed consent for participation was not required for this study in accordance with the national legislation and the institutional requirements.

Author contributions

All authors participated in the interpretation of the study results and in the drafting, critical revision, and approval of the final version of the manuscript. HZ, WG and YZ were involved in the study design, writing-original draft and writing-review and editing. HZ, MW and ZW were involved in the data curation and data analyses. XX, YD and DL were involved in sample and data collection. HZ, MW and XX were responsible for sample processing and methodology.

Funding

This work was supported by the grants from the Jiangsu Provincial Medical Youth Talent Project (No. QNRC2016013), National Natural Science Foundation of China (No. NSFC81803634 and No. NSFC81302849) and Nanjing Medical Science and Technique Development Foundation (No. YKK17075).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fphar.2022.1045843/full#supplementary-material

References

Agus, A., Planchais, J., and Sokol, H. (2018). Gut microbiota regulation of tryptophan metabolism in health and disease. Cell Host Microbe 23, 716–724. doi:10.1016/j.chom.2018.05.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Arora, A., and Shukla, Y. (2003). Modulation of vinca-alkaloid induced P-glycoprotein expression by indole-3-carbinol. Cancer Lett. 189, 167–173. doi:10.1016/S0304-3835(02)00550-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Birdwell, K. A., Decker, B., Barbarino, J. M., Peterson, J. F., Stein, C. M., Sadee, W., et al. (2015). Clinical Pharmacogenetics Implementation Consortium (CPIC) guidelines for CYP3A5 genotype and tacrolimus dosing. Clin. Pharmacol. Ther. 98, 19–24. doi:10.1002/cpt.113

PubMed Abstract | CrossRef Full Text | Google Scholar

Brunet, M., Van Gelder, T., Åsberg, A., Haufroid, V., Hesselink, D. A., Langman, L., et al. (2019). Therapeutic drug monitoring of tacrolimus-personalized therapy: Second consensus report. Ther. Drug Monit. 41, 261–307. doi:10.1097/FTD.0000000000000640

PubMed Abstract | CrossRef Full Text | Google Scholar

Cai, X., Li, R., Sheng, C., Tao, Y., Zhang, Q., Zhang, X., et al. (2020). Systematic external evaluation of published population pharmacokinetic models for tacrolimus in adult liver transplant recipients. Eur. J. Pharm. Sci. 145, 105237. doi:10.1016/j.ejps.2020.105237

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, J., Zhao, K. N., and Chen, C. (2014). The role of CYP3A4 in the biotransformation of bile acids and therapeutic implication for cholestasis. Ann. Transl. Med. 2, 7. doi:10.3978/j.issn.2305-5839.2013.03.02

PubMed Abstract | CrossRef Full Text | Google Scholar

Del Bello, A., Congy-Jolivet, N., Danjoux, M., Muscari, F., Lavayssière, L., Esposito, L., et al. (2018). High tacrolimus intra-patient variability is associated with graft rejection, and de novo donor-specific antibodies occurrence after liver transplantation. World J. Gastroenterol. 24, 1795–1802. doi:10.3748/wjg.v24.i16.1795

PubMed Abstract | CrossRef Full Text | Google Scholar

Fukudo, M., Yano, I., Fukatsu, S., Saito, H., Uemoto, S., Kiuchi, T., et al. (2003). Forecasting of blood tacrolimus concentrations based on the bayesian method in adult patients receiving living-donor liver transplantation. Clin. Pharmacokinet. 42, 1161–1178. doi:10.2165/00003088-200342130-00006

PubMed Abstract | CrossRef Full Text | Google Scholar

Ganesh, S., Almazroo, O. A., Tevar, A., Humar, A., and Venkataramanan, R. (2017). Drug metabolism, drug interactions, and drug-induced liver injury in living donor liver transplant patients. Clin. Liver Dis. 21, 181–196. doi:10.1016/j.cld.2016.08.011

PubMed Abstract | CrossRef Full Text | Google Scholar

Gonzales, H. M., Chandler, J. L., Mcgillicuddy, J. W., Rohan, V., Nadig, S. N., Dubay, D. A., et al. (2020). A comprehensive review of the impact of tacrolimus intrapatient variability on clinical outcomes in kidney transplantation. Am. J. Transpl. 20, 1969–1983. doi:10.1111/ajt.16002

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, Z. H., Murakami, T., Okochi, A., Yumoto, R., Nagai, J., and Takano, M. (2000). Expression and function of P-glycoprotein in rats with glycerol-induced acute renal failure. Eur. J. Pharmacol. 406, 453–460. doi:10.1016/S0014-2999(00)00699-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Hubbard, T. D., Murray, I. A., and Perdew, G. H. (2015). Indole and tryptophan metabolism: Endogenous and dietary routes to ah receptor activation. Drug Metab. Dispos. 43, 1522–1535. doi:10.1124/dmd.115.064246

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, B., Lee, J., Shin, K. H., Lee, S. H., Yu, K. S., Jang, I. J., et al. (2018). Identification of ω- or (ω-1)-hydroxylated medium-chain acylcarnitines as novel urinary biomarkers for CYP3A activity. Clin. Pharmacol. Ther. 103, 879–887. doi:10.1002/cpt.856

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, S., Park, Y., Duraisingham, S., Strobel, F. H., Khan, N., Soltow, Q. A., et al. (2013). Predicting network activity from high throughput metabolomics. PLoS Comput. Biol. 9, e1003123. doi:10.1371/journal.pcbi.1003123

PubMed Abstract | CrossRef Full Text | Google Scholar

Magliocco, G., Thomas, A., Desmeules, J., and Daali, Y. (2019). Phenotyping of human CYP450 enzymes by endobiotics: Current knowledge and methodological approaches. Clin. Pharmacokinet. 58, 1373–1391. doi:10.1007/s40262-019-00783-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Muduma, G., Saunders, R., Odeyemi, I., and Pollock, R. F. (2016). Systematic review and meta-analysis of tacrolimus versus ciclosporin as primary immunosuppression after liver transplant. PLoS One 11, e0160421. doi:10.1371/journal.pone.0160421

PubMed Abstract | CrossRef Full Text | Google Scholar

Niewczas, M. A., Mathew, A. V., Croall, S., Byun, J., Major, M., Sabisetti, V. S., et al. (2017). Circulating modified metabolites and a risk of ESRD in patients with type 1 diabetes and chronic kidney disease. Diabetes Care 40, 383–390. doi:10.2337/dc16-0173

PubMed Abstract | CrossRef Full Text | Google Scholar

Nigam, S. K., and Bush, K. T. (2019). Uraemic syndrome of chronic kidney disease: Altered remote sensing and signalling. Nat. Rev. Nephrol. 15, 301–316. doi:10.1038/s41581-019-0111-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Okabe, H., Hashimoto, Y., and Inui, K. I. (2000). Pharmacokinetics and bioavailability of tacrolimus in rats with experimental renal dysfunction. J. Pharm. Pharmacol. 52, 1467–1472. doi:10.1211/0022357001777676

PubMed Abstract | CrossRef Full Text | Google Scholar

Okabe, H., Yano, I., Hashimoto, Y., Saito, H., and Inui, K. (2002). Evaluation of increased bioavailability of tacrolimus in rats with experimental renal dysfunction. J. Pharm. Pharmacol. 54, 65–70. doi:10.1211/0022357021771931

PubMed Abstract | CrossRef Full Text | Google Scholar

Oweira, H., Lahdou, I., Opelz, G., Daniel, V., Terness, P., Schmidt, J., et al. (2018). Association of pre- and early post-transplant serum amino acids and metabolites of amino acids and liver transplant outcome. Transpl. Immunol. 46, 42–48. doi:10.1016/j.trim.2017.12.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Pang, H., Jia, W., and Hu, Z. (2019). Emerging applications of metabolomics in clinical Pharmacology. Clin. Pharmacol. Ther. 106, 544–556. doi:10.1002/cpt.1538

PubMed Abstract | CrossRef Full Text | Google Scholar

Patel, C. A. W. V. C., Shawcross, A. S. D. L., SingAnAyAgAm, A., and Shawcross, D. L. (2018). Review article : The gut microbiome as a therapeutic target in the pathogenesis and treatment of chronic liver disease. Aliment. Pharmacol. Ther. 47, 192–202. doi:10.1111/apt.14397

PubMed Abstract | CrossRef Full Text | Google Scholar

Phapale, P. B., Kim, S. D., Lee, H. W., Lim, M., Kale, D. D., Kim, Y. L., et al. (2010). An integrative approach for identifying a metabolic phenotype predictive of individualized pharmacokinetics of tacrolimus. Clin. Pharmacol. Ther. 87, 426–436. doi:10.1038/clpt.2009.296

PubMed Abstract | CrossRef Full Text | Google Scholar

Prokopienko, A. J., and Nolin, T. D. (2018). Microbiota-derived uremic retention solutes: Perpetrators of altered nonrenal drug clearance in kidney disease. Expert Rev. Clin. Pharmacol. 11, 71–82. doi:10.1080/17512433.2018.1378095

PubMed Abstract | CrossRef Full Text | Google Scholar

Qin, X. L., Chen, X., Wang, Y., Xue, X. P., Wang, Y., Li, J. L., et al. (2014). In vivo to in vitro effects of six bioactive lignans of wuzhi tablet (schisandra sphenanthera extract) on the CYP3A/pglycoprotein-mediated absorption and metabolism of tacrolimus. Drug Metab. Dispos. 42, 193–199. doi:10.1124/dmd.113.053892

PubMed Abstract | CrossRef Full Text | Google Scholar

Riva, N., Woillard, J. B., Distefano, M., Moragas, M., Dip, M., Halac, E., et al. (2019). Identification of factors affecting tacrolimus trough levels in Latin American pediatric liver transplant patients. Liver Transpl. 25, 1397–1407. doi:10.1002/lt.25495

PubMed Abstract | CrossRef Full Text | Google Scholar

Rohr, M. W., Narasimhulu, C. A., Rudeski-Rohr, T. A., and Parthasarathy, S. (2020). Negative effects of a high-fat diet on intestinal permeability: A review. Adv. Nutr. 11, 77–91. doi:10.1093/advances/nmz061

PubMed Abstract | CrossRef Full Text | Google Scholar

Sam, W. J., Lai, S. T., Holmes, M. J., Aw, M., Seng, H. Q., Kang, H. L., et al. (2006). Population pharmacokinetics of tacrolimus in whole blood and plasma in Asian liver transplant patients. Clin. Pharmacokinet. 45, 59–75. doi:10.2165/00003088-200645010-00004

PubMed Abstract | CrossRef Full Text | Google Scholar

Schymanski, E. L., Jeon, J., Gulde, R., Fenner, K., Ruff, M., Singer, H. P., et al. (2014). Identifying small molecules via high resolution mass spectrometry: Communicating confidence. Environ. Sci. Technol. 48, 2097–2098. doi:10.1021/es5002105

PubMed Abstract | CrossRef Full Text | Google Scholar

Shin, K. H., Choi, M. H., Lim, K. S., Yu, K. S., Jang, I. J., and Cho, J. Y. (2013). Evaluation of endogenous metabolic markers of hepatic CYP3A activity using metabolic profiling and midazolam clearance. Clin. Pharmacol. Ther. 94, 601–609. doi:10.1038/clpt.2013.128

PubMed Abstract | CrossRef Full Text | Google Scholar

Shuker, N., Van Gelder, T., and Hesselink, D. A. (2015). Intra-patient variability in tacrolimus exposure: Causes, consequences for clinical management. Transpl. Rev. 29, 78–84. doi:10.1016/j.trre.2015.01.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Tolstikov, V., James Moser, A., Sarangarajan, R., Narain, N. R., and Kiebish, M. A. (2020). Current status of metabolomic biomarker discovery: Impact of study design and demographic characteristics. Metabolites 10, 224. doi:10.3390/metabo10060224

PubMed Abstract | CrossRef Full Text | Google Scholar

Tomita, M., Doi, N., and Hayashi, M. (2010). Effects of acylcarnitines on efflux transporting system in Caco-2 cell monolayers. Eur. J. Drug Metab. Pharmacokinet. 35, 1–7. doi:10.1007/s13318-010-0001-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Watt, K. D. S., and McCashland, T. M. (2008). “Cholestasis post liver transplantation,” in Cholestatic liver disease. Editors K. D. Lindor, and J. A. Talwalkar (Totowa, NJ: Humana Press), 171–181. doi:10.1007/978-1-59745-118-5_10

CrossRef Full Text | Google Scholar

Woillard, J. B., Saint-Marcoux, F., Debord, J., and Åsberg, A. (2018). Pharmacokinetic models to assist the prescriber in choosing the best tacrolimus dose. Pharmacol. Res. 130, 316–321. doi:10.1016/j.phrs.2018.02.016

PubMed Abstract | CrossRef Full Text | Google Scholar

Ye, L., and Mao, W. (2016). Metabonomic biomarkers for risk factors of chronic kidney disease. Int. Urol. Nephrol. 48, 547–552. doi:10.1007/s11255-016-1239-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, M., Liu, M., Zhang, W., and Ming, Y. (2018). Pharmacokinetics, pharmacodynamics and pharmacogenetics of tacrolimus in kidney transplantation. Curr. Drug Metab. 19, 513–522. doi:10.2174/1389200219666180129151948

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, J., Gao, L. Z., Chen, Y. J., Zhu, P. P., Yin, S. S., Su, M. M., et al. (2019). Continuum of host-gut microbial co-metabolism: Host cyp3a4/3a7 are responsible for tertiary oxidations of deoxycholate species. Drug Metab. Dispos. 47, 283–294. doi:10.1124/dmd.118.085670

PubMed Abstract | CrossRef Full Text | Google Scholar