IOV modeling and simulation


I’ve been working on a project trying to model the PK of ertapenem in intermittent hemodialysis patients. The two issues I ran into was: 1) Simulating inter-occasional variability and 2) simulating a population where IOV is expected to occur at certain time points within the data. I’ve included a snapshot of my PUMAS code and the dataset I’ve been working with

Brief background: I have data for ertapenem concentrations during dialysis and off of dialysis. Each participant goes through 4-5 cycles of dialysis (every 48 hours for 6 hours) and we wanted to incorporate IOV on CL for each dialysis cycle.

Thanks in advance!

This is a snapshot of the dataset I’ve been working with.

Hello clee,

Please take a look at Defining NLME models in Pumas · Pumas. There are several examples demonstrating between/inter occasion variability. Hopefully that will help you. Please reply here if it requires further explanation.

Andreas Noack

I was just drafting an answer that included the link that @andreasnoack posted.

Additionally, I would also recommend the Simulation tutorial at that will show you how to perform simulations.

Hi @clee,

Please also see this post for more information about IOV in Pumas:

For your data set, you will want to have a single OCC column with 1-4/5 for each cycle of dialysis.

Regarding your 2), if I understand you correctly, you may want to do something like the choose_covariates function here: . So, you’ll determine each subject’s OCC column based on whatever criteria you choose for that given subject, and then pass OCC as a covariate.

Thanks everyone for your advice. I’ve been looking at all the resources you have provided and tried a couple of methods.

  1. I changed my OCC to a single column to have 1 -6 different occasions
  2. CL = tvcl * exp(tvdialDIAL) * exp(η[1]+DIALη_occ[OCC]) ; this seems to work for me where DIAL is a covariate, and when it is present (i.e. =1), add IOV.
    But, I found two different methods to define η_occ in the @random block and both gives different answers.
    η_occ ~ MvNormal(Diagonal(fill(Ω_occ,6))) vs.
    η_occ ~ MvNormal([Ω_occ, Ω_occ, Ω_occ, Ω_occ, Ω_occ,Ω_occ]).
    I was wonder if anyone why this is the case?

The difference is a little subtle. In the first version, you pass a matrix, which the MvNornal constructor interprets as a covariance matrix. In the second version, you pass a vector which the MvNormal constructor interprets as a vector of standard deviations. Compare

julia> cov(MvNormal(Diagonal([0.1, 0.1])))
2×2 Matrix{Float64}:
 0.1  0.0
 0.0  0.1

julia> cov(MvNormal([0.1, 0.1]))
2×2 Matrix{Float64}:
 0.01  0.0
 0.0   0.01

For inter/between occasion variability, we usually recommend a third constructor which is

julia> cov(MvNormal(2, 0.1))
2×2 Matrix{Float64}:
 0.01  0.0
 0.0   0.01

i.e. the first argument is the dimension and the second is the standard deviation. Personally, I also prefer to use lower case greek letters, i.e. ω, for standard deviations, squared versions, i.e. ω², for variances and only use the upper case version for the whole covariance matrix. I think that makes it easier to read.

1 Like

Thanks @donaldlee3 ,
Thanks for sharing the reference. Based on the provided example, you would only be assigning a specific covariate for each simulated subjected (same covariate value at all timepoints for the subject). I was wondering if there was a way to specify OCC to be present at specific time points for each individual (i.e. subject 1: OCC = 0 at t = 0-18h, OCC = 1 at t = 19-24h, OCC = 0 at t = 25-42h, and OCC =2 at t = 43-48h ).
Many thanks!

One way to do this is to convert your Pumas population into a data frame, and then add your OCC column that you made manually. You can then use read_pumas to convert the dataframe back to a Pumas population.

For example:

pop_df = DataFrame(Population(map(j -> Subject(id=j, events=ev), 1:Nsubj)))

pop_df[!, :OCC] .= pre_made_occs

pop = read_pumas(pop_df,
                 id                      =   :id,
                 time                    =   :time,
                 amt                     =   :amt,
                 cmt                     =   :cmt,
                 observations            =   [:dv],
                 covariates              =   [:OCC],
                 covariates_direction    =   :left,
                 evid                    =   :evid)