using Random
using Pumas
using PumasUtilities
using CairoMakie
using AlgebraOfGraphics
using DataFramesMeta
using CSV
using Dates
PK37 - Invitro/Invivo Extrapolation I
1 Learning Outcome
In this exercise, you will learn to assess the relationship between rate
and concentration
using a non-linear clearance model. The Clearance is impacted by two parallel non-linear processes, and both of them are taken into account for calculating rates. The use of a single enzyme system was not able to account for the deviations in the system.
2 Libraries
Call the necessary libraries to get started
3 Model
= @model begin
pk_37 @metadata begin
= "Invivo/Invitro Extrapolation Model"
desc end
@param begin
"""
Maximum Rate of System 1 (μmol/min/g)
"""
∈ RealDomain(lower = 0)
tvvmax1 """
Michaelis-Menten constant for system 1 (μmol/L)
"""
∈ RealDomain(lower = 0)
tvkm1 """
Maximum Rate of System 2 (μmol/min/g)
"""
∈ RealDomain(lower = 0)
tvvmax2 """
Michaelis-Menten constant for system 2 (μmol/L)
"""
∈ RealDomain(lower = 0)
tvkm2 ∈ PDiagDomain(4)
Ω """
Additive RUV
"""
∈ RealDomain(lower = 0)
σ_add end
@random begin
~ MvNormal(Ω)
η end
@covariates Conc
@pre begin
= tvvmax1 * exp(η[1])
Vmax1 = tvkm1 * exp(η[2])
Km1 = tvvmax2 * exp(η[3])
Vmax2 = tvkm2 * exp(η[4])
Km2 = Conc
_Conc end
@derived begin
= @. ((Vmax1 * _Conc / (Km1 + _Conc)) + (Vmax2 * _Conc / (Km2 + _Conc)))
rate """
Observed Concentrations (μmol/L)
"""
~ @. Normal(rate, σ_add)
dv_rate end
end
PumasModel
Parameters: tvvmax1, tvkm1, tvvmax2, tvkm2, Ω, σ_add
Random effects: η
Covariates: Conc
Dynamical system variables:
Dynamical system type: No dynamical model
Derived: rate, dv_rate
Observed: rate, dv_rate
4 Parameters
The parameters are as given below. Note that tv
represents the typical value for parameters.
tvvmax1
- Maximum Rate of System 1 (μmol/min/g)tvkm1
- Michaelis-Menten constant for system 1 (μmol/L)tvvmax2
- Maximum Rate of System 2 (μmol/min/g)tvkm2
- Michaelis-Menten constant for system 2 (μmol/L)
= (
param = 2.43386,
tvvmax1 = 256.384,
tvkm1 = 0.22523,
tvvmax2 = 2.23625,
tvkm2 = Diagonal([0.04, 0.04, 0.04, 0.04]),
Ω = 0.0441543,
σ_add )
5 Dosage Regimen
Since there is no dosage regimen, we will create a DataFrame
with no events and read the file using read_pumas
which will be used for simulation.
= DataFrame(id = 1, time = 0:1:40000, Conc = 0.00:0.05:2000, dv_rate = missing); df_sub1
=
sub1 read_pumas(df_sub1, observations = [:dv_rate], covariates = [:Conc], event_data = false)
Population
Subjects: 1
Covariates: Conc
Observations: dv_rate
6 Simulations
We will now simulate the rate using the simobs
function
Random.seed!(123)
The random effects are zero’ed out since we are simulating a single subject
= zero_randeffs(pk_37, sub1, param) zfx
1-element Vector{@NamedTuple{η::Vector{Float64}}}:
(η = [0.0, 0.0, 0.0, 0.0],)
= simobs(pk_37, sub1, param, zfx) sim_sub1
Simulated population (Vector{<:Subject})
Simulated subjects: 1
Simulated variables: rate, dv_rate
= DataFrame(sim_sub1)
df37 first(df37, 5)
Row | id | time | rate | dv_rate | evid | Conc | Vmax1 | Km1 | Vmax2 | Km2 | _Conc | η_1 | η_2 | η_3 | η_4 |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
String | Float64 | Float64? | Float64? | Int64 | Float64? | Float64? | Float64? | Float64? | Float64? | Float64? | Float64 | Float64 | Float64 | Float64 | |
1 | 1 | 0.0 | 0.0 | 0.100922 | 0 | 0.0 | 2.43386 | 256.384 | 0.22523 | 2.23625 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 |
2 | 1 | 1.0 | 0.00540031 | 0.0495946 | 0 | 0.05 | 2.43386 | 256.384 | 0.22523 | 2.23625 | 0.05 | 0.0 | 0.0 | 0.0 | 0.0 |
3 | 1 | 2.0 | 0.0105896 | 0.0203634 | 0 | 0.1 | 2.43386 | 256.384 | 0.22523 | 2.23625 | 0.1 | 0.0 | 0.0 | 0.0 | 0.0 |
4 | 1 | 3.0 | 0.0155811 | 0.0397764 | 0 | 0.15 | 2.43386 | 256.384 | 0.22523 | 2.23625 | 0.15 | 0.0 | 0.0 | 0.0 | 0.0 |
5 | 1 | 4.0 | 0.020387 | 0.0676771 | 0 | 0.2 | 2.43386 | 256.384 | 0.22523 | 2.23625 | 0.2 | 0.0 | 0.0 | 0.0 | 0.0 |
7 Visualization
Use the DataFrame
for plotting
@chain df37 begin
@rsubset :Conc ∈
0.02, 0.05, 0.1, 0.2, 0.5, 1, 10, 20, 50, 100, 180, 200, 250, 500, 1000, 2000]
[data(_) *
mapping(:Conc => "Concentration (μM)", :rate => "Rate (μmol/g/min)") *
visual(ScatterLines; linewidth = 4, markersize = 12, markercolor = :blue)
draw(
= (;
axis = "Metabolic Rate vs Substrate Concentration",
title = 0:200:2000,
xticks = 0:0.5:2.5,
yticks
),
)end
@chain df37 begin
@rsubset :Conc ∈
0.02, 0.05, 0.1, 0.2, 0.5, 1, 10, 20, 50, 100, 180, 200, 250, 500, 1000, 2000]
[data(_) *
mapping(:Conc => "Concentration (μM)", :rate => "Rate (μmol/g/min)") *
visual(ScatterLines; linewidth = 4, markersize = 12, markercolor = :blue)
draw(
= (;
axis = "Metabolic Rate vs Substrate Concentration",
title = 0:0.5:2.5,
yticks = log10,
xscale
),
)end
8 Population simulation
= (
par = 2.43386,
tvvmax1 = 256.384,
tvkm1 = 0.22523,
tvvmax2 = 2.23625,
tvkm2 = Diagonal([0.0426, 0.0124, 0.0256, 0.0371]),
Ω = 0.0441543,
σ_add
)
= map(
df1_pop -> DataFrame(
i = i,
id = (1:1:16),
time = [
Conc 0.02,
0.05,
0.1,
0.2,
0.5,
1,
10,
20,
50,
100,
180,
200,
250,
500,
1000,
2000,
],= missing,
dv_rate
),1:55,
)= vcat(DataFrame.(df1_pop)...)
df2_pop =
pop read_pumas(df2_pop, observations = [:dv_rate], covariates = [:Conc], event_data = false)
Random.seed!(1234)
= simobs(pk_37, pop, par)
sim_pop
= DataFrame(sim_pop);
df_sim
#CSV.write("pk_37.csv", df_sim)