A Digital Twin of Canagliflozin Pharmacokinetics and Pharmacodynamics in Type 2 Diabetes Mellitus

1. Introduction

Type 2 diabetes mellitus (T2DM) is a chronic metabolic disorder and a major global health challenge. Persistent hyperglycemia can lead to severe complications, and many patients do not achieve sustained glycemic control despite the availability of multiple glucose-lowering therapies [1,2]. These challenges highlight the continued need for effective and well-characterized therapeutic strategies.

The kidney is central to glucose homeostasis by reabsorbing nearly all filtered glucose within the proximal tubule [3]. This process is mediated mainly by the sodium–glucose co-transporter 2 (SGLT2, 90%) and to a lesser extent by the sodium–glucose co-transporter 1 (SGLT1, 10%) [4,5]. Inhibition of SGLT2 lowers the renal threshold for glucose (RTG), thereby increasing urinary glucose excretion (UGE) and reducing plasma glucose concentrations [5]. This insulin-independent mechanism has been successfully exploited by SGLT2 inhibitors [6].

Canagliflozin is an orally active SGLT2 inhibitor approved for the treatment of adults with T2DM. Clinically, canagliflozin is administered once daily at doses of 100–300 mg, with dosing guided primarily by renal function and therapeutic response [7]. It is characterized by oral efficacy and sustained glucose-lowering effects.

The pharmacokinetics of canagliflozin are well characterized. Following oral administration, it is rapidly absorbed, reaching peak plasma concentrations within 1–2 h, and exhibits a terminal half-life of approximately 11–15 h [5]. Canagliflozin is predominantly metabolized by uridine 5^′-diphospho-glucuronosyltransferases UGT1A9 and UGT2B4 to the inactive O-glucuronide metabolites M7 and M5, respectively, with a minor oxidative pathway via cytochrome P450 3A4 (CYP3A4) forming the hydroxylated metabolite M9 [5,6]. Several other metabolites have been identified, but these are generally considered to be of low significance and not major human metabolites. Most of the administered dose is recovered in feces as unchanged canagliflozin and M9, whereas urinary excretion occurs primarily in the form of the metabolites M5 and M7 [5,6]. Pharmacodynamically, canagliflozin produces dose-dependent increases in UGE, reflecting inhibition of renal glucose reabsorption.

Patients with T2DM exhibit substantial inter-individual variability in response to canagliflozin. Clinical studies have shown that renal and hepatic function can significantly influence the pharmacokinetics and pharmacodynamics of canagliflozin. Renal impairment increases systemic exposure to canagliflozin and its metabolites, while simultaneously reducing filtered glucose and UGE, leading to diminished glycemic efficacy [8]. Hepatic impairment can alter the formation and elimination of glucuronide metabolites, although effects on parent drug exposure are generally considered moderate in patients with mild to moderate impairment [9]. However, data on patients with severe hepatic impairment remain limited, and the impact on systemic exposure and therapeutic response is not well characterized. Consequently, dose recommendations for canagliflozin in the setting of hepatic dysfunction are less well defined. Together, these pathophysiological factors complicate dose selection in patients with impaired organ function. The impact of different dosing strategies in these populations cannot be comprehensively evaluated in clinical trials alone, highlighting the need for mechanistic approaches to assess drug exposure and response across diverse patient groups [6].

Several computational models of canagliflozin have been published, including population PK and PK/PD analyses, mechanistic PBPK and quantitative systems pharmacology models, and model-based meta-analyses addressing exposure–response relationships and class-wide glycemic effects [10,11,12,13,14,15,16]. While these studies provided valuable insights, they generally addressed specific sub-questions in isolation and did not offer an integrated, physiologically consistent description of canagliflozin pharmacokinetics and pharmacodynamics across populations.

A major limitation of existing canagliflozin models is their limited reusability and reproducibility. Most models are not openly available: the underlying software is closed, executable code is not provided, curated clinical datasets are inaccessible, and models and workflows do not follow FAIR principles. As a consequence, independent reproduction, long–term preservation, and systematic extension of these models are not ensured. In addition, many existing models rely on restricted or narrowly defined datasets, limiting their applicability across patient populations and clinical scenarios. Together, these limitations severely hamper reuse, cross–validation, and cumulative model development. A systematic overview of existing canagliflozin models, including their scope, availability, and reproducibility, is provided in Supplementary Table S1.

Reproducibility is a broader challenge in systems biology [17] and PBPK modeling more generally [18], with many published models lacking accessible equations, data, or executable workflows. To address these challenges, we developed the present model in accordance with FAIR principles [19] and Open Science practices, ensuring that all model files, simulation code, and curated datasets are openly available and reusable.

Physiologically based pharmacokinetic/pharmacodynamic (PBPK/PD) modeling provides a quantitative framework to integrate physiological, biochemical, and drug-specific data. By representing the processes of absorption, distribution, metabolism, and excretion (ADME) in mechanistic detail, such models can simulate the impact of organ dysfunction, and other patient-specific factors on exposure and response, offering insights beyond what can be obtained from clinical trials alone [20,21,22].

Here, we present a fully reproducible digital twin of canagliflozin developed using published clinical data from 22 studies (Table 1), including healthy individuals, patients with T2DM, and individuals with renal or hepatic impairment. The model is a physiologically based whole–body PBPK/PD framework that explicitly represents canagliflozin and its major metabolites M5, M7, and M9, as well as the key organs involved in absorption, metabolism, and elimination, in particular intestine, liver, and kidney, which are linked through the systemic circulation and renal glucose handling.

All clinical data are stored in a curated database, and the model is implemented in SBML [23] to ensure transparency and reproducible simulations. The digital twin was calibrated and evaluated across multiple datasets and applied to systematically investigate the effects of dose, renal impairment, and hepatic impairment on canagliflozin pharmacokinetics and pharmacodynamics. The model reproduces canagliflozin exposure and renal glucose handling across a wide range of clinical scenarios and quantitatively links parent drug and metabolites to pharmacodynamic response under normal and impaired organ function. This framework provides a mechanistic basis for analyzing inter-individual variability in drug response and supports future model–informed dose optimization.

2. Materials and Methods

The digital twin of canagliflozin was developed as a PBPK/PD model through a systematic workflow combining clinical data curation, mechanistic modeling, and in silico simulations. This involved a structured literature search, development an SBML-based model, parameter optimization with selected data subset, and simulation experiments reflecting clinical trial conditions. Pharmacokinetic and pharmacodynamic outcomes were then analyzed to characterize drug disposition and variability across physiological, pathophysiological and prandial states.

2.1. Systematic Literature Research and Data Curation

A systematic literature search was conducted to compile studies that reported pharmacokinetic and pharmacodynamic data of canagliflozin. PubMed was queried on 2024-07-02 using the terms canagliflozin AND pharmacokinetics, and the PKPDAI database was screened in parallel [24]. The clinical studies encompassed a diverse cohort of participants, including healthy volunteers, patients diagnosed with type 1 or type 2 diabetes, and studies investigating renal or hepatic impairment. Pediatric studies, animal studies, and reports lacking sufficient data were excluded from the analysis. The Supplementary Materials Figure S1 provides an overview of the literature review process. Data from eligible studies were curated in the open pharmacokinetics database PK-DB [25]. Information on demographics, disease status, dosing protocols, plasma and urine concentration-time profiles of canagliflozin and its metabolites, and pharmacodynamic outcomes such as UGE were extracted according to established protocols [25]. Figure-based data were digitized using WebPlotDigitizer [26], while tabular and textual data were reformatted into standardized PK-DB formats [25]. The curated dataset, including patient and study characteristics, dosing regimens, and PK/PD outcomes, provided the basis for PBPK/PD model development and is publicly available via PK-DB and within the model files. An overview of the curated studies is available in Table 1.

2.2. Computational Model

The PBPK/PD model was developed in the Systems Biology Markup Language (SBML) [23,27]. Programmatic model construction and visualization were performed using the sbmlutils [28] and cy3sbml [29,30] libraries. Numerical solutions of the underlying ordinary differential equations (ODEs) were obtained with sbmlsim [31], powered by the high-performance SBML simulation engine libRoadRunner [32,33]. The complete model, including simulation scripts and documentation, is available in SBML format under a CC-BY 4.0 license via GitHub (https://github.com/matthiaskoenig/canagliflozin-model) and archived on Zenodo 0.7.0 [34].

The model comprises a whole–body framework with submodels for the intestine, liver, and kidney, linked through the systemic circulation (Figure 1) to characterize canagliflozin absorption, distribution, metabolism, and excretion. The model has a hierarchical structure, with the whole–body model linking the individual organ submodels.

Key ADME processes include intestinal absorption, hepatic glucuronidation by UGT2B4 and UGT1A9 to form the major metabolites M5 and M7, respectively, and minor oxidative metabolism via CYP3A4 to form M9. The kidney submodel captures both metabolic and excretory processes, including additional renal glucuronidation by UGT1A9 and urinary excretion of canagliflozin and its metabolites M5 and M7. The pharmacodynamic component links canagliflozin plasma concentrations to UGE through inhibition of renal glucose reabsorption in the proximal tubule, parameterized by fasting plasma glucose and the RTG. Mathematical descriptions of the submodels and corresponding ordinary differential equations are provided in the Supplementary Materials Section S4.

Fractional organ volumes and blood flows were taken from literature sources [20]. The fractional compartment volumes were set to $F{V}_{\mathrm{gu}}=1.71\%$ for the gut, $F{V}_{\mathrm{ki}}=0.44\%$ for the kidneys, $F{V}_{\mathrm{li}}=2.10\%$ for the liver, and $F{V}_{\mathrm{lu}}=0.76\%$ for the lungs. Fractional blood flows were defined as $F{Q}_{\mathrm{gu}}=18.00\%$ for the gut, $F{Q}_{\mathrm{ki}}=19.00\%$ for the kidneys, $F{Q}_{\mathrm{h}}=21.50\%$ for the hepatic venous outflow, and $F{Q}_{\mathrm{lu}}=100\%$ for the lungs. Absolute organ volumes and blood flows were calculated by scaling the corresponding fractional values with body weight.

Tissue-to-plasma partition coefficients were assumed to be identical across all tissues for canagliflozin and its glucuronide metabolites. Fixed values of $K{p}^{\mathrm{CAN}}=1.00$, $K{p}^{\mathrm{M}5}=11.59$, and $K{p}^{\mathrm{M}7}=217.43$ were used. No tissue-to-plasma partition coefficient was defined for the minor metabolite M9. Transport processes of canagliflozin and its metabolites in the liver and kidneys were modeled explicitly.

Several key factors influencing pharmacokinetic and pharmacodynamic variability were implemented as scaling parameters.

Renal impairment was modeled as a progressive decline in renal function using the factor f_renal, which was applied to glomerular filtration rate (GFR) and renal clearance of parent drug and metabolites. Scaling values were mapped from KDIGO categories: normal (eGFR ≥90 ml/min, $f_{\mathrm{renal}}=1.00$), mild (GFR 50–89 ml/min, $f_{\mathrm{renal}}=0.69$), moderate (GFR 30–49 ml/min, $f_{\mathrm{renal}}=0.32$), severe (GFR ≤30 ml/min, $f_{\mathrm{renal}}=0.19$), derived from KDIGO guidelines and related modeling studies [35,36].

Hepatic impairment was modeled using the scaling factor f_cirrhosis, representing reduced functional liver parenchyma and blood shunting. Scaling values were mapped to Child–Turcotte–Pugh (CPT) classes: A (mild: 5–6 points, f_cirrhosis = 0.40), B (moderate: 7–9 points, f_cirrhosis = 0.70), and C (severe: 10–15 points, f_cirrhosis = 0.80) [37,38,39,40].

Subject- and study-specific physiological and clinical parameters were incorporated when available, including body weight, glomerular filtration rate to adjust renal function, and fasting plasma glucose for the pharmacodynamic component.

Multiple-dose regimens were implemented by stepwise numerical integration between dosing intervals, with dosing events applied according to the study-specific protocols. Oral and intravenous doses were specified using the parameters PODOSE_can and IVDOSE_can, respectively. Simulation time horizons and post-dose sampling windows were selected to match the corresponding clinical study designs.

All simulations were performed deterministically using the optimized parameter set representing the typical (mean) individual. Inter-individual or between-subject variability was not included, as the objective was to evaluate typical pharmacokinetic and pharmacodynamic behavior across studies rather than to perform population-based variability analyses.

2.3. Model Assumptions

Additional details on model equations and assumptions are provided in Supplementary Material S4. The key model assumptions and simplifications are summarized below.

• Canagliflozin absorption was modeled as a first-order process.
• Diurnal variation in plasma glucose concentrations was not modeled explicitly. Instead, a constant fasting plasma glucose (FPG) concentration was assumed and used for the calculation of UGE. When reported, study-specific FPG values were used. Otherwise, FPG values of 5 mM for healthy subjects, 8.0 mM for subjects with T1DM, and 8.0 mM for subjects with T2DM were assumed.
• The RTG was parameterized using parameter optimization, with optimized values reported in Supplementary Table S3.
• Renal filtration and tubular glucose reabsorption were not modelled explicitly. Renal elimination of canagliflozin (CAN) and its metabolites was instead described using first-order processes, depending on kidney volume, renal function (GFR), and compound-specific excretion rate constants. The parameters KI__CANEX_k, KI__M5EX_k and KI__M7EX_k were estimated via parameter optimization.
• The conversion of CAN to M5 by UGT2B4 and to M7 by UGT1A9 in the liver and kidneys was modelled using irreversible Michaelis–Menten kinetics.

2.4. Model Parameterization

A subset of curated data from healthy subjects, patients with T2DM, and individuals with renal impairment following single-dose administration was used to optimize parameters using a local standard optimizer, with the optimal parameters used in the final model. The final optimized parameter set was applied consistently across all subsequent simulations without further study-specific refitting. Optimized parameter values and optimization diagnostics are reported in the Supplementary Materials Section S5, including fitted pharmacokinetic and pharmacodynamic parameter sets as well as convergence behavior and goodness-of-fit assessments for both model components.

The cost function, defined as a function of the parameter vector $\overrightarrow{p}$, minimized the sum of squared, weighted residuals $r_{i,k}$ across all time courses k and data points i. Time courses were weighted by the number of participants in each study $n_k$, and individual time points were weighted by the inverse of the associated measurement uncertainty, represented by the standard deviation ${\sigma }_{i,k}$. This resulted in weights $w_{i,k}=n_k/{\sigma }_{i,k}$.

$$F\left(\overrightarrow{p}\right)=0.5\sum _{i,k}{(w_{i,k}·r_{i,k}\left(\overrightarrow{p}\right))}^2,$$

The weighting by uncertainty weights data points with smaller uncertainty higher, the weighting by participants weights data points depending on the number of subjects.

Multiple optimization runs ($n=100$) were performed using different initial parameter values. Optimization was conducted sequentially: pharmacokinetic parameters were estimated first, followed by pharmacodynamic parameters.

2.5. Simulations

For each curated clinical study (Table 1), a corresponding in silico experiment was implemented to reproduce the reported dosing regimen and study conditions. Parameters for oral and intravenous dosing, prandial state, bodyweight, fasting plasma glucose, and renal or hepatic function were adjusted according to study-specific information. Multiple dosing protocols were incorporated where applicable.

To further explore sources of variability, simulation experiments and parameter scans were performed across physiologically relevant ranges for renal function, hepatic function, and dose to enable systematic evaluation of the influence of key physiological and pathophysiological parameters on PK/PD outcomes.

2.6. Pharmacokinetic and Pharmacodynamic Parameters

Pharmacokinetic parameters of canagliflozin and its metabolites were calculated from plasma concentration-time curves and urinary excretion using standard non-compartmental methods. Pharmacodynamic outcomes were evaluated in terms of UGE, calculated from the simulated plasma concentration-time courses in combination with fasting plasma glucose and the renal threshold for glucose. Simulated profiles and derived pharmacokinetic and pharmacodynamic parameters were compared against the curated clinical data.

2.7. Sensitivity Analysis

The influence of model parameters on pharmacokinetic (PK) and pharmacodynamic (PD) outcomes was assessed using sensitivity analysis. A reference simulation corresponding to a single oral dose of 200 mg canagliflozin was used. PK readouts comprised the area under the concentration-time curve (AUC), maximum concentration (C_max), half-life, volume of distribution (V_d), clearance (CL), and elimination rate constant ($k_{\mathrm{el}}$) for canagliflozin, as well as AUC, C_max, and half-life for its metabolites M5 and M7. PD readout was urinary glucose excretion (UGE) at 24 hours. Parameters representing physical constants, unit conversion factors, and dosing parameters were excluded from the analysis. Results are reported in Supplementary Materials, Section S7.

2.7.1. Sampling-Based Sensitivity Analysis

To quantify uncertainty in PK and PD outcomes arising from variability in model parameters, a sampling-based uncertainty analysis was performed. Model parameters were sampled uniformly within their predefined bounds using Latin hypercube sampling (LHS). A total of $n=1000$ simulations were conducted, and the resulting distributions of PK and PD readouts were used to characterize parameter-induced uncertainty.

2.7.2. Local Sensitivity Analysis

Local sensitivities were computed by perturbing each parameter $p_i$ individually by $\pm 1\%$ relative to its reference value $p_{i,0}$. Sensitivities were calculated using a symmetric midpoint approximation,

$$S(q_k,p_i)={\displaystyle \frac{q_k\left(p_i^{+}\right)-q_k\left(p_i^{-}\right)}{p_i^{+}-p_i^{-}}},$$

where $p_i^{\pm }=p_{i,0}·(1\pm 0.01)$. Sensitivities were subsequently normalized to obtain dimensionless measures,

$$S_{\mathrm{norm}}(q_k,p_i)={\displaystyle \frac{q_k\left(p_i^{+}\right)-q_k\left(p_i^{-}\right)}{p_i^{+}-p_i^{-}}}·{\displaystyle \frac{p_{i,0}}{q_k\left(p_{i,0}\right)}},$$

representing the relative change in the model output per relative change in the parameter. Normalized sensitivities with absolute values below 0.1, as well as parameters without measurable effects on any PK or PD readout, were omitted from the heatmap visualization. For clarity, sensitivity matrices were hierarchically clustered by model parameters using single-linkage clustering.

2.7.3. Global Sensitivity Analysis

Global sensitivity analysis was performed using variance-based Sobol indices [41,42,43], as implemented in SALib [44,45]. First-order (S₁) and total-effect (S_T) indices were computed. Parameters with S₁ and S_T values below 0.05, as well as parameters without measurable effects on any readout, were omitted from the heatmap visualization. Samples were generated using Saltelli’s extension of the Sobol’ sequence with $n=4096$ base samples.

3. Results

3.1. Canagliflozin Database

An extensive database of canagliflozin’s pharmacokinetics and pharmacodynamics was created to develop and validate the model. A systematic literature search initially yielded 152 records. After screening according to predefined inclusion and exclusion criteria, 22 studies were selected for detailed curation and formed the core dataset used to evaluate the PBPK/PD model. See Supplementary Materials Figure S1 for the study selection process. All curated pharmacokinetic and pharmacodynamic data of the 22 studies are publicly available in the model files and in the PK-DB database with unique study identifiers as referenced in the manuscript (Table 1).

3.2. Computational Model

Using the curated dataset, a digital twin in the form of a PBPK/PD model was developed to describe the ADME of canagliflozin and its primary metabolites, as well as its pharmacodynamic effect on UGE (Figure 1). The complete model, including simulation scripts and documentation, is available in SBML format under a CC-BY 4.0 license via GitHub (https://github.com/matthiaskoenig/canagliflozin-model) and archived on Zenodo 0.7.0 [34].

Simulations followed the respective clinical study designs, accounting for single- and multiple-dose regimens, and renal or hepatic impairment. All simulations are provided in Supplementary Materials Section S6 (Figures S7–S67).

The final parameter set is summarized in the Supplementary Materials Tables S2 and S3. Parameter convergence and goodness-of-fit metrics are shown in the Supplementary Materials Figure S5 and S6. All submodel visualizations are provided in Supplementary Materials Section S3 (Figures S2-S4), and all model equations and ODEs are given in Supplementary Materials Section S4.

3.3. Dose Dependency

The pharmacokinetic and pharmacodynamic effects of canagliflozin over an oral dose range of 0–800 mg are shown in Figure 2 and Figure 3. With increasing doses, plasma concentrations of canagliflozin and its primary metabolites increased, as did the amounts excreted in urine and feces. Canagliflozin reached peak plasma levels after approximately 3.5 h and has a half-life of approximately 8 h, whereas the metabolites showed a slightly delayed peak and elimination. The parameter scan showed a clear dose-dependent rise in exposure metrics ($AU{C}_{0-inf}$ and $C_{\max }$), while $T_{\max }$ and half-lives remained largely unchanged. Higher doses of canagliflozin were associated with lower RTG and a nonlinear increase in urinary glucose excretion.

Time-course simulations were performed for all curated clinical dose-dependency studies (Chen2015 [46], Devineni2012 [47], Devineni2013 [48], Devineni2015a [49], Devineni2016 [50], Iijima2015 [51], Sha2011 [52], Sha2014 [53], Tamborlane2018 [54]). Simulated canagliflozin, M5 and M7 plasma concentrations, urinary excretion, and UGE are shown for both single dose and multiple dose regimens.

3.4. Renal Impairment

The impact of renal impairment on canagliflozin disposition and pharmacodynamics is shown in Figure 4. Simulations were performed for four renal function groups (normal, mild, moderate, severe). Plasma concentrations of canagliflozin were minimally affected by renal dysfunction, whereas exposure to its main metabolites increased with declining renal function. Urinary excretion of canagliflozin, M5 and M7 decreased with impairment, while fecal excretion remained unchanged. Across renal function states, parent $AU{C}_{0-inf}$ and $C_{max}$ changed modestly, while metabolite exposure increased with impairment.

The pharmacodynamic response showed a reduction in UGE with declining renal function, while RTG changed slightly. In severe impairment, simulated UGE was clearly lower than under normal function, whereas RTG was only minimally affected.

Simulated parent and metabolite plasma time courses and urinary excretion profiles across renal function groups are shown together with clinical data available from two studies (Devineni2015c [9] and Inagaki2014 [55]). Fecal excretion data were not available for direct comparison.

3.5. Hepatic Impairment

The effect of hepatic impairment on canagliflozin pharmacokinetics and pharmacodynamics is summarized in Figure 5. Simulations were performed for normal liver function and for mild, moderate, and severe cirrhosis. With increasing hepatic impairment, canagliflozin plasma exposure increased, reflected by higher $C_{\max }$ values, increased $AU{C}_{0-inf}$, and a prolonged plasma time course. Urinary excretion of canagliflozin increased with higher plasma exposure, while fecal excretion of the parent drug remained unchanged across hepatic function states.

For the metabolite M5, plasma concentrations decreased with increasing hepatic impairment, with reductions in both $C_{\max }$ and AUC. In contrast, M7 showed a reduced peak concentration but a prolonged exposure, resulting in lower $C_{\max }$ values and increased AUC with increasing severity of cirrhosis. Time to maximum concentration ($T_{\max }$) was stable for canagliflozin, M5, and M7 across all hepatic impairment categories. Parameter scans showed an increase in half-life for the parent drug and both metabolites with increasing hepatic impairment, without a linear relationship to cirrhosis severity.

Pharmacodynamic effects were minimal. UGE remained largely unchanged across hepatic function states. Simulated plasma concentration time courses and urinary excretion profiles for canagliflozin, M5, and M7 under different degrees of hepatic impairment are shown together with clinical data from Devineni2015c [9].

Table 1. Summary of studies for modeling. Overview of study identifiers, PK-DB IDs, PMID, route, dosing, and subject characteristics, including health status, renal impairment, hepatic impairment and T2DM.

Study	PK-DB	PMID	Route	Dosing	Dose [mg]	Healthy	Renal impairment	Hepatic impairment	T2DM
Chen2015 [46]	PKDB00858	26048186	PO	single	100, 300	✓
Devineni2012 [47]	PKDB00874	22226086	PO	single, multi	100, 300				✓
Devineni2013 [48]	PKDB00875	23670707	PO	single, multi	50, 100, 300				✓
Devineni2014 [56]	PKDB00879	24726680	PO	multi	300	✓
Devineni2015 [57]	PKDB00880	27140803	PO	single	300	✓
Devineni2015a [49]	PKDB00876	25500487	PO	single, multi	50, 100, 300	✓
Devineni2015b [58]	PKDB00877	27136910	PO, IV	single	0.010, 300	✓
Devineni2015c [9]	PKDB00881	25659911	PO	single	300	✓	✓	✓
Devineni2015d [59]	PKDB00882	25407255	PO	single, multi	300	✓
Devineni2015e [60]	PKDB00883	27136908	PO	single	50, 100, 300	✓
Devineni2016 [50]	PKDB00884	26687552	PO	single	200, 300	✓
Iijima2015 [51]	PKDB00885	26280756	PO	single, multi	25, 100, 200, 400				✓
Inagaki2014 [55]	PKDB00878	25200141	PO	single	100, 200		✓		✓
Kinoshita2015 [61]	PKDB00886	25424014	PO	single, multi	200	✓
Mamidi2014 [62]	PKDB00887	24568888	PO	single	188	✓
Mohamed2019 [63]	PKDB00888	31384829	PO	single	50	✓
Murphy2015 [64]	PKDB00857	25546166	PO	single	150	✓
Sha2011 [52]	PKDB00889	21457428	PO	single, multi	10, 30, 100, 200, 400, 600, 800	✓
Sha2014 [53]	PKDB00890	25166023	PO	single, multi	30, 100, 200, 300, 400				✓
Sha2015 [65]	PKDB00891	25421015	PO	multi	300	✓
Tamborlane2018 [54]	PKDB00892	29271103	PO	multi	100, 300				✓
Wattamwar2020 [66]	PKDB01127	32763846	PO	multi	150	✓

Table 1. Summary of studies for modeling. Overview of study identifiers, PK-DB IDs, PMID, route, dosing, and subject characteristics, including health status, renal impairment, hepatic impairment and T2DM.

Study	PK-DB	PMID	Route	Dosing	Dose [mg]	Healthy	Renal impairment	Hepatic impairment	T2DM
Chen2015 [46]	PKDB00858	26048186	PO	single	100, 300	✓
Devineni2012 [47]	PKDB00874	22226086	PO	single, multi	100, 300				✓
Devineni2013 [48]	PKDB00875	23670707	PO	single, multi	50, 100, 300				✓
Devineni2014 [56]	PKDB00879	24726680	PO	multi	300	✓
Devineni2015 [57]	PKDB00880	27140803	PO	single	300	✓
Devineni2015a [49]	PKDB00876	25500487	PO	single, multi	50, 100, 300	✓
Devineni2015b [58]	PKDB00877	27136910	PO, IV	single	0.010, 300	✓
Devineni2015c [9]	PKDB00881	25659911	PO	single	300	✓	✓	✓
Devineni2015d [59]	PKDB00882	25407255	PO	single, multi	300	✓
Devineni2015e [60]	PKDB00883	27136908	PO	single	50, 100, 300	✓
Devineni2016 [50]	PKDB00884	26687552	PO	single	200, 300	✓
Iijima2015 [51]	PKDB00885	26280756	PO	single, multi	25, 100, 200, 400				✓
Inagaki2014 [55]	PKDB00878	25200141	PO	single	100, 200		✓		✓
Kinoshita2015 [61]	PKDB00886	25424014	PO	single, multi	200	✓
Mamidi2014 [62]	PKDB00887	24568888	PO	single	188	✓
Mohamed2019 [63]	PKDB00888	31384829	PO	single	50	✓
Murphy2015 [64]	PKDB00857	25546166	PO	single	150	✓
Sha2011 [52]	PKDB00889	21457428	PO	single, multi	10, 30, 100, 200, 400, 600, 800	✓
Sha2014 [53]	PKDB00890	25166023	PO	single, multi	30, 100, 200, 300, 400				✓
Sha2015 [65]	PKDB00891	25421015	PO	multi	300	✓
Tamborlane2018 [54]	PKDB00892	29271103	PO	multi	100, 300				✓
Wattamwar2020 [66]	PKDB01127	32763846	PO	multi	150	✓

4. Discussion

In this study, we established a comprehensive clinical dataset of canagliflozin pharmacokinetics and pharmacodynamics and used it to develop a mechanistic PBPK/PD digital twin. In total, 22 clinical trials were curated, covering a broad spectrum of dosing regimens and study populations, including healthy individuals, patients with type 2 diabetes, and cohorts with renal and hepatic impairment. Overall, data availability was sufficient for model development and evaluation, with consistent reporting of plasma concentration time courses, urinary excretion, and UGE across both single dose and multiple dose designs. Fecal excretion data and RTG related endpoints were sufficient to support model calibration and validation.

The PBPK/PD framework integrates systemic distribution, hepatic metabolism, and renal elimination into a coherent representation of canagliflozin disposition and pharmacodynamic effect. Robust parameter optimization and good agreement between simulations and observed data across a wide range of studies support the reliability of the model. A key strength of this mechanistic approach is its ability to investigate scenarios that are difficult to address systematically in clinical trials, such as direct comparisons between renal and hepatic impairment or controlled exploration of dose dependent behavior across physiological states. Pharmacodynamic effects are driven by plasma canagliflozin concentrations through inhibition of SGLT2, resulting in a modulation of RTG and a mechanistic, physiologically interpretable link between exposure and UGE.

The simulations reproduced the expected dose dependent pharmacokinetic and pharmacodynamic behavior of canagliflozin, with increasing doses leading to higher plasma exposure and greater UGE. Across the investigated dose range, time to maximum concentration and elimination half life remained largely unchanged, indicating dose proportional kinetics within the therapeutic window. The modeled dose dependent increase in UGE reflects progressive inhibition of renal glucose reabsorption and is consistent with the observed nonlinear relationship between exposure and pharmacodynamic response reported in clinical studies.

The pharmacokinetics and pharmacodynamics of canagliflozin under conditions of hepatic impairment are less extensively characterized than in other clinical populations. Evaluation in this setting was constrained by the availability of clinical data, as only one study reported plasma concentration–time profiles. Corresponding pharmacodynamic measurements, including UGE and RTG, were not reported, which limited direct assessment of pharmacodynamic effects. In addition, the available plasma concentration–time profiles exhibited substantial variability, complicating quantitative comparisons between simulations and observations.

Clinical observations under hepatic impairment showed little to no change in parent canagliflozin plasma exposure across degrees of liver dysfunction. In contrast, metabolite behavior differed, with increased plasma exposure of M7, while M5 profiles were noisy and did not exhibit a consistent systematic trend. Model simulations predicted increased parent canagliflozin exposure with worsening hepatic impairment, reflecting reduced hepatic metabolic capacity. In parallel, simulated M5 exposure decreased due to diminished hepatic formation, whereas M7 showed prolonged exposure driven by reduced clearance and the additional formation pathway in the kidney. As a result, the clinical parent drug exposure was not fully reproduced by the model, and agreement for metabolite profiles remained qualitative and limited by data variability.

Renal metabolism was represented as part of the mechanistic model structure. Published studies report renal formation of M7, and this pathway was therefore implemented in the kidney submodel. Inclusion of renal M7 metabolism improved agreement between simulations and clinical data, particularly for urinary excretion and, to a lesser extent, for plasma concentration–time profiles. In contrast, renal metabolism of M5 was not implemented, consistent with reports indicating negligible renal conversion despite expression of the corresponding UGT enzymes in kidney tissue [67].

Despite these considerations, discrepancies between simulations and clinical observations persisted under hepatic impairment. These findings suggest that additional processes may contribute to canagliflozin disposition, for example altered renal handling, extrahepatic metabolism, or compensatory clearance pathways that are not captured by the current model.

Renal impairment had a pronounced influence on the pharmacodynamics of canagliflozin and represented an important source of variability in treatment response. Evaluation in this population was supported by two clinical studies reporting both pharmacokinetic and pharmacodynamic measurements. Across levels of renal impairment, parent canagliflozin plasma exposure changed only modestly, whereas exposure to the major metabolites increased with declining renal function. This pattern is consistent with reduced renal elimination of the metabolites, indicating that systemic exposure to the parent drug is relatively resistant to moderate reductions in renal function. The PBPK/PD model reproduced these trends, supporting the mechanistic representation of renal clearance pathways.

Unlike the limited pharmacokinetic changes, renal impairment strongly affected the pharmacodynamic effects. Both clinical data and simulations revealed a significant decrease in UGE as renal function declined, despite relatively stable parent drug concentrations. These results highlight the dependence of the pharmacodynamic response on renal glucose filtration and tubular handling, emphasizing the critical role of kidney function in determining the efficacy of SGLT2 inhibition.

Agreement between simulations and clinical observations under renal impairment was generally stronger than for hepatic impairment. Although the available studies were restricted to single–dose regimes, the presence of plasma, urinary, and pharmacodynamic endpoints allowed evaluation of both exposure and response. The consistent reproduction of these trends across renal function categories supports the validity of the renal impairment implementation within the PBPK/PD framework.

Future work should focus on improving the breadth and consistency of clinical data available for evaluating canagliflozin under organ impairment. In particular, more comprehensive pharmacokinetic and pharmacodynamic reporting under hepatic impairment is needed to enable mechanistic interpretation beyond plasma exposure alone. Standardized measurement of UGE and RTG would substantially strengthen assessment of pharmacodynamic behavior in this setting.

Although the model performed more robustly under renal impairment, the available evidence remains limited. A wider range of renal function datasets, including multiple-dose regimens and aligned pharmacokinetic and pharmacodynamic endpoints, would further support the refinement and validation of the model. Systematic reporting across impairment studies would enhance the utility of mechanistic PBPK/PD models and improve the general understanding of drug behavior and therapeutic response in patients with altered organ function.

Beyond the pharmacological findings, this study also addresses a broader challenge in computational pharmacology and systems biology related to model transparency and reproducibility. A substantial fraction of published PBPK models cannot be independently reproduced because essential components, such as model equations, executable code, and curated calibration datasets, are not publicly available [17,18]. This limits independent verification and substantially restricts reuse, extension, and cumulative model development within the community.

To address these limitations, reproducibility and accessibility were treated as core design principles in the present work. The complete modeling framework, including the SBML model, simulation scripts, and curated clinical datasets, is openly available and structured in accordance with FAIR principles [19,68]. This enables independent reproduction, transparent evaluation, and systematic reuse of the model across different contexts. All resources are released under permissive MIT and CC-BY licenses, lowering barriers to reuse in both academic and industrial settings and allowing straightforward integration into existing workflows, including commercial applications.

In summary, this PBPK/PD digital twin of canagliflozin integrates diverse clinical data into a mechanistic framework that captures key pharmacokinetic and pharmacodynamic behavior across dosing regimens and patient populations. The model provides quantitative insight into dose dependency and the effects of organ impairment, supporting analysis of variability in drug response under clinically relevant conditions. By providing full open access to the model, simulation code, and curated datasets, this work establishes a transparent and reproducible reference framework that supports independent validation, reuse, and future PBPK/PD model development.