From: Mats Karlsson <*mats.karlsson*>

Date: Sat, 26 May 2012 11:59:37 +0200

Hi Martin,

If you write P(3) =1-P(1)-P(2) you reduce the risk of writing something

incorrect J

Mats Karlsson

Mats Karlsson, PhD

Professor of Pharmacometrics

FIRST WORLD CONFERENCE ON PHARMACOMETRICS, 5-7 September 2012, Seoul (

<http://www.go-wcop.org/> www.go-wcop.org)

Dept of Pharmaceutical Biosciences

Faculty of Pharmacy

Uppsala University

Box 591

75124 Uppsala

Phone: +46 18 4714105

Fax + 46 18 4714003

From: owner-nmusers

Behalf Of Martin Bergstrand

Sent: 26 May 2012 09:06

To: 'Andy Stein'; yapingz2011

Cc: nmusers

Subject: RE: [NMusers] Question about handling BLOQ data with mixture model

Dear Andy and Yaping,

I am sorry for my embarrassingly bad advice. I can now see that it is my

clumsy code that does not sum up to 1. That I should have written and how

normally code it is this:

$MIX

NSPOP=3

P(1) = THETA(8)/100

P(2) = (1-THETA(8)/100)*THETA(9)/100

P(3) = 1-THETA(8)/100 -((1-THETA(8)/100)*THETA(9)/100)

$THETA (0, 3.78, 100) ; PMIX1

$THETA (80, 91, 100) ; PMIX2/(1-PMIX1)

This is as the sober me can see mathematically identical to what you did. I

should get an alcohol lock for outlook so that I can't do NMusers postings

after more than 1 beer.

For what it is worth coming from me you could try to add a very small number

to the CUMD (e.g. 10^-6). If you have a BQL observation that NONMEM predicts

a very low probability NONMEM will in my experience crash. When you get a

run working with this arbitrary imputation you can investigate what

observations that has this extremely low probability and investigate if

there seem to be something fishy about them. In general it can be said that

the Laplacian estimation method in NONMEM is notoriously instable and

slightly prone to local minima. It can for that reason be useful to by

default try running with perturbed initial estimates.

Good luck,

Martin

From: owner-nmusers

Behalf Of Andy Stein

Sent: den 26 maj 2012 01:45

To: LGibiansky

Cc: nmusers

Subject: Re: [NMusers] Question about handling BLOQ data with mixture model

I wanted to follow up on the comments to Yaping's email. First, the three

probabilities below from the original code do in fact sum to 1.

P(1)=THETA(8)/100

P(2)=(1-THETA(8)/100)*THETA(9)/1000

P(3)=(1-THETA(8)/100)*(1-THETA(9)/1000)

Note that: P(2) + P(3) = 1-THETA(8)/100

And thus P(1) + P(2) + P(3) = 1

Also, the model worked completely fine when the BLOQ part of the code was

left out and only the mixture was modeled. That is what led us to think

that the combination of BLOQ with a Mixture was causing the problem.

Andy

On Fri, May 25, 2012 at 1:07 PM, LGibiansky

<LGibiansky

Sum of probabilities should sum to 1. More standard way would be to use

P(1) = 1/(1+THETA(8)+THETA(9))

P(2) = THETA(8)/(1+THETA(8)+THETA(9))

P(3) = THETA(9)/(1+THETA(8)+THETA(9))

where THETA(8) and THETA(9) are any positive numbers

Regards

Leonid

Original Message:

-----------------

From: Martin Bergstrand martin.bergstrand

Date: Fri, 25 May 2012 22:59:45 +0700

To: yapingz2011

Subject: RE: [NMusers] Question about handling BLOQ data with mixture model

Dear Yaping,

I can see that you need to make any particular consideration because you

are applying a mixture model. CUMD is dependent on IPRED that in its turn is

dependent on the assigned mixture. That should be enough.

However, I spot what must be an error in your way of defining your mixture

probabilities. As it is now you total probability does not sum up to 1. Why

don't you parameterize it his way:

$MIX

NSPOP=3

P(1) = THETA(8)/100

P(2) = (1-THETA(8)/100)*THETA(9)/100

P(3) = 1-THETA(8)/100 -THETA(9)/100

$THETA (0, 3.78, 100) ; PMIX1

$THETA (80, 91, 100) ; PMIX2/(1-PMIX1)

I have kept the division of THETAs by 100 since I assume that you want

estimates in %. However the division with 1000 did not make any sense to me

despite the correctly assigned THETA boundaries?

Finally a word of caution, be careful with the use of NONMEM reserved

variables such as T (integrated time in $DES) and F (default model

prediction with some ADVANS). In my experience things can go wrong when you

use them outside the way it was intended, conflicts can occur.

Kind regards,

Martin Bergstrand, PhD

Pharmacometrics Research Group

Dept of Pharmaceutical Biosciences

Uppsala University

Sweden

martin.bergstrand

Visiting scientist:

Mahidol-Oxford Tropical Medicine Research Unit,

Bangkok, Thailand

From: owner-nmusers

Behalf Of Yaping Zhang

Sent: den 25 maj 2012 02:04

To: nmusers

Subject: [NMusers] Question about handling BLOQ data with mixture model

Hello NMUsers,

I am trying to analyze BLOQ data using the M3 method (Stuart Beal, Ways to

Fit a PK Model with Some Data Below the Quantification Limit, 2001). The

model I have is a mixture model to describe PD response of three

subpopulations. The complete control stream is pasted below.

I have implemented just the mixture model (no BLOQ) and just the BLOQ error

model (no mixture) and it works fine with nonmem 6. But the run crashed

immediately using nonmem 6 if including both BLOQ and the mixture model.

I am wondering if I need to account for the mixture in the section of the

code below when putting a mixture model together with the BLOQ error model

IF (BLOQ.EQ.1) THEN

F_FLAG=1

Y =CUMD + something based on MIXNUM?

ENDIF

Any ideas are gratefully received!

Many thanks,

Yaping

$PROB AMN107A2303

$INPUT NUM=DROP STUD=DROP SUBJ=DROP ID AGE0=DROP SEX=DROP RACE=DROP

DART=DROP ARM ACTT=DROP POP PPK=DROP COUN=DROP

SOK OTIM=DROP TIME TVIS=DROP

AMT=DROP DOSE=DROP SCHD=DROP AUC=DROP CMIN=DROP

ODV=DROP MDV DV BLOQ STY=DROP EVDT=DROP

$DATA ../data/AMN2303_ENEST.csv ; currently modified with matlab and saved

with oocalc

IGNORE=

IGNORE=(MDV.EQ.1)

IGNORE=(ID.EQ.66, ID.EQ.92, ID.EQ.335, ID.EQ.346, ID.EQ.348, ID.EQ.416,

ID.EQ.496, ID.EQ.527, ID.EQ.762, ID.EQ.790)

$PRED

mu =THETA(1)*EXP(ETA(1))

AA =THETA(2)*EXP(ETA(2))

alpha =THETA(3)*EXP(ETA(3))

BB =THETA(4)*EXP(ETA(4))

beta =THETA(5)*EXP(ETA(5))

InCC =THETA(6)+ETA(6)

gamma =THETA(7)*EXP(ETA(7))

T = TIME

IF (TIME.LE.0) T = 0

EST=MIXEST

IF (MIXNUM.EQ.1) F =AA*EXP(mu*T/8766)

IF (MIXNUM.EQ.2) F =AA*EXP(alpha*T/8766)+BB*EXP(beta*T/8766)

IF (MIXNUM.EQ.3) F

=AA*EXP(alpha*T/8766)+BB*EXP(beta*T/8766)+EXP(InCC)*EXP(gamma*T/8766)

PROP=THETA(10)

W=SQRT(PROP*PROP)

IPRED=-2.8

IF(F.GT.0)IPRED =LOG10(F)

LLOQ=-2.5

DUM=(LLOQ-IPRED)/W

CUMD=PHI(DUM)

IF (BLOQ.EQ.0) THEN

F_FLAG=0

Y =IPRED+W*ERR(1)

ENDIF

IF (BLOQ.EQ.1) THEN

F_FLAG=1

Y =CUMD

ENDIF

$MIX

NSPOP=3

P(1)=THETA(8)/100

P(2)=(1-THETA(8)/100)*THETA(9)/1000

P(3)=(1-THETA(8)/100)*(1-THETA(9)/1000)

$THETA (-10,-0.57,0) ; mu THETA(1)

$THETA (0.0001,50.7,100) ; AA THETA(2)

$THETA (-100,-14.2,0) ; alpha THETA(3)

$THETA (0.0001,0.196,10) ; BB THETA(4)

$THETA (-10,-0.678,0) ; beta THETA(5)

$THETA (-100,-6.92,0) ; InCC THETA(6)

$THETA (0, 2.15) ; gamma THETA(7)

$THETA (0, 3.78, 100) ; THETA(8)

$THETA (800, 910, 1000) ; THETA(9)

$THETA (0.001,0.104, 5) ; ERR THETA(10)

$OMEGA 0.001 FIX ; mu

$OMEGA 0.767 ; AA

$OMEGA 0.236 ; alpha

$OMEGA 3.58 ; BB

$OMEGA 0.148 ; beta

$OMEGA 5 ; CC

$OMEGA 0.858 ; gamma

$SIGMA 1 FIX

$EST METHOD=1 LAPLACIAN NOABORT MAXEVAL=9990 PRINT=1 MSFO=run012.nmmsf SIG=2

$COV

$TABLE ID TIME DV MDV IPRED EST

FILE=run012.nmfit NOPRINT ONEHEADER

$TABLE ID ETA1 ETA2 ETA3 ETA4 ETA5 ETA6 ETA7 mu AA alpha BB beta InCC gamma

NOAPPEND FIRSTONLY FILE=run012.nmpar NOPRINT ONEHEADER

--------------------------------------------------------------------

mail2web LIVE - Free email based on MicrosoftR Exchange technology -

http://link.mail2web.com/LIVE

Received on Sat May 26 2012 - 05:59:37 EDT

Date: Sat, 26 May 2012 11:59:37 +0200

Hi Martin,

If you write P(3) =1-P(1)-P(2) you reduce the risk of writing something

incorrect J

Mats Karlsson

Mats Karlsson, PhD

Professor of Pharmacometrics

FIRST WORLD CONFERENCE ON PHARMACOMETRICS, 5-7 September 2012, Seoul (

<http://www.go-wcop.org/> www.go-wcop.org)

Dept of Pharmaceutical Biosciences

Faculty of Pharmacy

Uppsala University

Box 591

75124 Uppsala

Phone: +46 18 4714105

Fax + 46 18 4714003

From: owner-nmusers

Behalf Of Martin Bergstrand

Sent: 26 May 2012 09:06

To: 'Andy Stein'; yapingz2011

Cc: nmusers

Subject: RE: [NMusers] Question about handling BLOQ data with mixture model

Dear Andy and Yaping,

I am sorry for my embarrassingly bad advice. I can now see that it is my

clumsy code that does not sum up to 1. That I should have written and how

normally code it is this:

$MIX

NSPOP=3

P(1) = THETA(8)/100

P(2) = (1-THETA(8)/100)*THETA(9)/100

P(3) = 1-THETA(8)/100 -((1-THETA(8)/100)*THETA(9)/100)

$THETA (0, 3.78, 100) ; PMIX1

$THETA (80, 91, 100) ; PMIX2/(1-PMIX1)

This is as the sober me can see mathematically identical to what you did. I

should get an alcohol lock for outlook so that I can't do NMusers postings

after more than 1 beer.

For what it is worth coming from me you could try to add a very small number

to the CUMD (e.g. 10^-6). If you have a BQL observation that NONMEM predicts

a very low probability NONMEM will in my experience crash. When you get a

run working with this arbitrary imputation you can investigate what

observations that has this extremely low probability and investigate if

there seem to be something fishy about them. In general it can be said that

the Laplacian estimation method in NONMEM is notoriously instable and

slightly prone to local minima. It can for that reason be useful to by

default try running with perturbed initial estimates.

Good luck,

Martin

From: owner-nmusers

Behalf Of Andy Stein

Sent: den 26 maj 2012 01:45

To: LGibiansky

Cc: nmusers

Subject: Re: [NMusers] Question about handling BLOQ data with mixture model

I wanted to follow up on the comments to Yaping's email. First, the three

probabilities below from the original code do in fact sum to 1.

P(1)=THETA(8)/100

P(2)=(1-THETA(8)/100)*THETA(9)/1000

P(3)=(1-THETA(8)/100)*(1-THETA(9)/1000)

Note that: P(2) + P(3) = 1-THETA(8)/100

And thus P(1) + P(2) + P(3) = 1

Also, the model worked completely fine when the BLOQ part of the code was

left out and only the mixture was modeled. That is what led us to think

that the combination of BLOQ with a Mixture was causing the problem.

Andy

On Fri, May 25, 2012 at 1:07 PM, LGibiansky

<LGibiansky

Sum of probabilities should sum to 1. More standard way would be to use

P(1) = 1/(1+THETA(8)+THETA(9))

P(2) = THETA(8)/(1+THETA(8)+THETA(9))

P(3) = THETA(9)/(1+THETA(8)+THETA(9))

where THETA(8) and THETA(9) are any positive numbers

Regards

Leonid

Original Message:

-----------------

From: Martin Bergstrand martin.bergstrand

Date: Fri, 25 May 2012 22:59:45 +0700

To: yapingz2011

Subject: RE: [NMusers] Question about handling BLOQ data with mixture model

Dear Yaping,

I can see that you need to make any particular consideration because you

are applying a mixture model. CUMD is dependent on IPRED that in its turn is

dependent on the assigned mixture. That should be enough.

However, I spot what must be an error in your way of defining your mixture

probabilities. As it is now you total probability does not sum up to 1. Why

don't you parameterize it his way:

$MIX

NSPOP=3

P(1) = THETA(8)/100

P(2) = (1-THETA(8)/100)*THETA(9)/100

P(3) = 1-THETA(8)/100 -THETA(9)/100

$THETA (0, 3.78, 100) ; PMIX1

$THETA (80, 91, 100) ; PMIX2/(1-PMIX1)

I have kept the division of THETAs by 100 since I assume that you want

estimates in %. However the division with 1000 did not make any sense to me

despite the correctly assigned THETA boundaries?

Finally a word of caution, be careful with the use of NONMEM reserved

variables such as T (integrated time in $DES) and F (default model

prediction with some ADVANS). In my experience things can go wrong when you

use them outside the way it was intended, conflicts can occur.

Kind regards,

Martin Bergstrand, PhD

Pharmacometrics Research Group

Dept of Pharmaceutical Biosciences

Uppsala University

Sweden

martin.bergstrand

Visiting scientist:

Mahidol-Oxford Tropical Medicine Research Unit,

Bangkok, Thailand

From: owner-nmusers

Behalf Of Yaping Zhang

Sent: den 25 maj 2012 02:04

To: nmusers

Subject: [NMusers] Question about handling BLOQ data with mixture model

Hello NMUsers,

I am trying to analyze BLOQ data using the M3 method (Stuart Beal, Ways to

Fit a PK Model with Some Data Below the Quantification Limit, 2001). The

model I have is a mixture model to describe PD response of three

subpopulations. The complete control stream is pasted below.

I have implemented just the mixture model (no BLOQ) and just the BLOQ error

model (no mixture) and it works fine with nonmem 6. But the run crashed

immediately using nonmem 6 if including both BLOQ and the mixture model.

I am wondering if I need to account for the mixture in the section of the

code below when putting a mixture model together with the BLOQ error model

IF (BLOQ.EQ.1) THEN

F_FLAG=1

Y =CUMD + something based on MIXNUM?

ENDIF

Any ideas are gratefully received!

Many thanks,

Yaping

$PROB AMN107A2303

$INPUT NUM=DROP STUD=DROP SUBJ=DROP ID AGE0=DROP SEX=DROP RACE=DROP

DART=DROP ARM ACTT=DROP POP PPK=DROP COUN=DROP

SOK OTIM=DROP TIME TVIS=DROP

AMT=DROP DOSE=DROP SCHD=DROP AUC=DROP CMIN=DROP

ODV=DROP MDV DV BLOQ STY=DROP EVDT=DROP

$DATA ../data/AMN2303_ENEST.csv ; currently modified with matlab and saved

with oocalc

IGNORE=

IGNORE=(MDV.EQ.1)

IGNORE=(ID.EQ.66, ID.EQ.92, ID.EQ.335, ID.EQ.346, ID.EQ.348, ID.EQ.416,

ID.EQ.496, ID.EQ.527, ID.EQ.762, ID.EQ.790)

$PRED

mu =THETA(1)*EXP(ETA(1))

AA =THETA(2)*EXP(ETA(2))

alpha =THETA(3)*EXP(ETA(3))

BB =THETA(4)*EXP(ETA(4))

beta =THETA(5)*EXP(ETA(5))

InCC =THETA(6)+ETA(6)

gamma =THETA(7)*EXP(ETA(7))

T = TIME

IF (TIME.LE.0) T = 0

EST=MIXEST

IF (MIXNUM.EQ.1) F =AA*EXP(mu*T/8766)

IF (MIXNUM.EQ.2) F =AA*EXP(alpha*T/8766)+BB*EXP(beta*T/8766)

IF (MIXNUM.EQ.3) F

=AA*EXP(alpha*T/8766)+BB*EXP(beta*T/8766)+EXP(InCC)*EXP(gamma*T/8766)

PROP=THETA(10)

W=SQRT(PROP*PROP)

IPRED=-2.8

IF(F.GT.0)IPRED =LOG10(F)

LLOQ=-2.5

DUM=(LLOQ-IPRED)/W

CUMD=PHI(DUM)

IF (BLOQ.EQ.0) THEN

F_FLAG=0

Y =IPRED+W*ERR(1)

ENDIF

IF (BLOQ.EQ.1) THEN

F_FLAG=1

Y =CUMD

ENDIF

$MIX

NSPOP=3

P(1)=THETA(8)/100

P(2)=(1-THETA(8)/100)*THETA(9)/1000

P(3)=(1-THETA(8)/100)*(1-THETA(9)/1000)

$THETA (-10,-0.57,0) ; mu THETA(1)

$THETA (0.0001,50.7,100) ; AA THETA(2)

$THETA (-100,-14.2,0) ; alpha THETA(3)

$THETA (0.0001,0.196,10) ; BB THETA(4)

$THETA (-10,-0.678,0) ; beta THETA(5)

$THETA (-100,-6.92,0) ; InCC THETA(6)

$THETA (0, 2.15) ; gamma THETA(7)

$THETA (0, 3.78, 100) ; THETA(8)

$THETA (800, 910, 1000) ; THETA(9)

$THETA (0.001,0.104, 5) ; ERR THETA(10)

$OMEGA 0.001 FIX ; mu

$OMEGA 0.767 ; AA

$OMEGA 0.236 ; alpha

$OMEGA 3.58 ; BB

$OMEGA 0.148 ; beta

$OMEGA 5 ; CC

$OMEGA 0.858 ; gamma

$SIGMA 1 FIX

$EST METHOD=1 LAPLACIAN NOABORT MAXEVAL=9990 PRINT=1 MSFO=run012.nmmsf SIG=2

$COV

$TABLE ID TIME DV MDV IPRED EST

FILE=run012.nmfit NOPRINT ONEHEADER

$TABLE ID ETA1 ETA2 ETA3 ETA4 ETA5 ETA6 ETA7 mu AA alpha BB beta InCC gamma

NOAPPEND FIRSTONLY FILE=run012.nmpar NOPRINT ONEHEADER

--------------------------------------------------------------------

mail2web LIVE - Free email based on MicrosoftR Exchange technology -

http://link.mail2web.com/LIVE

Received on Sat May 26 2012 - 05:59:37 EDT