using PumasUtilities
using Random
using Pumas
using CairoMakie
using AlgebraOfGraphics
using DataFramesMeta
using CSV
using Dates
PK44 (Part 1) - Estimation of inhibitory constant Ki
1 Background
- Structural model - Estimation of inhibitory rate constant in competitive enzyme inhibition model
- Number of subjects - 1
- Number of compounds - 1
2 Learning Outcome
This exercise demonstrates how to build a competitive inhibitory model and understand the relationship between rate of metabolite formation and concentration.
3 Objectives
To build a competitive inhibitory model and simulate the model for a single subject.
4 Libraries
Load the necessary libraries
5 Model definition
Note the expression of the model parameters with helpful comments. The model is expressed with differential equations. Residual variability is an additive error model.
pk_44_cim = @model begin
@metadata begin
desc = "Competitive Inhibitory Model"
timeu = u"minute"
end
@param begin
"""
Maximum metabolic rate (μM*gm/min)
"""
tvvmax ∈ RealDomain(lower = 0)
"""
Michaelis-Mentons constant (μmol/L)
"""
tvkm ∈ RealDomain(lower = 0)
"""
Inhibitory constant (μmol/L)
"""
tvki ∈ RealDomain(lower = 0)
Ω ∈ PDiagDomain(2)
"""
Additive RUV
"""
σ_add ∈ RealDomain(lower = 0)
end
@random begin
η ~ MvNormal(Ω)
end
@covariates conc I
@pre begin
Vmax = tvvmax * exp(η[1])
Km = tvkm * exp(η[2])
Ki = tvki
_conc = conc
_I = I
end
@derived begin
## Competitive Inhibiton Model
rate_cim = @. Vmax * _conc / ((Km * (1 + (_I / Ki)) + _conc))
"""
Metabolic rate (nmol/min/mg)
"""
dv_rate_cim ~ @. Normal(rate_cim, σ_add)
end
endPumasModel
Parameters: tvvmax, tvkm, tvki, Ω, σ_add
Random effects: η
Covariates: conc, I
Dynamical system variables:
Dynamical system type: No dynamical model
Derived: rate_cim, dv_rate_cim
Observed: rate_cim, dv_rate_cim
6 Initial Estimates of Model Parameters
The model parameters for simulation are the following. Note that tv represents the typical value for parameters.
Vmax- Maximum metabolic rate (μM*gm_protein/min)Km- Michaelis-Mentons constant (μmol/L)Ki- Inhibitory constant (μmol/L)I- Inhibitor concentration/Exposure
param = (
tvvmax = 99.8039,
tvkm = 11.3192,
tvki = 5.05605,
Ω = Diagonal([0.01, 0.01]),
σ_add = 1.256,
)7 Creating a Dataset
The time, concentration data and exposure (I) will be used to estimate the rate of metabolite concentration.
df_sub1 = map(
i -> DataFrame(id = i, time = 1:1:1000, dv_rate_cim = missing, conc = 1:1:1000, I = 0),
1:6,
)
df = vcat(DataFrame.(df_sub1)...)
df[!, :I] = ifelse.(df.id .== 2, 10, df.I)
df[!, :I] = ifelse.(df.id .== 3, 25, df.I)
df[!, :I] = ifelse.(df.id .== 4, 50, df.I)
df[!, :I] = ifelse.(df.id .== 5, 75, df.I)
df[!, :I] = ifelse.(df.id .== 6, 100, df.I)
df_sub1 = df
sub1 = read_pumas(
df_sub1,
observations = [:dv_rate_cim],
covariates = [:conc, :I],
event_data = false,
)Population
Subjects: 6
Covariates: conc, I
Observations: dv_rate_cim
8 Single-Subject Simulation
Simulate the rate of metabolite formation
Initialize the random number generator with a seed for reproducibility of the simulation.
Random.seed!(123)sim_sub1 = simobs(pk_44_cim, sub1, param)Simulated population (Vector{<:Subject})
Simulated subjects: 6
Simulated variables: rate_cim, dv_rate_cim
pkdata_44_1 = DataFrame(sim_sub1)
first(pkdata_44_1, 5)| Row | id | time | rate_cim | dv_rate_cim | evid | conc | I | η_1 | η_2 | Vmax | Km | Ki | _conc | _I |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| String? | Float64 | Float64? | Float64? | Int64? | Int64? | Int64? | Float64? | Float64? | Float64? | Float64? | Float64? | Int64? | Int64? | |
| 1 | 1 | 1.0 | 9.02751 | 9.43209 | 0 | 1 | 0 | -0.0453924 | -0.1684 | 95.3748 | 9.56491 | 5.05605 | 1 | 0 |
| 2 | 1 | 2.0 | 16.4938 | 17.0399 | 0 | 2 | 0 | -0.0453924 | -0.1684 | 95.3748 | 9.56491 | 5.05605 | 2 | 0 |
| 3 | 1 | 3.0 | 22.7717 | 23.5659 | 0 | 3 | 0 | -0.0453924 | -0.1684 | 95.3748 | 9.56491 | 5.05605 | 3 | 0 |
| 4 | 1 | 4.0 | 28.124 | 28.7451 | 0 | 4 | 0 | -0.0453924 | -0.1684 | 95.3748 | 9.56491 | 5.05605 | 4 | 0 |
| 5 | 1 | 5.0 | 32.7413 | 32.4165 | 0 | 5 | 0 | -0.0453924 | -0.1684 | 95.3748 | 9.56491 | 5.05605 | 5 | 0 |
9 Visualize Results
@chain pkdata_44_1 begin
@rsubset :conc ∈ [1, 5, 10, 15, 26, 104, 251, 502, 1000]
data(_) *
mapping(
:conc => "Concentration (μM)",
:rate_cim => "Metabolic rate (nmol/min/mg protein)",
color = :I => nonnumeric => "Exposure",
) *
visual(ScatterLines, linewidth = 4, markersize = 12)
draw(; figure = (; fontsize = 22), axis = (; xticks = 0:100:1000, yticks = 0:10:100))
end@chain pkdata_44_1 begin
@rsubset :conc ∈ [1, 5, 10, 15, 26, 104, 251, 502, 1000]
data(_) *
mapping(
:conc => "log Concentration (μM)",
:rate_cim => "Metabolic rate (nmol/min/mg protein)",
color = :I => nonnumeric => "Exposure",
) *
visual(ScatterLines, linewidth = 4, markersize = 12)
draw(;
figure = (; fontsize = 22),
axis = (;
xscale = log10,
xtickformat = i -> (@. string(round(i; digits = 1))),
yticks = 0:10:100,
),
)
end10 Conclusion
This tutorial showed how to build a competitive inhibitory model and perform a rate of metabolite formation versus concentration simulation.