From: Sebastien Bihorel <*Sebastien.Bihorel*>

Date: Fri, 14 May 2010 10:01:11 -0400

Thanks Martin,

Running the model in NM7 using your parameterization returns non-NaN

values for DV, PRED and RES. I agree with you that my parametrization

would have force the residues to be positive but I'm still puzzled by

NM7 results: W is within 0.001 and 1, which should not create numerical

problems... As a side note, I am a bit concerned with the use of a

single sigma parameter for simulation purposes. I believe this will

induce a correlation between the additive and the pseudo-proportional

parts of the RV. Shouldn't we use the following parameterization for

simulation?

DFLG=0

IF(AMT.GT.0)DFLG=1

IPRED=LOG(F+DFLG)

W1=THETA(5)*EPS(1)

W2=THETA(6)*EPS(2)*(1-F/(THETA(7)+F))

Y=IPRED+W1+W2

This should theoretically give the same variance for Y as your

parameterization and enable the two components of the RV model to be

independent during the simulation.

Sebastien

Martin Bergstrand wrote:

*> Dear Sebastien,
*

*>
*

*> I can't answer your primary question but I want to make a correction on the
*

*> way that I specified the residual error model when using it (see updated
*

*> code below). Your implementation will in my opinion cause bias since it will
*

*> make all residuals positive (-x*-x = x^2).
*

*>
*

*> The reason that I used this model for simulations was to avoid the problems
*

*> with the estimation error model, i.e. [W = SQRT(W1**2+(W2/F)**2]. These
*

*> problems have been mentioned before at NMusers
*

*> (http://www.mail-archive.com/nmusers *

*>
*

*> The error model that you are referring wasn't used for estimation since it
*

*> depends on tree parameters (W1,W2 and WH) and these are in most cases not
*

*> all identifiable. I didn't get this error model out of the literature
*

*> (haven't seen it described), it was simply my own solution to the problems I
*

*> was experiencing.
*

*>
*

*> ;; --- Altered error model code ---------
*

*> IPRED=LOG(F+DFLG)
*

*> W1=THETA(5)
*

*> W2=THETA(6)
*

*> WH=THETA(7) ; Equal to 0.5 in examples out of the publication)
*

*>
*

*> W=SQRT(W1*W1+(W2*(1-F/(WH+F)))**2)
*

*>
*

*> Y=IPRED+W*EPS(1)
*

*>
*

*> $SIMGA 1 FIX
*

*> ;; --------------------------------------
*

*>
*

*> Best regards,
*

*>
*

*> Martin Bergstrand, MSc, PhD student
*

*> -----------------------------------------------
*

*> Pharmacometrics Research Group,
*

*> Department of Pharmaceutical Biosciences,
*

*> Uppsala University
*

*> -----------------------------------------------
*

*> martin.bergstrand *

*> -----------------------------------------------
*

*> Work: +46 18 471 4639
*

*> Mobile: +46 709 994 396
*

*>
*

*>
*

*> -----Original Message-----
*

*> From: owner-nmusers *

*> Behalf Of Sebastien Bihorel
*

*> Sent: Thursday, May 13, 2010 6:06 PM
*

*> To: nmusers; Martin Bergstrand
*

*> Subject: [NMusers] Unexpected differences in predictions between NM 6.2.0
*

*> and NM 7.1.2
*

*>
*

*> Dear R-users,
*

*>
*

*> I have recently observed a very puzzling difference of behavior between NM
*

*> 6.2.0 and NM 7.1.2, while trying to reproduce the simulation examples
*

*> described in Martin Bergstrand's paper 'Handling data below the limit of
*

*> quantification in mixed effect models' (AAPS Journal 2009, 11-2). My model
*

*> corresponds to Martins' model B (2-cmt linear model with single IV
*

*> dosing) and implement my interpretation of the residual variability
*

*> described in his equation 3 (see code below). When I run the code in NONMEM
*

*> 7.1.2,, all the PRED values are reported as NaN and RES values either as NaN
*

*> or 0. NM6 returns non-NaN values for PRED and RES, but there does not seem
*

*> to be any differences in simulated DV.
*

*>
*

*> Did anybody experienced the same differences? In my case, could they be
*

*> explained by an improper implementation of the RV model?
*

*>
*

*> As a side note, I would also be interested to know if a reference in the
*

*> literature would describe the properties of this particular RUV model.
*

*>
*

*> Thank you
*

*>
*

*> Sebastien Bihorel
*

*>
*

*> --------
*

*>
*

*> $PROBLEM base-2cmt-sim
*

*>
*

*> $DATA basedata.csv IGNORE= *

*>
*

*> $INPUT ID TIME AMT RATE CMT EVID DV MDV STDY
*

*>
*

*> $THETA
*

*> (0.,5.) ;1- clearance
*

*> (0.,20.) ;2- central compartment volume
*

*> (0.,5.) ;3- distribution clearance
*

*> (0.,100.) ;4- peripheral compartment volume
*

*> (0.,0.1) ;5- 'additive' log RV
*

*> (0.,0.1) ;6- 'second' log RV term
*

*>
*

*> $OMEGA
*

*> 0.3 ;1- IIV in clearance
*

*> 0.3 ;2- IIV in central compartment volume
*

*> 0.3 ;3- IIV in distribution clearance
*

*> 0.3 ;4- IIV in peripheral compartment volume
*

*>
*

*> $SIGMA
*

*> 1 FIX ;1- 'additive' log RV
*

*> 1 FIX ;2- 'second' log RV term
*

*>
*

*> $SUBROUTINES ADVAN3 TRANS4
*

*>
*

*> $PK
*

*>
*

*> ; Model parameter assignment
*

*>
*

*> TVCL=THETA(1)
*

*> TVV1=THETA(2)
*

*> TVQ =THETA(3)
*

*> TVV2=THETA(4)
*

*>
*

*> ECL=EXP(ETA(1))
*

*> EV1=EXP(ETA(2))
*

*> EQ =EXP(ETA(3))
*

*> EV2=EXP(ETA(4))
*

*>
*

*> ; PREDPP required variables
*

*>
*

*> F1=1.0 ; bioavailability in central compartment
*

*> CL=TVCL*ECL ; elimination clearance
*

*> V1=TVV1*EV1 ; volume of the central compartment
*

*> Q =TVQ*EQ ; inter-compartment distribution clearance
*

*> V2=TVV2*EV2 ; volume of the peripheral compartment
*

*>
*

*> S1=V1
*

*>
*

*> $ERROR
*

*>
*

*> ;set up Dose flag
*

*> DFLG=0
*

*> IF(AMT.GT.0)DFLG=1
*

*>
*

*> IPRED=LOG(F+DFLG)
*

*> W1=THETA(5)*EPS(1)
*

*> W2=THETA(6)*EPS(2)
*

*> W=SQRT(W1*W1+(W2*(1-F/(0.5+F)))**2)
*

*>
*

*> Y=IPRED+W
*

*>
*

*>
*

*> $SIM (123456) ONLYSIM
*

*>
*

*> $TABLE ID TIME AMT RATE CMT EVID DV MDV STDY W NOPRINT ONEHEADER
*

*> FILE=base-2cmt-sim.tbl
*

*>
*

*>
*

*> *

Received on Fri May 14 2010 - 10:01:11 EDT

Date: Fri, 14 May 2010 10:01:11 -0400

Thanks Martin,

Running the model in NM7 using your parameterization returns non-NaN

values for DV, PRED and RES. I agree with you that my parametrization

would have force the residues to be positive but I'm still puzzled by

NM7 results: W is within 0.001 and 1, which should not create numerical

problems... As a side note, I am a bit concerned with the use of a

single sigma parameter for simulation purposes. I believe this will

induce a correlation between the additive and the pseudo-proportional

parts of the RV. Shouldn't we use the following parameterization for

simulation?

DFLG=0

IF(AMT.GT.0)DFLG=1

IPRED=LOG(F+DFLG)

W1=THETA(5)*EPS(1)

W2=THETA(6)*EPS(2)*(1-F/(THETA(7)+F))

Y=IPRED+W1+W2

This should theoretically give the same variance for Y as your

parameterization and enable the two components of the RV model to be

independent during the simulation.

Sebastien

Martin Bergstrand wrote:

Received on Fri May 14 2010 - 10:01:11 EDT