From: Bob Leary <*bleary*>

Date: Mon, 10 Aug 2009 09:05:47 -0400

I suspect that the reason that log(1-R) works but log(R) does not is the =

real*4 type

problem, not the fact that R may take on the value 0 (I do not know of =

any pseudorandom

uniform number generators where this can happen). The 1 may be treated =

as a double precision

constant, in which case fortran casts the mixed quantity 1-R (real*8 - =

real*4) as a real*8

quantity, which agrees in type with what dlog is expecting.

Robert H. Leary, PhD

Fellow

Pharsight - A Certara(tm) Company

5625 Dillard Dr., Suite 205

Cary, NC 27511

Phone/Voice Mail: (919) 852-4625, Fax: (919) 859-6871

Email: bleary

*> This email message (including any attachments) is for the sole use of =
*

the intended recipient and may contain confidential and proprietary =

information. Any disclosure or distribution to third parties that is =

not specifically authorized by the sender is prohibited. If you are not =

the intended recipient, please contact the sender by reply email and =

destroy all copies of the original message.

-----Original Message-----

From: owner-nmusers

[mailto:owner-nmusers

Sent: Monday, August 10, 2009 1:51 AM

To: 'nmusers'

Subject: Re: [NMusers] Re: count simulations

Elodie,

I've tried using LOG(R) in NONMEM 6 with the Intel V 10 and G77 =

compilers.

T=0

NEVENT=0

DO WHILE (T.LT.1)

CALL RANDOM (2,R)

T=T-LOG(R)/COUNT

IF (T.LT.1) NEVENT=NEVENT+1

ENDDO

DV=NEVENT

The Intel compiler gave this error:

FSUBS.obj : error LNK2019: unresolved external symbol _DLOG referenced

in function _ERROR

count_sim_R.exe : fatal error LNK1120: 1 unresolved externals

The G77 compiler gave this error:

FSUBS.for: In subroutine `error':

FSUBS.for:221:

B00009=DLOG(R)

^

Reference to intrinsic `DLOG' at (^) invalid -- one or more arguments

have incorrect type

With NONMEM 7 the Intel V 10 and Gfortran compilers had no problems.

This seems to be a problem caused by the code generated by NONMEM 6. R

is declared as a REAL by NONMEM 6 but the DLOG function wants a DOUBLE

and the compilers won't play the game.

NONMEM 7 declares REAL(KIND=DPSIZE) :: R which seems to let =

DLOG

work OK.

In any event its safer (to avoid LOG(0)) to use LOG(1-R) even in NONMEM =

7.

Nick

Elodie Plan wrote:

*> Dear Nick,
*

*>
*

*> I agree, it must be "remembered" somewhere that log(R) is forbidden, =
*

because

*> the simple following trick works (but won't once every 2 billion =
*

times, when

*> 0 is picked):
*

*> IF (ICALL.EQ.4) THEN
*

*> T=0
*

*> N=0
*

*> DO WHILE (T.LT.1)
*

*> CALL RANDOM (2,R)
*

*> RA=R
*

*> T=T-LOG(RA)/LAMB
*

*> IF (T.LT.1) N=N+1
*

*> END DO
*

*> DV=N
*

*> ENDIF
*

*>
*

*> Best regards,
*

*> Elodie
*

*>
*

*> Elodie L. Plan, PharmD, MSc,
*

*> PhDstudent
*

*> *******************************************
*

*> Department of Pharmaceutical Biosciences,
*

*> Faculty of Pharmacy, Uppsala University
*

*> Box 591 - 751 24 Uppsala - SWEDEN
*

*> Office +46 18 4714385 - Fax +46 18 4714003
*

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

*>
*

*>
*

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

*> From: owner-nmusers *

[mailto:owner-nmusers

*> Behalf Of Nick Holford
*

*> Sent: Friday, August 07, 2009 1:43 PM
*

*> To: nmusers
*

*> Subject: Re: [NMusers] Re: count simulations
*

*>
*

*> Mats,
*

*>
*

*> A uniform random deviate from the interval 0-1 is by convention =
*

allowed

*> to include 0 but not include 1. So in theory R could be 0 and LOG(0)
*

*> would be invalid, however 1-R if would be OK if R was 0.
*

*>
*

*> I don't understand what triggered the error message that Elodie shows
*

*> below. It implies that the language processor (is it NM-TRAN or
*

*> FORTRAN?) 'knows' that R could be 0 and therefore does not allow =
*

LOG(R).

*> But surely both NM-TRAN and FORTRAN allow the LOG of any variable and =
*

it

*> is only at run-time that a variable <=0 would throw an exception?
*

*>
*

*> In practice the probability of having R=0 is approaching 0 so I =
*

cannot

*> see that there is a real risk of LOG(R) being a problem. On the other
*

*> hand there is certainly a performance overhead in computing 1-R
*

*> everytime when R would do just as well.
*

*>
*

*> Nick
*

*>
*

*> Mats Karlsson wrote:
*

*>
*

*>> Hi Nick,
*

*>>
*

*>> Elodie came up with one reason why it may have been 1-R rather than R =
*

(see

*>> below). I rationalize the code as simulating one event after another =
*

on a

*>> time interval standardized by lambda. You sample a survival =
*

probability,

*>> translate that into a time, check if it is beyond the standardized
*

*>>
*

*> interval,
*

*>
*

*>> if not increase N and add the event time to the elapsed time in the
*

*>> interval.
*

*>> SUR = EXP(-LAMB*T)
*

*>> R = EXP(-LAMB*T)
*

*>> LOG(R) = -LAMB*T
*

*>> T = -LOG(R)/LAMB
*

*>>
*

*>> Here is elodies mail
*

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

*>> Really in theory having (R) or (1-R) cannot give simulated datasets =
*

more

*>> different to each other than if they were simulated with 2 different =
*

seed

*>> number.
*

*>>
*

*>> But there was a reason in practice in nonmem because having only (R)
*

*>>
*

*> gives:
*

*>
*

*>> FSUBS.f: In subroutine `pred':
*

*>> FSUBS.f:34:
*

*>> B00006=DLOG(R)
*

*>> ^
*

*>> Reference to intrinsic `DLOG' at (^) invalid -- one or more arguments =
*

have

*>> incorrect type
*

*>> No nonmem execution.
*

*>>
*

*>> Is it that nonmem can give random number 0 and not 1? I can see only =
*

the

*>> numerical issue of log(0).
*

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

*>>
*

*>> 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: owner-nmusers *

[mailto:owner-nmusers

*>>
*

*> On
*

*>
*

*>> Behalf Of Nick Holford
*

*>> Sent: Thursday, August 06, 2009 9:48 PM
*

*>> To: nmusers
*

*>> Subject: [NMusers] Re: count simulations
*

*>>
*

*>> Mats,
*

*>>
*

*>> Thanks for pointing out that R and 1-R are equivalent when R is a
*

*>> uniform 0-1 random deviate.
*

*>>
*

*>> There is an NM-TRAN example using 1-R in this paper:
*

*>>
*

*>> Frame B, Miller R, Lalonde RL. Evaluation of Mixture Modeling with =
*

Count

*>> Data
*

*>> using NONMEM. Journal of Pharmacokinetics and Pharmacodynamics.
*

*>> 2003;30(3):167-83.
*

*>>
*

*>> I have to admit to having cut and pasted this example and used it to
*

*>> show others how to simulate count data so it may have propogated that =
*

*>> way too.
*

*>>
*

*>> Do you know of a clear explanation of why this simple algorithm =
*

produces

*>> Poisson distribution samples?
*

*>>
*

*>>
*

*>> Nick
*

*>>
*

*>>
*

*>> Mats Karlsson wrote:
*

*>>
*

*>>
*

*>>> Dear both,
*

*>>>
*

*>>> You have both simulated count data using the code below (or very
*

*>>> similar). My question is why do you use LOG(1-R) rather than the
*

*>>> simpler LOG(R)? If you've done it because you inherited the code,
*

*>>> where did you get the code.
*

*>>>
*

*>>> *IF (ICALL.EQ.4) THEN*
*

*>>>
*

*>>> * T=0*
*

*>>>
*

*>>> * N=0*
*

*>>>
*

*>>> * DO WHILE (T.LT.1)*
*

*>>>
*

*>>> * CALL RANDOM (2,R)*
*

*>>>
*

*>>> * T=T-LOG(1-R)/LAMB*
*

*>>>
*

*>>> * IF (T.LT.1) N=N+1*
*

*>>>
*

*>>> * END DO*
*

*>>>
*

*>>> * DV=N*
*

*>>>
*

*>>> *ENDIF*
*

*>>>
*

*>>> 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
*

*>>>
*

*>>>
*

*>>>
*

*>>
*

*>>
*

*>
*

*>
*

--

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

Received on Mon Aug 10 2009 - 09:05:47 EDT

Date: Mon, 10 Aug 2009 09:05:47 -0400

I suspect that the reason that log(1-R) works but log(R) does not is the =

real*4 type

problem, not the fact that R may take on the value 0 (I do not know of =

any pseudorandom

uniform number generators where this can happen). The 1 may be treated =

as a double precision

constant, in which case fortran casts the mixed quantity 1-R (real*8 - =

real*4) as a real*8

quantity, which agrees in type with what dlog is expecting.

Robert H. Leary, PhD

Fellow

Pharsight - A Certara(tm) Company

5625 Dillard Dr., Suite 205

Cary, NC 27511

Phone/Voice Mail: (919) 852-4625, Fax: (919) 859-6871

Email: bleary

the intended recipient and may contain confidential and proprietary =

information. Any disclosure or distribution to third parties that is =

not specifically authorized by the sender is prohibited. If you are not =

the intended recipient, please contact the sender by reply email and =

destroy all copies of the original message.

-----Original Message-----

From: owner-nmusers

[mailto:owner-nmusers

Sent: Monday, August 10, 2009 1:51 AM

To: 'nmusers'

Subject: Re: [NMusers] Re: count simulations

Elodie,

I've tried using LOG(R) in NONMEM 6 with the Intel V 10 and G77 =

compilers.

T=0

NEVENT=0

DO WHILE (T.LT.1)

CALL RANDOM (2,R)

T=T-LOG(R)/COUNT

IF (T.LT.1) NEVENT=NEVENT+1

ENDDO

DV=NEVENT

The Intel compiler gave this error:

FSUBS.obj : error LNK2019: unresolved external symbol _DLOG referenced

in function _ERROR

count_sim_R.exe : fatal error LNK1120: 1 unresolved externals

The G77 compiler gave this error:

FSUBS.for: In subroutine `error':

FSUBS.for:221:

B00009=DLOG(R)

^

Reference to intrinsic `DLOG' at (^) invalid -- one or more arguments

have incorrect type

With NONMEM 7 the Intel V 10 and Gfortran compilers had no problems.

This seems to be a problem caused by the code generated by NONMEM 6. R

is declared as a REAL by NONMEM 6 but the DLOG function wants a DOUBLE

and the compilers won't play the game.

NONMEM 7 declares REAL(KIND=DPSIZE) :: R which seems to let =

DLOG

work OK.

In any event its safer (to avoid LOG(0)) to use LOG(1-R) even in NONMEM =

7.

Nick

Elodie Plan wrote:

because

times, when

[mailto:owner-nmusers

allowed

LOG(R).

it

cannot

(see

on a

probability,

more

seed

have

the

[mailto:owner-nmusers

Count

produces

--

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

Received on Mon Aug 10 2009 - 09:05:47 EDT