Integrated Population Pharmacokinetic Analysis of Rivaroxaban Across Multiple Patient Populations.

The population pharmacokinetics (PK) of rivaroxaban have been evaluated in several population-specific models. We developed an integrated population PK model using pooled data from 4,918 patients in 7 clinical trials across all approved indications. Effects of gender, age, and weight on apparent clearance (CL/F) and apparent volume of distribution (V/F), renal function, and comedication on CL/F, and relative bioavailability as a function of dose (F) were analyzed. Virtual subpopulations for exposure simulations were defined by age, creatinine clearance (CrCL) and body mass index (BMI). Rivaroxaban PK were adequately described by a one-compartment disposition model with a first-order absorption rate constant. Significant effects of CrCL, use of comedications, and study population on CL/F, age, weight, and gender on V/F, and dose on F were identified. CrCL had a modest influence on exposure, whereas age and BMI had a minor influence. The model was suitable to predict rivaroxaban exposure in patient subgroups of special interest.

WHAT IS THE CURRENT KNOWLEDGE ON THE TOPIC?

☑ Population PK models in healthy individuals and several disease indications have been developed. However, it is not known if the covariate relationships identified can be applied to the whole population.

WHAT QUESTION DID THIS STUDY ADDRESS?

☑ Do covariate relationships identified in previous rivaroxaban population PK models apply to the whole population in the global clinical trials?

WHAT DOES THIS STUDY ADD TO OUR KNOWLEDGE?

☑ Our model uses covariates consistently across indications and accounts for PK differences between populations. The PK of rivaroxaban was demonstrated to be predictable by the developed integrated population PK model. Renal function had a greater effect on rivaroxaban exposure than age, body weight, and comedication use.

HOW MIGHT THIS CHANGE DRUG DISCOVERY, DEVELOPMENT, AND/OR THERAPEUTICS?

☑ Exposure to rivaroxaban is stable in response to age, body weight, renal function, and comedication use. Individual rivaroxaban exposure in approved indications can be estimated based on dose and patient characteristics.

Rivaroxaban, a direct oral anticoagulant that reversibly inhibits factor Xa, is approved in the United States and Europe for the following indications: prevention of venous thromboembolism (VTE) in adults undergoing elective hip or knee replacement surgery; treatment of deep vein thrombosis (DVT), pulmonary embolism (PE), and prevention of recurrent DVT and PE in adults; and prevention of stroke and systemic embolism in adults with nonvalvular atrial fibrillation (AF). The European Medicines Agency has also approved rivaroxaban for the prevention of atherothrombotic events in adults with acute coronary syndrome (ACS), co‐administered with aspirin alone or in combination with clopidogrel or ticlopidine. Rivaroxaban is classified under the Biopharmaceutics Classification System1 as a low solubility, high permeability compound (class 2). Approximately one‐third of the dose is eliminated renally as unchanged drug and the remainder is subject to metabolic degradation via cytochrome P450 (CYP)3A4, CYP2J2 and CYP‐independent biotransformation processes.2

The extensive global clinical development program for rivaroxaban was based on the concept that fixed dose regimens of rivaroxaban can provide reliable anticoagulation without the need for routine coagulation assays or measurement of drug levels to guide dose adjustment. In healthy volunteers, rivaroxaban exhibited a highly predictable and dose‐proportional pharmacokinetic/pharmacodynamic response with few drug–drug and food–drug interactions.2, 3 Phase II dose‐ranging studies in patients undergoing elective hip or knee replacement surgery and receiving rivaroxaban for the prevention of VTE,4, 5, 6 patients receiving rivaroxaban for the treatment of DVT,7, 8 and patients with ACS receiving rivaroxaban for the prevention of atherothrombotic events9 demonstrated that rivaroxaban has a wide therapeutic window. Dosing regimens were subsequently confirmed in large‐scale phase III studies.10, 11, 12, 13, 14, 15, 16, 17

Population pharmacokinetic (PK) models are routinely used to identify sources of variability between patients or subgroups of patients, and the knowledge gained from population PK modeling and simulation may be used to optimize dose and ultimately patient response. Previously, population PK models have been developed separately to describe the PK of rivaroxaban in healthy individuals,18 for VTE prevention,19, 20 VTE treatment,21 and for patients with AF21, 22 and ACS.23 A two‐compartment model was used to describe the PK of rivaroxaban in healthy volunteers; however, as is commonly the case, sparse sampling was used to collect PK samples in the global clinical trials. A one‐compartment model that is simpler than the two‐compartment model was consequently developed for use across the patient populations. These models have identified the factors that influence rivaroxaban PK and together have shown that they are similar across indications. For example, age and renal function have been shown to influence rivaroxaban clearance, which is expected, given that approximately one‐third of the rivaroxaban dose is eliminated renally2, 24 and renal function decreases with advancing age,25 whereas volume of distribution is affected by age and a size parameter, for example, lean body weight.18, 19, 20, 21, 22, 23 In addition, the oral bioavailability of rivaroxaban has been demonstrated to be dose‐dependent.19, 20, 21, 23 However, the impact of these covariates on rivaroxaban PK has not been simultaneously determined across the patient populations.

The aims of this study were: (1) to develop an integrated and comprehensive population PK model for rivaroxaban across all approved indications using pooled PK data obtained from global clinical trials; (2) to harmonize the relationships between relevant covariates and rivaroxaban exposure across all indications; and (3) to use the model to perform exposure simulations in patient subgroups of special interest, including elderly patients, patients with renal impairment, and patients with obesity. As exposure‐response studies across all rivaroxaban indications are planned, such a population PK model is a prerequisite for providing reliable exposure estimates owing to the sparse PK data obtained during the global clinical studies.

MATERIALS AND METHODS

Studies included in the analyses

Available rivaroxaban concentration–time data from 4,918 patients from 7 global clinical trials in different indications were pooled (Table 1). Briefly, these studies were three phase II studies for the prevention of VTE in patients undergoing elective hip or knee replacement surgery (ODIXa‐Hip2 [dose‐ranging],4 ODIXa‐OD‐Hip [dose‐ranging],5 ODIXa Knee6), two phase II studies for the treatment of acute symptomatic VTE (ODIXa‐DVT,7 EINSTEIN‐DVT [dose ranging]8), one phase II study for the prevention of cardiovascular events in patients with ACS (ACS TIMI‐469), and one phase III study of stroke prevention in patients with AF (ROCKET AF14). Data from the immediate postsurgical phase of the VTE prevention studies was excluded on the basis of previous results that demonstrated that patients can be categorized as slow or fast absorbers of rivaroxaban during the first day after surgery.20

Table 1

Description of studies and demographic and baseline characteristics of patients included in the integrated population pharmacokinetic model

Indication	VTE‐P			VTE‐T		AF	ACS	Total
Study	ODIXa‐Hip2	ODIXa‐OD‐Hip	ODIXa‐Knee	ODIXa‐DVT	EINSTEIN DVT	ROCKET AF	ACS TIMI‐46
Study number	10944	11527	10945	11223	11528	3001	2001	NA
Phase	II	II	II	II	II	III	II	NA
No.	553	666	501	470	400	161	2,290	5,041
Dose, mg and regimen	2.5, 5, 10, 20, 30 b.i.d. dose‐ranging	5, 10, 20, 30, 40 o.d. dose‐ranging	2.5, 5, 10, 20, 30 b.i.d.	10, 20, 30 b.i.d. 40 o.d.	20, 30, 40 o.d. dose‐ranging	20 (15 in patients with CrCL 30–49 mL/min) o.d.	2.5, 5, 7.5, 10, 15, 20 o.d. or b.i.d.	NA
Treatment duration	9 ± 2 days	9 ± 2 days	8 ± 2 days	12 weeks	12 weeks	Median 590 days	180 days	NA
No. of patients with PK observations	512	661	463	470	400	161	2,251	4,918
Total number of PK observations	2,197	3,995	1,841	3,524	1,110	800	9,376	22,843
Median number of PK observations per patient (min, max)	4 (1, 10)	6 (1, 8)	4 (1, 10)	8 (1, 9)	3 (1, 5)	5 (2, 15)	4 (1, 13)	NA
Age, years	65.01 (9.92)	64.68 (10.47)	66.64 (9.64)	59.29 (16.04)	58.18 (17.00)	65.46 (9.51)	57.23 (9.51)	60.53 (11.82)
Weight, kg	77.68 (14.17)	76.75 (14.05)	88.31 (18.79)	80.07 (17.87)	80.50 (17.23)	85.40 (18.60)	84.67 (16.50)	82.48 (16.87)
CrCL, mL/min	99.00 (30.89)	85.45 (28.67)	110.5 (39.16)	93.37 (37.35)	96.63 (38.65)	81.76 (32.06)	100.4 (31.78)	97.74 (33.97)
CrCLtruncated, mL/min	98.25 (29.03)	85.27 (27.99)	108.4 (34.65)	92.68 (35.65)	96.16 (37.62)	81.76 (32.06)	100.1 (30.73)	97.17 (32.34)
Lean body mass, kg	53.46 (9.40)	53.50 (9.39)	54.72 (10.32)	56.48 (10.45)	56.03 (10.41)	57.53 (9.86)	59.72 (9.30)	57.05 (10.00)
Serum creatinine, mg/dL	0.78 (0.19)	0.90 (0.21)	0.79 (0.22)	0.97 (0.25)	0.94 (0.24)	1.09 (0.29)	0.98 (0.24)	0.93 (0.25)
Body surface area, m2	1.86 (0.20)	1.86 (0.20)	1.96 (0.22)	1.91 (0.24)	1.92 (0.23)	1.96 (0.22)	1.97 (0.21)	1.93 (0.22)
Female, %	60.2	58.0	62.7	38.7	49.5	37.9	22.2	39.3

ACS, acute coronary syndromes; AF, atrial fibrillation; b.i.d., twice daily; CrCL, creatinine clearance; CrCLtruncated, Tietz‐truncated creatinine clearance; DVT, deep vein thrombosis; NA, not applicable; o.d., once daily; PK, pharmacokinetic; VTE‐P, venous thromboembolism prevention; VTE‐T, venous thromboembolism treatment

Data are presented as mean (SD) unless stated otherwise.

All studies were conducted in accordance with the Declaration of Helsinki and the principles of Good Clinical Practice. Study protocols and amendments were approved by independent ethics committees. All participants provided written informed consent prior to study enrollment.

Sampling and investigations

Details of the studies in which PK sampling was performed are shown in Table 1. Blood samples were collected using a sparse sampling approach, as described previously.20, 21, 22, 23 Rivaroxaban plasma concentrations were determined using a selective chromatographic assay combined with tandem mass spectrometry, as described previously.20, 21, 22 The lower limit of quantification was 0.5 μg/L. Any data below the lower limit of quantification for the assay were excluded.

Population pharmacokinetic modeling

Analysis was conducted using a nonlinear mixed‐effects modeling approach in NONMEM (ICON Development Solutions, version 7.3) on Windows Server 2012 R2. This modeling technique quantifies random effects, such as unexplained interindividual variability and residual variability, as well as the influence of measured patient characteristics or covariates (fixed effects) on model parameters and allows population means to be estimated. The first‐order conditional estimation algorithm with η‐ɛ interaction was used for all analyses. The R (The R Foundation for Statistical Computing, versions 2.31 and 3.2.2) and PsN (version 4.2.0 with ActiveState Perl version 5.16) software packages were used for model evaluation.

Structural pharmacokinetic model

Previous population PK modeling studies for the phase II and III populations19, 20, 21, 23 using a one‐compartment model suggested that the model structure and PK parameters were consistent across different patient populations. Thus, the same model structure was used for the current integrated population PK study. Model refinement was based on the change in objective function value (OFV) and/or the model qualification assessments described below. A comparison of log‐transformed vs. nontransformed data was also made.

Covariate model development

The difference in OFV was used for selection between models. Structural model components (e.g., covariate influence) were incorporated into the model if the likelihood ratio test (i.e., difference in the OFV, which is an approximately χ2 distribution) showed significance at a critical level of P ≤ 0.01 (change in OFV >6.63 with 1 degree of freedom). Components remained in the final model when, following backward elimination, the difference in OFV showed significance at a critical level of P ≤ 0.001 (a change in OFV of >10.8 with 1 degree of freedom).

The following covariates identified from prior analyses18, 19, 20, 21, 22, 23 were used in the base model as a starting point: effects of age at baseline on apparent clearance (CL/F); effects of renal function at baseline as estimated by creatinine clearance (CrCL) on CL/F; effects of weight on apparent volume of distribution (V/F); and relative bioavailability as a function of dose (F). The influence of selected comedications (weak, moderate, or strong CYP3A4 inhibitors, CYP3A4 inducers, and P‐glycoprotein [P‐gp] inhibitors) on CL/F was also tested.

The covariates and model properties tested before declaration of the final covariate model are shown in Supplementary Table S1. For CrCL, a modified Cockroft–Gault equation was used, as suggested by Tietz et al.26 to avoid high and physiologically implausible CrCL values. The Tietz‐truncated CrCL is set to a value of 140 mL/min × body surface area/(1.73 m2), if the CrCL calculated by the Cockroft–Gault equation exceeds the value of the Tietz‐truncated CrCL.

Continuous covariates were included in the model in power form according to the following equation where is the continuous covariate, is the population mean or median, is the model parameter, P is the typical value for that parameter, and is the power estimate for the covariate effect. Categorical covariates (i.e., gender and comedications) with M categories were included according to the following equation where represents the fractional change in P for category m (m in 2…M) of covariate X, and is an indicator variable with value 1 if patient i falls into category m.

Model qualification

Diagnostic plots of observed data vs. population prediction (PRED) and individual predictions (IPRED) were examined for adequate fit. Plots of conditional weighted residuals (CWRES) vs. PRED and vs. time (times after last and first doses) were inspected for evidence of systematic lack of fit, and to confirm the absence of bias in the error distributions. Individual deviations from the population mean were expected to be normally distributed with a mean of zero and variance σ2. The distribution of CWRES was checked using the quantile–quantile normal plot to assess the assumption of normality. To verify absence of bias, between‐subject random effects were graphed in scatter plots vs. key continuous model covariates with potential trends. Prediction‐corrected visual predictive checks were also used for model diagnostics to allow for comparison across dose levels.

Exposure simulation

The final integrated population PK model was used to simulate rivaroxaban exposure for subpopulations of special interest based on the large pool of patients from the treatment and control groups in the phase II and phase III studies for whom complete covariate information was available. Virtual subpopulations were defined by age (18 to <65 years, 65 to ≤75 years, and >75 years), renal function (>80 mL/min, 50 to ≤80 mL/min, 30 to <50 mL/min, and <30 mL/min), and body mass index (BMI) as a marker of body weight (<18.5 kg/m2, 18.5 to <25 kg/m2, 25 to <30 kg/m2, 30 to <40 kg/m2, and ≥40 kg/m2). The classifications yielded 60 possible subgroups. The target number of individuals per subgroup was 1,000.

Any subgroups of fewer than 60 patients were omitted from the exposure predictions; for subgroups of 60–1,000 patients, sampling with replacement was performed to obtain 1,000 individuals, for a subgroup of more than 1,000 patients, and sampling without replacement was performed to obtain 1,000 individuals. Individual steady state exposure estimates of area under the plasma concentration–time curve (AUC) from time 0 to 24 hours (AUC0–24), maximum plasma concentration (Cmax), and trough plasma concentration (Ctrough) at approved doses for each indication were then calculated for each subgroup and for each indication using the developed model.

RESULTS

Patient demographics

A total of 22,843 PK observations from 4,918 patients who received rivaroxaban across six phase II studies and one phase III study contributed to development of the integrated population PK model (Table 1). Rivaroxaban doses ranged from 2.5 mg once daily (o.d.) to 30 mg twice daily (b.i.d.). Demographic data and baseline characteristics for the patients included in the model are presented in Table 1. With the exception of gender, these parameters were broadly consistent across studies. The overall proportion of female patients across all studies was 39.3%, but ranged from 22.2–62.7% in individual studies.

Pharmacokinetic model development

The final model structure is shown in Supplementary Figure S1 and is mathematically described using Eq. 1 in the Supplementary Materials.

As discussed earlier, a one‐compartment disposition model with first order absorption and elimination, parameterized in terms of CL/F, V/F, and a first‐order absorption rate constant (ka) was chosen as the structural model.

The key model‐building steps are shown in Table 2. Dose had a highly significant effect on F; a nonlinear function of the dose was used to capture the saturation in F (Figure 1 a). Inclusion of a study‐specific effect on CL/F resulted in a distinct improvement in the fit of the model to the data, with factors of 0.85 for AF, 1.14 for ACS, 1.04 for VTE prevention (≤72 hours), and 1.29 for VTE prevention (>72 hours), relative to VTE treatment (study‐specific factor = 1). A proportion residual variance approach using nontransformed data led to an improvement in the residual patterns compared with the use of log‐transformed data.

Table 2

Key modeling steps

Run	Based on	Description	ΔOFV	Condition number
0	–	1 comp IIV on ka+CL/F+V/F prop. res. Error nonlog transform	–	13
1	0	F estimated per dose level in reference to 10 mg (F = 1 fixed)a	−724.1	32
2	1	Covariance between CL/F and V/F	−556.2	58
3	2	Continuous function parameterized by two fixed effects parameters to describe dose‐dependent Fb	11.7	223
4	3	Tietz‐truncated CrCL on CL/F	−206.6	220
5	4	Weight on V/F	−212.2	239
6	5	Age on CL/F	0	245
7	6	Age on V/F	−55.0	320
8	7	Gender on V/F	−78.1	NAc
9	8	Weight on CL/F	−72.8	215
10	9	Removing age on CL/F	1.3	212
11	10	Comedication on CL/F	−80.8	207
12d	11	Study population on CL/F	−550.3	184

CL/F, apparent clearance; CrCL, creatinine clearance; ΔOFV, delta objective function value; F, relative oral bioavailability; IIV, interindividual variability; ka, first‐order absorption rate constant; NA, not applicable; V/F, apparent volume of distribution.

Seven mixed effects parameters.

Degrees of freedom reduced by 5.

Run 8 did not render a successful covariance step.

Final population pharmacokinetic model.

Figure 1

(a) Visualization of estimated relative bioavailability function used in the integrated population pharmacokinetic model compared with the 10 mg dose (F = 1). Dots represent the doses used in the studies. Values indicate the number of patients who received each dose. (b) Visualization of the estimated effects of patient characteristic covariates on apparent clearance (CL/F) and apparent volume of distribution (V/F). Estimated effects are given as fold‐change compared with reference category, symbols represent point estimates in fold‐change, and bars represent the 95% confidence interval of the point estimate obtained via bootstrapping. For continuous covariates (CrCLtruncated, weight and age), low and high are defined as the 5th and 95th percentile of the covariate distribution of the analysis population, respectively. The effects in these groups are shown relative to the median of the respective distribution. AF, atrial fibrillation; CrCLtruncated, Tietz‐truncated creatinine clearance; CYP3A4, cytochrome P450 3A4; PGP, p‐glycoprotein; VTE‐P, venous thromboembolism prevention; VTE‐T, venous thromboembolism treatment.

The covariate analysis demonstrated a significant effect of CrCL, comedications (proportion of patients taking weak, moderate, or strong CYP3A4 inhibitors, CYP3A4 inducers, or P‐gp inhibitors) and study population on CL/F and of age, weight, and gender on V/F (Table 3).The estimated effects of CrCL, weight, comedication, and study population on CL/F, and of age, weight, and gender on V/F are also visualized in a forest plot (Figure 1 b). The resulting model, after applying the effects of covariates, was taken as the final population PK model (run 12, Table 2) for the integrated patient population. The parameter estimates based on the final model are given in Table 3.

Table 3

Parameter estimates in the final population pharmacokinetic model

Parameter	Unit	Estimate	Relative SE (%)a	Lower confidence interval (2.5%)b	Upper confidence interval (97.5%)b
ka	1/h	0.821	2.36	0.780	0.860
CL/F	L/h	6.58	2.33	6.29	6.86
V/F	L	62.5	2.04	59.6	64.4
Fmin	Proportion	0.590	5.99	0.51	0.653
Fmax	Proportion	1.25	– (fixed)c	–	–
D50	mg	14.4	14.8	10.7	19.7
θCL/F, CrCL	Power	0.406	6.03	0.351	0.453
θCL/F, weight	Power	−0.278	15.4	−0.359	−0.187
θV/F, weight	Power	0.216	17.4	0.143	0.278
θV/F, Age	Power	−0.189	16.3	−0.246	−0.127
θV/F, Sex	Proportion	0.889	1.46	0.867	0.914
θCL/F, PGP	Proportion	0.966	1.73	0.933	1.00
θCL/F, Strong CYP3A4 inhibitor	Proportion	0.978	5.40	0.902	1.97
θCL/F, Medium CYP3A4 inhibitor	Proportion	0.863	3.79	0.793	0.920
θCL/F, Weak CYP3A4 inhibitor	Proportion	0.939	2.17	0.900	0.975
θCL/F, CYP3A4 inducer	Proportion	1.30	6.30	1.16	1.48
θCL/F, AF	Proportion	0.849	3.48	0.793	0.900
θCL/F, ACS	Proportion	1.14	1.93	1.10	1.18
θCL/F, VTE ≤72 h	Proportion	1.04	1.96	0.992	1.08
θCL/F, VTE >72 h	Proportion	1.29	1.76	1.24	1.34
ω2 ka	Variance	0.628 (32.8)d	5.39	0.544	0.703
ω2 CL/F	Variance	0.167 (10.3)d	2.96	0.157	0.176
ω2 CL/F, V/F	Covariance	0.0674	6.84	0.0586	0.0782
ω2 V/F	Variance	0.0391 (25.7)d	10.6	0.0318	0.048
σ2 prop	Variance	0.203 (10.3)d	1.54	0.197	0.21

ACS, acute coronary syndrome; AF, atrial fibrillation; CL/F, apparent clearance; CrCL, creatinine clearance; CYP3A4, cytochrome P450 3A4; D, dose; F, relative bioavailability as a function of dose; ka, first‐order absorption rate constant; V/F, apparent volume of distribution; PGP, P‐glycoprotein; VTE, venous thromboembolism.

Relative standard errors (SE) were obtained from the estimated covariance matrix.

The 95% confidence intervals of the estimates were obtained from a nonparametric bootstrap evaluation with 200 runs.

Relative bioavailability was fixed at 1.0 for a 10 mg dose.

For interindividual variability terms (ω2) and residual variability (σ2), the respective shrinkage estimates are displayed in parentheses.

Pharmacokinetic model qualification

The final model characterized the data well (Figure 2 a,b, and Supplementary Figure S2). CWRES were randomly scattered around zero across predicted range and time (Figure 3 , and Supplementary Figure S3). The prediction‐corrected visual predictive checks showed that there was good overall agreement between the final model and observed data, and the general trend was well described overall (Figure 2 b). The prediction‐corrected visual predictive checks also indicated a slight trend toward underestimation of Cmax. The impact of covariates on key parameters was qualitatively similar across all indications, as demonstrated by prediction‐corrected visual predictive checks stratified by age, CrCL, weight, gender, rivaroxaban dose, and study‐specific effects (Supplementary Figure S4). For CL/F and V/F, shrinkage was within the range of 20–30%, which is commonly accepted27; for ka, the numerical shrinkage value was marginally above 30% (Table 3). Thus, qualification of the model demonstrated that it is suitable for prediction of rivaroxaban exposure in patient subgroups of special interest.

Figure 2

(a) Observations vs. population predictions (PRED) and individual predictions (IPRED), for all studies. DV, dependent variable. (b) Prediction‐corrected visual predictive checks (VPCs) of the pooled dataset. The orange circles represent the prediction‐corrected observations; the black horizontal bars show the 5th percentiles, medians, and 95th percentiles of the prediction‐corrected observations within the binning interval; the widths of the horizontal bars represent the range between the 5th–95th percentiles of the sampling time points of the observations within the binning interval; the vertical black lines indicate the medians of the sampling times within the binning interval; the blue shaded area represents the 95% confidence band around the 5th and 95th percentiles of the prediction‐corrected simulations; the red shaded area represents the 95% confidence band around the median of the prediction‐corrected simulations; the values indicate the numbers of observations at each binning interval.

Figure 3

Conditionally weighted residuals (CWRES) vs. population predictions (PRED) and time after dose for all studies.

Using the final integrated population PK model, exposure at steady state was predicted for each virtual subpopulation (defined by age, renal function, and BMI) in each indication. Figure 4 shows exposure simulation results for AUC according to renal function, age, and weight for the AF indication; results for Cmax and Ctrough are shown in Supplementary Figure S5. The median AUC0–24 was 53% higher and Ctrough was 2.1‐fold higher at steady state in individuals with severe renal impairment compared with individuals with normal kidney function. The Cmax was increased by 35% in patients with severe renal impairment (Supplementary Figure S5). Apart from these trends in the population medians, the 90% prediction intervals for AUC0–24, Cmax, and Ctrough were largely overlapping. The influence of age and body weight on rivaroxaban PK was minor. Population medians varied by 15%, 10%, and 23% for AUC0–24, Cmax, and Ctrough, respectively, among the predefined age groups, and 7.5%, 7.6%, and 44%, respectively, among the BMI categories compared with the reference groups.

Figure 4

Simulation results (area under the plasma concentration–time curve from time 0–24 hours (AUC0–24) at steady state following rivaroxaban 20 mg (left‐hand graphs) and 15 mg once daily (right‐hand graphs)) for the subgroups* according to renal function (top row), age (middle row), and weight (bottom row) for the atrial fibrillation indication. Calculations were performed using analytical equations. Boxes show the 25th–75th percentiles; the horizontal line in the box indicates the median; the open rhombus in the box indicates the mean; the whiskers represent 1.5 × interquartile range. *Virtual subgroups were defined by age (adults, 18 to <65 years; elderly, 65 to ≤75 years; and very elderly, >75 years), renal function (normal, creatinine clearance (CrCL) >80 mL/min; mild renal impairment (RI), CrCL 50 to ≤ 80 mL/min; moderate RI, CrCL 30 to <50 mL/min; and severe RI, CrCL <30 mL/min) and body mass index (underweight, <18.5 kg/m2; normal weight; 18.5 to <25 kg/m2; overweight, 25 to <30 kg/m2; obese, 30 to <40 kg/m2; and morbidly obese ≥40 kg/m2).

DISCUSSION

We have developed the first integrated population PK model for rivaroxaban across all approved indications. The model was based on 22,843 PK observations from 4,918 patients who received rivaroxaban across six phase II studies4, 5, 6, 7, 8, 9 and one phase III study,14 which evaluated the efficacy and safety of rivaroxaban in patients with VTE, AF, and ACS. The demographic and clinical characteristics can, therefore, be considered to be representative of the range of patients in clinical practice. The results demonstrate that rivaroxaban has a predictable PK profile, confirming the results from early clinical studies in healthy volunteers3 and from previous indication‐specific models.19, 20, 21

The data presented here show that renal function is the main driver of rivaroxaban exposure, which is consistent with the observation that approximately one‐third of the rivaroxaban dose is excreted unchanged by the kidneys.2, 24 In comparison, the influence of two other covariates, age and body weight, on rivaroxaban PK is minor. These findings are in agreement with subgroup analyses of several phase III clinical trials, which demonstrated that the efficacy and safety profile of rivaroxaban was maintained across prespecified subgroups classified according to age, weight and BMI.14, 15, 16, 17, 28, 29

Population PK modeling has been used to predict exposure to several direct oral anticoagulants30, 31, 32, 33, 34, 35 and the population PK modeling approach has been used throughout the development of rivaroxaban, resulting in a lineage of modeling efforts18, 19, 20, 21, 22, 23; however, this is the first population PK model of rivaroxaban, and also of a direct oral anticoagulant, across all approved indications. The model uses covariates consistently across indications and accounts for PK differences between studies and/or populations. Prediction‐corrected visual predictive checks showed that there was good overall agreement between the final model and observed data, and that overall, the general trend was well described. Importantly, qualification of the model demonstrates that it is suitable to be used to predict rivaroxaban exposure in patient subgroups of special interest and will provide reliable exposure estimates for exposure–response analyses across all approved indications.

The effects of age, CrCL, and comedication on CL/F, and weight and sex on V/F were included as covariates, as well as F. As shown in Table 3, the effects of age and CrCL on CL/F and weight on V/F were tested initially; the latter two resulted in significant improvement in the OFV, whereas the effect of age on CL/F did not. This relation was kept in the model while evaluating the effects of age and gender on V/F and of weight on CL/F. These relationships rendered substantial drops in OFV. A backward deletion test showed that the effect of age on CL/F did not contribute to the predictive performance of the model after the addition of age and gender on V/F and weight on CL/F; this relationship was, therefore, removed. It should be noted that the effect of weight on CL/F did not render any improvement in OFV when tested univariately without the addition of any other parameter‐covariate relationship. A likely explanation for the strong improvement in OFV is that this relationship refines the weight dependence introduced through CrCL via the Cockcroft–Gault equation. The effects of age and gender on V/F were likewise tested univariately and showed an improvement in OFV. A ΔOFV of −80.8 was observed with the addition of the comedication effect on CL/F. The trend of the comedication effects is consistent with findings in phase I drug–drug interaction studies, which recommend that coadministration of rivaroxaban with strong inhibitors of CYP3A4 and P‐gp (e.g., ketoconazole) should be avoided owing to an increase in exposure and increased risk of bleeding complications.2

In agreement with findings from previous indication‐specific rivaroxaban population PK models,19, 20, 21, 23 and as expected for a Biopharmaceutics Classification System class 2 substance, dose had a highly significant influence on F2; this dependence was captured as a nonlinear function. Indeed, previous studies have demonstrated that, at doses of 15 mg or 20 mg rivaroxaban in the absence of food, bioavailability and absorption rate were less than dose proportional.18, 36 The decreased absorption rate of rivaroxaban at these doses is likely a result of limited aqueous solubility at high vs. low doses.37

The main covariate effects have been consistently characterized in the model. We have shown that renal function has the most significant effect on exposure differences in the subgroups; however, the prediction intervals for AUC0–24, Cmax, and Ctrough are largely overlapping, owing to the variability within one subpopulation. The prediction‐corrected visual predictive checks indicate a slight trend toward an underestimation of Cmax, which is not unexpected and most likely originates from the selection of the structural model, which, in turn, is a consequence of the sparsely sampled PK data and differences between PK sampling schemes in the studies that were used to build this integrated population PK model.18, 19, 21, 22, 23 The CrCL in the current model is truncated according to body surface area to reduce the impact of high, nonphysiological, values of CrCL and extremes of body surface area. In addition, weight has been included in the expression to describe individual rivaroxaban CL/F values. Inclusion of the weight effect accounts for the intrinsic confounding between the two variables because weight is used in the modified Cockcroft–Gault calculation of CrCL.

The inclusion of study/population‐specific factors on CL/F resulted in a distinct improvement of the fit of the model to the data. These factors were 0.85 for AF and 1.14 for ACS, indicating a moderate variation of approximately +/‐ 15% relative to the indication VTE treatment (in which the CL/F factor was set to 1). Most notably, clearance of rivaroxaban increases in the first 3 days after elective hip or knee replacement surgery, up to a study/population specific factor on CL/F of 1.29. The sources of these observed study/population specificities of CL/F cannot be fully elucidated, but they are likely multifactorial. Differences in the overall health status of patients with different indications for rivaroxaban treatment can plausibly affect PK, but study‐specific differences, for example, in the time windows for peak and trough sampling, can also emphasize different parts of the PK profiles and, in turn, influence the estimation of PK parameters, such as CL/F per study or indication.

Strengths and limitations of the study

One of the strengths of our population PK model is the large dataset (almost 5,000 patients for whom the PK data and the complete covariate information was available) from the rivaroxaban phase II/III study program. The covariate relationship is consistent across all patient populations, indicating that the integrated model can be used across indications to predict individual exposure based on covariates and, if available, sparse PK data.

During the model refinement, a qualitative shift of CL/F predictions with concomitant weak, moderate, and strong CYP3A4 inhibitors was seen (Table 2). However, this was smaller than expected based on results from drug‐drug interaction studies (Supplementary Table S2).2 For example, ketoconazole increases the AUC for rivaroxaban 2.6‐fold in drug–drug interaction studies, whereas, in our model, the AUC was increased 1.02‐fold in the presence of strong CYP3A4 inhibitors. This may reflect the small number of patients who received CYP3A4 inhibitors in the phase II/III clinical trials from which data were obtained for the integrated population PK model. Indeed, strong CYP3A4 inhibitors were specifically excluded by the phase II/III program protocols and, accordingly, the proportion of patients who used concomitant strong CYP3A4 inhibitors across the phase II/III program was very low (6 of 5,041 patients). In addition, most patients who took these medications used them for only a short period, although duration of use varied considerably. Furthermore, several comedications may have been administered topically and are, therefore, unlikely to have resulted in significant changes in rivaroxaban clearance. These factors could contribute to the difference in AUC changes seen in this analysis and dedicated drug–drug interaction studies with CYP3A4 inhibitors. The lack of an apparent effect with the strong inhibitors is, therefore, consistent with the conduct of the studies; the small proportion of patients who used concomitant strong CYP3A4 inhibitors translated to a high degree of uncertainty of the point estimate for the effect on CL/F, as shown in Figure 1 b. Coadministration of rivaroxaban with strong inhibitors of CYP3A4 and P‐gp (e.g., ketoconazole) should be avoided,38, 39 owing to increased rivaroxaban exposure and an increased risk of bleeding complications.2

CONCLUSIONS

We have developed an integrated population PK model for rivaroxaban based on the phase II/III PK pool in four indications that can be applied to multiple patient populations. This model has been used for simulations to examine the relationship between covariates and exposure to rivaroxaban and has shown that, as expected, renal function has the most significant effect on exposure. The influence of age and body weight on rivaroxaban PK was minor, demonstrating that fixed doses of rivaroxaban can be prescribed in adult patients without adjustment for age or body weight. The model will be used to estimate individual patient exposure based on covariate information and provides a foundation for further exposure–response analyses based on clinical outcomes in large phase III studies.

Source of Funding

This study was funded by Bayer AG and Janssen Research and Development LLC.

Conflicts of Interest

M.F., D.G., D.K., W.M., A.S., and S.W. are employees of Bayer AG and may own limited stock of Bayer AG. This manuscript was developed within the scope of their employment and no additional payment was received. L.Z. is an employee of Janssen Research and Development and owns stock in Johnson and Johnson. X.Y. was an employee of Janssen Research and Development at the time that the study was carried out. This manuscript was developed within the scope of their employment and no additional payment was received. X.Y. is now an employee of Regeneron Pharmaceuticals, Tarrytown, NY, USA and owns stock in Johnson and Johnson. S.S. is a paid consultant for Bayer AG.

Author Responsibility

The corresponding author, Stefan Willmann, confirms that he had full access to all data in the study and takes final responsibility for the manuscript.

Author Contributions

S.W., L.Z., M.F., D.K., W.M., S.S., A.S., X.Y., and D.G. wrote the manuscript. D.G., D.K., W.M., A.S., S.W., X.Y., and L.Z. designed the research. D.G., A.S., S.W., X.Y., and L.Z. performed the research. S.W., L.Z., M.F., D.K., W.M., S.S., A.S., X.Y., and D.G. analyzed the data.

Supplementary Table S1 Covariates and model properties tested in the integrated population PK model

Click here for additional data file.

Supplementary Table S2 Changes in area under the curve in response to CYP3A4 and P‐gp inhibitors and inducers from the drug–drug interaction studies in the rivaroxaban phase I program analysis and in the integrated population pharmacokinetic model

Click here for additional data file.

Supplementary Figure S1 Graphical representation of the final integrated population pharmacokinetic model

Ka, first–order absorption rate constant; K10, elimination rate constant from the central compartment

Click here for additional data file.

Equation 1. The integrated population pharmacokinetic model is described by a one compartment disposition model with first order absorption and elimination.

Click here for additional data file.

Supplementary Figure S2 Observations vs. population (PRED) and, individual predictions (IPRED), for (a) ODIXa‐Hip2, (b) ODIXa‐OD‐Hip, (c) ODIXa‐Knee, (d) ODIXa‐DVT, (e) EINSTEIN DVT, (f) ROCKET AF, and (g) ACS TIMI‐46.

Click here for additional data file.

Supplementary Figure S3 Conditionally weighted residuals (CWRES) vs. population predictions (PRED) and time after dose (TAD) for (a) ODIXa‐Hip2, (b) ODIXa‐OD‐Hip, (c) ODIXa‐Knee, (d) ODIXa‐DVT, (e) EINSTEIN DVT, (f) ROCKET AF, and (g) ACS TIMI‐46.

Click here for additional data file.

Supplementary Figure S4 Prediction‐corrected visual predictive checks (VPCs) stratified by (a) age, (b) creatinine clearance (CrCL), (c) weight, (d) gender, (e) rivaroxaban dose, and (f) study specific effects.

Click here for additional data file.

Supplementary Figure S5 Simulation results for (a) Cmax and (b) Ctrough at steady state following rivaroxaban 20 mg [left‐hand graphs] and 15 mg once daily [right‐hand graphs]) for the subgroups* according to renal function (top row), age (middle row), and weight (bottom row) for the atrial fibrillation indication.

Click here for additional data file.

Model Code File

Click here for additional data file.