Parent-metabolite + dual absorption: fitting error

chrisv · November 16, 2023, 11:46am

Hello. I am trying to fit a parent-metabolite model with dual absorption (ka1 and ka2). In the database I have two administrations at time t=0 with the total dose at CMT=1 (ev) and CMT=2 (ev). I get the error that the gradient is 0 and therefore the parameter is not identified:

First I get the error for tvlag. Then I set FIXED the tvlag parameter.
I get the error for tvbio. Then I set FIXED the tvio parameter.
I get the error for ETA(5). etc.

And so on but the model does not achieve the fit.

I will be very grateful for any possible help.

This is the code:
mymodel = @model begin
@metadata begin
desc = “Two parallel first-order absorptions - Two Compartment Model with Metabolite Compartment”
timeu = u"hr"
end

@param begin
  "Volume of Central Compartment (L)"
  tvvc            ∈ RealDomain(lower=0)
  "Clearance of parent drug (L/hr)"
  tvclp          ∈ RealDomain(lower=0)
  "Inter Compartmental Clearance (L/hr)"
  tvq           ∈ RealDomain(lower=0)
  "Volume of Peripheral Compartment (L)"
  tvvp            ∈ RealDomain(lower=0)
  "Volume of Compartment of Metabolite (L)"
  tvvm          ∈ RealDomain(lower=0)
  "Clearance of metabolite (L/hr)"
  tvclm_M           ∈ RealDomain(lower=0)
  tvclm_F           ∈ RealDomain(lower=0)
  "Lag time (hr)"
  tvlag           ∈ RealDomain(lower=0)
  "Absorption rate constant 1 (hr⁻¹)"
  tvka1            ∈ RealDomain(lower=0)
  "Absorption rate constant 2 (hr⁻¹)"
  tvka2            ∈ RealDomain(lower=0)
  "fraction metabolized"
  tvfmet            ∈ RealDomain(lower=0)
  "Fraction"
  tvbio ∈ RealDomain(lower=0)
  # variance
   Ω              ∈ PDiagDomain(5)
  "Proportional RUV - Plasma"
  σ²_prop_cp      ∈ RealDomain(lower=0.001)
  "Proportional RUV - Metabolite"
  σ²_prop_met     ∈ RealDomain(lower=0.001)
end

@random begin
  η               ~ MvNormal(Ω)
end

@covariates begin
  sex
end

@pre begin
  tvclm = ifelse(sex == "M", tvclm_M, tvclm_F)
  Vc              = tvvc * exp(η[1])
  Clp              = tvclp * exp(η[2])
  Q            = tvq
  Vp              = tvvp
  Vm             = tvvm
  Clm             = tvclm * exp(η[3])
  Ka1              = tvka1 * exp(η[4])
  Ka2              = tvka2
  lags          = (SR = tvlag,)
  bioav = (IR = tvbio*exp(η[5]), SR = 1-tvbio*exp(η[5]))
  fmet              = tvfmet
end

@dynamics begin
  IR'      = -Ka1*IR
  SR'      = -Ka2*SR
  Central'        =  Ka1*IR + Ka2*SR - (1-fmet)*(Clp/Vc)*Central -(fmet)*(Clm/Vc)*Central - (Q/Vc)*Central + (Q/Vp)*Peripheral
  Peripheral'     = (Q/Vc)*Central - (Q/Vp)*Peripheral
  Metabolite'     = (fmet)*(Clp/Vc)*Central - (Clm/Vm)*Metabolite
end

@derived begin
  cp              = @. Central/Vc
  dv_cp           ~ @. Normal(cp, sqrt(cp^2*σ²_prop_cp))
  met             = @. Metabolite/Vm
  dv_met          ~ @. Normal(met, sqrt(met^2*σ²_prop_met))
end

end

iparams = (
tvka1 = 1,
tvka2 = 10,
tvclp = 40,
tvvc = 50,
tvq = 20,
tvvp = 120,
tvbio= 0.1,
tvlag = 0.14,
tvfmet = 0.7,
tvclm_M = 100,
tvclm_F = 100,
tvvm = 5,
Ω = Diagonal([0.1, 0.1, 0.1,0.1,0.0001]),
σ²_prop_cp = 0.1,
σ²_prop_met = 0.1
)

chrisv · November 16, 2023, 12:32pm

With this fit: fit_results = fit(model20, pop, iparams, Pumas.FOCE(),constantcoef = (tvfmet = 0.7, ))

I get this error:
ERROR: gradient is exactly zero in tvlag. This indicates that tvlag isn’t identified.

So I fix tvlag. With this fit: fit_results = fit(model20, pop, iparams, Pumas.FOCE(),constantcoef = (tvfmet = 0.7,tvlag=0.14, ))

I get this error: ERROR: gradient is exactly zero in tvbio. This indicates that tvbio isn’t identified.

And so on parameter by parameter. I can never get the fit to work.

parssh.mehta · November 16, 2023, 9:47pm

Hi Chris,

There are a few parameters that need to come under the @dosecontrol block and not under @pre

@dosecontrol begin
  lags        = (SR = tvlag,)
  bioav       = (;IR = tvbio*exp(η[5]), SR = 1-tvbio*exp(η[5]))
end

This should help fix your problem. You can also read more about the @dosecontrol block Docs-Dose Control Parameters

Best,

Parssh

chrisv · November 27, 2023, 11:42pm

Thank you very much for your help. I have made that change and I have also simplified the model a little, the code is now a little simpler, here it can be seen:

model24 = @model begin

@param begin
  tvka1            ∈ RealDomain(lower=0)
  tvka2            ∈ RealDomain(lower=0)
  tvclp          ∈ RealDomain(lower=0) #Clearance of parent drug
  tvvc            ∈ RealDomain(lower=0) #Volume of central compartment (parent drug)
  tvq           ∈ RealDomain(lower=0) 
  tvvp            ∈ RealDomain(lower=0) #Volume of peripheral compartment (parent drug)
  tvvm          ∈ RealDomain(lower=0) #Volume of Compartment of Metabolite
  tvclm           ∈ RealDomain(lower=0)
  tvlag           ∈ RealDomain(lower=0)
  tvbio ∈ RealDomain(lower=0)
   Ω              ∈ PDiagDomain(5)
   σ_cp ∈ RealDomain(lower = 0.0001)
   σ_met ∈ RealDomain(lower = 0.0001)
  
end

@random begin
  η               ~ MvNormal(Ω)
end

@pre begin
Vc = tvvc * exp(η[1])
Clp = tvclp * exp(η[2])
Q = tvq
Vp = tvvp
Vm = tvvm
Clm = tvclm * exp(η[3])
Ka1 = tvka1
Ka2 = tvka2 * exp(η[4])
end

@dosecontrol begin
  lags          = (SR = tvlag,)
  bioav = (IR = tvbio*exp(η[5]), SR = 1-tvbio*exp(η[5]))
end

@dynamics begin
  IR'      = -Ka1*IR
  SR'      = -Ka2*SR
  Central'        =  Ka1*IR + Ka2*SR - (Clp/Vc)*Central - (Q/Vc)*Central + (Q/Vp)*Peripheral
  Peripheral'     = (Q/Vc)*Central - (Q/Vp)*Peripheral
  Metabolite'     = (Clp/Vc)*Central - (Clm/Vm)*Metabolite
end

@derived begin
  cp              = @. Central/Vc
  dv_cp ~ @. Normal(cp, σ_cp*abs(cp))
  met             = @. Metabolite/Vm
  dv_met ~ @. Normal(met, σ_met*abs(met))
end

end

iparams = (
tvka1 = 0.5,
tvka2 = 10,
tvclp = 40,
tvvc = 40,
tvq = 20,
tvvp = 160,
tvbio= 0.5,
tvlag = 0.15,
tvclm = 100,
tvvm = 5,
Ω = Diagonal([0.1, 0.1, 0.1,0.1,0.1]),
σ_cp = 0.1,
σ_met = 0.1
)

fitted_model = fit(model24, pop, iparams, FOCE())

But in this time I get a stranger error: “TaskedFailedException”. And no additional explanation.

chrisv · November 28, 2023, 12:04am

I found the error explanation in the terminal:
nested task error: DomainError with -Inf:
Duration must be non-negative. Got a duration of -Inf and a bioavailability of -0.0036060595751163937 for compartment 2.

But the initial parameters are okay. And I understand that the differential equations are okay also.

I don’t understand why I get this domain error.

Thank you for help !

storopoli · November 28, 2023, 12:22pm

Probably because

SR = 1-tvbio*exp(η[5]))

can be negative if tvbio*exp(η[5])) > 1.

Why don’t you try a NaivePool fit first to see if the error is with the η[5] interaction?
Or even better, why not a FOCE fit without an eta on tvbio? Reduce the Ω ∈ PDiagDomain(4) and see if the error still happens.

Topic		Replies	Views
Dual absorption 2 CMT model fitting error	3	125	November 28, 2023
PK model with transit compartments for absorption	4	365	January 31, 2022
Error while running working example	8	701	November 12, 2019
Parent-metab modeling	3	159	September 1, 2023
Model fitting giving error: gradient is exactly zero in θtd50 Simulation	4	277	June 29, 2022

Parent-metabolite + dual absorption: fitting error

Related Topics