Dear NONMEM Users,
I'm developing a model with NONMEM6|7. In this model, I need to generate a random number at every time step. The problem is that the use of the function "CALL RANDOM (2,R)" is supposed to be only for simulation, isn't it?
Thus, is it possible to generate a random number at every time step with NONMEM? Does anybody know how?
This is part of the code that does not work (it's not working because the condition ICALL.EQ.4 never happens and for that it's not getting into the "if", but on the other hand, in order to use the function RANDOM(2,R), such function requires ICALL.EQ.4):
...
$SUBS ADVAN6 TOL=5
$MODEL
...
$PK
...
IF (ICALL.EQ.4) THEN
CALL RANDOM (2,R) ;Rand number in[0,1[
T=1-R
...
ENDIF
$DES
...
$ERROR
...
$ESTIMATION MAXEVAL=0 NUMERICAL METHOD=COND LAPLACE LIKE CENTERING PRINT=2 MSFO=msfo3
Thank you!
Nieves
--
--------------------------------
Nieves Velez de Mendizabal, Ph.D
Departamento de Farmacia y Tecnología Farmacéutica
Facultad de Farmacia
Universidad de Navarra
Phone: (+34) 658 732 851
Phone: (+34) 948 255 400 ext. 5827
[email protected]
How to generate a random number with $EST
Nieves,
I have no idea why you want to do this but if you want to have the same random number at every time then the simplest and most efficient way to do this is to create the random number as a data item in your data set. You can use the rand() function in Excel to generate a pseudo-random uniform sample at each time. It is a waste of computer cycles to have to keep calling the random number generator for every record if you want the same random number value.
On the other hand if you want a different random number every time you reach a particular record then that is a more complicated matter. I expect it can be done with a verbatim code call but I do not know how. Perhaps you can use?
" CALL RANDOM(2,R)
It would be interesting to know why you want the random number when you are not doing simulation.
Best wishes,
Nick
Quoted reply history
On 14/03/2011 10:32 a.m., Nieves Vélez de Mendizabal wrote:
> Dear NONMEM Users,
>
> I'm developing a model with NONMEM6|7. In this model, I need to generate a random number at every time step. The problem is that the use of the function "CALL RANDOM (2,R)" is supposed to be only for simulation, isn't it?
>
> Thus, is it possible to generate a random number at every time step with NONMEM? Does anybody know how?
>
> This is part of the code that does not work (it's not working because the condition ICALL.EQ.4 never happens and for that it's not getting into the "if", but on the other hand, in order to use the function RANDOM(2,R), such function requires ICALL.EQ.4):
>
> ...
>
> $SUBS ADVAN6 TOL=5
>
> $MODEL
>
> ...
>
> $PK
>
> ...
>
> IF (ICALL.EQ.4) THEN
>
> CALL RANDOM (2,R) ;Rand number in[0,1[
>
> T=1-R
>
> ...
>
> ENDIF
>
> $DES
>
> ...
>
> $ERROR
>
> ...
>
> $ESTIMATION MAXEVAL=0 NUMERICAL METHOD=COND LAPLACE LIKE CENTERING PRINT=2 MSFO=msfo3
>
> Thank you!
>
> Nieves
>
> --
> --------------------------------
> Nieves Velez de Mendizabal, Ph.D
> Departamento de Farmacia y Tecnología Farmacéutica
> Facultad de Farmacia
> Universidad de Navarra
> Phone: (+34) 658 732 851
> Phone: (+34) 948 255 400 ext. 5827
> [email protected]
> --------------------------------
--
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
Nieves,
I'm also curious why you wish to generate random numbers. However, FORTRAN has a built-in uniform random number generator:
http://gcc.gnu.org/onlinedocs/gfortran/RAND.html
I'm assuming you wish to use a uniform distribution, so this should work just fine.
However, if you wish to generate Normal deviates, the procedure is slightly more complicated. You can transform uniform deviates to normal deviates using the Box-Muller method. The basic idea is to take two uniform deviates from -1 to 1, until they are in the unit circle. When they are, you can transform them to two normal deviates; The second deviate is saved for the next call.
Fortran 90 code is available in Numeric Recipes in Fortran on page 279-280, see:
http://apps.nrbook.com/fortran/index.html
Matt.
Quoted reply history
From: owner-nmusers
Behalf Of Nick Holford
Sent: Monday, March 14, 2011 4:09 AM
To: nmusers
Subject: Re: [NMusers] How to generate a random number with $EST
Nieves,
I have no idea why you want to do this but if you want to have the same random number at every time then the simplest and most efficient way to do this is to create the random number as a data item in your data set. You can use the rand() function in Excel to generate a pseudo-random uniform sample at each time. It is a waste of computer cycles to have to keep calling the random number generator for every record if you want the same random number value.
On the other hand if you want a different random number every time you reach a particular record then that is a more complicated matter. I expect it can be done with a verbatim code call but I do not know how. Perhaps you can use?
" CALL RANDOM(2,R)
It would be interesting to know why you want the random number when you are not doing simulation.
Best wishes,
Nick
On 14/03/2011 10:32 a.m., Nieves Vlez de Mendizabal wrote:
Dear NONMEM Users,
I'm developing a model with NONMEM6|7. In this model, I need to generate a random number at every time step. The problem is that the use of the function "CALL RANDOM (2,R)" is supposed to be only for simulation, isn't it?
Thus, is it possible to generate a random number at every time step with NONMEM? Does anybody know how?
This is part of the code that does not work (it's not working because the condition ICALL.EQ.4 never happens and for that it's not getting into the "if", but on the other hand, in order to use the function RANDOM(2,R), such function requires ICALL.EQ.4):
...
$SUBS ADVAN6 TOL=5
$MODEL
...
$PK
...
IF (ICALL.EQ.4) THEN
CALL RANDOM (2,R) ;Rand number in[0,1[
T=1-R
...
ENDIF
$DES
...
$ERROR
...
$ESTIMATION MAXEVAL=0 NUMERICAL METHOD=COND LAPLACE LIKE CENTERING PRINT=2 MSFO=msfo3
Thank you!
Nieves
--
--------------------------------
Nieves Velez de Mendizabal, Ph.D
Departamento de Farmacia y Tecnologa Farmacutica
Facultad de Farmacia
Universidad de Navarra
Phone: (+34) 658 732 851
Phone: (+34) 948 255 400 ext. 5827
nvelez
--------------------------------
--
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: n.holford
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
________________________________
This e-mail (including any attachments) is confidential and may be legally privileged. If you are not an intended recipient or an authorized representative of an intended recipient, you are prohibited from using, copying or distributing the information in this e-mail or its attachments. If you have received this e-mail in error, please notify the sender immediately by return e-mail and delete all copies of this message and any attachments.
Thank you.
Nieves,
I'm also curious why you wish to generate random numbers. However, FORTRAN
has a built-in uniform random number generator:
http://gcc.gnu.org/onlinedocs/gfortran/RAND.html
I'm assuming you wish to use a uniform distribution, so this should work just
fine.
However, if you wish to generate Normal deviates, the procedure is slightly
more complicated. You can transform uniform deviates to normal deviates using
the Box-Muller method. The basic idea is to take two uniform deviates from -1
to 1, until they are in the unit circle. When they are, you can transform them
to two normal deviates; The second deviate is saved for the next call.
Fortran 90 code is available in Numeric Recipes in Fortran on page 279-280, see:
http://apps.nrbook.com/fortran/index.html
Matt.
Quoted reply history
From: [email protected] [mailto:[email protected]] On
Behalf Of Nick Holford
Sent: Monday, March 14, 2011 4:09 AM
To: nmusers
Subject: Re: [NMusers] How to generate a random number with $EST
Nieves,
I have no idea why you want to do this but if you want to have the same random
number at every time then the simplest and most efficient way to do this is to
create the random number as a data item in your data set. You can use the
rand() function in Excel to generate a pseudo-random uniform sample at each
time. It is a waste of computer cycles to have to keep calling the random
number generator for every record if you want the same random number value.
On the other hand if you want a different random number every time you reach a
particular record then that is a more complicated matter. I expect it can be
done with a verbatim code call but I do not know how. Perhaps you can use?
" CALL RANDOM(2,R)
It would be interesting to know why you want the random number when you are not
doing simulation.
Best wishes,
Nick
On 14/03/2011 10:32 a.m., Nieves Vélez de Mendizabal wrote:
Dear NONMEM Users,
I'm developing a model with NONMEM6|7. In this model, I need to generate a
random number at every time step. The problem is that the use of the function
"CALL RANDOM (2,R)" is supposed to be only for simulation, isn't it?
Thus, is it possible to generate a random number at every time step with
NONMEM? Does anybody know how?
This is part of the code that does not work (it's not working because the
condition ICALL.EQ.4 never happens and for that it's not getting into the "if",
but on the other hand, in order to use the function RANDOM(2,R), such function
requires ICALL.EQ.4):
...
$SUBS ADVAN6 TOL=5
$MODEL
...
$PK
...
IF (ICALL.EQ.4) THEN
CALL RANDOM (2,R) ;Rand number in[0,1[
T=1-R
...
ENDIF
$DES
...
$ERROR
...
$ESTIMATION MAXEVAL=0 NUMERICAL METHOD=COND LAPLACE LIKE CENTERING PRINT=2
MSFO=msfo3
Thank you!
Nieves
--
--------------------------------
Nieves Velez de Mendizabal, Ph.D
Departamento de Farmacia y Tecnología Farmacéutica
Facultad de Farmacia
Universidad de Navarra
Phone: (+34) 658 732 851
Phone: (+34) 948 255 400 ext. 5827
[email protected]<mailto:[email protected]>
--------------------------------
--
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]<mailto:[email protected]>
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
________________________________
This e-mail (including any attachments) is confidential and may be legally
privileged. If you are not an intended recipient or an authorized
representative of an intended recipient, you are prohibited from using, copying
or distributing the information in this e-mail or its attachments. If you have
received this e-mail in error, please notify the sender immediately by return
e-mail and delete all copies of this message and any attachments.
Thank you.
Nieves,
Two different thoughts using typical code. I'm not sure how you would provide a seed for the CALL RANDOM using verbatim code.
(1) Generate the random numbers (one for all obs. or one for each time value) in another software and import them as part of the data (as Nick suggested).
(2) If it is important that NONMEM generate the random numbers then you could try:
a) Run your control stream once using SIMONLY to obtain the random numbers. Also have a line that saves the input DV to a new variable name (i.e. CP =DV). In the table file use the NOAPPEND option to prevent adding (DV,PRED,RES,WRES), output all the variables in the input datasets + the random number variable + CP.
---------------
$PK
CP=DV
IF(ICALL.EQ.4)THEN
CALL RANDOM (2,R)
RN=1-R
ENDIF
....
$SIMULATION (1234123) (12343 UNIFORM) ONLYSIM
$TABLE ID TIME DAY TSLD AMT RATE DUR DOSE N CP CVLQ CMT MDV EVID WTKG SAMP NUM RN FORMAT=sF8.4 NOAPPEND NOPRINT NOHEADER FILE=xxx.tbl
-------------------
b) Use the table file from (a) as the input dataset. However, you will need to use CP as DV (since the original DV values were replaced in the simulation run).
-------------------
$DATA trash-1.tbl
$INPUT ID TIME DAY=DROP TSLD AMT RATE DUR DOSE N DV CVLQ
CMT MDV EVID WTKG SAMP NUM RN
$PK
Normal model statements
Use RN as your random number
Note: TIME is now the time that NONMEM computes (elapsed time from first observation), so I dropped DAY. DV is really CP (the original DV). Also, note that you may want to specify the format (option in NM7) on the output table to maintain the precision of the numbers in your original input file.
---------------------
Best regards,
Luann Phillips
Director, PK/PD
Cognigen Corporation
Nieves Vélez de Mendizabal wrote:
> Dear NONMEM Users,
>
> I’m developing a model with NONMEM6|7. In this model, I need to generate a random number at every time step. The problem is that the use of the function “CALL RANDOM (2,R)” is supposed to be only for simulation, isn’t it?
>
> Thus, is it possible to generate a random number at every time step with NONMEM? Does anybody know how?
>
> This is part of the code that does not work (it’s not working because the condition ICALL.EQ.4 never happens and for that it’s not getting into the “if”, but on the other hand, in order to use the function RANDOM(2,R), such function requires ICALL.EQ.4):
>
> …
>
> $SUBS ADVAN6 TOL=5
>
> $MODEL
>
> …
>
> $PK
>
> …
>
> IF (ICALL.EQ.4) THEN
>
> CALL RANDOM (2,R) ;Rand number in[0,1[
>
> T=1-R
>
> …
>
> ENDIF
>
> $DES
>
> …
>
> $ERROR
>
> …
>
> $ESTIMATION MAXEVAL=0 NUMERICAL METHOD=COND LAPLACE LIKE CENTERING PRINT=2 MSFO=msfo3
>
> Thank you!
>
> Nieves
>
> --
> --------------------------------
> Nieves Velez de Mendizabal, Ph.D
> Departamento de Farmacia y Tecnología Farmacéutica
> Facultad de Farmacia
> Universidad de Navarra
> Phone: (+34) 658 732 851
> Phone: (+34) 948 255 400 ext. 5827
> [email protected]
>
>
Dear all,
A single control stream with $SIM followed by $EST could possibly do the trick
(similar to what Luann suggests below, using two control streams). However,
before we air additional suggestions on how we may or may not achieve
random-number generation during estimation; maybe it is better to hold off
until we hear back from Nieves? - As several persons have already indicated
they can not see much use of this feature and we do not know specifically how
Nieves meant to use it.
Best regards
Jakob
Quoted reply history
-----Original Message-----
From: [email protected] [mailto:[email protected]] On
Behalf Of Luann Phillips
Sent: 14 March 2011 15:00
To: Nieves Vélez de Mendizabal
Cc: [email protected]
Subject: Re: [NMusers] How to generate a random number with $EST
Nieves,
Two different thoughts using typical code. I'm not sure how you would
provide a seed for the CALL RANDOM using verbatim code.
(1) Generate the random numbers (one for all obs. or one for each time
value) in another software and import them as part of the data (as Nick
suggested).
(2) If it is important that NONMEM generate the random numbers then you
could try:
a) Run your control stream once using SIMONLY to obtain the random
numbers. Also have a line that saves the input DV to a new variable name
(i.e. CP =DV). In the table file use the NOAPPEND option to prevent
adding (DV,PRED,RES,WRES), output all the variables in the input
datasets + the random number variable + CP.
---------------
$PK
CP=DV
IF(ICALL.EQ.4)THEN
CALL RANDOM (2,R)
RN=1-R
ENDIF
....
$SIMULATION (1234123) (12343 UNIFORM) ONLYSIM
$TABLE ID TIME DAY TSLD AMT RATE DUR DOSE N CP CVLQ CMT MDV EVID WTKG
SAMP NUM RN FORMAT=sF8.4 NOAPPEND NOPRINT NOHEADER FILE=xxx.tbl
-------------------
b) Use the table file from (a) as the input dataset. However, you will
need to use CP as DV (since the original DV values were replaced in the
simulation run).
-------------------
$DATA trash-1.tbl
$INPUT ID TIME DAY=DROP TSLD AMT RATE DUR DOSE N DV CVLQ
CMT MDV EVID WTKG SAMP NUM RN
$PK
Normal model statements
Use RN as your random number
Note: TIME is now the time that NONMEM computes (elapsed time from first
observation), so I dropped DAY. DV is really CP (the original DV). Also,
note that you may want to specify the format (option in NM7) on the
output table to maintain the precision of the numbers in your original
input file.
---------------------
Best regards,
Luann Phillips
Director, PK/PD
Cognigen Corporation
Nieves Vélez de Mendizabal wrote:
> Dear NONMEM Users,
>
>
> I'm developing a model with NONMEM6|7. In this model, I need to generate
> a random number at every time step. The problem is that the use of the
> function "CALL RANDOM (2,R)" is supposed to be only for simulation,
> isn't it?
>
> Thus, is it possible to generate a random number at every time step with
> NONMEM? Does anybody know how?
>
> This is part of the code that does not work (it's not working because
> the condition ICALL.EQ.4 never happens and for that it's not getting
> into the "if", but on the other hand, in order to use the function
> RANDOM(2,R), such function requires ICALL.EQ.4):
>
> ...
>
> $SUBS ADVAN6 TOL=5
>
> $MODEL
>
> ...
>
> $PK
>
> ...
>
> IF (ICALL.EQ.4) THEN
>
> CALL RANDOM (2,R) ;Rand number in[0,1[
>
> T=1-R
>
> ...
>
> ENDIF
>
> $DES
>
> ...
>
> $ERROR
>
> ...
>
> $ESTIMATION MAXEVAL=0 NUMERICAL METHOD=COND LAPLACE LIKE CENTERING
> PRINT=2 MSFO=msfo3
>
>
>
> Thank you!
>
> Nieves
>
> --
> --------------------------------
> Nieves Velez de Mendizabal, Ph.D
> Departamento de Farmacia y Tecnología Farmacéutica
> Facultad de Farmacia
> Universidad de Navarra
> Phone: (+34) 658 732 851
> Phone: (+34) 948 255 400 ext. 5827
> [email protected]
> --------------------------------
>
Dear nmusers,
Thanks everybody. The generation of a random number from a uniform distribution finally works. Below you can see the code. You also you can see the table that NONMEM7 generated. The generated random values can be seen at the 6th column. The aim of this model is to estimate a parameter (PERT) that I use to predict a probability (PROB). Such probability represents the probability of suffering one event. There is individual variability in such probability. During the estimation process I have to produce a train of events for every individual. I don't know yet if it makes sense to do that in the way I'm trying to, but I thought it could be worth it.
$PK
.....
PERT=THETA(3)*EXP(ETA(1))
PROB=EXP(PERT)/(1+EXP(PERT))
...
$DES
PER=0
KK=0
"KK=RAND(0)
IF (KK.LT.PROB) THEN ;perturbation
PER=1
ENDIF
....
$ESTIMATION MAXEVAL=9990 NUMERICAL METHOD=COND LAPLACE LIKE CENTERING PRINT=1 NOABORT MSFO=msfo5
$TABLE ID TIME DV PER PROB KK PER NOAPPEND ONEHEADER NOPRINT FILE=sdtab5
TABLE NO. 1
ID TIME DV PER PROB KK
1.0000E+00 1.0000E+00 0.0000E+00 1.0000E+00 5.2498E-01 2.1428E-01
1.0000E+00 2.0000E+00 1.0000E+00 0.0000E+00 5.2498E-01 6.6192E-01
1.0000E+00 3.0000E+00 1.0000E+00 0.0000E+00 5.2498E-01 7.7577E-01
1.0000E+00 4.0000E+00 0.0000E+00 0.0000E+00 5.2498E-01 6.1920E-01
1.0000E+00 5.0000E+00 0.0000E+00 1.0000E+00 5.2498E-01 1.6852E-01
1.0000E+00 6.0000E+00 1.0000E+00 1.0000E+00 5.2498E-01 2.1436E-01
1.0000E+00 7.0000E+00 0.0000E+00 0.0000E+00 5.2498E-01 6.9345E-01
1.0000E+00 8.0000E+00 0.0000E+00 0.0000E+00 5.2498E-01 5.6443E-01
1.0000E+00 9.0000E+00 1.0000E+00 0.0000E+00 5.2498E-01 6.0616E-01
1.0000E+00 1.0000E+01 0.0000E+00 0.0000E+00 5.2498E-01 8.2914E-01
1.0000E+00 1.1000E+01 0.0000E+00 0.0000E+00 5.2498E-01 5.9900E-01
1.0000E+00 1.2000E+01 1.0000E+00 0.0000E+00 5.2498E-01 7.6669E-01
1.0000E+00 1.3000E+01 0.0000E+00 0.0000E+00 5.2498E-01 8.5406E-01
1.0000E+00 1.4000E+01 0.0000E+00 0.0000E+00 5.2498E-01 8.0210E-01
--
--------------------------------
Nieves Velez de Mendizabal, Ph.D
Departamento de Farmacia y Tecnología Farmacéutica
Facultad de Farmacia
Universidad de Navarra
Phone: (+34) 658 732 851
Phone: (+34) 948 255 400 ext. 5827
[email protected]