Building a Compartmental PK Model

Author

Vijay Ivaturi

1 Introduction

This tutorial walks through the process of building a pharmacokinetic (PK) model of the combined PKPD model for Warfarin that was introduced in the previous tutorial. It is assumed that the one-compartment model with first-order absorption and absorption lag, inter-individual variability on clearance (\(CL)\), central volume of distribution (\(V_c\)), and absorption half-life (\(t_{abs}\)), allometric scaling on PK parameters, and a combined proportional and additive residual error model has already been established. Conceptually, the structural model can be summarized as follows:

  1. Depot (absorption compartment)
  2. Central (systemic)

2 Breakdown of the Model

2.1 Model Structure

A schematic representation might look like this:

G Depot Depot Central Central Depot->Central Ka Nothing Central->Nothing CL

  • Depot is where the drug is administered (e.g., via oral dosing).
  • Central is the primary systemic compartment with clearance \(CL\).
  • \(Ka\) is the first-order absorption rate constant from the Depot to the Central compartment.

Allometric (size-based) scaling is often applied to parameters such as clearance (\(CL\)) and volume (\(V\)) to account for patient size differences or other covariate information. In this tutorial, scaling factors provided by covariates (e.g., FSZCL for clearance and FSZV for volume) will serve that role.

2.2 Writing the Differential Equations

Given the two-compartment PK structure and first-order absorption, the ordinary differential equations (ODEs) can be written as follows:

\[ \begin{aligned} \frac{d(\text{Depot})}{dt} &= - K_a \times \text{Depot}, \\[6pt] \frac{d(\text{Central})}{dt} &= K_a \times \text{Depot} \;-\; \frac{CL}{V_c} \times \text{Central} \end{aligned} \]

Here:

  • \(K_a = \log(2)/t_{abs}\) is derived from an absorption half-life \(t_{abs}\).
  • \(CL\) (clearance), \(V_c\) (central volume of distribution) can be made functions of covariates and random effects.
  • Allometric scaling is handled by multiplying each base population parameter by the appropriate factor (e.g., FSZCL or FSZV) and exponentiating random effects if desired.

2.3 Converting to a @dynamics Block

In Pumas, the above system of ODEs can be written directly in the @dynamics block. For example:

@dynamics begin
    Depot' = -Ka * Depot
    Central' = Ka * Depot - (CL / Vc) * Central
end

This block automatically interprets each line as one of the state equations (e.g. \(\frac{d(\text{state})}{dt} = \dots\)).

2.4 Using the Coefficients in the @pre Block

Before these equations are solved, the model needs to specify how each parameter is computed from the population values, covariates, and random effects. This is typically performed in a @pre block. Each parameter in the ODE system (e.g., \(CL\), \(V_c\), \(K_a\)) is expressed in terms of:

  • Population parameters (e.g., pop_CL, pop_Vc, pop_tabs)
  • Covariates (e.g., FSZCL, FSZV)
  • Random effects (e.g., exp(pk_η[i]))

An example set of algebraic expressions is shown below:

@pre begin
    # Example of allometric scaling using covariates:
    CL = pop_CL * FSZCL * exp(pk_η[1])
    Vc = pop_Vc * FSZV * exp(pk_η[2])

    # Absorption half-life and rate constant
    tabs = pop_tabs * exp(pk_η[3])
    Ka = log(2) / tabs
end

In this illustration:

  • FSZCL and FSZV are the covariate-based scaling factors.
  • pop_CL, pop_Vc, and pop_tabs are the typical population values.
  • pk_η[1], pk_η[2], and pk_η[3] are the individual random effects for \(CL\), \(V_c\), and \(t_{abs}\), respectively.

2.5 Specifying the Lag Time in the @dosecontrol Block

The lag time is specified in the @dosecontrol block. For example:

@dosecontrol begin
    lags = (Depot = pop_lag,)
end

2.6 Specifying Population Parameters in the @param Block

These population parameters must be declared so that the modeling environment recognizes them, along with any initial guesses, bounds, or descriptions. Below is a @param block that focuses on the PK parameters:

@param begin
    """
    Clearance (L/h/70kg)
    """
    pop_CL ∈ RealDomain(lower = 0.0, init = 0.134)

    """
    Central Volume (L/70kg)
    """
    pop_Vc ∈ RealDomain(lower = 0.0, init = 8.11)
    """
    Absorption time (h)
    """
    pop_tabs ∈ RealDomain(lower = 0.0, init = 0.523)

    """
    Lag time (h)
    """
    pop_lag ∈ RealDomain(lower = 0.0, init = 0.1)

    # Inter-individual variability (diagonal covariance for CL, Vc, tabs)
    """
      - ΩCL
      - ΩVc
      - ΩTabs
    """
    pk_Ω ∈ PDiagDomain([0.01, 0.01, 0.01])

    # Residual variability for concentration
    """
    Proportional error
    """
    σ_prop ∈ RealDomain(lower = 0.0, init = 0.00752)

    """
    Additive error
    """
    σ_add ∈ RealDomain(lower = 0.0, init = 0.0661)
end
  • Each parameter has a lower bound (e.g., 0.0) and an initial estimate (init) used during model fitting.
  • pk_Ω captures the variances for the three random effects on \(CL\), \(V_c\), and \(t_{abs}\).

2.7 Specifying the Random Effects in the @random Block

The random effects are drawn from specified distributions in the @random block. For instance:

@random begin
    # mean = 0, variances = pk_Ω
    pk_η ~ MvNormal(pk_Ω)
end
  • pk_η is a vector of random effects for the PK parameters (CL, Vc, tabs, etc.).

During model simulation or estimation, each individual in the data set will have its own sample from these distributions, thereby introducing inter-individual variability.

2.8 Constructing the @derived Block for Predicted Concentration

The final step is to specify the model’s predictions for observed data—in this case, the drug concentration in plasma. The predicted concentration (cp) is calculated by dividing the mass in the central compartment by the central volume of distribution. An associated @vars block might look like this:

@vars begin
    cp := Central / Vc
end

This tells the model framework how to calculate cp from the states and parameters.

A residual error model is then applied to link cp to the observed concentration. For example:

@derived begin
    """
    Warfarin Concentration (mg/L)
    """
    conc ~ @. Normal(cp, sqrt((σ_prop * cp)^2 + σ_add^2))
end
  • cp is the predicted concentration in the central compartment.
  • σ_prop and σ_add capture proportional and additive measurement error, respectively.

2.9 Putting It All Together

Below is a concise example model that incorporates the steps above. All PD components from the original warfarin model have been removed, leaving only the one-compartment PK system (with a depot for absorption), IIV on \(CL\), \(V_c\), and \(t_{abs}\), allometric scaling factors (FSZCL, FSZV), and a combined proportional and additive residual variability for concentrations:

warfarin_pk_model = @model begin

    @metadata begin
        desc = "Warfarin 1-compartment PK model (PD removed)"
        timeu = u"hr"
    end

    @param begin
        # PK parameters
        pop_CL ∈ RealDomain(lower = 0.0, init = 0.134)
        pop_Vc ∈ RealDomain(lower = 0.0, init = 8.11)
        pop_tabs ∈ RealDomain(lower = 0.0, init = 0.523)
        pop_lag ∈ RealDomain(lower = 0.0, init = 0.1)

        # Inter-individual variability
        pk_Ω ∈ PDiagDomain([0.01, 0.01, 0.01])

        # Residual variability
        σ_prop ∈ RealDomain(lower = 0.0, init = 0.00752)
        σ_add ∈ RealDomain(lower = 0.0, init = 0.0661)
    end

    @random begin
        # mean = 0, covariance = pk_Ω
        pk_η ~ MvNormal(pk_Ω)
    end

    @covariates begin
        # Scaling factors
        FSZCL
        FSZV
    end

    @pre begin
        CL = FSZCL * pop_CL * exp(pk_η[1])
        Vc = FSZV * pop_Vc * exp(pk_η[2])

        tabs = pop_tabs * exp(pk_η[3])
        Ka = log(2) / tabs
    end

    @dosecontrol begin
        lags = (Depot = pop_lag,)
    end

    @vars begin
        cp := Central / Vc
    end

    @dynamics begin
        Depot' = -Ka * Depot
        Central' = Ka * Depot - (CL / Vc) * Central
    end

    @derived begin
        conc ~ @. Normal(cp, sqrt((σ_prop * cp)^2 + σ_add^2))
    end
end
  1. The @param block lists all population-level parameters.
  2. The @random block defines how random effects are distributed.
  3. The @covariates block provides covariate inputs for allometric scaling.
  4. The @pre block transforms population parameters and random effects into individual-level parameters (\(CL, V_c, K_a\)).
  5. The @dosecontrol block defines the lag time for the absorption compartment.
  6. The @vars block defines the predicted concentration.
  7. The @dynamics block contains the ODEs for the one-compartment PK model.
  8. The @derived block generates the predicted concentration conc with a specified residual error model.

Through these steps, the tutorial demonstrates how to move from known compartmental equations to a fully specified model that can be used for fitting or simulation in a population PK context.