Effect of ABCB1 SNP polymorphisms on the plasma concentrations and clinical outcomes of rivaroxaban in Chinese NVAF patients: a population pharmacokinetic-based study

Background: This study utilized a population pharmacokinetic (PPK) approach to assess the influence of ABCB1 genetic polymorphisms on the plasma concentrations and clinical outcomes of rivaroxaban.

Methods: The PPK model for rivaroxaban was developed using the nonlinear mixed-effects modelling approach and Monte Carlo simulations were employed to derive peak concentration (C_max) and trough concentration (C_trough). ABCB1 genetic variants were analyzed for their impact on the plasma concentrations and clinical outcomes.

Results: Analysis of 287 rivaroxaban plasma concentrations from 228 non-valvular atrial fibrillation (NVAF) patients revealed significant associations between AST (aspartate aminotransferase)/ALT (alanine aminotransferase) ratios and the apparent clearance (CL/F), the apparent volume of distribution (V/F). ABCB1 1236C>T TT and ABCB1 c.2482-2236C>T CC genotypes exhibited higher dose-adjusted C_max (C_max/D) compared to other relevant genotypes. Additionally, the ABCB1 3435C>T TT genotype showed lower dose-adjusted C_trough (C_trough/D) compared to CC or CT genotypes. For clinical outcomes, the ABCB1 c.2482–2236C>T CC genotype had a higher bleeding risk compared to TT (RR = 1.99, 95% CI 1.08–3.69) or CT genotypes (RR = 1.42, 95% CI 1.04–1.92), and ABCB1 3435C>T TT genotype showed a higher thromboembolic risk compared to CC genotype (RR = 3.48, 95% CI 1.02–11.85).

Conclusion: The PPK model incorporated CL/F and V/F with the covariate AST/ALT. Model-based simulations revealed that ABCB1 1236C>T, ABCB1 c.2482–2236C>T, and ABCB1 3435C>T genotypes had significant impacts on the plasma concentrations of rivaroxaban. Specifically, ABCB1 c.2482-2236C>T and ABCB1 3435C>T genotypes were associated with bleeding events and thromboembolic events, respectively.

Introduction

Atrial fibrillation (AF) stands as the most prevalent sustained cardiac arrhythmia observed in clinical settings, impacting millions globally (Chugh et al., 2014). AF poses substantial risks to cardiovascular health and overall wellness. The irregular ventricular response associated with AF can lead to hemodynamic instability, significantly elevating the likelihood of thromboembolic events, particularly ischemic stroke (Chugh et al., 2014). Anticoagulant therapy, commonly using oral anticoagulants such as warfarin or non-vitamin K antagonist oral anticoagulants (NOACs), like rivaroxaban, dabigatran, or apixaban, plays a crucial role in mitigating thromboembolic risks among AF patients (January et al., 2019).

In recent years, rivaroxaban, a direct inhibitor of factor Xa, has gained considerable attention for its effectiveness and safety profile in managing thromboembolic disorders among patients with non-valvular atrial fibrillation (NVAF) (Patel et al., 2011). Despite the generally predictable nature of rivaroxaban’s PK, individual variability in drug exposure has been observed, leading to challenges in optimizing dosage regimens. This variability can impact both the efficacy and safety of rivaroxaban therapy, potentially increasing the risk of bleeding or reducing its antithrombotic effects (Testa et al., 2016). Specifically, genetic polymorphisms in drug-metabolizing enzymes and drug transporters may contribute to this variability.

Cytochrome P450 (CYP) isoforms 3A4, 3A5, and 2J2 primarily metabolize rivaroxaban (Daly and King, 2003). Investigations have indicated that the activity of CYP3A4 affects the peak concentration (C_max) and trough concentration (C_trough) of rivaroxaban (Sychev et al., 2018), suggesting that genetic variations in CYP genes may affect the pharmacokinetic (PK) of rivaroxaban. A recent study highlighted that individuals with mutant genotypes of CYP3A4 (rs2242480, rs2246709, and rs3735451) and CYP3A5 (rs776746) displayed higher C_trough of rivaroxaban compared to those with wild-type genotypes, and the minor allele (C) carriers on rs3735451 and the minor allele (A) carriers on rs2246709 were correlated with the minor bleeding events (Li X. Y. et al., 2023). Furthermore, rivaroxaban is transported by P-glycoprotein, encoded by the ATP Binding Cassette Subfamily B Member 1 (ABCB1) gene (Gnoth et al., 2011). Variations in this gene may influence P-glycoprotein activity and expression levels, thereby impacting rivaroxaban’s absorption, distribution, and elimination processes, finally leading to variations in plasma concentrations and clinical outcomes. Therefore, we aim to investigate the impact of ABCB1 gene polymorphisms on the plasma concentrations and clinical outcomes of rivaroxaban in this study.

However, accurately determining the C_max and C_trough can be challenging due to imprecise experimental operations in clinical experimental procedures. Meanwhile, Population pharmacokinetics (PPK) has emerged as a robust analytical tool for elucidating variability in drug exposure across diverse individuals within a population (Mould and Upton, 2012). PPK models can identify sources of variability and facilitate accurate prediction of PK parameters, even with limited sampling instances or variable dosing histories (Mould and Upton, 2012). In light of these considerations, the present study aims to investigate the impact of ABCB1 gene polymorphisms on the plasma concentrations and clinical outcomes of rivaroxaban in Chinese NVAF patients using a PPK-based approach. Specifically, we sought to determine whether genetic variations in ABCB1 influence the dose-adjusted C_max (C_max/D) and C_trough (C_trough/D) of rivaroxaban, and whether these variations are associated with bleeding or thromboembolic events. By addressing these questions, our study provides valuable insights into the role of ABCB1 genetic polymorphisms in rivaroxaban therapy, paving the way for personalized dosing strategies in NVAF patients.

Methods

Subjects and therapeutic interventions

This was prospective research conducted at the Fujian Provincial Hospital, which was approved by the institutional ethics committee (No. k2022-09-014). Inclusion criteria encompassed adult individuals diagnosed with NVAF receiving rivaroxaban treatment, while notable exclusion criteria included patients with poorly controlled hypertension (systolic blood pressure ≥ 160 mmHg), a history of cerebral hemorrhage or arterial dissection, or abnormal coagulation function. Specifically, poorly controlled hypertension is associated with an increased risk of bleeding events, such as intracerebral hemorrhage, which could complicate the interpretation of rivaroxaban’s safety profile (Lane and Lip, 2012). Additionally, poorly controlled hypertension may lead to secondary complications, such as renal impairment, which could indirectly influence rivaroxaban’s PK by reducing its elimination (Kvasnicka et al., 2017).

Patients received rivaroxaban doses of 5, 7.5, 10, 15, or 20 mg once daily. Blood samples were collected from residual blood after routine biochemical tests during patient clinical care, with inclusion limited to samples with validated sampling information. Plasma samples were centrifuged at 2,500×g at 4°C for 10 min and stored at −70°C. Additionally, the clinical characteristics and outcomes of the patients were extracted from medical records and the estimated glomerular filtration rate (eGFR) was calculated using the Chronic Kidney Disease Epidemiology Collaboration equation (Kong et al., 2013). Bleeding events comprised minor bleeding (epistaxis, gingival bleeding, skin ecchymosis), urinary and gastrointestinal bleeding (hematemesis, melena, stool routine positive, hematuria), and intracranial hemorrhage. Thrombotic events included arterial thrombosis (acute myocardial infarction, ischemic stroke), and venous thrombosis (lower extremity deep vein thrombosis, pulmonary embolism).

Measurement of the collected rivaroxaban plasma concentrations

The rivaroxaban plasma concentrations were analyzed using a highly sensitive and specific high-performance liquid chromatography-tandem mass spectrometry (HPLC-MS/MS) method. This method was adapted from a previously published protocol (Derogis et al., 2017), with modifications to chromatographic conditions and the linear calibration range to suit the specific requirements of our study. Validation of the analytical method was performed in accordance with the guidelines outlined by the U.S. Food and Drug Administration (FDA) and the European Medicines Agency (EMA) for bioanalytical method validation. The calibration curve was linear over 3 ∼ 1,600 ng/mL, with an LLOQ of 3 ng/mL and r² > 0.99. Both intra-day and inter-day precision and accuracy were assessed. Precision (coefficient of variation, CV) was < 5%, and accuracy ranged from 93% to 100%. Recovery was consistent and reproducible at low, medium, and high-quality control levels. No matrix interference was observed, and dilution integrity (1:2, 1:4) met the acceptance criteria. Stability studies under various conditions (short-term, freeze-thaw, long-term) were satisfactory. No carry-over was detected. Standard samples of rivaroxaban and apixaban (used as an internal standard) were obtained from Sigma-Aldrich (St. Louis, MO, United States). Methanol and formic acid of chromatography grade were sourced from Merck Company. These standards were used throughout the validation process to ensure consistency and reproducibility.

Assessment of gene polymorphisms

Gene polymorphisms at four gene loci, specifically ABCB1 3435C>T, ABCB1 1236C>T, ABCB1 2677G>T/A, and ABCB1 c.2482-2236C>T, were investigated in this research. Massarray SNP typing technology, conducted by BGI Tech Solutions (Beijing Liuhe) Co., Limited, was employed for the detection of gene polymorphisms.

Pharmacokinetic model development

The compartment model and nonlinear mixed-effects modelling strategy were employed in the development of the PPK model. The Phoenix^® NLME™ 7.0 software from Certara (St. Louis, MO) was utilized, employing the FOCE-ELS (first-order conditional estimation and extended least squares) method. It is important to note that a log-normal distribution was assumed for the inter-individual variability (IIV) of the PK parameters. Specifically, the random effect term (η_i) associated with each individual’s parameter was assumed to follow a normal distribution with a mean of zero and a variance of ω² (Equation 1).

$P_i=P\times e^{{\eta }_i}(1)$

where P is the typical value of a PK parameter and P_i represents the ith patient’s individual PK parameter.

Furthermore, the residual error is described by the proportional, additive, or combined error model (Equations 2–4).

$C_i=C\times \left(1+{\epsilon }_1\right)(2)$

$C_i=C+{\epsilon }_2(3)$

$C_i=C\times \left(1+{\epsilon }_1\right)+{\epsilon }_2(4)$

where Ci and C represent the individual plasma concentrations and predictions of plasma concentrations, respectively. ε₁ and ε₂ represent the proportional errors and additive errors of predictions for drug concentrations, respectively, which are normally distributed with a mean of zero and a variance of σ².

Interindividual variability in rivaroxaban PK was estimated, and covariate analyses were performed to identify factors that could explain or reduce this variability. The dataset included both categorical covariates and continuous covariates. The impact of categorical covariates on each parameter was assessed using a scale model (Equation 5) (Liu et al., 2022), whereas continuous covariates were evaluated using an exponential function model (Equation 6) (Li et al., 2016).

$P_i=\left\{\begin{array}{c}P\times {\theta }_1\left(\text{Category\hspace{0.17em}}1\right)\\ P\times {\theta }_2\left(\text{Category\hspace{0.17em}}2\right)\\ P\times {\theta }_3\left(\text{Category\hspace{0.17em}}3\right)\end{array}\right.(5)$

$P_i=P\times {\left(\frac{Cov}{{Cov}_{median}}\right)}^{\theta }\times e^{{\eta }_i}(6)$

where Cov and Cov_median represent the individual and median values of a covariate, respectively, while θ represents the estimated value of the continuous covariate effect. θ₁, θ₂ and θ₃ represent the estimated values of the different categorical covariate effects, respectively.

The correlation between covariates and the parameters of the base model was used to construct the initial population model. The final population model was obtained through the forward inclusion-backward elimination approach. A covariate was deemed significant if its inclusion resulted in a reduction of more than 6.635 in the objective function value (OFV) (P < 0.01), and its exclusion led to an increase of more than 10.828 in the OFV (P < 0.001).

Model validation

The final model’s reliability was evaluated using visual assessment methods, specifically employing goodness-of-fit (GOF) plots (Li Z. et al., 2023). These primary GOF plots included four scatter plots: conditional weighted residual errors (CWRES) versus time after last dose (TAD), CWRES versus population-predicted concentration (PRED), observations versus PRED, and observations versus individual predicted concentration (IPRED). To assess the final model’s robustness, a bootstrap resampling technique was employed. A total of 2000 datasets were generated by randomly selecting different patient combinations, and parameters were re-estimated using the final population model (Li et al., 2016). Median parameter values and their 95% confidence intervals (CI), derived from the 2.5th and 97.5th percentiles of the 2000 bootstrap-estimated parameters, were compared with the estimates of the final model. Additionally, visual predictive checks (VPCs) were performed using 1,000 simulations to assess the predictive performance of the final model (Li et al., 2015). In these VPCs, the observed drug concentrations were plotted against the 5th, 50th, and 95th percentiles of the simulated concentrations, which represent the 90% prediction intervals. Furthermore, the 90% CIs for the 5th, 50th, and 95th percentiles were calculated across the 1,000 simulated datasets to evaluate the robustness of the model predictions. Lastly, for each patient, observed plasma concentrations were plotted against time, and the corresponding predicted concentration-time profiles were generated. The individual fits were visually evaluated by assessing the alignment between observed data points and predicted concentration-time profiles, with the plots used to identify potential outliers or systematic deviations in model predictions at the individual level.

Monte Carlo simulations

Monte Carlo simulations were performed using the final PPK model to estimate individual C_max and C_trough for each patient. A total of 1,000 simulations were conducted per patient to account for IIV and residual error. The primary objectives of the Monte Carlo simulation were: (1) to evaluate the effects of ABCB1 genetic polymorphisms on C_max/D and C_trough/D; and (2) to explore potential correlations between rivaroxaban plasma concentrations and clinical outcomes, including bleeding and thromboembolic events.

Statistical analysis

The Hardy-Weinberg equilibrium test was conducted using the χ² test. To evaluate the impact of ABCB1 SNPs on rivaroxaban plasma concentrations, continuous variables such as C_max/D and C_trough/D were compared among genotypes using the Kruskal-Wallis test. Pairwise comparisons were performed using the Mann-Whitney U test with Bonferroni correction to adjust for multiple testing. For clinical outcomes, categorical variables such as bleeding events and thromboembolic events were analyzed using the χ² test or Fisher’s Exact test. Relative risks (RR) with 95% CI were calculated to assess the risk of clinical events associated with specific genotypes. All analyses were performed using Stata 17.0 (StataCorp, College Station, TX, United States). Data are presented as median and range or interquartile range, or mean ± standard deviation (SD). Statistical significance was defined as P < 0.05. However, for the univariate analysis, Bonferroni correction-adjusted P-values were considered statistically significant to address multiple testing.

Results

Participants

A cohort of 228 patients was successfully recruited for participation in this study, in accordance with the predetermined eligibility criteria. Throughout the study duration, a total of 287 plasma concentration data points were collected and documented. The median number of blood samples per participant was 1 (range: 1 ∼ 3; mean ± SD: 1.26 ± 0.54). Approximately 78% of participants contributed one sample, 18% contributed two samples, and 4% contributed three samples. Detailed information regarding the demographic and clinical characteristics of the enrolled patients has been summarized in Table 1. Genotyping for four SNPs was performed, and no observed frequencies deviated significantly from Hardy-Weinberg equilibrium, as shown in Table 2. Analyzing and summarizing confounding factors potentially impacting drug metabolism or clinical outcomes were conducted and presented in Table 3, where no significant differences were observed among different ABCB1 genotypes.

Table 1

Table 1. Demographic and clinical characteristics of participants.

Table 2

Table 2. Variation in genotypes and allele frequencies.

Table 3

Table 3. Confounding factors among different ABCB1 genotypes.

Population pharmacokinetic model

In this study, we evaluated six different combinations of compartment models and error models to determine the most appropriate structure for our PPK model. Specifically, we tested one-compartment and two-compartment models with additive, proportional, and combined error models. The optimal model was selected based on the OFV and diagnostic plots, with the combination yielding the lowest OFV and most consistent diagnostic plots being chosen as the final model. The OFV values for each combination are presented in Table 4. Through this comparative analysis, a one-compartment model with a proportional error model was ultimately identified as the best-fitting model for describing the PK of rivaroxaban in the studied population. Additionally, the absorption rate constant (k_a) was fixed at 0.617 h^-1 based on a prior PPK study involving Japanese patients (Kaneko et al., 2013). This decision was made due to the limited sampling during the absorption phase in our dataset, which could lead to high uncertainty in k_a estimation. Fixing k_a ensured stability in parameter estimation while maintaining consistency with prior knowledge.

Table 4

Table 4. Objective function values for different combinations of compartment models and error models.

A comprehensive set of covariates was evaluated for their potential influence on the apparent clearance (CL/F) and apparent volume of distribution (V/F) of rivaroxaban. The covariates analyzed included age, sex, body weight (BW), body mass index (BMI), albumin (ALB), bilirubin (BIL), alanine aminotransferase (ALT), aspartate aminotransferase (AST), the AST/ALT, ratio, serum creatinine (S_CR), and estimated glomerular filtration rate (eGFR). Notably, the k_a was fixed at 0.617 h⁻¹ based on findings from prior studies, thereby precluding the need for covariate analysis related to k_a.

Covariate analysis revealed significant associations between several variables and the PK parameters of rivaroxaban. Specifically, age, AST, ALT, eGFR, and the AST/ALT ratio demonstrated notable effects on CL/F or V/F, as evidenced by diagnostic covariate plots. Consequently, these covariates were incorporated into the initial PPK model. A stepwise approach, involving forward inclusion and backward elimination, was utilized to refine the model. This iterative process identified AST/ALT ratios as the sole covariate significantly influencing CL/F and V/F of rivaroxaban. The refinement resulted in a decrease in the OFV value from 2,661 to 2,506. Consequently, the final PPK model, which incorporates the AST/ALT covariate, is represented mathematically by Equations 7, 8.

$CL/F\left(L/h\right)=5.64\times {\left(\frac{AST/ALT}{1.188}\right)}^{-0.074}(7)$

$V/F\left(L\right)=41.7\times {\left(\frac{AST/ALT}{1.188}\right)}^{0.213}(8)$

In Equations 7, 8, the typical values for CL/F and V/F are 5.64 L/h and 41.7 L, respectively, while the median value of AST/ALT is 1.188. The coefficients (f_CL/F-AST/ALT) and (f_V/F-AST/ALT) are −0.074 and 0.213, respectively, indicating the relationship between AST/ALT and these two PK parameters. The equations indicate a clear decrease in CL/F and an increase in V/F as AST/ALT ratios rise. Initially, a diagonal OMEGA matrix was used to estimate the variances of the PK parameters (CL/F and V/F). To assess potential correlations between parameters, a full OMEGA variance-covariance matrix was subsequently evaluated. However, the correlation between CL/F and V/F was found to be negligible, and the inclusion of covariances did not significantly improve the model fit. Thus, a diagonal OMEGA matrix was retained in the final model to maintain parsimony.

Table 5 displays the estimates, relative standard errors (RSE), IIV, and shrinkage of the final PPK model, alongside the residual errors. It is noteworthy that the estimates’ precision is acceptable, given the RSE values ranging from 5.49% to 15.87%. The inter-individual variability (ETA) shrinkage values for CL/F and V/F were 15.8% and 22.4%, respectively, while the residual error (EPSILON) shrinkage value for the proportional error was 18.7%. These values are below the commonly accepted threshold of 30%, indicating that the model is well-supported by the data and the parameter estimates are robust.

Table 5

Table 5. Parameter estimates and bootstrap analysis of rivaroxaban population pharmacokinetic model in the NVAF population.

Goodness-of-fit and model evaluation

The final PPK model underwent a favorable GOF assessment. When comparing observed concentrations (DV) with both PRED and IPRED, data points were generally distributed symmetrically around the Y = X axis, indicating good agreement between model predictions and actual observations (Figures 1A,B). However, it is worth noting that for observed concentrations above approximately 200 μg/L, the model tends to underpredict the observed values. This systematic bias may be attributed to the limited number of high-concentration data points in the dataset, which could affect the model’s ability to accurately predict extreme values. Moreover, CWRES distribution versus PRED or TAD exhibited symmetry, with the majority of values falling within the −2 to +2 range (Figures 1C,D). The reliability and robustness of the final PPK model were further confirmed through bootstrap analysis. Median parameters derived from bootstrap samples closely resembled original parameter estimates, with the 95% CI encompassing these estimates as well (Table 5). Additionally, the VPCs demonstrated that the 90% CI of predictions from the final PPK model aligned well with observed data (Figure 2). Lastly, the individual fits demonstrated strong agreement between observed and predicted concentrations across most patients, indicating that the model adequately captured the PK profiles at the individual level. Representative individual fit plots are shown in Supplementary Figure S1.

Figure 1

Figure 1. The goodness-of-fit plots of the final population pharmacokinetic model. (A) DV versus PRED; (B) DV versus IPRED; (C) CWRES versus PRED; (D) CWRES versus TAD. DV, dependent value (observations); PRED, population predicted concentration; IPRED, individual predicted concentration; CWRES, conditional weighted residuals; TAD, time after last dose. The red lines represent theoretical reference lines (y = −2, 0, +2), while the blue lines indicate the observed trends in the data.

Figure 2

Figure 2. The visual predictive check plots of the final population pharmacokinetic model. The blue dots are the actual original observations; the red line is the 50% quantile from the population predicted concentrations, and the upper and lower blue dashed lines are the 95% and 5% quantiles from the population predicted concentrations, respectively; the shade areas are the 90% confidence intervals of the relevant quantiles. TAD, time after last dose.

Correlation analyses between rivaroxaban plasma concentrations and ABCB1 SNP polymorphisms

Using Monte Carlo simulations based on the final PPK model, we simulated C_max at 2 ∼ 4 h post-dose and C_trough at 24 h post-dose for each patient. Subsequently, we assessed the impact of ABCB1 SNPs on the C_max/D and C_trough/D of rivaroxaban. The results are summarized in Table 6, where ABCB1 1236C>T (P = 0.0001) and ABCB1 c.24822236C>T demonstrated (P = 0.0001) a significant impact on C_max/D. Furthermore, ABCB1 3435C>T showed a significant effect on C_trough/D (P = 0.0001). Further analyses in pairwise comparisons among ABCB1 SNPs revealed that ABCB1 1236C>T TT carriers exhibited higher C_max/D values compared to CC (P < 0.000001) or CT (P < 0.000001) carriers. Similarly, ABCB1 c.2482–2236C>T CC carriers displayed higher C_max/D values than TT (P = 0.000042) or CT (P = 0.000103) carriers. Regarding C_trough/D, ABCB1 3435C>T TT carriers showed lower C_trough/D compared to CC (P < 0.000001) or CT (P = 0.001744) carriers. Detailed results are provided in Table 7.

Table 6

Table 6. Impacts of ABCB1 SNPs on the plasma concentrations of rivaroxaban.

Table 7

Table 7. Pair comparisons of plasma concentrations among ABCB1 SNPs.

Impacts of ABCB1 SNP polymorphisms on the clinical outcomes

The effects of ABCB1 SNPs on clinical outcomes in NVAF patients receiving rivaroxaban are summarized in Table 8. We observed that ABCB1 c.2482-2236C>T and ABCB1 3435C>T were significantly associated with impacts on bleeding events and thromboembolic events, respectively. However, no difference was observed for the other ABCB1 genotypes. Further pairwise comparisons of RR were conducted for the significant SNPs and are depicted in Figure 3. The analysis revealed that ABCB1 c.24822236C>T CC carriers were at a higher risk of bleeding events compared to TT carriers (RR = 1.99, 95% CI 1.08–3.69) or CT carriers (RR = 1.42, 95% CI 1.04–1.92). Additionally, ABCB1 3435C>T TT carriers showed a higher risk of thromboembolic events compared to CC carriers (RR = 3.48, 95% CI 1.02–11.85).

Table 8

Table 8. Impacts of ABCB1 SNPs on the bleeding events and thromboembolic events.

Figure 3

Figure 3. Pair comparisons of risk ratio in significant SNPs for bleeding events and thromboembolic events. RR, relative risk; CI, confidence interval; SNP, Single nucleotide polymorphism.

Discussion

In this study, successful characterization of the PK of rivaroxaban was achieved using a one-compartment model. The final PPK model yielded typical values of 5.64 L/h for CL/F and 41.7 L for V/F. Through an extensive covariate search, AST/ALT ratios were identified as a significant factor affecting CL/F and V/F in the final PPK model. The relationship between CL/F, V/F and AST/ALT was best described by the exponential model, with a scaling factor (f_CL/F-AST/ALT) of −0.074 and (f_V/F-AST/ALT) of 0.213, indicating the decrease in CL/F and increase in V/F as AST/ALT ratios rise. Collectively, evidence from GOF analysis, bootstrap evaluation, VPCs and individual fits supports the conclusion that the final PPK model possesses the adequate predictive ability to accurately describe the PK behavior of rivaroxaban in these patients.

Prior rivaroxaban PPK studies typically excluded patients with body weight (BW) below 45 kg and severe renal or hepatic impairment (Liu et al., 2022; Li Z. et al., 2023). In contrast, our prospective study encompassed patients across a spectrum of BW (33.5 kg–99 kg), eGFR (13.3 mL/min–130.4 mL/min), and AST/ALT ratio (0.37–6.5). In our study, The AST/ALT ratio was identified as a significant covariate influencing both CL/F and V/F in the final PPK model. Beyond its statistical significance, the AST/ALT ratio has a well-established physiological relevance as an indicator of liver function (Cohen and Kaplan, 1979; Gitlin, 1982; Williams and Hoofnagle, 1988). Elevated AST/ALT ratios are associated with progressive liver impairment (Giannini et al., 1999), which can directly affect the activity of metabolic enzymes (e.g., CYP3A4, CYP3A5, and CYP2J2) and transporters (e.g., P-glycoprotein) involved in rivaroxaban disposition. Our findings suggest that higher AST/ALT ratios are associated with reduced CL/F and increased V/F, consistent with the expected effects of liver dysfunction on drug metabolism and protein binding. Incorporating the AST/ALT ratio into the PPK model thus provides valuable insights into the individualization of rivaroxaban dosing, particularly in patients with varying degrees of liver function.

Previous studies have emphasized the significant impact of genetic polymorphisms on both plasma concentration and clinical outcomes, such as bleeding or thromboembolic events, in NVAF patients undergoing rivaroxaban treatment. Xiang et al. suggested an association between SUSD3 rs76292544 and 12-month bleeding events, along with potential links of genetic variants from 52 SNPs in 36 genes (including GOT2 rs14221 and MMP13 rs640198) with the peak anti-FXa level related to bleeding events (Xiang et al., 2023). Zhang et al. found that rivaroxaban C_max/D and ABCB1 2677G variation correlated with a higher incidence of bleeding events (Zhang et al., 2022). Additionally, Lähteenmäki et al. emphasized associations of ABCB1 3435C>T SNP and 1236T-2677T-3435T (rs1128503-rs2032582-rs1045642) haplotype with reduced thromboembolic risks in rivaroxaban users, while ABCB1 c.2482-2236C>T SNP was associated with lower bleeding risk in apixaban users (Lähteenmäki et al., 2021).

Our investigation extends these findings, highlighting the significant impacts of ABCB1 1236C>T, ABCB1 c.2482–2236C>T, and ABCB1 3435C>T on the plasma concentrations of rivaroxaban. Specifically, ABCB1 1236C>T TT and ABCB1 c.2482-2236C>T CC genotypes were associated with higher C_max/D compared to other genotypes. Furthermore, ABCB1 3435C>T significantly influenced C_trough/D, with ABCB1 3435C>T TT demonstrating lower C_trough/D compared to CC or CT genotypes. Importantly, different ABCB1 genotypes showed significant differences in clinical outcomes. For instance, ABCB1 c.2482-2236C>T CC genotype had a higher bleeding risk than CT or TT genotypes, whereas ABCB1 3435C>T TT genotype showed a higher thromboembolic risk than CC genotypes.

In addition to genetic polymorphisms, various confounding factors may contribute to the variability in the plasma concentrations or clinical outcomes of rivaroxaban. Our study identified and analyzed several types of confounding factors. Firstly, rivaroxaban undergoes elimination through both renal excretion and hepatic metabolism (Kvasnicka et al., 2017), thereby age, BW, liver and kidney function were considered confounding factors. Additionally, our study identified AST/ALT ratios as a significant covariate influencing CL/F and V/F, it is important to consider the potential for confounding when interpreting the relationship between ABCB1 genotypes and rivaroxaban exposure. Specifically, differences in AST/ALT ratios across ABCB1 genotype groups could contribute to the observed genotype-exposure correlations. Secondly, CHA₂DS₂-VASc and HAS-BLED scores are clinical tools used to assess stroke/thromboembolism or bleeding risk in AF patients, respectively (Lane and Lip, 2012; Pisters et al., 2010). Variations in these scores among different genetic polymorphisms may lead to diverse clinical outcomes; thus, they were also included as confounding factors. In addition, polypharmacy is common among AF patients and is associated with increased mortality and bleeding risk (Proietti et al., 2016) due to drug-drug interactions (DDIs) (Ferri et al., 2022). DDIs can manifest as both pharmacokinetic and pharmacodynamic interactions. For instance, amiodarone can alter rivaroxaban plasma concentrations, while drugs like clopidogrel, aspirin, cilostazol, nonsteroidal anti-inflammatory drugs (NSAIDs), and blood-activating traditional Chinese medicines (BATCMs) can increase bleeding risk in AF patients taking rivaroxaban (Ferri et al., 2022). Therefore, these concomitant medications were considered confounding factors in our analysis. In general, the results of confounding factors analyses indicated no significant differences among different ABCB1 genotypes in the included patients (Table 3), suggesting that the observed effects of ABCB1 polymorphisms on rivaroxaban exposure are unlikely to be confounded by these factors. However, given the inherent limitations of observational studies, residual confounding cannot be entirely ruled out (Euser et al., 2009). Future studies incorporating larger sample sizes may help further disentangle the independent contributions of genetic factors.

In this study, we observed simultaneous impacts of ABCB1 genetic polymorphisms on both plasma concentrations and clinical outcomes. Specifically, carriers of the ABCB1 c.2482-2236C>T CC genotype exhibited higher C_max/D and bleeding risk compared to those with CT or TT genotypes. On the other hand, carriers of the ABCB1 3435C>T TT genotype showed lower C_trough/D and higher thromboembolic risk than CC carriers. These findings suggest a hypothesis that ABCB1 genotypes may influence clinical outcomes by altering the plasma concentrations of rivaroxaban. However, it is important to note that this hypothesis is based on trend observations from model simulation results. To validate this hypothesis more conclusively, a well-designed prospective study with a larger number of participants is necessary. Such a study would provide stronger evidence regarding the impact of ABCB1 genetic variants on both plasma concentrations and clinical outcomes of rivaroxaban.

Limitation

It is important to recognize several potential limitations that could affect the applicability and precision of our findings.

Firstly, our study lacked adequate samples during the absorption phase, leading to the adoption of a fixed k_a value of 0.617 h⁻¹ based on a previous rivaroxaban PPK study involving Japanese patients (Kaneko et al., 2013). Given potential differences between Chinese and Japanese populations (Johnson, 1997), this fixed parameter may introduce inaccuracies into the PPK model, consequently affecting the outcomes of model-based simulations.

Secondly, While the one-compartment model demonstrated superior performance in describing the PK of rivaroxaban in Chinese NVAF patients, it is important to recognize that the choice of model structure depends on the specific dataset and population characteristics. The observed biphasic decline in plasma concentrations (Figure 2) suggests that a two-compartment model might be more appropriate in certain scenarios. However, in our study, the one-compartment model provided a better fit based on OFV, while maintaining parsimony and avoiding overparameterization. Despite rigorous model development and validation procedures being followed, the potential for structural misspecification remains a limitation of the current study. Thus, the findings and conclusions of this study should be interpreted with caution.

Thirdly, our study examined individual ABCB1 SNPs’ effects on rivaroxaban plasma levels and outcomes, but it is important to note that SNPs often occur in linkage disequilibrium. Haplotypes, such as 1236T-2677T-3435T, may better explain rivaroxaban variability than single SNPs (Lähteenmäki et al., 2021). The absence of haplotype analysis is a limitation, as haplotypes capture more comprehensive genetic interactions. Future studies should explore haplotypes to clarify their combined impact on rivaroxaban exposure and risks, aiding personalized dosing strategies.

Fourthly, while we analyzed CHA₂DS₂-VASc and HAS-BLED scores across ABCB1 genotypes, finding no significant differences, this approach has limitations. These composite indices do not fully capture individual patient characteristics or confounders like socioeconomic status, medication adherence, or lifestyle factors, which may influence genotype-outcome associations. Our study used cross-sectional comparisons instead of time-to-event analyses (e.g., Cox models), risking bias by ignoring temporal variations in exposure or outcomes. Additionally, limited adjustments beyond composite scores may leave residual confounding unaddressed. Future studies should employ advanced methods, including time-to-event analysis and multivariable adjustments, to better assess ABCB1 polymorphisms’ effects on rivaroxaban outcomes and reduce biases.

Fifthly, Our study links ABCB1 polymorphisms to rivaroxaban plasma levels and clinical outcomes like bleeding or thromboembolic events. However, these correlations are inferred indirectly via model-based simulations, not direct pharmacokinetic-pharmacodynamic (PK-PD) relationships. A key limitation is the lack of a formal PK-PD model linking rivaroxaban exposure to outcomes. While ABCB1 polymorphisms may modulate outcomes via altered PK, results should be interpreted cautiously. Future studies with PK-PD models, integrating anti-Factor Xa activity or prothrombin time, are needed for clearer mechanistic insights.

Lastly, our study population comprised exclusively Chinese individuals, necessitating further investigations involving diverse racial backgrounds to validate the generalizability of our results. In general, given these limitations, caution should be exercised in interpreting and generalizing the conclusions drawn from our findings. Future larger and well-designed studies are warranted to corroborate and expand upon our observations.

Conclusion

We have successfully developed a PPK model of rivaroxaban in Chinese NVAF patients, incorporating AST/ALT with CL/F and V/F. The model-based simulations indicated that ABCB1 1236C>T TT or ABCB1 c.2482-2236C>T CC genotypes exhibited higher C_max/D compared to other genotypes, respectively. Additionally, the ABCB1 3435C>T TT genotype showed lower C_trough/D than CC or CT genotypes. For clinical outcomes, ABCB1 c.2482-2236C>T CC genotype exhibited a higher bleeding risk than CT or TT genotypes, while ABCB1 3435C>T TT genotype demonstrated a higher thromboembolic risk than CC genotype. However, it is important to note that our conclusions should be interpreted cautiously due to inherent limitations in our study design. These findings warrant validation through larger well-designed prospective studies in the future.

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 humans were approved by Fujian Provincial Hospital institutional ethics committee (No. k2022-09-014). The studies were conducted in accordance with the local legislation and institutional requirements. The participants provided their written informed consent to participate in this study.

Author contributions

FW: Data curation, Formal Analysis, Methodology, Validation, Writing – original draft. ZL: Data curation, Formal Analysis, Methodology, Software, Validation, Visualization, Writing – original draft, Writing – review and editing. YH: Data curation, Formal Analysis, Investigation, Resources, Validation, Writing – review and editing. QL: Data curation, Resources, Validation, Writing – review and editing. LZ: Data curation, Resources, Validation, Writing – review and editing. HW: Data curation, Methodology, Resources, Validation, Writing – review and editing. HG: Resources, Validation, Writing – review and editing. MC: Data curation, Validation, Writing – review and editing. YL: Data curation, Validation, Writing – review and editing. XL: Conceptualization, Formal Analysis, Methodology, Supervision, Writing – review and editing, Software, Validation, Visualization. MC: Conceptualization, Data curation, Formal Analysis, Funding acquisition, Investigation, Methodology, Project administration, Resources, Supervision, Writing – review and editing.

Funding

The author(s) declare that financial support was received for the research and/or publication of this article. This study was financed by Guiding Projects issued by Fujian Provincial Department of Science and Technology, China (Project NO: 2023Y0045), the Department of Finance, Fujian Province, China (Project. No (2024)881-2024058) and Fujian Provincial Joint Fund for Scientific and Technological Innovation Project (2024Y9009).

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.

Generative AI statement

The author(s) declare that no Generative AI was used in the creation of this manuscript.

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.2025.1574949/full#supplementary-material

References

Chugh, S. S., Havmoeller, R., Narayanan, K., Singh, D., Rienstra, M., Benjamin, E. J., et al. (2014). Worldwide Epidemiology of atrial fibrillation A global burden of disease 2010 study. Circulation 129 (8), 837–847. doi:10.1161/CIRCULATIONAHA.113.005119

PubMed Abstract | CrossRef Full Text | Google Scholar

Cohen, J. A., and Kaplan, M. M. (1979). The SGOT/SGPT ratio—an indicator of alcoholic liver disease. Dig. Diseases&Sciences 24 (11), 835–838. doi:10.1007/BF01324898

PubMed Abstract | CrossRef Full Text | Google Scholar

Daly, A. K., and King, B. P. (2003). Pharmacogenetics of oral anticoagulants. Pharmacogenetics 13 (5), 247–252. doi:10.1097/00008571-200305000-00002

PubMed Abstract | CrossRef Full Text | Google Scholar

Derogis, P. B. M., Sanches, L. R., de Aranda, V. F., Colombini, M. P., Mangueira, C. L. P., Katz, M., et al. (2017). Determination of rivaroxaban in patient's plasma samples by anti-Xa chromogenic test associated to High Performance Liquid Chromatography tandem Mass Spectrometry (HPLC-MS/MS). Plos One 12 (2), e0171272. doi:10.1371/journal.pone.0171272

PubMed Abstract | CrossRef Full Text | Google Scholar

Euser, A. M., Zoccali, C., Jager, K. J., and Dekker, F. W. (2009). Cohort studies: prospective versus retrospective. Nephron Clin. Pract. 113 (3), C214–C217. doi:10.1159/000235241

PubMed Abstract | CrossRef Full Text | Google Scholar

Ferri, N., Colombo, E., Tenconi, M., Baldessin, L., and Corsini, A. (2022). Drug-drug interactions of direct oral anticoagulants (DOACs): from pharmacological to clinical practice. Pharmaceutics 14 (6), 1120. doi:10.3390/pharmaceutics14061120

PubMed Abstract | CrossRef Full Text | Google Scholar

Giannini, E., Botta, F., Fasoli, A., Ceppa, P., Risso, D., Lantieri, P. B., et al. (1999). Progressive liver functional impairment is associated with an increase in AST/ALT ratio. Dig. Dis. and Sci. 44 (6), 1249–1253. doi:10.1023/a:1026609231094

PubMed Abstract | CrossRef Full Text | Google Scholar

Gitlin, N. (1982). The serum glutamic oxaloacetic transaminase/serum glutamic pyruvic transaminase ratio as a prognostic index in severe acute viral hepatitis. American J. Gastroenterology 77 (1), 2–4.

PubMed Abstract | Google Scholar

Gnoth, M. J., Buetehorn, U., Muenster, U., Schwarz, T., and Sandmann, S. (2011). In vitro and in vivo P-glycoprotein transport characteristics of rivaroxaban. J. Pharmacol. Exp. Ther. 338 (1), 372–380. doi:10.1124/jpet.111.180240

PubMed Abstract | CrossRef Full Text | Google Scholar

January, C. T., Wann, L. S., Calkins, H., Chen, L. Y., Cigarroa, J. E., Cleveland, J. C., et al. (2019). 2019 AHA/ACC/HRS focused update of the 2014 AHA/ACC/HRS guideline for the management of patients with atrial fibrillation: a report of the American College of cardiology/American heart association task force on clinical practice guidelines and the heart rhythm society in collaboration with the society of thoracic surgeons. Circulation 140 (2), E125–E151. doi:10.1161/CIR.0000000000000665

PubMed Abstract | CrossRef Full Text | Google Scholar

Johnson, J. A. (1997). Influence of race or ethnicity on pharmacokinetics of drugs. J. Pharm. Sci. 86 (12), 1328–1333. doi:10.1021/js9702168

PubMed Abstract | CrossRef Full Text | Google Scholar

Kaneko, M., Tanigawa, T., Hashizume, K., Kajikawa, M., Tajiri, M., and Mueck, W. (2013). Confirmation of model-based dose selection for a Japanese phase III study of rivaroxaban in non-valvular atrial fibrillation patients. Drug Metabolism Pharmacokinet. 28 (4), 321–331. doi:10.2133/dmpk.dmpk-12-rg-109

PubMed Abstract | CrossRef Full Text | Google Scholar

Kong, X. L., Ma, Y. C., Chen, J. H., Luo, Q., Yu, X., Li, Y., et al. (2013). Evaluation of the Chronic Kidney Disease Epidemiology Collaboration equation for estimating glomerular filtration rate in the Chinese population. Nephrol. Dial. Transplant. 28 (3), 641–651. doi:10.1093/ndt/gfs491

PubMed Abstract | CrossRef Full Text | Google Scholar

Kvasnicka, T., Malikova, I., Zenahlikova, Z., Kettnerova, K., Brzezkova, R., Zima, T., et al. (2017). Rivaroxaban - metabolism, pharmacologic properties and drug interactions. Curr. Drug Metab. 18 (7), 636–642. doi:10.2174/1389200218666170518165443

PubMed Abstract | CrossRef Full Text | Google Scholar

Lähteenmäki, J., Vuorinen, A. L., Pajula, J., Harno, K., Lehto, M., Niemi, M., et al. (2021). Pharmacogenetics of bleeding and thromboembolic events in direct oral anticoagulant users. Clin. Pharmacol. and Ther. 110 (3), 768–776. doi:10.1002/cpt.2316

PubMed Abstract | CrossRef Full Text | Google Scholar

Lane, D. A., and Lip, G. Y. H. (2012). Use of the CHA2DS2-VASc and HAS-BLED scores to aid decision making for thromboprophylaxis in nonvalvular atrial fibrillation. Circulation 126 (7), 860–865. doi:10.1161/CIRCULATIONAHA.111.060061

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, X. G., Wu, Y. X., Sun, S. S., Mei, S., Wang, J., Wang, Q., et al. (2015). Population pharmacokinetics of vancomycin in postoperative neurosurgical patients. J. Pharm. Sci. 104 (11), 3960–3967. doi:10.1002/jps.24604

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, X. G., Wu, Y. X., Sun, S. S., Zhao, Z., and Wang, Q. (2016). Population pharmacokinetics of vancomycin in postoperative neurosurgical patients and the application in dosing recommendation. J. Pharm. Sci. 105 (11), 3425–3431. doi:10.1016/j.xphs.2016.08.012

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, X. Y., Gu, Z. C., Wang, Z., Xu, Q., Ma, C.L., and Lv, Q.Z. (2024). Mutant CYP3A4/5 correlated with clinical outcomes by affecting rivaroxaban pharmacokinetics and pharmacodynamics in patients with atrial fibrillation, in Cardiovasc. Drugs. Ther. 38 (6), 1315–1325. doi:10.1007/s10557-023-07495-4

CrossRef Full Text | Google Scholar

Li, Z., Yang, S. Y., Hua, Z. X., Lu, Y., and Li, X. (2023). Population pharmacokinetics of rivaroxaban in Chinese deep vein thrombosis patients and the exposure simulation for dosing recommendation. Naunyn-Schmiedebergs Archives Pharmacol. 397, 3351–3362. doi:10.1007/s00210-023-02798-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, X. Q., Zhang, Y. F., Ding, H. Y., Yan, M. M., Jiao, Z., Zhong, M. K., et al. (2022). Population pharmacokinetic and pharmacodynamic analysis of rivaroxaban in Chinese patients with non-valvular atrial fibrillation. Acta Pharmacol. Sin. 43 (10), 2723–2734. doi:10.1038/s41401-022-00892-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Mould, D. R., and Upton, R. N. (2012). Basic concepts in population modeling, simulation, and model-based drug development. Cpt Pharmacometrics and Syst. Pharmacol. 1 (9), e6. doi:10.1038/psp.2012.4

PubMed Abstract | CrossRef Full Text | Google Scholar

Patel, M. R., Mahaffey, K. W., Garg, J., Pan, G., Singer, D. E., Hacke, W., et al. (2011). Rivaroxaban versus warfarin in nonvalvular atrial fibrillation. New England J. Med. 365 (10), 883–891. doi:10.1056/NEJMoa1009638

PubMed Abstract | CrossRef Full Text | Google Scholar

Pisters, R., Lane, D. A., Nieuwlaat, R., de Vos, C. B., Crijns, H. J. G. M., and Lip, G. Y. H. (2010). A novel user-friendly score (HAS-BLED) to assess 1-year risk of major bleeding in patients with atrial fibrillation: the Euro Heart Survey. Chest 138 (5), 1093–1100. doi:10.1378/chest.10-0134

PubMed Abstract | CrossRef Full Text | Google Scholar

Proietti, M., Raparelli, V., Olshansky, B., and Lip, G. Y. H. (2016). Polypharmacy and major adverse events in atrial fibrillation: observations from the AFFIRM trial. Clin. Res. Cardiol. 105 (5), 412–420. doi:10.1007/s00392-015-0936-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Sychev, D. A., Vardanyan, A., Rozhkov, A., Hachatryan, E., Badanyan, A., Smirnov, V., et al. (2018). CYP3A activity and rivaroxaban serum concentrations in Russian patients with deep vein thrombosis. Genet. Test. Mol. Biomarkers 22 (1), 51–54. doi:10.1089/gtmb.2017.0152

PubMed Abstract | CrossRef Full Text | Google Scholar

Testa, S., Tripodi, A., Legnani, C., Pengo, V., Abbate, R., Dellanoce, C., et al. (2016). Plasma levels of direct oral anticoagulants in real life patients with atrial fibrillation: results observed in four anticoagulation clinics. Thrombosis Res. 137, 178–183. doi:10.1016/j.thromres.2015.12.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Williams, A. L. B., and Hoofnagle, J. H. (1988). Ratio of serum aspartate to alanine aminotransferase in chronic hepatitis. Relationship to cirrhosis. Gastroenterology 95 (3), 734–739. doi:10.1016/s0016-5085(88)80022-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Xiang, Q., Wang, Z., Mu, G. Y., Xie, Q., Liu, Z., Zhou, S., et al. (2023). Genetic variants influenced the risk of bleeding and pharmacodynamics of rivaroxaban in patients with nonvalvular atrial fibrillation: a multicentre prospective cohort study. Clin. Transl. Med. 13 (5), e1263. doi:10.1002/ctm2.1263

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, D., Qin, W., Du, W. W., Wang, X., Chen, W., and Li, P. (2022). Effect of ABCB1 gene variants on rivaroxaban pharmacokinetic and hemorrhage events occurring in patients with non-valvular atrial fibrillation. Biopharm. and Drug Dispos. 43 (4), 163–171. doi:10.1002/bdd.2327

PubMed Abstract | CrossRef Full Text | Google Scholar