Hello,

I have a question - once the dug is absorbed into circulation, the fate should be same - same CL & Volume of distribution despite the dose route, normal PK or flip-flop PK. Is this true? If IM displays flip-flop PK and PO displays normal/conventional PK, is it good idea to simultaneously fit PO/IM PK to avoid switch of kel and ka? Here is the code I used to simulate IM PK profile using fitted parameters. Depot is for PO route here (Set oral bioavailability Fpo=1, solve FIM for fitting). IM PK (IR and SR) is described as parallel absorption - slow process with transit compartment delay. The code works for simulation but the simulation output doesn’t look right (pkdata_im = DataFrame(obs), some parameters show NA in the output csv file) . Could you help check the code? Many Thanks!

```
two_parallel_foabs = @model begin
@param begin
tvcl ∈ RealDomain(lower=0)
tvvc ∈ RealDomain(lower=0)
tvka ∈ RealDomain(lower=0)
tvka1 ∈ RealDomain(lower=0)
tvka2 ∈ RealDomain(lower=0)
tvlag ∈ RealDomain(lower=0)
tvFracIM ∈ RealDomain(lower=0)
tvvp ∈ RealDomain(lower=0)
tvq ∈ RealDomain(lower=0)
tvvp2 ∈ RealDomain(lower=0)
tvq2 ∈ RealDomain(lower=0)
tvktr ∈ RealDomain(lower=0)
tvfIM ∈ RealDomain(lower=0)
Ω ∈ PDiagDomain(8)
σ_prop ∈ RealDomain(lower=0)
end
@random begin
η ~ MvNormal(Ω)
end
@pre begin
CL = tvcl*exp(η[1])
Vc = tvvc
Vp = tvvp
Q = tvq
Vp2 = tvvp2
Q2 = tvq2*exp(η[2])
Ka = tvka*exp(η[3])
Ka1 = tvka1*exp(η[4])
Ka2 = tvka2*exp(η[5])
Ktr = tvktr*exp(η[7])
end
@dosecontrol begin
lags = (Depot = tvlag,)
bioav = ( Depot=1, IR = tvfIM*exp(η[8])*(1-tvFracIM)*exp(η[6]), SR = tvfIM*exp(η[8])*tvFracIM*exp(η[6]))
end
@dynamics begin
Depot' = -Ka*Depot
IR' = -Ka1*IR
SR' = -Ktr*SR
Transit1' = Ktr*SR - Ktr*Transit1
Transit2' = Ktr*Transit1 - Ktr*Transit2
Transit3' = Ktr*Transit2 - Ktr*Transit3
Transit4' = Ktr*Transit3 - Ktr*Transit4
Transit5' = Ktr*Transit4 - Ktr*Transit5
Transit6' = Ktr*Transit5 - Ktr*Transit6
Transit7' = Ktr*Transit6 - Ka2*Transit7
Central' = Ka*Depot + Ka1*IR + Ka2*Transit7 - (CL/Vc)*Central - (Q/Vc)*Central + (Q/Vp)*Peripheral1 - (Q2/Vc)*Central + (Q2/Vp2)*Peripheral2
Peripheral1' = (Q/Vc)*Central - (Q/Vp)*Peripheral1
Peripheral2' = (Q2/Vc)*Central - (Q2/Vp2)*Peripheral2
end
@derived begin
cp = @. (Central/Vc)
dv ~ @. Normal(abs(cp), abs(cp)*σ_prop)
end
end
param = ( tvcl = 14, tvvc = 43, tvq = 1,
tvvp = 8.7,
tvq2 = 0.7,
tvvp2 = 206,
tvka = 1.28, tvka1 = 0.0024, tvka2 = 0.001, tvlag = 0.18, tvFracIM = 0.49, tvktr=0.0166, tvfIM=3,
Ω = Diagonal([0, 0, 0, 0, 0, 0,0,0]),σ_prop =0.00001)
## simulate a profile with this type of absorption:
dose_fo1 = DosageRegimen(112500, cmt = 2, time = 0)
dose_fo2 = DosageRegimen(112500, cmt = 3, time = 0)
dose = DosageRegimen(dose_fo1, dose_fo2) # Actual dose is made up of 2 virtual doses
s3_IM = Subject(id=3, events=dose,observations = (conc = nothing,))
subj_with_covariates = map(1:10) do i
Subject(id = i,
events = dose,
covariates = choose_covariates(),
observations = (conc = nothing,))
end
obs = simobs(two_parallel_foabs, s3_IM, param, obstimes = 0:1:4800)
## plot the results
sim_plot(two_parallel_foabs, obs, observations = :cp)
pkdata_im = DataFrame(obs)
CSV.write("../data/pk_im.csv", pkdata_im)
```