Dear NMuser
I am trying to analyze time to repeated event data using NONMEM.
The response were obtained till 24 hours after drug administration.
Inhibitory Emax model was implemented.
I am getting unreasonable parameter estimates which is far beyond what data
say.
If some body can point out what i am doing wrong, it would be very helpful.
Thank you in advance.
Hyewon
Data set (# of observations =62, # of patients=50 )
C ID TIME CS MDV AUC
. 101 0 . 1 1.111
. 101 0.05 0 0 1.111
. 101 2 0 0 1.111
. 102 0 . 1 0
. 102 24 1 0 0
. 103 0 . 1 0.999
. 103 0.75 . 0 0.999
....
Model File
$PROB RUN# 101
$INPUT C ID TIME CS MDV AUC
;CS:0=having event,1=censored
$DATA .data.csv IGNORE=C
$SUBROUTINE ADVAN=6 TOL=6
$MODEL
COMP=(HAZARD)
$PK
LAMD=THETA(1)*EXP(ETA(1)) ;scale factor
ALPH=THETA(2) ;shape factor
EC=THETA(3) ;AUC when effect is half of
its max
EFF=1-AUC/(EC+AUC) ;drug effect
$DES
DEL=1E-6
DADT(1)=LAMD*ALPH*(TIME+DEL)**(ALPH-1)*EXP(-LAMD*(TIME+DEL)**ALPH)*EFF
$ERROR
DEL2=1E-6
IF(NEWIND.NE.2) OLDCHZ=0
CHZ=A(1)-OLDCHZ
OLDCHZ=A(1)
SUR=EXP(-CHZ)
HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*EXP(-LAMD*(TIME+DEL2)**ALPH)*EFF
IF(CS.EQ.1) Y=SUR ;survival prob of censored
IF(CS.EQ.0) Y=SUR*HAZNOW ;pdf of event
$THETA
(0,10) ;[scale factor]
(0,1.0) ;[shaph factor]
(0,0.01) ;[AUC at half of Emax]
$OMEGA
0.01 ;[p] omega(1,1)
$EST METHOD=COND LIKE LAPLACIAN PRINT=5 SIG=3 MAX=9999 MSFO=101.msf
$COV PRINT=E
$TABLE ID TIME SUR AUC ONEHEADER NOPRINT FILE=101.tab
Time-to-event analysis with $DES
Hi Hyewon
My 2 cents of advice:
1. Note that you need to estimate the baseline hazard. You can then add
AUC as a covariate on the baseline hazard, or perhaps the Weibull shape
factor.
2. Because of my limited experience, I am not sure of the type of
survival distribution you are trying to model, but if it is a Weibull
distribution, I would write the code as follows:
$PK
IF(NEWIND.LE.1) OT=0 ; OT is zero at start for each individual
HAZBASE=THETA(X)*EXP(ETA(X)) ; baseline hazard for RTTE
ALPH=THETA(Y) ; weibull shape parameter
EFF=HAZBASE*(1-AUC/(EC+AUC) ; Imax model on hazbase (or on ALPH if you
want)
etc
$DES
DADT(1) = EFF*ALPH*(EFF*(T-OT+DEL))**(ALPH-1) ; HAZARD CHANGES WITH TIME
$ERROR
OT=TIME ; so that NONMEM remembers this value to subtract at the next
record in $DES
etc
Hope that helps.
Regards,
Emmanuel
>>> Hyewon Kim 03/24/11 1:18 PM >>>
Dear NMuser
I am trying to analyze time to repeated event data using NONMEM.
The response were obtained till 24 hours after drug administration.
Inhibitory Emax model was implemented.
I am getting unreasonable parameter estimates which is far beyond what
data say.
If some body can point out what i am doing wrong, it would be very
helpful.
Thank you in advance.
Hyewon
Data set (# of observations =62, # of patients=50 )
C ID TIME CS MDV AUC
. 101 0 . 1 1.111
. 101 0.05 0 0 1.111
. 101 2 0 0 1.111
. 102 0 . 1 0
. 102 24 1 0 0
. 103 0 . 1 0.999
. 103 0.75 . 0 0.999
....
Model File
$PROB RUN# 101
$INPUT C ID TIME CS MDV AUC
;CS:0=having event,1=censored
$DATA .data.csv IGNORE=C
$SUBROUTINE ADVAN=6 TOL=6
$MODEL
COMP=(HAZARD)
$PK
LAMD=THETA(1)*EXP(ETA(1)) ;scale factor
ALPH=THETA(2) ;shape factor
EC=THETA(3) ;AUC when effect is half of its max
EFF=1-AUC/(EC+AUC) ;drug effect
$DES
DEL=1E-6
DADT(1)=LAMD*ALPH*(TIME+DEL)**(ALPH-1)*EXP(-LAMD*(TIME+DEL)**ALPH)*EFF
$ERROR
DEL2=1E-6
IF(NEWIND.NE.2) OLDCHZ=0
CHZ=A(1)-OLDCHZ
OLDCHZ=A(1)
SUR=EXP(-CHZ)
HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*EXP(-LAMD*(TIME+DEL2)**ALPH)*EFF
IF(CS.EQ.1) Y=SUR ;survival prob of censored
IF(CS.EQ.0) Y=SUR*HAZNOW ;pdf of event
$THETA
(0,10) ;[scale factor]
(0,1.0) ;[shaph factor]
(0,0.01) ;[AUC at half of Emax]
$OMEGA
0.01 ;[p] omega(1,1)
$EST METHOD=COND LIKE LAPLACIAN PRINT=5 SIG=3 MAX=9999 MSFO=101.msf
$COV PRINT=E
$TABLE ID TIME SUR AUC ONEHEADER NOPRINT FILE=101.tab
###
UNIVERSITY OF CAPE TOWN
This e-mail is subject to the UCT ICT policies and e-mail disclaimer
published on our website at
http://www.uct.ac.za/about/policies/emaildisclaimer/ or obtainable from
+27 21 650 9111. This e-mail is intended only for the person(s) to whom
it is addressed. If the e-mail has reached you in error, please notify
the author. If you are not the intended recipient of the e-mail you may
not use, disclose, copy, redirect or print the content. If this e-mail
is not related to the business of UCT it is sent by the sender in the
sender's individual capacity.
###
Hyewon,
The most obvious problem with your code is in $DES. You must use the variable T not the variable TIME when referring to a time varying hazard. TIME is the time on each data record. It does not change in $DES. T is the time from the last data record upto the current data record and changes within $DES.
You are predicting the likelihood of an event at the exact time of the event observation record with multiple events per subject. As you seem to realize you need tocompute the cumulative hazard (CUMHAZ) either from TIME=0 or from the TIME of the last non-censored event for each subject.
In my opinion the following code is clearer and less dependent on your data structure than the method you are using which works only if your data has just event records after the TIME=0 record.
IF (MDV.EQ.0.AND.CS.EQ.0) THEN
OLDCHZ=A(1) ; cum haz upto time of this event
ELSE
OLDCHZ=OLDCHZ ; need to do this if OLDCHZ is a random variable
ENDIF
Your hazard model looks rather complicated. It seems to be based on the product of a Weibull baseline hazard
LAMD*ALPH*(TIME+DEL2)**(ALPH-1)
then something odd involving the Weibull parameters
*EXP(-LAMD*(TIME+DEL2)**ALPH)
and then forces the hazard to be zero if AUC is zero.
*EFF
Is this really what you want? Do you know the hazard of event is zero if AUC is zero? Or is the last right parenthesis in the wrong place?
I would suggest something like this where LAMD and ALPH are the two parameters of the Weibull baseline hazard (when AUC is zero) and BETAE is a parameter describing the effect of the drug on the overall hazard.
HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*EXP(BETAE*EFF)
You may also find it easier to develop the model if you do not try to estimate a random effect on LAMD until you have got reasonable estimates for the other parameters.
You may find it helpful to look at this presentation showing how to code time to event models in NM-TRAN:
http://pkpdrx.com/holford/docs/time-to-event-analysis.pdf
Best wishes,
Nick
Quoted reply history
On 24/03/2011 3:25 a.m., Hyewon Kim wrote:
> Dear NMuser
>
> I am trying to analyze time to repeated event data using NONMEM.
> The response were obtained till 24 hours after drug administration.
> Inhibitory Emax model was implemented.
>
> I am getting unreasonable parameter estimates which is far beyond what data say. If some body can point out what i am doing wrong, it would be very helpful.
>
> Thank you in advance.
>
> Hyewon
>
> Data set (# of observations =62, # of patients=50 )
> C ID TIME CS MDV AUC
> . 101 0 . 1 1.111
> . 101 0.05 0 0 1.111
> . 101 2 0 0 1.111
> . 102 0 . 1 0
> . 102 24 1 0 0
> . 103 0 . 1 0.999
> . 103 0.75 . 0 0.999
> ....
>
> Model File
> $PROB RUN# 101
> $INPUT C ID TIME CS MDV AUC
> ;CS:0=having event,1=censored
> $DATA .data.csv IGNORE=C
> $SUBROUTINE ADVAN=6 TOL=6
> $MODEL
> COMP=(HAZARD)
> $PK
> LAMD=THETA(1)*EXP(ETA(1)) ;scale factor
> ALPH=THETA(2) ;shape factor
>
> EC=THETA(3) ;AUC when effect is half of its max
>
> EFF=1-AUC/(EC+AUC) ;drug effect
>
> $DES
> DEL=1E-6
> DADT(1)=LAMD*ALPH*(TIME+DEL)**(ALPH-1)*EXP(-LAMD*(TIME+DEL)**ALPH)*EFF
>
> $ERROR
> DEL2=1E-6
> IF(NEWIND.NE.2) OLDCHZ=0
> CHZ=A(1)-OLDCHZ
> OLDCHZ=A(1)
> SUR=EXP(-CHZ)
> HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*EXP(-LAMD*(TIME+DEL2)**ALPH)*EFF
>
> IF(CS.EQ.1) Y=SUR ;survival prob of censored
> IF(CS.EQ.0) Y=SUR*HAZNOW ;pdf of event
>
> $THETA
> (0,10) ;[scale factor]
> (0,1.0) ;[shaph factor]
> (0,0.01) ;[AUC at half of Emax]
>
> $OMEGA
> 0.01 ;[p] omega(1,1)
>
> $EST METHOD=COND LIKE LAPLACIAN PRINT=5 SIG=3 MAX=9999 MSFO=101.msf
>
> $COV PRINT=E
>
> $TABLE ID TIME SUR AUC ONEHEADER NOPRINT FILE=101.tab
--
Nick Holford, Professor Clinical Pharmacology
Dept Pharmacology& Clinical Pharmacology
University of Auckland,85 Park Rd,Private Bag 92019,Auckland,New Zealand
tel:+64(9)923-6730 fax:+64(9)373-7090 mobile:+64(21)46 23 53
email: [email protected]
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
Hi All
I would just like to point out that in the Weibull hazard that Nick
proposes:
HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)
The units of LAMD would have to be TIME**(-ALPH). That is, as soon as ALPH
changes then the units of LAMD change. This seems strange and unnecessary
to include this type of correlation in the parameters. If we instead
reparameterize the model to:
HAZNOW=LAMD*ALPH*(LAMD*(TIME+DEL2))**(ALPH-1)
Then LAMD has units of 1/TIME no matter the value of ALPH. I have tested
the two parameterizations in one simple simulation example and found that
the relative standard errors were significantly lower with this second
parameterization.
Best regards
Andy
Andrew Hooker, Ph.D.
Associate Professor of Pharmacometrics
Dept. of Pharmaceutical Biosciences
Uppsala University
Box 591, 751 24, Uppsala, Sweden
Phone: +46 18 471 4355
http://www.farmbio.uu.se/research/Research+groups/Pharmacometrics/
www.farmbio.uu.se/research/Research+groups/Pharmacometrics/
Quoted reply history
From: [email protected] [mailto:[email protected]] On
Behalf Of Nick Holford
Sent: den 24 mars 2011 07:54
To: nmusers
Subject: Re: [NMusers] Time-to-event analysis with $DES
Hyewon,
Sorry -- I just realized that EFF is not the drug effect but 1-drug effect
so that it has a value of 1 when AUC is zero. Please ignore my remarks about
EFF being zero when AUC is zero in my comments below and consider this code
to describe the effect of AUC on hazard:
HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*(1-BETAE*AUC/(EC+AUC))
I think it is a very strong assumption that the drug will reduce the hazard
to zero at infinite AUC which you can test by estimating a parameter (such
as BETAE). If BETAE is significantly less than 1 then this means the
assumption is not supported.
Nick
On 24/03/2011 8:17 a.m., Nick Holford wrote:
Hyewon,
The most obvious problem with your code is in $DES. You must use the
variable T not the variable TIME when referring to a time varying hazard.
TIME is the time on each data record. It does not change in $DES. T is the
time from the last data record upto the current data record and changes
within $DES.
You are predicting the likelihood of an event at the exact time of the event
observation record with multiple events per subject. As you seem to realize
you need to compute the cumulative hazard (CUMHAZ) either from TIME=0 or
from the TIME of the last non-censored event for each subject.
In my opinion the following code is clearer and less dependent on your data
structure than the method you are using which works only if your data has
just event records after the TIME=0 record.
IF (MDV.EQ.0.AND.CS.EQ.0) THEN
OLDCHZ=A(1) ; cum haz upto time of this event
ELSE
OLDCHZ=OLDCHZ ; need to do this if OLDCHZ is a random variable
ENDIF
Your hazard model looks rather complicated. It seems to be based on the
product of a Weibull baseline hazard
LAMD*ALPH*(TIME+DEL2)**(ALPH-1)
then something odd involving the Weibull parameters
*EXP(-LAMD*(TIME+DEL2)**ALPH)
and then forces the hazard to be zero if AUC is zero.
*EFF
Is this really what you want? Do you know the hazard of event is zero if AUC
is zero? Or is the last right parenthesis in the wrong place?
I would suggest something like this where LAMD and ALPH are the two
parameters of the Weibull baseline hazard (when AUC is zero) and BETAE is a
parameter describing the effect of the drug on the overall hazard.
HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*EXP(BETAE*EFF)
You may also find it easier to develop the model if you do not try to
estimate a random effect on LAMD until you have got reasonable estimates for
the other parameters.
You may find it helpful to look at this presentation showing how to code
time to event models in NM-TRAN:
http://pkpdrx.com/holford/docs/time-to-event-analysis.pdf
Best wishes,
Nick
On 24/03/2011 3:25 a.m., Hyewon Kim wrote:
Dear NMuser
I am trying to analyze time to repeated event data using NONMEM.
The response were obtained till 24 hours after drug administration.
Inhibitory Emax model was implemented.
I am getting unreasonable parameter estimates which is far beyond what data
say.
If some body can point out what i am doing wrong, it would be very helpful.
Thank you in advance.
Hyewon
Data set (# of observations =62, # of patients=50 )
C ID TIME CS MDV AUC
. 101 0 . 1 1.111
. 101 0.05 0 0 1.111
. 101 2 0 0 1.111
. 102 0 . 1 0
. 102 24 1 0 0
. 103 0 . 1 0.999
. 103 0.75 . 0 0.999
....
Model File
$PROB RUN# 101
$INPUT C ID TIME CS MDV AUC
;CS:0=having event,1=censored
$DATA .data.csv IGNORE=C
$SUBROUTINE ADVAN=6 TOL=6
$MODEL
COMP=(HAZARD)
$PK
LAMD=THETA(1)*EXP(ETA(1)) ;scale factor
ALPH=THETA(2) ;shape factor
EC=THETA(3) ;AUC when effect is half of
its max
EFF=1-AUC/(EC+AUC) ;drug effect
$DES
DEL=1E-6
DADT(1)=LAMD*ALPH*(TIME+DEL)**(ALPH-1)*EXP(-LAMD*(TIME+DEL)**ALPH)*EFF
$ERROR
DEL2=1E-6
IF(NEWIND.NE.2) OLDCHZ=0
CHZ=A(1)-OLDCHZ
OLDCHZ=A(1)
SUR=EXP(-CHZ)
HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*EXP(-LAMD*(TIME+DEL2)**ALPH)*EFF
IF(CS.EQ.1) Y=SUR ;survival prob of censored
IF(CS.EQ.0) Y=SUR*HAZNOW ;pdf of event
$THETA
(0,10) ;[scale factor]
(0,1.0) ;[shaph factor]
(0,0.01) ;[AUC at half of Emax]
$OMEGA
0.01 ;[p] omega(1,1)
$EST METHOD=COND LIKE LAPLACIAN PRINT=5 SIG=3 MAX=9999 MSFO=101.msf
$COV PRINT=E
$TABLE ID TIME SUR AUC ONEHEADER NOPRINT FILE=101.tab
--
Nick Holford, Professor Clinical Pharmacology
Dept Pharmacology & Clinical Pharmacology
University of Auckland,85 Park Rd,Private Bag 92019,Auckland,New Zealand
tel:+64(9)923-6730 fax:+64(9)373-7090 mobile:+64(21)46 23 53
email: [email protected]
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
--
Nick Holford, Professor Clinical Pharmacology
Dept Pharmacology & Clinical Pharmacology
University of Auckland,85 Park Rd,Private Bag 92019,Auckland,New Zealand
tel:+64(9)923-6730 fax:+64(9)373-7090 mobile:+64(21)46 23 53
email: [email protected]
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
Andy,
Thanks for this response. In fact I did not propose that particular Weibull parameterization but simply pointed out that it seemed to be what the original questioner seemed to be trying to use.
My own preference for parameterizing the Weibull is
HAZARD=LAMDA0*EXP(BETA*LOG(TIME))
This parameterization makes is easier to understand that the Weibull (ln(time)) is just a transformation of the Gompertz (time):
HAZARD=LAMDA0*EXP(BETA*TIME)
Its then easy to extend hazards to more interesting and possibly more useful hazards such as the Gombull (aka the Weipertz):
HAZARD=LAMDA0*EXP(BETAG*TIME + BETAW*LOG(TIME))
One of the features of the Weibull that I dont like is that the hazard is predicted to be zero when time is zero. This seems a particularly unlikely biological occurrence and also requires extra coding to deal with this special case when time is zero to avoid NaNs.
I'd be interested in knowing if you find any difference in the standard errors with the ln(time) parameterization. But whatever you find it won't change my behaviour because as you know standard errors are not much use for making predictions :-)
Best wishes,
Nick
Quoted reply history
On 28/03/2011 3:05 p.m., Andrew Hooker wrote:
> Hi All
>
> I would just like to point out that in the Weibull hazard that Nick proposes:
>
> HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)
>
> The units of LAMD would have to be TIME**(-ALPH). That is, as soon as ALPH changes then the units of LAMD change. This seems strange and unnecessary to include this type of correlation in the parameters. If we instead reparameterize the model to:
>
> HAZNOW=LAMD*ALPH*(LAMD*(TIME+DEL2))**(ALPH-1)
>
> Then LAMD has units of 1/TIME no matter the value of ALPH. I have tested the two parameterizations in one simple simulation example and found that the relative standard errors were significantly lower with this second parameterization.
>
> Best regards
>
> Andy
>
> Andrew Hooker, Ph.D.
>
> Associate Professor of Pharmacometrics
>
> Dept. of Pharmaceutical Biosciences
>
> Uppsala University
>
> Box 591, 751 24, Uppsala, Sweden
>
> Phone: +46 18 471 4355
>
> www.farmbio.uu.se/research/Research+groups/Pharmacometrics/ < http://www.farmbio.uu.se/research/Research+groups/Pharmacometrics/ >
>
> *From:* [email protected] [ mailto: [email protected] ] *On Behalf Of *Nick Holford
>
> *Sent:* den 24 mars 2011 07:54
> *To:* nmusers
> *Subject:* Re: [NMusers] Time-to-event analysis with $DES
>
> Hyewon,
>
> Sorry -- I just realized that EFF is not the drug effect but 1-drug effect so that it has a value of 1 when AUC is zero. Please ignore my remarks about EFF being zero when AUC is zero in my comments below and consider this code to describe the effect of AUC on hazard:
>
> HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*(1-BETAE*AUC/(EC+AUC))
>
> I think it is a very strong assumption that the drug will reduce the hazard to zero at infinite AUC which you can test by estimating a parameter (such as BETAE). If BETAE is significantly less than 1 then this means the assumption is not supported.
>
> Nick
>
> On 24/03/2011 8:17 a.m., Nick Holford wrote:
>
> Hyewon,
>
> The most obvious problem with your code is in $DES. You must use the variable T not the variable TIME when referring to a time varying hazard. TIME is the time on each data record. It does not change in $DES. T is the time from the last data record upto the current data record and changes within $DES.
>
> You are predicting the likelihood of an event at the exact time of the event observation record with multiple events per subject. As you seem to realize you need to compute the cumulative hazard (CUMHAZ) either from TIME=0 or from the TIME of the last non-censored event for each subject.
>
> In my opinion the following code is clearer and less dependent on your data structure than the method you are using which works only if your data has just event records after the TIME=0 record.
>
> IF (MDV.EQ.0.AND.CS.EQ.0) THEN
> OLDCHZ=A(1) ; cum haz upto time of this event
> ELSE
> OLDCHZ=OLDCHZ ; need to do this if OLDCHZ is a random variable
> ENDIF
>
> Your hazard model looks rather complicated. It seems to be based on the product of a Weibull baseline hazard
>
> LAMD*ALPH*(TIME+DEL2)**(ALPH-1)
> then something odd involving the Weibull parameters
> *EXP(-LAMD*(TIME+DEL2)**ALPH)
> and then forces the hazard to be zero if AUC is zero.
> *EFF
>
> Is this really what you want? Do you know the hazard of event is zero if AUC is zero? Or is the last right parenthesis in the wrong place?
>
> I would suggest something like this where LAMD and ALPH are the two parameters of the Weibull baseline hazard (when AUC is zero) and BETAE is a parameter describing the effect of the drug on the overall hazard.
>
> HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*EXP(BETAE*EFF)
>
> You may also find it easier to develop the model if you do not try to estimate a random effect on LAMD until you have got reasonable estimates for the other parameters.
>
> You may find it helpful to look at this presentation showing how to code time to event models in NM-TRAN:
>
> http://pkpdrx.com/holford/docs/time-to-event-analysis.pdf
>
> Best wishes,
>
> Nick
>
> On 24/03/2011 3:25 a.m., Hyewon Kim wrote:
>
> Dear NMuser
>
> I am trying to analyze time to repeated event data using NONMEM.
> The response were obtained till 24 hours after drug administration.
> Inhibitory Emax model was implemented.
>
> I am getting unreasonable parameter estimates which is far beyond what data say. If some body can point out what i am doing wrong, it would be very helpful.
>
> Thank you in advance.
>
> Hyewon
>
> Data set (# of observations =62, # of patients=50 )
> C ID TIME CS MDV AUC
> . 101 0 . 1 1.111
> . 101 0.05 0 0 1.111
> . 101 2 0 0 1.111
> . 102 0 . 1 0
> . 102 24 1 0 0
> . 103 0 . 1 0.999
> . 103 0.75 . 0 0.999
> ....
>
> Model File
> $PROB RUN# 101
> $INPUT C ID TIME CS MDV AUC
> ;CS:0=having event,1=censored
> $DATA .data.csv IGNORE=C
> $SUBROUTINE ADVAN=6 TOL=6
> $MODEL
> COMP=(HAZARD)
> $PK
> LAMD=THETA(1)*EXP(ETA(1)) ;scale factor
> ALPH=THETA(2) ;shape factor
>
> EC=THETA(3) ;AUC when effect is half of its max
>
> EFF=1-AUC/(EC+AUC) ;drug effect
>
> $DES
> DEL=1E-6
> DADT(1)=LAMD*ALPH*(TIME+DEL)**(ALPH-1)*EXP(-LAMD*(TIME+DEL)**ALPH)*EFF
>
> $ERROR
> DEL2=1E-6
> IF(NEWIND.NE.2) OLDCHZ=0
> CHZ=A(1)-OLDCHZ
> OLDCHZ=A(1)
> SUR=EXP(-CHZ)
> HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*EXP(-LAMD*(TIME+DEL2)**ALPH)*EFF
>
> IF(CS.EQ.1) Y=SUR ;survival prob of censored
> IF(CS.EQ.0) Y=SUR*HAZNOW ;pdf of event
>
> $THETA
> (0,10) ;[scale factor]
> (0,1.0) ;[shaph factor]
> (0,0.01) ;[AUC at half of Emax]
>
> $OMEGA
> 0.01 ;[p] omega(1,1)
>
> $EST METHOD=COND LIKE LAPLACIAN PRINT=5 SIG=3 MAX=9999 MSFO=101.msf
>
> $COV PRINT=E
>
> $TABLE ID TIME SUR AUC ONEHEADER NOPRINT FILE=101.tab
>
> --
> Nick Holford, Professor Clinical Pharmacology
> Dept Pharmacology& Clinical Pharmacology
> University of Auckland,85 Park Rd,Private Bag 92019,Auckland,New Zealand
> tel:+64(9)923-6730 <tel:+64%289%29923-6730> fax:+64(9)373-7090
> mobile:+64(21)46 23 53
> email:[email protected] <mailto:[email protected]>
> http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
>
> --
> Nick Holford, Professor Clinical Pharmacology
> Dept Pharmacology& Clinical Pharmacology
> University of Auckland,85 Park Rd,Private Bag 92019,Auckland,New Zealand
> tel:+64(9)923-6730 <tel:+64%289%29923-6730> fax:+64(9)373-7090
> mobile:+64(21)46 23 53
> email:[email protected] <mailto:[email protected]>
> http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
--
Nick Holford, Professor Clinical Pharmacology
Dept Pharmacology& Clinical Pharmacology
University of Auckland,85 Park Rd,Private Bag 92019,Auckland,New Zealand
tel:+64(9)923-6730 fax:+64(9)373-7090 mobile:+64(21)46 23 53
email: [email protected]
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
Hi Nick, Andy,
The Weibull hazard is not necessarily my first choice either. However, I
should point out that it is not a much of an issue that the Weibull predicts
a zero hazard at time 0 (only defined for ALP>1). One typically assumes the
survival fraction is = 1 at time = 0, so it doesn't matter what the hazard
is for any model at time =0. Also, If Beta = 0 (in yours or ALP=1 in
Andy's) it is well known that the Weibull = Exponential model, which has a
constant hazard throughout. Note that if 0<ALP<1 or -1<Beta<0 the hazard
can be rather large (not near 0) for small times and is decreasing, again
indicating the value at time = 0 is not much of an issue. Perhaps I am
missing something, but for the Gombull or Weipertz, isn't the hazard 0 at
time = 0 as well (rewriting HAZARD = LAMBD0 * EXP(BETAG*TIME)*(TIME)^BETAW
[for BETAW>0])?
Best regards,
Matt
Quoted reply history
From: [email protected] [mailto:[email protected]] On
Behalf Of Nick Holford
Sent: Monday, March 28, 2011 9:35 AM
To: Andrew Hooker
Cc: 'nmusers'
Subject: Re: [NMusers] Time-to-event analysis with $DES
Andy,
Thanks for this response. In fact I did not propose that particular Weibull
parameterization but simply pointed out that it seemed to be what the
original questioner seemed to be trying to use.
My own preference for parameterizing the Weibull is
HAZARD=LAMDA0*EXP(BETA*LOG(TIME))
This parameterization makes is easier to understand that the Weibull
(ln(time)) is just a transformation of the Gompertz (time):
HAZARD=LAMDA0*EXP(BETA*TIME)
Its then easy to extend hazards to more interesting and possibly more useful
hazards such as the Gombull (aka the Weipertz):
HAZARD=LAMDA0*EXP(BETAG*TIME + BETAW*LOG(TIME))
One of the features of the Weibull that I dont like is that the hazard is
predicted to be zero when time is zero. This seems a particularly unlikely
biological occurrence and also requires extra coding to deal with this
special case when time is zero to avoid NaNs.
I'd be interested in knowing if you find any difference in the standard
errors with the ln(time) parameterization. But whatever you find it won't
change my behaviour because as you know standard errors are not much use for
making predictions :-)
Best wishes,
Nick
On 28/03/2011 3:05 p.m., Andrew Hooker wrote:
Hi All
I would just like to point out that in the Weibull hazard that Nick
proposes:
HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)
The units of LAMD would have to be TIME**(-ALPH). That is, as soon as ALPH
changes then the units of LAMD change. This seems strange and unnecessary
to include this type of correlation in the parameters. If we instead
reparameterize the model to:
HAZNOW=LAMD*ALPH*(LAMD*(TIME+DEL2))**(ALPH-1)
Then LAMD has units of 1/TIME no matter the value of ALPH. I have tested
the two parameterizations in one simple simulation example and found that
the relative standard errors were significantly lower with this second
parameterization.
Best regards
Andy
Andrew Hooker, Ph.D.
Associate Professor of Pharmacometrics
Dept. of Pharmaceutical Biosciences
Uppsala University
Box 591, 751 24, Uppsala, Sweden
Phone: +46 18 471 4355
www.farmbio.uu.se/research/Research+groups/Pharmacometrics/
From: [email protected] [mailto:[email protected]] On
Behalf Of Nick Holford
Sent: den 24 mars 2011 07:54
To: nmusers
Subject: Re: [NMusers] Time-to-event analysis with $DES
Hyewon,
Sorry -- I just realized that EFF is not the drug effect but 1-drug effect
so that it has a value of 1 when AUC is zero. Please ignore my remarks about
EFF being zero when AUC is zero in my comments below and consider this code
to describe the effect of AUC on hazard:
HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*(1-BETAE*AUC/(EC+AUC))
I think it is a very strong assumption that the drug will reduce the hazard
to zero at infinite AUC which you can test by estimating a parameter (such
as BETAE). If BETAE is significantly less than 1 then this means the
assumption is not supported.
Nick
On 24/03/2011 8:17 a.m., Nick Holford wrote:
Hyewon,
The most obvious problem with your code is in $DES. You must use the
variable T not the variable TIME when referring to a time varying hazard.
TIME is the time on each data record. It does not change in $DES. T is the
time from the last data record upto the current data record and changes
within $DES.
You are predicting the likelihood of an event at the exact time of the event
observation record with multiple events per subject. As you seem to realize
you need to compute the cumulative hazard (CUMHAZ) either from TIME=0 or
from the TIME of the last non-censored event for each subject.
In my opinion the following code is clearer and less dependent on your data
structure than the method you are using which works only if your data has
just event records after the TIME=0 record.
IF (MDV.EQ.0.AND.CS.EQ.0) THEN
OLDCHZ=A(1) ; cum haz upto time of this event
ELSE
OLDCHZ=OLDCHZ ; need to do this if OLDCHZ is a random variable
ENDIF
Your hazard model looks rather complicated. It seems to be based on the
product of a Weibull baseline hazard
LAMD*ALPH*(TIME+DEL2)**(ALPH-1)
then something odd involving the Weibull parameters
*EXP(-LAMD*(TIME+DEL2)**ALPH)
and then forces the hazard to be zero if AUC is zero.
*EFF
Is this really what you want? Do you know the hazard of event is zero if AUC
is zero? Or is the last right parenthesis in the wrong place?
I would suggest something like this where LAMD and ALPH are the two
parameters of the Weibull baseline hazard (when AUC is zero) and BETAE is a
parameter describing the effect of the drug on the overall hazard.
HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*EXP(BETAE*EFF)
You may also find it easier to develop the model if you do not try to
estimate a random effect on LAMD until you have got reasonable estimates for
the other parameters.
You may find it helpful to look at this presentation showing how to code
time to event models in NM-TRAN:
http://pkpdrx.com/holford/docs/time-to-event-analysis.pdf
Best wishes,
Nick
On 24/03/2011 3:25 a.m., Hyewon Kim wrote:
Dear NMuser
I am trying to analyze time to repeated event data using NONMEM.
The response were obtained till 24 hours after drug administration.
Inhibitory Emax model was implemented.
I am getting unreasonable parameter estimates which is far beyond what data
say.
If some body can point out what i am doing wrong, it would be very helpful.
Thank you in advance.
Hyewon
Data set (# of observations =62, # of patients=50 )
C ID TIME CS MDV AUC
. 101 0 . 1 1.111
. 101 0.05 0 0 1.111
. 101 2 0 0 1.111
. 102 0 . 1 0
. 102 24 1 0 0
. 103 0 . 1 0.999
. 103 0.75 . 0 0.999
....
Model File
$PROB RUN# 101
$INPUT C ID TIME CS MDV AUC
;CS:0=having event,1=censored
$DATA .data.csv IGNORE=C
$SUBROUTINE ADVAN=6 TOL=6
$MODEL
COMP=(HAZARD)
$PK
LAMD=THETA(1)*EXP(ETA(1)) ;scale factor
ALPH=THETA(2) ;shape factor
EC=THETA(3) ;AUC when effect is half of
its max
EFF=1-AUC/(EC+AUC) ;drug effect
$DES
DEL=1E-6
DADT(1)=LAMD*ALPH*(TIME+DEL)**(ALPH-1)*EXP(-LAMD*(TIME+DEL)**ALPH)*EFF
$ERROR
DEL2=1E-6
IF(NEWIND.NE.2) OLDCHZ=0
CHZ=A(1)-OLDCHZ
OLDCHZ=A(1)
SUR=EXP(-CHZ)
HAZNOW=LAMD*ALPH*(TIME+DEL2)**(ALPH-1)*EXP(-LAMD*(TIME+DEL2)**ALPH)*EFF
IF(CS.EQ.1) Y=SUR ;survival prob of censored
IF(CS.EQ.0) Y=SUR*HAZNOW ;pdf of event
$THETA
(0,10) ;[scale factor]
(0,1.0) ;[shaph factor]
(0,0.01) ;[AUC at half of Emax]
$OMEGA
0.01 ;[p] omega(1,1)
$EST METHOD=COND LIKE LAPLACIAN PRINT=5 SIG=3 MAX=9999 MSFO=101.msf
$COV PRINT=E
$TABLE ID TIME SUR AUC ONEHEADER NOPRINT FILE=101.tab
--
Nick Holford, Professor Clinical Pharmacology
Dept Pharmacology & Clinical Pharmacology
University of Auckland,85 Park Rd,Private Bag 92019,Auckland,New Zealand
tel:+64(9)923-6730 <tel:+64%289%29923-6730> fax:+64(9)373-7090
mobile:+64(21)46 23 53
email: [email protected]
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
--
Nick Holford, Professor Clinical Pharmacology
Dept Pharmacology & Clinical Pharmacology
University of Auckland,85 Park Rd,Private Bag 92019,Auckland,New Zealand
tel:+64(9)923-6730 <tel:+64%289%29923-6730> fax:+64(9)373-7090
mobile:+64(21)46 23 53
email: [email protected]
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
--
Nick Holford, Professor Clinical Pharmacology
Dept Pharmacology & Clinical Pharmacology
University of Auckland,85 Park Rd,Private Bag 92019,Auckland,New Zealand
tel:+64(9)923-6730 fax:+64(9)373-7090 mobile:+64(21)46 23 53
email: [email protected]
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford