Hi Mats, Nick and all.
Two questions on the proposed model for simulation using the uncertainty =
in the covariance matrix that was discussed last summer.
1. I would like to perform simulation with uncertainty in parameter =
estimates followed by estimation with alternative models. Would it be =
possible to use the model Mats proposed as the simulation model using =
the SSE module in PsN?
2. I did run the model proposed by Mats.
When I look at the output files 52 (estimated thetas) and 53 (similated =
thetas), the mean of the estimated thetas are only around half of the =
simulated values. Also there is no correlation between simulated and =
estimated thetas. How can this be? Is there something I am =
misunderstanding.
I used a datafile with 20 subjects with two observations/subject as =
explained below.
The model I used was:
$PROB
$DATA data1.csv IGNORE=
$INPUT ID DV TIME TVCL TVV AMT
;data1=
; 1 0 0 0 0 100
; 1 0 1 0 0 0
; 1 0 21 0 0 0
; 2 0 0 0 0 100
; etc
$SIM (20090726 NEW) SUBPROBLEMS=100
$THETA ; estimates of THETA and OMEGA from previous run
(0,1) ; POP_CL theta1
(0,10) ; POP_V theta2
$OMEGA
0.5 ; PPV_CL eta1
0.5 ; PPV_V eta2
;variancecovariance matrix of the THETA estimates from previous run
$OMEGA BLOCK(2)
0.2 ; UNC_POP_CL eta3
0.1 3 FIX; UNC_POP_V eta 4
$SIGMA .001
$SUB ADVAN1
$PK
CL=THETA(1)*EXP(ETA(1)) ; with uncertainty for CL
V =THETA(2)*EXP(ETA(2)) ; with uncertainty for V
TCLE=THETA(1); estimated typical value
TVE=THETA(2); estimated typical value
K=CL/V
IF (ICALL.EQ.4.AND.NEWIND.EQ.0) THEN ; do this just once per =
subproblem
DO WHILE (ETA(3).LT..99.OR.ETA(4).LT.9.9) CALL SIMETA(ETA) ;resample =
if "bad" values
END DO
TVCL=THETA(1)+ETA(3)
TVV =THETA(2)+ETA(4)
ENDIF
IF (ICALL.EQ.4) THEN ; do this for all simulated values
CL=TVCL*EXP(ETA(1)) ; with uncertainty for CL
V =TVV*EXP(ETA(2)) ; with uncertainty for V
ENDIF
$ERROR
Y =F*EXP(EPS(1))
$INFN
IF (NEWIND.EQ.0) THEN
WRITE (52,*) THETA ; Output estimated THETAs to file 52
WRITE (53,*) TVCL, TVV ; Output simulated THETAs to file 53
ENDIF
$EST MAX=9999 METH=1 INTER
$TABLE ID DV TIME TVCL TVV TCLE TVE CL V ETA1 ETA2 ETA3 ETA4
ONEHEADER NOPRINT FILE=patab37
=

Confidentiality Notice: This message is private and may contain =
confidential and proprietary information. If you have received this =
message in error, please notify us and remove it from your system and =
note that you must not copy, distribute or take any action in reliance =
on it. Any unauthorized use or disclosure of the contents of this =
message is not permitted and may be unlawful.
Original Message
From: ownernmusers
On Behalf Of Nick Holford
Sent: den 27 juli 2009 09:42
To: nmusers
Subject: Re: AW: [NMusers] Simulations with/without residual error
Mats,
Thanks for pointing out the error in my original code fragement  it
was, of course, an intentional model misspecification :)
The code as you write it below won't run for me because of PREDPP errors =
so I dont know how you get any output at all but I do accept that it
does need to have ICALL.EQ.4 in the block that calculates the initial
values of UNCCL and UNCV.
IF (ICALL.EQ.4.AND.NEWIND.EQ.0) THEN ; do this just once per subproblem
UNCCL=THETA(1)*EXP(ETA(3))
UNCV=THETA(2)*EXP(ETA(4))
ENDIF
I don't understand why ICALL.EQ.4 is needed because this is a $SIM
ONLYSIM problem. I expected the code to run as if it was all in an
ICALL.EQ.4 block but it seems that once again NONMEM does the unexpected =
thing  no doubt consist with some deeply hidden rule in the =
documentation.
When ICALL.EQ.4 is included and ETA is resampled to get reasonable
values for CL and V then this is the output with the same value of UNCCL =
and UNCV for all records in each subproblem.
TABLE NO. 1
REP ID UNCCL CL UNCV V
1 1 1.3889 0.50415 11.591 6.1054
1 1 1.3889 0.50415 11.591 6.1054
1 1 1.3889 0.50415 11.591 6.1054
1 2 1.3889 3.1992 11.591 10.021
1 2 1.3889 3.1992 11.591 10.021
1 2 1.3889 3.1992 11.591 10.021
TABLE NO. 1
REP ID UNCCL CL UNCV V
2 1 1.5026 2.6093 11.656 2.4064
2 1 1.5026 2.6093 11.656 2.4064
2 1 1.5026 2.6093 11.656 2.4064
2 2 1.5026 0.66622 11.656 4.6768
2 2 1.5026 0.66622 11.656 4.6768
2 2 1.5026 0.66622 11.656 4.6768
Best wishes,
Nick
Mats Karlsson wrote:
> Nick,
>
> I said the code was incorrect because it actually simulates values =
dependent
> on the uncertainty for the very first record in each problem only. =
Here is
> output from your code:
>
>
>
> LINE NO. ID TIME DV UNCC UNCV CL =
V
> DV PRED RES WRES
>
> 1
> + 1.00E+00 0.00E+00 0.00E+00 1.36E+00 9.38E+00 2.36E+00
> 3.23E+00 0.00E+00 1.00E+02 0.00E+00 0.00E+00
>
> 2
> + 1.00E+00 1.00E+00 1.#RE+00 0.00E+00 0.00E+00 0.00E+00
> 0.00E+00 1.#RE+00 1.#RE+00 1.#RE+00 0.00E+00
>
> 3
> + 1.00E+00 2.10E+01 1.#RE+00 0.00E+00 0.00E+00 0.00E+00
> 0.00E+00 1.#RE+00 1.#RE+00 1.#RE+00 0.00E+00
>
> 4
> + 2.00E+00 0.00E+00 0.00E+00 0.00E+00 0.00E+00 0.00E+00
> 0.00E+00 0.00E+00 1.00E+02 0.00E+00 0.00E+00
>
>
> And here is that code:
> $PROB
> $DATA data2.csv
> $INPUT ID DV TIME DROP DROP AMT
>
> $SIM (20090709) ONLYSIM SUBPROBLEMS=100
> ; estimates of THETA and OMEGA from previous run
> $THETA
> 1 ; POP_CL theta1
> 10 ; POP_V theta2
> $OMEGA
> 0.5 ; PPV_CL eta1
> 0.5 ; PPV_V eta2
> ;variancecovariance matrix of the THETA estimates from previous run
> $OMEGA BLOCK(2)
> 0.2 ; UNC_POP_CL eta3
> 0.1 3 ; UNC_POP_V eta 4
> $SIGMA .1
> $SUB ADVAN1
> $PK
> ; get CL and V uncertainties
> IF (NEWIND.EQ.0) THEN ; do this just once per subproblem
> UNCCL=THETA(1)+ETA(3)
> UNCV=THETA(2)+ETA(4)
> ENDIF
> CL=UNCCL*EXP(ETA(1)) ; with uncertainty for CL
> V =UNCV*EXP(ETA(2)) ; with uncertainty for V
>
> K=CL/V
> $ERROR
> Y =F*EXP(EPS(1))
> $TABLE ID TIME DV UNCCL UNCV CL V
>
> As for the error model this was of course intentional  checking =
impact of
> model misspecification is one of the reasons we do these kinds of =
simulation
> studies. :)
>
> Best regards,
> Mats
>
> Mats Karlsson, PhD
> Professor of Pharmacometrics
> Dept of Pharmaceutical Biosciences
> Uppsala University
> Box 591
> 751 24 Uppsala Sweden
> phone: +46 18 4714105
> fax: +46 18 471 4003
>
>
> Original Message
> From: ownernmusers
[mailto:ownernmusers
> Behalf Of Nick Holford
> Sent: Sunday, July 26, 2009 10:42 PM
> To: nmusers
> Subject: Re: AW: [NMusers] Simulations with/without residual error
>
> Mats,
>
> Thanks for extending the code fragment I proposed earlier to add some
> extra useful features.
>
> I don't know what was incorrect about the code I suggested but
> resampling the uncertainty ETA values is a pragmatic way to avoid =
PREDPP
> errors caused by randomly generated nonpositive CL or V values. The
> code below produces the same results as your code but it shows how to
> check on reasonable values for CL and V that do not depend on the =
value
> of CL or V defined by THETA or other things such as the effect of
> covariates on the group parameter value.
>
> I note that your example using simulation and estimation does not use
> the same residual error model for simulation and estimation. The
> simulation residual error is exp(EPS(1)) but METHOD=CONDITIONAL
> INTERACTION would approximate this by (1+EPS(1)). You may remember the =
> commentary on this point that Stuart Beal wrote in 2002 :)
>
> Best wishes,
>
> Nick
>
> Beal SL. Commentary on Significance Levels for Covariate Effects in
> NONMEM. Journal of Pharmacokinetics & Pharmacodynamics. =
2002;29(4):40310.
>
> $SIM (20090726 NEW) ONLYSIM SUBPROBLEMS=100
> ; estimates of THETA and OMEGA from previous run
> $THETA
> 1 ; POP_CL theta1
> 10 ; POP_V theta2
> $OMEGA
> 0.5 ; PPV_CL eta1
> 0.5 ; PPV_V eta2
> ;variancecovariance matrix of the THETA estimates from previous run
> $OMEGA BLOCK(2)
> 0.2 ; UNC_POP_CL eta3
> 0.1 3 ; UNC_POP_V eta 4
> $SIGMA .01
>
> $SUB ADVAN1
>
> $PK
> ; get CL and V uncertainties
> IF (ICALL.EQ.4.AND.NEWIND.EQ.0) THEN ; do this just once per =
subproblem
> UNCCL=THETA(1)+ETA(3)
> UNCV=THETA(2)+ETA(4)
> DO WHILE (UNCCL.LE.0.OR.UNCV.LE.0) ;resample if values unreasonable
> CALL SIMETA(ETA)
> UNCCL=THETA(1)+ETA(3)
> UNCV=THETA(2)+ETA(4)
> END DO
> ENDIF
>
> CL=UNCCL*EXP(ETA(1)) ; with uncertainty for CL
> V =UNCV*EXP(ETA(2)) ; with uncertainty for V
> K=CL/V
>
> $ERROR
> Y =F*EXP(EPS(1))
>
>
> Mats Karlsson wrote:
>
>> Dear Nick, Marc and all,
>>
>> An even later addition to this thread. Nick pointed out the
>> possibility to do simulations in NONMEM that incorporated error in
>> simulation parameters. The code provided gave a good hint at how to =
do
>> that but would work incorrectly unless some additional NONMEM =
switches
>> were flicked. Also, as pointed out by Marc and others, one might like =
>> to truncate the parameter distribution for simulation (avoid negative =
>> CL etc). Also it may be convenient to do simulation followed by
>> reestimation in a single run. For appropriate postprocessing this
>> requires that simulation parameters are output somewhere. A model
>> file, extending Nick's example, that does these things are provided
>> below.
>>
>> $PROB
>>
>> $DATA data1
>>
>> $INPUT ID DV TIME TVCL TVV AMT
>>
>> ;data1=
>>
>> ; 1 0 0 0 0 100
>>
>> ; 1 0 1 0 0 0
>>
>> ; 1 0 21 0 0 0
>>
>> ; 2 0 0 0 0 100
>>
>> ; etc
>>
>> $SIM (20090726 NEW) SUBPROBLEMS=100
>>
>> $THETA ; estimates of THETA and OMEGA from previous run
>>
>> (0,1) ; POP_CL theta1
>>
>> (0,10) ; POP_V theta2
>>
>> $OMEGA
>>
>> 0.5 ; PPV_CL eta1
>>
>> 0.5 ; PPV_V eta2
>>
>> ;variancecovariance matrix of the THETA estimates from previous run
>>
>> $OMEGA BLOCK(2)
>>
>> 0.2 ; UNC_POP_CL eta3
>>
>> 0.1 3 FIX; UNC_POP_V eta 4
>>
>> $SIGMA .01
>>
>> $SUB ADVAN1
>>
>> $PK
>>
>> CL=THETA(1)*EXP(ETA(1)) ; with uncertainty for CL
>>
>> V =THETA(2)*EXP(ETA(2)) ; with uncertainty for V
>>
>> K=CL/V
>>
>> IF (ICALL.EQ.4.AND.NEWIND.EQ.0) THEN ; do this just once per =
subproblem
>>
>> DO WHILE (ETA(3).LT..99.OR.ETA(4).LT.9.9) CALL SIMETA(ETA) =
;resample
>> if "bad" values
>>
>> END DO
>>
>> TVCL=THETA(1)+ETA(3)
>>
>> TVV =THETA(2)+ETA(4)
>>
>> ENDIF
>>
>> IF (ICALL.EQ.4) THEN ; do this for all simulated values
>>
>> CL=TVCL*EXP(ETA(1)) ; with uncertainty for CL
>>
>> V =TVV*EXP(ETA(2)) ; with uncertainty for V
>>
>> ENDIF
>>
>> $ERROR
>>
>> Y =F*EXP(EPS(1))
>>
>> $INFN
>>
>> IF (NEWIND.EQ.0) THEN
>>
>> WRITE (52,*) THETA ; Output estimated THETAs to file 52
>>
>> WRITE (53,*) TVCL, TVV ; Output simulated THETAs to file 53
>>
>> ENDIF
>>
>> $EST MAX=9999 METH=1 INTER
>>
>> Mats Karlsson, PhD
>>
>> Professor of Pharmacometrics
>>
>> Dept of Pharmaceutical Biosciences
>>
>> Uppsala University
>>
>> Box 591
>>
>> 751 24 Uppsala Sweden
>>
>> phone: +46 18 4714105
>>
>> fax: +46 18 471 4003
>>
>>
> Nick Holford wrote:
>
>> $SIM (20090709) ONLYSIM SUBPROBLEMS=100
>> ; estimates of THETA and OMEGA from previous run
>> $THETA
>> 1 ; POP_CL theta1
>> 10 ; POP_V theta2
>> $OMEGA
>> 0.5 ; PPV_CL eta1
>> 0.5 ; PPV_V eta2
>> ;variancecovariance matrix of the THETA estimates from previous run
>> $OMEGA BLOCK(2)
>> 0.2 ; UNC_POP_CL eta3
>> 0.1 3 ; UNC_POP_V eta 4
>>
>> $PK
>> ; get CL and V uncertainties
>> IF (NEWIND.EQ.0) THEN ; do this just once per subproblem
>> UNCCL=THETA(1)+ETA(3)
>> UNCV=THETA(2)+ETA(4)
>> ENDIF
>> CL=UNCCL*EXP(ETA(1)) ; with uncertainty for CL
>> V =UNCV*EXP(ETA(2)) ; with uncertainty for V
>> ...
>>
>>
>
>

Nick Holford, Professor Clinical Pharmacology
Dept Pharmacology & Clinical Pharmacology
University of Auckland, 85 Park Rd, Private Bag 92019, Auckland, New =
Zealand
n.holford
mobile: +64 21 46 23 53
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
 application/octetstream attachment: fort52
 application/octetstream attachment: fort53
Received on Thu Mar 11 2010  09:12:17 EST