Jump to content


Photo

self-define model structure with if, else

time nlme if else regressor

  • Please log in to reply
8 replies to this topic

#1 linda

linda

    Member

  • Members
  • PipPip
  • 11 posts

Posted 11 January 2024 - 03:58 PM

Hi Everyone,

 

I am trying to code a self-define model structure with ifelse statements in Phoenix nlme for a pop Pk model.

 

I would like to have a statement like this:

 

t_0=10
Ac_0=y0
 
if t < Tin
  ddt_Ac = k * Ac - kel * Ac
else 
  ddt_Ac = - kel * Ac
 
Q1: Ac is target PK concentration, starting at day 10, and initial value should be fixed as a regression variable with subject ID which could be read from inputs.
How should I code to define start time and its subject dependent initial value?
 
k, kel and Tin are parameters needed to be fitted.
 
Q2: How to code if else statement with a uncertain condition (Tin needs to be fitted)?
 
I wrote a partial code for this model structure, please give me some suggestions for these two questions.
 

test(){

covariate(Time)

 

stparm (Ac = Time=10? Ac=y0)

stparm(Dev = Time<Tin ? k*Ac-kel*Ac:  -kel*Ac)

deriv(Ac = Dev)

 

error(AcEps = 1)

observe(AcObs(Time) = Ac+ AcEps )

 

fixef(tvTin= c(,1,))

fixef(tvk = c(, 1, ))

fixef(tvkel = c(, 1, ))

 

stparm(Tin = tvTin * exp(nTin))

stparm(kel = tvkel * exp(nkel))

stparm(k = tvk* exp(nk))

 

ranef(diag(nTin, nk, nkel) = c(1, 1, 1))

}

 

 

Best,

Linda

 

 


  • DanielGon, Williambus, LavillKaH and 12 others like this

#2 bwendt@certara.com

bwendt@certara.com

    Advanced Member

  • Administrators
  • 282 posts

Posted 12 January 2024 - 03:47 PM

Hi Linda,

 

there is something a bit unclear about your model:

 

you are defining Dev as structural parameter:

 

stparm(Dev = Time<Tin ? k*Ac-kel*Ac:  -kel*Ac)

 

This should actually be a deriv() statement, since you are coupling concentrations (Ac) with micro constants (k, kel)

 

Should Ac be fixed until Time=10 ?

When you have a first order input (k*Ac) what is the dose?

Perhaps, you can provide more details?

 

Thanks,

Bernd


  • Thomasgaks likes this

#3 linda

linda

    Member

  • Members
  • PipPip
  • 11 posts

Posted 12 January 2024 - 06:51 PM

Hi Bernd,

 

Thank you for your reply. I agree that I should use deriv(), instead of stparm()

The model is just a kinetics model, where is no dose. The model is used to fitting the data with two phases, production phase with k and decay phase with kel. I can solve the Q1 by chaning the initial time points. But I still need help for Q2, which is how to sturcture a two phase model.

 

I rewrote a code, but still couldn't run successfully.

 

test(){

covariate(Time)

 

deriv(Ac = Time<Tin ? k*Ac-kel*Ac:  -kel*Ac)

 

error(AcEps = 1)

observe(AcObs(Time) = Ac+ AcEps )

 

fixef(tvTin= c(,1,))

fixef(tvk = c(, 1, ))

fixef(tvkel = c(, 1, ))

 

stparm(Tin = tvTin * exp(nTin))

stparm(kel = tvkel * exp(nkel))

stparm(k = tvk* exp(nk))

 

ranef(diag(nTin, nk, nkel) = c(1, 1, 1))

}

 

 

Thanks,

Linda


Edited by linda, 12 January 2024 - 06:51 PM.


#4 bwendt@certara.com

bwendt@certara.com

    Advanced Member

  • Administrators
  • 282 posts

Posted 12 January 2024 - 07:45 PM

Hi Linda,

 

please try this as a first draft:

 

test(){
deriv(Ac = k*Ac-kel*Ac)
error(AcEps = 1)
observe(AcObs(t) = Ac+ AcEps )
sequence{Ac=10}
fixef(tvTin= c(,10,))
fixef(tvk = c(, 1.5, ))
fixef(tvkel = c(, 1, ))
stparm(Tin = tvTin * exp(nTin))
stparm(kel = tvkel * exp(nkel))
stparm(k = t<Tin? tvk* exp(nk):0)
ranef(diag(nTin, nk, nkel) = c(1, 1, 1))
}
 
It will produce a two-phase curve:
 
Ac.png
 
Bernd


#5 linda

linda

    Member

  • Members
  • PipPip
  • 11 posts

Posted 12 January 2024 - 08:02 PM

Hi Berna,

 

Thank you for your quick reply. The code can execute. What is sequence{Ac=10}, is that for set an initial value? I would like to set initial value at time 0 from data. Is that possible? I attached input dataset to make you more clear.

 

Thanks,

Linda


Edited by linda, 12 January 2024 - 08:02 PM.


#6 bwendt@certara.com

bwendt@certara.com

    Advanced Member

  • Administrators
  • 282 posts

Posted 12 January 2024 - 09:00 PM

Sequence{AC=0} sets the initial value for Ac, you are right.

 

You can read in values from a file with:

 

covariate(Ac0)

 

then you can read those in through:

 

sequence{Ac=Ac0}

 

Here is the output:

 

ac0.png



#7 linda

linda

    Member

  • Members
  • PipPip
  • 11 posts

Posted 15 January 2024 - 06:21 PM

Hi Berna,

 

Thank you so much. The initial question Q1 is solved. 

 

However, when I check the fitted population parameters, I found a problem:

I used Monolix for this fitting before. They can support if else statement as

if t < Tin
  ddt_Ac = k * Ac - kel * Ac
else 
  ddt_Ac = - kel * Ac
 
And I found our approach gets a smaller value of k. I am thinking whether because we have some timepoints setting k as 0. Whether we have another way for coding this?
 
Bests,
Linda


#8 bwendt@certara.com

bwendt@certara.com

    Advanced Member

  • Administrators
  • 282 posts

Posted 15 January 2024 - 06:34 PM

Hi Linda,

 

the implementation of if-else condition in Phoenix is just using a different syntax:

 

stparm(k = t<Tin? tvk* exp(nk):0)

 

means, that while time is smaller than Tin, the value of k is estimated, e.g. tvk. When time is greater or equal to Tin the value of k =0

 

When you look at the differential equation, the first term, k*Ac will be 0 when k=0:

 

deriv(Ac = k*Ac-kel*Ac)

 

Hope, this helps.

 

Bernd



#9 linda

linda

    Member

  • Members
  • PipPip
  • 11 posts

Posted 19 January 2024 - 07:02 PM

Hi Bernd,

 

Thank you for your answer. I think if we set k =0 for phase 1, I will effect the population estimation value of k because theotically there should be no k value on phase1. 

 

Thanks,

Linda







Also tagged with one or more of these keywords: time, nlme, if else, regressor

0 user(s) are reading this topic

0 members, 0 guests, 0 anonymous users