From: "Eleveld, DJ" d.j.eleveld@anest.umcg.nl
Subject: [NMusers] Modelling twitch potentiation in NONMEM, problems with $DES
Date: Thu, March 10, 2005 6:22 am
Hello everyone,
I am wondering if anyone can give me advice or insight into using
NONMEM for modelling twitch potentiation along with a pharmacodynamic
model. I am trying to use NONMEM to model the method described in "Model
to describe the degree of twitch potentiation during neuromuscular
monitoring, D. J. Eleveld, A. F. Kopman, J. H. Proost, J. M. K. H. Wierda,
British journal of Anaesthesia, Br J Anaesth 2004; 92: 373-80".
The basic idea is that the degree of twitch potentiation increases as a
result of a twitch and decreases exponentially in the inter-twitch period.
I have been using a compartment to model the degree of potentiation. Because
I need to increase the degree of potentiation step-wise immediately after
each twitch observation, I have used another compartment to act as a
counter with the twitch observation count.
The fitting is related to the standard LN-PD model fitting for muscle
relaxants. It is slightly more complex because the twitch height calculation
has an additional factor varies dependant on an 'internal' variable, the NMB.
What I mean by 'internal' here is that NMB is dependent on Keo, which I am
also fitting. It seems to be logical to model the degree of potentiation
with a compartment, but then I have to make sure this compartment gets a
'dose' at every twitch observation. I can't calculate the 'dose' size
before hand because it depends on the NMB, which depends on the Keo, which
I am also fitting. So I think I have to incorporate the 'dose' into the
potentiation compartment within the $DES code. This doesn't seem to be the
usual way of writing $DES code, but I can't think of any other way to get it to work.
There is the slight added complication that there are some obviously incorrect
twitches removed from the data file, but these twitches cannot simply be
ignored as they do include a twitch and therefore influence the time course
of potentiation. These have been marked as 999 in the data file, a number
that otherwise does not occur.
My control file looks like:
$PROB Potentiation LN-PD fitting
$DATA data2.prn
IGNORE=C
$INPUT ID TIME DV PSTP PCON
$SUBROUTINES ADVAN6 TOL=6
$MODEL COMP(EFFECT) ; NMB effect compartment
COMP(POTENT) ; Potentiation compartment
COMP(HACK) ; Ugly (but necessary) twitch counter
$PK
KEO=THETA(1)*EXP(ETA(1)) ; Link Keo
EC50=THETA(2)*EXP(ETA(2)) ; PD Ec50
GAMM=THETA(3)*EXP(ETA(3)) ; PD Gamma
SCAL=THETA(4)+ETA(4) ; Potentiation model Kscale
POTR=ABS(THETA(5)+ETA(5)) ; Potentiation model r must be +ve
POTK=THETA(6)*EXP(ETA(6)) ; Potentiation model k
$DES
DADT(1)=KEO*(PCON-A(1)) ; Effect site concentration decay
PDT1=A(1)**GAMM ; Calculate NMB
NMB=PDT1/(PDT1+(EC50**GAMM))
PINC=POTR*(1-NMB) ; Increase potentiation
IF (A(3).EQ.PSTP) PINC=0 ; only once per twitch
DADT(2)=-POTK*A(2)+PINC ; The degree of potentiation (P)
DADT(3)=-A(3)+PSTP ; Keep twicth count flag correct
$ERROR
EDT1=A(1)**GAMM ; NMB must be diff names than $DES, sigh....
ENMB=EDT1/(EDT1+(EC50**GAMM))
BASE=100*(1+A(2))
Y=SCAL*(BASE/100)*(1-ENMB)+ERR(1) ; Predicted twitch response
TWIT=Y
IF (DV.EQ.999) Y=999+ERR(1); Missing twitches are fit almost excatly
$THETA (0.05,,0.25)(1000.,,2500.)(2.,,6.)(90.,,120.)
(0.001,,0.02)(0.01,,0.5)
$OMEGA 0.5 0.5 0.5 10 0.01 0.5
$SIGMA 20
This way to model increases in potentiation seems logical, if a bit awkward. The amount
in compartment 2 is the degree of twitch potentiation and compartment 3 functions as a
'counter' along with the PSTP record (different for each observation) to make sure
potentiation only increases once per observation even if the $DES code is called more
than once per observation.
Possibly my assumptions about how the $DES code is called are incorrect. Does anyone
see any obvious problems with my control file or reasons why it might not work as expected?
Thank you,
Doug Eleveld
_______________________________________________________