Re: Time to event analysis
Gaurav,
1. Error Message
When I run your code I get this error from NM-TRAN.
460 $ESTIM: WITH LIKELIHOOD, $OMEGA MUST INCLUDE INITIAL ESTIMATE.
I get this with NONMEM VI, 7.1 and 7.2. I don't know how you got the message "laplace required with likelihood, Pop. data, unless omega=0". The error occurs because you do not have a $OMEGA record. This is because you have joined the $THETA and $OMEGA record like this:
$THETA (0,.025) ; BASE $OMEGA 0 FIX
If you separate the line as follows:
$THETA (0,.025) ; BASE
$OMEGA 0 FIX
then your code will run. The error message is certainly not very helpful in diagnosing this problem. This is something the NONMEM developers might want to look at.
2. Using NM-TRAN for human communication
One of the reasons for having an easy to use modelling language like NM-TRAN is the ability to use the language to communicate to others (and to yourself -- see below) what you are doing. While your code and data can be made to run it contains several 'features' that may make it harder for others to understand what you did.
The model does not actually use the time as the DV. The observation is either that an event occurred at a particular time or that it had not occurred at the last time and is therefore a censored event. So I would recommend you use the CS data item as the DV to make it clear that this shows the observation is an event or a censored event.
By convention, in the survival analysis literature censored events are given the value 0 and events are give the value 1. If you use this convention it makes it simpler for people who know about survival analysis to understand your data and your data could be used in another program that does survival analysis.
When you use the LIKE option in $EST it makes no difference which column you choose for DV because NONMEM does not use the value in this data item . The only thing that NONMEM uses is the likelihood returned in $ERROR when MDV=0. You may continue to use TIME as the DV but this is not the usual way to do this. If you use the DV to indicate the type of event then the code and data become easier to read for humans (IMHO - see below).
I don't know if it is an error in your data but all the values of CS in your code when MDV=0 have a value of 1 which means your model only computes the likelihood for censored events and therefore has little information about the hazard.
Your $MODEL record describes the hazard compartment as RISK. This is OK because this compartment is used to describe the cumulative hazard which is also known as risk. However in $ERROR you use
HZ = A(1) ;hazard
This is not correct. A(1) is the cumulative hazard not the hazard.
So it would be clearer to write:
RISK=A(1); cumulative hazard
The hazard in your model is BASE.
Your use of HZ as if it was the hazard then leads to an error in computing the likelihood of an event. The likelihood of an event at an exact time is calculated from the product of the survivor function and the hazard ie. SUR* BASE. You have used the product of SUR times the cumulative hazard which is incorrect.
Below you will find a suggestion for writing the code. I have also changed your data example so that DV is 0 for a censored event and 1 for an event. I have put some values of DV=1 and 0 to make this realistic. Note that when MDV=1 I have put a "." in the DV column to indicate to human readers that the DV is a missing value. Doing this also makes it easier to see what events are recorded in the data.
You can find more examples here:
http://www.pkpdrx.com/holford/docs/time-to-event-analysis.pdf
and here:
http://www.page-meeting.org/?abstract=2281
Note that the URL in your email associated with your reference to the nmusers archive was rather strange and caused my email client to flag your email as a possible scam. I have reformatted the link below in your message so that it no longer causes this behaviour.
Best wishes,
Nick
$PROB Time To Event data
$INPUT C STUD ID TIME DV DOSE MDV
$DATA hz1_nh.csv IGNORE=C
$SUBR ADVAN=6 TOL=6
$MODEL COMP=(RISK)
$PK
BASE = THETA(1)*EXP(ETA(1)) ;the ETA is a placeholder here
$DES
DADT(1)=BASE ;hazard
$ERROR
RISK = A(1) ; cumulative hazard
SUR = EXP(-RISK) ;survival probability
IF (DV.EQ.0) THEN ; censored event
Y=SUR
ELSE ; event at exact time
Y=SUR*BASE
ENDIF
$THETA (0,.025) ; BASE
$OMEGA 0 FIX
$ESTIM MAXEVAL=9990 METHOD = COND LAPLACE LIKE PRINT=1
;METHOD=ZERO works exactly the same way because there are no random effects in this model
;$ESTIM MAXEVAL=9990 METHOD = ZERO LIKE PRINT=1
MSFO=msfb01
$COV PRINT=E
------------- hz1_nh.csv ---------------------------------------------------------
C STUD ID TIME DV DOSE MDV
. 101 4 0 . 0.4 1
. 101 4 111 0 0.4 0
. 101 5 0 . 0.4 1
. 101 5 51 1 0.4 0
. 101 6 0 . 0.8 1
. 101 6 31 1 0.8 0
. 101 7 0 . 0.8 1
. 101 7 50 1 0.8 0
. 101 8 0 . 1.6 1
. 101 8 59 1 1.6 0
. 101 9 0 . 3.2 1
. 101 9 477 0 3.2 0
. 101 11 0 . 3.2 1
. 101 11 278 0 3.2 0
. 101 12 0 . 6.4 1
. 101 12 58 1 6.4 0
Quoted reply history
On 29/07/2011 3:10 a.m., Gaurav Bajaj wrote:
> Dear All,
>
> I am trying to run an example of time to single event data from NONMEM archive, ( http://www.mail-archive.com/ [email protected] /msg01743.html ) using NONMEM 7.2. It gives an error in control statement: $EST "laplace required with likelihood, Pop. data, unless omega=0". I am guessing that it might be a NONMEM version issue. Since I am new to NONMEM, it will be good to have your suggestions.
>
> *Dataset:*
>
> C STUD ID TIME=DV CS DOSE MDV
>
> . 101 4 0 1
> 0.4 1
> . 101 4 111 1
> 0.4 0
> . 101 5 0 1
> 0.4 1
> . 101 5 51 1
> 0.4 0
> . 101 6 0 1
> 0.8 1
> . 101 6 31 1
> 0.8 0
> . 101 7 0 1
> 0.8 1
> . 101 7 . 50
> 1
> 0.8 0
> . 101 8 0 0 1.6 1
> . 101 8 59 1 1.6 0
> . 101 9 0 0 3.2 1
> . 101 9 477 1
> 3.2 0
> . 101 11 0 1
> 3.2 1
> . 101 11 278 1 3.2 0
> . 101 12 0 1
> 6.4 1
> . 101 12 58 1
> 6.4 0
>
> CS=0 is the censored event
> Time is in days
> ---------------------------------------------------------------------------------------------
> *Model:*
>
> $PROB Time To Event data
>
> $INPUT C STUD ID TIME=DV CS DOSE MDV
>
> $DATA hz1.csv IGNORE=C
>
> $SUBR ADVAN=6 TOL=6
> $MODEL COMP=(RISK)
>
> $PK
>
> BASE = THETA(1)*EXP(ETA(1)) ;the ETA is a placeholder here
>
> $DES DADT(1)=BASE ;hazard
>
> $ERROR
> HZ = A(1) ;hazard
> SUR = EXP(-HZ) ;survival probability
> DENS=HZ*SUR
>
> Y = (1-CS)*DENS + CS*SUR
>
> $THETA (0,.025) ; BASE $OMEGA 0 FIX
>
> $ESTIM MAXEVAL=9990 METHOD = COND LAPLACE LIKE PRINT=1 MSFO=msfb01
> $COV PRINT=E
>
> ----
> Thanks,
>
> Gaurav
--
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