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 _______________________________________________________