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{(:η,), Tuple{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	η_1	η_2	η_3	η_4	Vmax1	Km1	Vmax2	Km2	_Conc
	String?	Float64	Float64?	Float64?	Int64?	Float64?	Float64?	Float64?	Float64?	Float64?	Float64?	Float64?	Float64?	Float64?	Float64?
1	1	0.0	0.0	-0.0151318	0	0.0	0.0	0.0	0.0	0.0	2.43386	256.384	0.22523	2.23625	0.0
2	1	1.0	0.00540031	0.0864815	0	0.05	0.0	0.0	0.0	0.0	2.43386	256.384	0.22523	2.23625	0.05
3	1	2.0	0.0105896	0.0922725	0	0.1	0.0	0.0	0.0	0.0	2.43386	256.384	0.22523	2.23625	0.1
4	1	3.0	0.0155811	0.0200618	0	0.15	0.0	0.0	0.0	0.0	2.43386	256.384	0.22523	2.23625	0.15
5	1	4.0	0.020387	-0.0343151	0	0.2	0.0	0.0	0.0	0.0	2.43386	256.384	0.22523	2.23625	0.2

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)