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

using Random
using Pumas
using PumasUtilities
using CairoMakie
using AlgebraOfGraphics
using DataFramesMeta
using CSV
using Dates

3 Model

pk_37 = @model begin
    @metadata begin
        desc = "Invivo/Invitro Extrapolation Model"
    end

    @param begin
        """
        Maximum Rate of System 1 (μmol/min/g)
        """
        tvvmax1  RealDomain(lower = 0)
        """
        Michaelis-Menten constant for system 1 (μmol/L)
        """
        tvkm1  RealDomain(lower = 0)
        """
        Maximum Rate of System 2 (μmol/min/g)
        """
        tvvmax2  RealDomain(lower = 0)
        """
        Michaelis-Menten constant for system 2 (μmol/L)
        """
        tvkm2  RealDomain(lower = 0)
        Ω  PDiagDomain(4)
        """
        Additive RUV
        """
        σ_add  RealDomain(lower = 0)
    end

    @random begin
        η ~ MvNormal(Ω)
    end

    @covariates Conc

    @pre begin
        Vmax1 = tvvmax1 * exp(η[1])
        Km1 = tvkm1 * exp(η[2])
        Vmax2 = tvvmax2 * exp(η[3])
        Km2 = tvkm2 * exp(η[4])
        _Conc = Conc
    end

    @derived begin
        rate = @. ((Vmax1 * _Conc / (Km1 + _Conc)) + (Vmax2 * _Conc / (Km2 + _Conc)))
        """
        Observed Concentrations (μmol/L)
        """
        dv_rate ~ @. Normal(rate, σ_add)
    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 = (
    tvvmax1 = 2.43386,
    tvkm1 = 256.384,
    tvvmax2 = 0.22523,
    tvkm2 = 2.23625,
    Ω = Diagonal([0.04, 0.04, 0.04, 0.04]),
    σ_add = 0.0441543,
)

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.

df_sub1 = DataFrame(id = 1, time = 0:1:40000, Conc = 0.00:0.05:2000, dv_rate = missing);
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

zfx = zero_randeffs(pk_37, sub1, param)
1-element Vector{@NamedTuple{η::Vector{Float64}}}:
 (η = [0.0, 0.0, 0.0, 0.0],)
sim_sub1 = simobs(pk_37, sub1, param, zfx)
Simulated population (Vector{<:Subject})
  Simulated subjects: 1
  Simulated variables: rate, dv_rate
df37 = DataFrame(sim_sub1)
first(df37, 5)
5×15 DataFrame
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 = (;
            title = "Metabolic Rate vs Substrate Concentration",
            xticks = 0:200:2000,
            yticks = 0:0.5:2.5,
        ),
    )
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 = (;
            title = "Metabolic Rate vs Substrate Concentration",
            yticks = 0:0.5:2.5,
            xscale = log10,
        ),
    )
end

8 Population simulation

par = (
    tvvmax1 = 2.43386,
    tvkm1 = 256.384,
    tvvmax2 = 0.22523,
    tvkm2 = 2.23625,
    Ω = Diagonal([0.0426, 0.0124, 0.0256, 0.0371]),
    σ_add = 0.0441543,
)

df1_pop = map(
    i -> DataFrame(
        id = i,
        time = (1:1:16),
        Conc = [
            0.02,
            0.05,
            0.1,
            0.2,
            0.5,
            1,
            10,
            20,
            50,
            100,
            180,
            200,
            250,
            500,
            1000,
            2000,
        ],
        dv_rate = missing,
    ),
    1:55,
)
df2_pop = vcat(DataFrame.(df1_pop)...)
pop =
    read_pumas(df2_pop, observations = [:dv_rate], covariates = [:Conc], event_data = false)

Random.seed!(1234)
sim_pop = simobs(pk_37, pop, par)

df_sim = DataFrame(sim_pop);

#CSV.write("pk_37.csv", df_sim)