TMDD Part 3: Model Approximations (QE, QSS, Michaelis-Menten)

Author

Vijay Ivaturi

1 Introduction

This is Part 3 of our comprehensive TMDD tutorial series. Here we explore the practical approximations of the full TMDD model that are commonly used when the full model is either overparameterized or when data are insufficient for estimating all parameters.

Tutorial Series Overview:

  1. Part 1: Understanding TMDD - Conceptual foundations and parameter meaning 1b. Part 1b: TMDD Dynamics Through the Four Phases - Deep dive into phase characterization
  2. Part 2: The Full TMDD Model - Complete mechanistic model implementation
  3. Part 3: TMDD Approximations (this tutorial) - QE, QSS, and Michaelis-Menten models
  4. Part 4: Special Cases and Extensions - Advanced TMDD scenarios
  5. Part 5: Practical Workflows - Fitting, diagnostics, and model selection

2 Learning Objectives

By the end of this tutorial, you will be able to:

  • Explain why TMDD model simplifications are necessary
  • Implement Quasi-Equilibrium (QE), Quasi-Steady-State (QSS), and Michaelis-Menten models in Pumas
  • Understand the mathematical basis and assumptions behind each approximation
  • Select the appropriate approximation based on your data and scientific context
  • Compare model predictions across different approximation levels

3 Libraries

using Pumas
using PumasUtilities
using AlgebraOfGraphics
using CairoMakie
using Random
using DataFrames

4 Why Simplify the Full TMDD Model?

The full TMDD model contains nine kinetic parameters: \(CL\), \(V_c\), \(Q\), \(V_p\), \(k_{syn}\), \(k_{deg}\), \(k_{on}\), \(k_{off}\), and \(k_{int}\). In practice, estimating all these parameters is often challenging because:

  1. Parameter non-identifiability: With typical PK data (free ligand concentrations only), several parameters are structurally or practically non-identifiable.

  2. Correlation between parameters: \(k_{on}\) and \(k_{off}\) are often highly correlated, making individual estimation unreliable.

  3. Data limitations: Receptor and complex concentrations are rarely measured, limiting the information available for parameter estimation.

  4. Numerical instability: The full model can have stiff ODEs that are computationally expensive and prone to convergence issues.

Gibiansky et al. (2008) provided critical insights into when simplifications are appropriate and how to select the right approximation.

5 The Approximation Hierarchy

The TMDD model approximations form a hierarchy from most mechanistic to most simplified:

Full TMDD → Quasi-Steady-State (QSS) → Quasi-Equilibrium (QE) → Michaelis-Menten (MM)
         ↑                           ↑                        ↑
     Most general               Assumes rapid             Assumes receptor
                              binding kinetics          concentration negligible

Each approximation makes additional assumptions that reduce model complexity but restrict applicability.

6 The Full TMDD Model (Reference)

For comparison, here is the one-compartment full TMDD model:

\[ \frac{dL}{dt} = -k_{el} \cdot L - k_{on} \cdot L \cdot R + k_{off} \cdot LR \]

\[ \frac{dR}{dt} = k_{syn} - k_{deg} \cdot R - k_{on} \cdot L \cdot R + k_{off} \cdot LR \]

\[ \frac{dLR}{dt} = k_{on} \cdot L \cdot R - k_{off} \cdot LR - k_{int} \cdot LR \]

Where:

  • \(L\) = Free ligand concentration
  • \(R\) = Free receptor concentration
  • \(LR\) = Ligand-receptor complex concentration
  • \(k_{el}\) = Linear elimination rate constant
  • \(k_{on}\), \(k_{off}\) = Association and dissociation rate constants
  • \(k_{syn}\), \(k_{deg}\) = Receptor synthesis and degradation rates
  • \(k_{int}\) = Complex internalization rate
tmdd_full_model = @model begin
    @param begin
        tvCl  RealDomain(lower = 0)
        tvVc  RealDomain(lower = 0)
        tvKsyn  RealDomain(lower = 0)
        tvKdeg  RealDomain(lower = 0)
        tvKon  RealDomain(lower = 0)
        tvKoff  RealDomain(lower = 0)
        tvKint  RealDomain(lower = 0)
    end

    @pre begin
        Cl = tvCl
        Vc = tvVc
        Ksyn = tvKsyn
        Kdeg = tvKdeg
        Kon = tvKon
        Koff = tvKoff
        Kint = tvKint
    end

    @dosecontrol begin
        bioav = (Ligand = 1 / tvVc,)
    end

    @init begin
        Receptor = Ksyn / Kdeg
    end

    @vars begin
        Kel := Cl / Vc
        R0 := Ksyn / Kdeg
        Ltot := Ligand + Complex
        Rtot := Receptor + Complex
    end

    @dynamics begin
        Ligand' = -Kel * Ligand - Kon * Ligand * Receptor + Koff * Complex
        Receptor' = Ksyn - Kdeg * Receptor - Kon * Ligand * Receptor + Koff * Complex
        Complex' = Kon * Ligand * Receptor - Koff * Complex - Kint * Complex
    end

    @derived begin
        Free_Ligand = @. Ligand
        Free_Receptor = @. Receptor
        Complex_Conc = @. Complex
        Total_Ligand = @. Ltot
        Total_Receptor = @. Rtot
    end
end
PumasModel
  Parameters: tvCl, tvVc, tvKsyn, tvKdeg, tvKon, tvKoff, tvKint
  Random effects:
  Covariates:
  Dynamical system variables: Ligand, Receptor, Complex
  Dynamical system type: Nonlinear ODE
  Derived: Free_Ligand, Free_Receptor, Complex_Conc, Total_Ligand, Total_Receptor
  Observed: Free_Ligand, Free_Receptor, Complex_Conc, Total_Ligand, Total_Receptor

7 Quasi-Equilibrium (QE) Approximation

7.1 Mathematical Basis

The QE approximation assumes that the binding and unbinding of ligand to receptor is much faster than all other processes (elimination, distribution, internalization). Mathematically, this means:

\[ k_{on} \cdot L \cdot R \approx k_{off} \cdot LR \]

This implies that the binding reaction is at equilibrium at all times, and we can define the equilibrium dissociation constant:

\[ K_D = \frac{k_{off}}{k_{on}} = \frac{L \cdot R}{LR} \]

7.2 When is QE Appropriate?

The QE approximation is valid when:

  • \(k_{on}\) and \(k_{off}\) are both very large compared to \(k_{int}\), \(k_{deg}\), and \(k_{el}\)
  • Specifically, \(k_{int} \ll k_{off}\) is a key requirement
Critical Assumption

QE assumes that the complex dissociates much faster than it is internalized. If \(k_{int}\) is comparable to or larger than \(k_{off}\), use QSS instead.

7.3 Deriving the QE Equations

From the equilibrium assumption and mass balance:

Total ligand: \(L_{tot} = L + LR\)

Total receptor: \(R_{tot} = R + LR\)

Equilibrium: \(K_D = \frac{L \cdot R}{LR}\)

Solving for free ligand \(L\) as a function of \(L_{tot}\) and \(R_{tot}\):

\[ L = \frac{1}{2}\left[(L_{tot} - R_{tot} - K_D) + \sqrt{(L_{tot} - R_{tot} - K_D)^2 + 4 \cdot K_D \cdot L_{tot}}\right] \]

The complex concentration is then:

\[ LR = L_{tot} - L \]

And free receptor:

\[ R = R_{tot} - LR \]

7.4 QE Model in Pumas

tmdd_qe_model = @model begin
    @param begin
        tvCl  RealDomain(lower = 0)
        tvVc  RealDomain(lower = 0)
        tvKsyn  RealDomain(lower = 0)
        tvKdeg  RealDomain(lower = 0)
        tvKD  RealDomain(lower = 0)       # KD replaces Kon/Koff
        tvKint  RealDomain(lower = 0)
    end

    @pre begin
        Cl = tvCl
        Vc = tvVc
        Ksyn = tvKsyn
        Kdeg = tvKdeg
        KD = tvKD
        Kint = tvKint
    end

    @dosecontrol begin
        bioav = (Ltot = 1 / tvVc,)
    end

    @init begin
        Rtot = Ksyn / Kdeg
    end

    @vars begin
        Kel := Cl / Vc
        R0 := Ksyn / Kdeg
        # Solve for free ligand from equilibrium
        L := 0.5 * ((Ltot - Rtot - KD) + sqrt((Ltot - Rtot - KD)^2 + 4 * KD * Ltot))
        # Complex and free receptor from mass balance
        LR := Ltot - L
        R := Rtot - LR
    end

    @dynamics begin
        Ltot' = -Kel * L - Kint * LR
        Rtot' = Ksyn - Kdeg * R - Kint * LR
    end

    @derived begin
        Free_Ligand = @. L
        Free_Receptor = @. R
        Complex_Conc = @. LR
        Total_Ligand = @. Ltot
        Total_Receptor = @. Rtot
    end
end
PumasModel
  Parameters: tvCl, tvVc, tvKsyn, tvKdeg, tvKD, tvKint
  Random effects:
  Covariates:
  Dynamical system variables: Ltot, Rtot
  Dynamical system type: Nonlinear ODE
  Derived: Free_Ligand, Free_Receptor, Complex_Conc, Total_Ligand, Total_Receptor
  Observed: Free_Ligand, Free_Receptor, Complex_Conc, Total_Ligand, Total_Receptor
Key Simplification

Notice that QE reduces the number of differential equations from 3 to 2, and replaces \(k_{on}\) and \(k_{off}\) with a single parameter \(K_D\). This improves identifiability when only ligand concentrations are measured.

8 Quasi-Steady-State (QSS) Approximation

8.1 Mathematical Basis

The QSS approximation, introduced by Gibiansky et al. (2008), is more general than QE. It assumes that the rate of complex formation equals the rate of complex removal:

\[ k_{on} \cdot L \cdot R = (k_{off} + k_{int}) \cdot LR \]

This leads to the steady-state constant:

\[ K_{SS} = \frac{k_{off} + k_{int}}{k_{on}} \]

Note that \(K_{SS} \geq K_D\) always, with equality only when \(k_{int} = 0\).

8.2 When is QSS Appropriate?

The QSS approximation is valid when:

  • Binding kinetics are fast relative to other processes
  • \(k_{int}\) may be significant (unlike QE)
  • Both ligand and receptor concentrations are measured, OR
  • There is a strong scientific rationale for the target biology
QSS vs QE

The key difference is in the steady-state constant:

  • QE: \(K_D = k_{off} / k_{on}\) (assumes \(k_{int} \ll k_{off}\))
  • QSS: \(K_{SS} = (k_{off} + k_{int}) / k_{on}\) (accounts for internalization)

QSS is more general and should be preferred when \(k_{int}\) is not negligible.

8.3 QSS Model in Pumas

tmdd_qss_model = @model begin
    @param begin
        tvCl  RealDomain(lower = 0)
        tvVc  RealDomain(lower = 0)
        tvKsyn  RealDomain(lower = 0)
        tvKdeg  RealDomain(lower = 0)
        tvKSS  RealDomain(lower = 0)       # KSS replaces Kon/Koff/Kint combination
        tvKint  RealDomain(lower = 0)
    end

    @pre begin
        Cl = tvCl
        Vc = tvVc
        Ksyn = tvKsyn
        Kdeg = tvKdeg
        KSS = tvKSS
        Kint = tvKint
    end

    @dosecontrol begin
        bioav = (Ltot = 1 / tvVc,)
    end

    @init begin
        Rtot = Ksyn / Kdeg
    end

    @vars begin
        Kel := Cl / Vc
        R0 := Ksyn / Kdeg
        # Solve for free ligand from quasi-steady-state
        L := 0.5 * ((Ltot - Rtot - KSS) + sqrt((Ltot - Rtot - KSS)^2 + 4 * KSS * Ltot))
        # Complex and free receptor from mass balance
        LR := Ltot - L
        R := Rtot - LR
    end

    @dynamics begin
        Ltot' = -Kel * L - Kint * LR
        Rtot' = Ksyn - Kdeg * R - Kint * LR
    end

    @derived begin
        Free_Ligand = @. L
        Free_Receptor = @. R
        Complex_Conc = @. LR
        Total_Ligand = @. Ltot
        Total_Receptor = @. Rtot
    end
end
PumasModel
  Parameters: tvCl, tvVc, tvKsyn, tvKdeg, tvKSS, tvKint
  Random effects:
  Covariates:
  Dynamical system variables: Ltot, Rtot
  Dynamical system type: Nonlinear ODE
  Derived: Free_Ligand, Free_Receptor, Complex_Conc, Total_Ligand, Total_Receptor
  Observed: Free_Ligand, Free_Receptor, Complex_Conc, Total_Ligand, Total_Receptor

8.4 Relationship Between Full Model and QSS

The full model parameters can be recovered from QSS parameters if both \(K_D\) and \(k_{int}\) are known:

\[ K_D = K_{SS} - \frac{k_{int}}{k_{on}} \]

Or equivalently:

\[ k_{on} = \frac{k_{int}}{K_{SS} - K_D} \]

In practice, QSS is often used when \(k_{on}\) and \(k_{off}\) cannot be separately identified.

9 Michaelis-Menten (MM) Approximation

9.1 Mathematical Basis

The Michaelis-Menten approximation is the most simplified TMDD model. It applies when the ligand concentration is much greater than the total receptor concentration:

\[ L \gg R_{tot} \]

Under this assumption, the target is essentially fully saturated at therapeutic concentrations, and the target-mediated elimination reduces to a capacity-limited (saturable) process:

\[ \frac{dL}{dt} = -k_{el} \cdot L - \frac{V_{max} \cdot L}{K_M + L} \]

Two Derivation Paths for MM

The MM model can be derived from different parent models, leading to different interpretations of \(K_M\):

1. From QSS (reversible binding): \[V_{max} = R_0 \cdot k_{int}, \quad K_M = K_{SS} = \frac{k_{off} + k_{int}}{k_{on}}\]

2. From Irreversible Binding (Gibiansky & Gibiansky, 2009): \[V_{max} = k_{syn}, \quad K_M = K_{IB} = \frac{k_{deg}}{k_{on}}\]

The irreversible binding derivation applies when \(k_{int} \gg k_{off}\) (complex degradation much faster than dissociation). In practice, both derivations yield similar MM kinetics, but the parameter relationships differ.

9.2 When is MM Appropriate?

The Michaelis-Menten approximation is valid when:

  • Free receptor and complex concentrations are not measured
  • Ligand concentrations are always above target concentrations
  • Target-mediated clearance follows saturable kinetics
  • Only ligand PK data are available
Limitation of MM

The MM approximation cannot describe receptor dynamics. If target engagement or receptor occupancy predictions are needed, use the QSS or full model.

9.3 MM Model in Pumas

tmdd_mm_model = @model begin
    @param begin
        tvCl  RealDomain(lower = 0)
        tvVc  RealDomain(lower = 0)
        tvVmax  RealDomain(lower = 0)     # Maximum elimination rate
        tvKM  RealDomain(lower = 0)       # Michaelis constant
    end

    @pre begin
        Cl = tvCl
        Vc = tvVc
        Vmax = tvVmax
        KM = tvKM
    end

    @dosecontrol begin
        bioav = (Ligand = 1 / tvVc,)
    end

    @vars begin
        Kel := Cl / Vc
    end

    @dynamics begin
        Ligand' = -Kel * Ligand - Vmax * Ligand / (KM + Ligand)
    end

    @derived begin
        Free_Ligand = @. Ligand
    end
end
PumasModel
  Parameters: tvCl, tvVc, tvVmax, tvKM
  Random effects:
  Covariates:
  Dynamical system variables: Ligand
  Dynamical system type: Nonlinear ODE
  Derived: Free_Ligand
  Observed: Free_Ligand
Fewest Parameters

The MM model has only 4 parameters (\(CL\), \(V_c\), \(V_{max}\), \(K_M\)), making it the most identifiable but least mechanistic approximation. It is often a good starting point for model development.

9.4 Terminal Slope Limitations of MM

A critical insight from Peletier & Gabrielsson (2012) is that the Michaelis-Menten model can overestimate the terminal elimination rate compared to the full TMDD model.

9.4.1 Why This Matters

At low ligand concentrations (terminal phase), the elimination behavior differs:

Full TMDD Model (internalization-limited):

\[ \lambda_z^{TMDD} \approx k_{int} \]

When receptor turnover is slow relative to internalization, the terminal slope is governed by the internalization rate of the drug-receptor complex.

Michaelis-Menten Model:

\[ \lambda_z^{MM} = k_{el} + \frac{V_{max}}{K_M \cdot V_c} = k_{el} + \frac{k_{int} \cdot R_0}{K_M} \]

The MM model includes an additional term that can significantly increase the apparent terminal slope.

Practical Implication

If you fit a Michaelis-Menten model to TMDD data and use the terminal slope to predict long-term drug exposure, you may underestimate the drug concentrations at late time points. This is because the MM approximation assumes \(L \gg R_{tot}\) throughout, which breaks down in the terminal phase.

9.4.2 When MM Terminal Slope is Accurate

The MM approximation for terminal slope is accurate when:

  1. \(k_{el} \gg k_{int}\) (linear elimination dominates)
  2. \(R_0 \ll K_M\) (low receptor relative to binding constant)
  3. Data only extend to the transition phase (phase C) and don’t capture the true terminal phase (phase D)

9.4.3 Quantifying the Difference

The discrepancy ratio is approximately:

\[ \frac{\lambda_z^{MM}}{\lambda_z^{TMDD}} \approx \frac{k_{el}}{k_{int}} + \frac{R_0}{K_M} \]

For a typical mAb with \(k_{el} \approx 0.002\) hr\(^{-1}\), \(k_{int} \approx 0.002\) hr\(^{-1}\), \(R_0 = 5\) nM, and \(K_M = 1\) nM:

\[ \text{Ratio} \approx 1 + 5 = 6 \]

This means the MM model might predict a terminal half-life 6 times shorter than reality!

Recommendation

When accurate prediction of terminal phase kinetics is important (e.g., for washout periods, immunogenicity assessment, or exposure at trough), prefer the QSS or full TMDD model over Michaelis-Menten.

10 Parameter Comparison Across Models

The following table shows how parameters relate across approximation levels:

Full Model QE Model QSS Model MM Model Relationship
\(k_{on}\), \(k_{off}\) \(K_D\) \(K_{SS}\) \(K_M\) \(K_D = k_{off}/k_{on}\); \(K_{SS} = (k_{off}+k_{int})/k_{on}\)
\(k_{syn}\), \(k_{deg}\) \(k_{syn}\), \(k_{deg}\) \(k_{syn}\), \(k_{deg}\) Implicit in \(V_{max}\) \(R_0 = k_{syn}/k_{deg}\)
\(k_{int}\) \(k_{int}\) \(k_{int}\) Implicit in \(V_{max}\) \(V_{max} = R_0 \cdot k_{int}\)
\(CL\), \(V_c\) \(CL\), \(V_c\) \(CL\), \(V_c\) \(CL\), \(V_c\) Same

11 Comparing Model Predictions

Let’s simulate all four models with equivalent parameters to visualize their predictions.

Show/Hide Code
# Define common parameters based on typical mAb values
R0 = 5.0         # Baseline receptor (nM)
Kdeg = 0.008     # Receptor degradation (hr⁻¹) ~0.2 day⁻¹
Ksyn = R0 * Kdeg # Receptor synthesis (nM/hr)
Kon = 0.33       # Association rate (nM⁻¹ hr⁻¹) ~8 day⁻¹
Koff = 0.33      # Dissociation rate (hr⁻¹) ~8 day⁻¹
Kint = 0.002     # Internalization rate (hr⁻¹) ~0.04 day⁻¹

KD = Koff / Kon  # = 1.0 nM
KSS = (Koff + Kint) / Kon  # ≈ 1.006 nM
Vmax = R0 * Kint  # = 0.01 nM/hr

# Parameter sets for each model
param_full = (
    tvCl = 0.006,
    tvVc = 3.0,
    tvKsyn = Ksyn,
    tvKdeg = Kdeg,
    tvKon = Kon,
    tvKoff = Koff,
    tvKint = Kint,
)

param_qe =
    (tvCl = 0.006, tvVc = 3.0, tvKsyn = Ksyn, tvKdeg = Kdeg, tvKD = KD, tvKint = Kint)

param_qss =
    (tvCl = 0.006, tvVc = 3.0, tvKsyn = Ksyn, tvKdeg = Kdeg, tvKSS = KSS, tvKint = Kint)

param_mm = (tvCl = 0.006, tvVc = 3.0, tvVmax = Vmax, tvKM = KSS)

# Create population with varying doses (wide range for clear separation)
doses = [5, 50, 500]
Random.seed!(42)

pop_full = [
    Subject(
        id = i,
        events = DosageRegimen(dose, time = 0, cmt = 1),
        observations = (
            Free_Ligand = nothing,
            Free_Receptor = nothing,
            Complex_Conc = nothing,
        ),
    ) for (i, dose) in enumerate(doses)
]

pop_qe = [
    Subject(
        id = i,
        events = DosageRegimen(dose, time = 0, cmt = 1),
        observations = (
            Free_Ligand = nothing,
            Free_Receptor = nothing,
            Complex_Conc = nothing,
        ),
    ) for (i, dose) in enumerate(doses)
]

pop_qss = [
    Subject(
        id = i,
        events = DosageRegimen(dose, time = 0, cmt = 1),
        observations = (
            Free_Ligand = nothing,
            Free_Receptor = nothing,
            Complex_Conc = nothing,
        ),
    ) for (i, dose) in enumerate(doses)
]

pop_mm = [
    Subject(
        id = i,
        events = DosageRegimen(dose, time = 0, cmt = 1),
        observations = (Free_Ligand = nothing,),
    ) for (i, dose) in enumerate(doses)
]

# Simulate all models (extended time for mAb kinetics)
obstimes = 0.1:2:672

sim_full = simobs(tmdd_full_model, pop_full, param_full, obstimes = obstimes)
sim_qe = simobs(tmdd_qe_model, pop_qe, param_qe, obstimes = obstimes)
sim_qss = simobs(tmdd_qss_model, pop_qss, param_qss, obstimes = obstimes)
sim_mm = simobs(tmdd_mm_model, pop_mm, param_mm, obstimes = obstimes)
Simulated population (Vector{<:Subject})
  Simulated subjects: 3
  Simulated variables: Free_Ligand

11.1 Free Ligand Comparison

Show/Hide Code
# Convert all simulations to DataFrames and combine
df_full = DataFrame(sim_full)
df_full.Model .= "Full"

df_qe = DataFrame(sim_qe)
df_qe.Model .= "QE"

df_qss = DataFrame(sim_qss)
df_qss.Model .= "QSS"

df_mm = DataFrame(sim_mm)
df_mm.Model .= "MM"

# Create dose labels
dose_label_df = DataFrame(id = string.(1:3), Dose = string.(doses) .* " units")

df_full = innerjoin(df_full, dose_label_df, on = :id)
df_qe = innerjoin(df_qe, dose_label_df, on = :id)
df_qss = innerjoin(df_qss, dose_label_df, on = :id)
df_mm = innerjoin(df_mm, dose_label_df, on = :id)

# Combine all models
df_models = vcat(
    select(df_full, :time, :Free_Ligand, :Model, :Dose),
    select(df_qe, :time, :Free_Ligand, :Model, :Dose),
    select(df_qss, :time, :Free_Ligand, :Model, :Dose),
    select(df_mm, :time, :Free_Ligand, :Model, :Dose),
)

plt =
    data(df_models) *
    mapping(
        :time => "Time (hr)",
        :Free_Ligand => "Free Ligand Concentration (nM)";
        color = :Model,
        linestyle = :Model,
        layout = :Dose,
    ) *
    visual(Lines; linewidth = 2)

draw(plt; axis = (yscale = Makie.pseudolog10,))

Comparison of free ligand predictions across TMDD model approximations for three dose levels.

11.2 Receptor and Complex Dynamics

The MM model cannot predict receptor or complex dynamics. Let’s compare the other three models:

Show/Hide Code
# Use middle dose for comparison
dose_idx = 2

# Prepare data for models that have receptor/complex
df_receptor = vcat(
    select(
        filter(row -> row.id == string(dose_idx), df_full),
        :time,
        :Free_Receptor,
        :Complex_Conc,
        :Model,
    ),
    select(
        filter(row -> row.id == string(dose_idx), df_qe),
        :time,
        :Free_Receptor,
        :Complex_Conc,
        :Model,
    ),
    select(
        filter(row -> row.id == string(dose_idx), df_qss),
        :time,
        :Free_Receptor,
        :Complex_Conc,
        :Model,
    ),
)

# Stack for faceted plot
df_long = stack(
    df_receptor,
    [:Free_Receptor, :Complex_Conc],
    variable_name = :Species,
    value_name = :Concentration,
)
df_long.Species =
    replace.(
        string.(df_long.Species),
        "Free_Receptor" => "Free Receptor",
        "Complex_Conc" => "Complex",
    )

plt =
    data(df_long) *
    mapping(
        :time => "Time (hr)",
        :Concentration => "Concentration (nM)";
        color = :Model,
        linestyle = :Model,
        layout = :Species,
    ) *
    visual(Lines; linewidth = 2)

draw(plt; figure = (title = "Dose = $(doses[dose_idx]) units",))

Comparison of receptor (left) and complex (right) dynamics across models. Note: MM model cannot predict these quantities.

12 Model Selection Decision Framework

Choosing the right TMDD approximation depends on several factors. Here is a practical decision framework:

12.1 Step 1: What Data Are Available?

Data Available Recommended Starting Model
Free ligand only Michaelis-Menten
Free ligand + total receptor QSS
Free ligand + total ligand QSS
Free ligand + free receptor + complex Full TMDD

12.2 Step 2: Consider Scientific Context

  • Do you need receptor dynamics? Use QSS or Full model
  • Is internalization significant? Use QSS over QE
  • Is target concentration known to be low relative to drug? MM may suffice
  • Do you need mechanistic interpretation? Use Full or QSS model

12.3 Step 3: Assess Model Fit and Identifiability

Start with the simplest appropriate model (often MM) and increase complexity if:

  • The model fails to capture key PK features
  • Additional data become available
  • Mechanistic predictions are needed
Practical Recommendation

Start simple, add complexity as needed.

  1. Fit Michaelis-Menten model first

  2. If MM fails, try QSS

  3. If QSS fails and you have receptor data, try Full TMDD

  4. Compare models using information criteria (AIC, BIC)

    :::

13 Identifiability Considerations

Following Gibiansky et al. (2008), here are key identifiability insights:

13.1 With Free Ligand Data Only

  • MM model: All 4 parameters identifiable
  • QSS model: \(CL\), \(V_c\) identifiable; \(k_{syn}\), \(k_{deg}\), \(K_{SS}\), \(k_{int}\) may require constraints
  • Full model: Generally not identifiable

13.2 With Total Ligand and Total Receptor Data

  • QSS model: All parameters typically identifiable
  • Full model: Still may have issues with \(k_{on}\) and \(k_{off}\) separately

13.3 Practical Strategies

  1. Fix known parameters: Use literature values for \(k_{deg}\), \(K_D\) from in vitro studies
  2. Reduce model: If parameters are not identifiable, use a simpler approximation
  3. Collect more data: Measure receptor concentrations when possible
  4. Profile likelihood: Assess practical identifiability of estimated parameters

14 Adding Peripheral Compartment (Two-Compartment Models)

All TMDD models can be extended to include peripheral distribution. Here is the two-compartment QSS model:

tmdd_qss_2cmt_model = @model begin
    @param begin
        tvCl  RealDomain(lower = 0)
        tvVc  RealDomain(lower = 0)
        tvQ  RealDomain(lower = 0)        # Intercompartmental clearance
        tvVp  RealDomain(lower = 0)       # Peripheral volume
        tvKsyn  RealDomain(lower = 0)
        tvKdeg  RealDomain(lower = 0)
        tvKSS  RealDomain(lower = 0)
        tvKint  RealDomain(lower = 0)
    end

    @pre begin
        Cl = tvCl
        Vc = tvVc
        Q = tvQ
        Vp = tvVp
        Ksyn = tvKsyn
        Kdeg = tvKdeg
        KSS = tvKSS
        Kint = tvKint
    end

    @dosecontrol begin
        bioav = (Ltot = 1 / tvVc,)
    end

    @init begin
        Rtot = Ksyn / Kdeg
    end

    @vars begin
        Kel := Cl / Vc
        K12 := Q / Vc
        K21 := Q / Vp
        R0 := Ksyn / Kdeg
        # Solve for free ligand from quasi-steady-state
        L := 0.5 * ((Ltot - Rtot - KSS) + sqrt((Ltot - Rtot - KSS)^2 + 4 * KSS * Ltot))
        # Complex and free receptor from mass balance
        LR := Ltot - L
        R := Rtot - LR
    end

    @dynamics begin
        Ltot' = -Kel * L - Kint * LR - K12 * L + K21 * (Peripheral / Vc)
        Peripheral' = K12 * L * Vc - K21 * Peripheral
        Rtot' = Ksyn - Kdeg * R - Kint * LR
    end

    @derived begin
        Free_Ligand = @. L
        Free_Receptor = @. R
        Complex_Conc = @. LR
    end
end
PumasModel
  Parameters: tvCl, tvVc, tvQ, tvVp, tvKsyn, tvKdeg, tvKSS, tvKint
  Random effects:
  Covariates:
  Dynamical system variables: Ltot, Peripheral, Rtot
  Dynamical system type: Nonlinear ODE
  Derived: Free_Ligand, Free_Receptor, Complex_Conc
  Observed: Free_Ligand, Free_Receptor, Complex_Conc

15 Summary

In this tutorial, we covered:

  • Why simplify: Identifiability, data limitations, numerical stability
  • QE approximation: Assumes binding at equilibrium; uses \(K_D\)
  • QSS approximation: More general; uses \(K_{SS} = (k_{off} + k_{int})/k_{on}\)
  • MM approximation: Simplest; assumes \(L \gg R_{tot}\); uses \(V_{max}\) and \(K_M\)
  • Model selection: Start simple, add complexity based on data and needs
  • Identifiability: Key consideration for parameter estimation

In Part 4, we will explore special cases and extensions including constant receptor models, irreversible binding, and multi-target scenarios.

16 References

References

Gibiansky, L., & Gibiansky, E. (2009). Target-mediated drug disposition model: Approximations, identifiability of model parameters and applications to the population pharmacokinetic-pharmacodynamic modeling of biologics. Expert Opinion on Drug Metabolism & Toxicology, 5(7), 803–812. https://doi.org/10.1517/17425250903030043
Gibiansky, L., Gibiansky, E., Kakkar, T., & Ma, P. (2008). Approximations of the target-mediated drug disposition model and identifiability of model parameters. Journal of Pharmacokinetics and Pharmacodynamics, 35(5), 573–591. https://doi.org/10.1007/s10928-008-9102-8
Peletier, L. A., & Gabrielsson, J. (2012). Dynamics of target-mediated drug disposition: Characteristic profiles and parameter identification. Journal of Pharmacokinetics and Pharmacodynamics, 39(5), 429–451. https://doi.org/10.1007/s10928-012-9260-6

Reuse