```
From: "Tsai, Max" max.tsai@spcorp.com
Subject: [NMusers] Simulating different populations
Date: Wed, 10 Aug 2005 10:35:42 -0400

Can NONMEM perform simulations to incorporate variability associated with
fixed effects across populations?  When one develops a NM model, NONMEM
output lists the final parameter estimates and the standard error estimates
for those parameters.  Can NONMEM use those standard errors to simulate
different populations so that the fixed effects are not so "fixed"?  For
instance, if SUBPROBS=2, I would like to see one population have one set of
fixed effects and a second population have a different set of fixed effects.
For a given population, the fixed effects would be determined by the
population parameter estimates plus some population variability whose
distribution is specified by the standard error.  This would be on the
population level between replicate studies and thus a step higher than BSV
and WSV.  Is this capability built into NONMEM?  If so, how would it be
done?  Thanks.

-Max

Max Tsai, Ph.D.
Senior Scientist (DM/PK)
Schering-Plough Corporation
K-15-2-2650
Kenilworth, NJ 07033-0530
*:  (908) 740-3911
*:  (908) 740-2916
*:  max.tsai@spcorp.com
_______________________________________________________

From: "Gastonguay, Marc"
Subject: Re: [NMusers] Simulating different populations
Date: Wed, 10 Aug 2005 11:56:36 -0400

Max,
If I understand you correctly, it sounds like you want to simulate from the posterior
(or uncertainty) distributions on the population parameters (THETA, OMEGA and SIGMA).
NONMEM is an approximate maximum likelihood method, so you don't get true posteriors but
you could use the var-cov matrix of the estimates as an approximation. Implementation is
not so straightforward. To do this correctly you need to simulate from a multivariate
normal dstribution for THETAs and an Inverse Wishart distribution for OMEGA and SIGMA.
There is an experimental, unsupported and undocumented subroutine in NONMEM, known as
the PRIOR subroutine, which could implement this sort of simulation, but this is not
publicly available and I have no idea if/when it may be included in future versions on
NONMEM. You could also do this by sampling THETAs OMEGAs and SIMGAs from a posterior
distribution generated from a parametric or non-parametric bootstrap - but you'd have
to write your own program to do so. We (and I believe others) have been working on our
own utilities to implement these types of simulations in NONMEM and will share results
with the user community in the near future.

Another option is to move to a full Bayesian MCMC method, such as WinBUGS, which generates
posterior distributions for any parameters and derived quantities you're interested in.
You could then simulate from these distributions.

In the meantime, you could try this bit of code (below) in the NMTRAN control stream, which
allows for a simple simulation from independent posterior distributions for any population
fixed effects parameters (could also be applied to variance terms if you use a log-normal
distribution - but this is a crude approximation). The result is a simulation with inter-trial
uncertainty in the population fixed effects parameters.

Hope this helps.

Marc
=========================
Marc R. Gastonguay, Ph.D.
Metrum Research Group LLC
www.metrumrg.com

;To define inter-trial uncertainty in a parameter,use code similar to this  in \$PRED or \$PK...
; This is from \$PRED example where uncertainty in relative bio fraction (FR1) is simulated

IF(IREP.EQ.1.AND.NEWIND.EQ.0) COUNT = 0
IF(ICALL.EQ.4) THEN
IF(IREP.NE.COUNT) THEN  ;SHOULD ONLY HAPPEN AT THE START OF EACH SIMULATION REPLICATE
FR1 = -1
DOWHILE(FR1.LE.0)
CALL RANDOM(2,R)
FR1 = THETA(2) + THETA(3)*R
ENDDO
ENDIF      COUNT = IREP
;Define model using FR1...
AUCR = THETA(1) + ETA(1)
AUCT = AUCR * FR1 + ETA(2)
Y = AUCR*(1-FLG) + AUCT*FLG + EPS(1)
ENDIF

TRL=IREP ; trial replicate number

\$THETA 100 ; TYPICAL REFERENCE AUC
0.8   ; MODE OF POSTERIOR DISTRIBUTION ( POINT ESTIMATE OF FR1)
0.075  ;THETA(3) IS SD OF UNCERTAINTY, OR STANDARD ERROR FROM NONMEM \$COV
\$OMEGA 100 100
\$SIGMA 10 1
\$SIMULATION (437565 NEW) (47003 NORMAL) ONLYSIM SUBPROBLEMS=100
\$TABLE ID DAY FLG AUCR AUCT F1 TRL DV FILE = SIMUNCERT.tab NOAPPEND NOPRINT NOHEADER

_______________________________________________________

From: Liping Zhang ZHANG_LIPING@lilly.com
Subject: Re: [NMusers] Simulating different populations
Date: Wed, 10 Aug 2005 11:09:16 -0500

Hi, Max,

what I typically do, is to sample (let's say 1000 replicates) from the
distribution of "fixed effects" parameters in Splus, given estimated
population parameter value and their variance-covariance matrix generated
by NONMEM. So now you have a 1000 different "population", each has
different "fixed effects", but overall they reflect the population
parameter estimates and their uncertainty. Then through some simple Unix
C-shell script I simulate the outcome from each "population" (with fixed
population parameters and BSV and WSV)  and do statistical summary
according to modeling needs. As Ken Kowalski pointed out not long ago,
this procedure is computer-labor intensive, you need to simulate large
number of population, and within each population, a sufficient number
trials. I would not go for this route unless I know for my modeling need
it requires this kind of "truthful representation".

best regards, Liping

Liping Zhang, PhD

Global PK/PD Modeling and Trial Simulation
ELi Lilly & Company
Work: 317-277-8687       Fax: 317-433-6661
DC 0734
Email: zhang_liping@lilly.com
_______________________________________________________

From: "Kowalski, Ken" Ken.Kowalski@pfizer.com
Subject: RE: [NMusers] Simulating different populations
Date: Wed, 10 Aug 2005 13:12:09 -0400

Liping,

Just to clarify.  I indicated that simulating different sets of population
estimates based on the var-cov matrix from NONMEM and then assuming
multivariate normal and possibly inverse-Wishart (for the variance
parameters in Omega and/or Sigma) distributions for the parameter estimates
is NOT computationally intensive.  However, it is laborious to use these
population estimates and feed them back in to NONMEM to perform simulations
that take into account this parameter uncertainty because it requires custom
coding.  Based on Marc's comments it sounds like he and others are working
on developing utilities that will automate this process in NONMEM which
would make it a lot easier to do this as a general practice.

If on the other hand we do not wish to make a parametric assumption using
the var-cov matrix and multivariate Normal/Inverse-Wishart distribution we
could perform parametric or nonparametric bootstrap simulations and re-fit
the model in NONMEM for each of say 1000 bootstrap datasets to obtain 1000
estimates of the population parameters from the posterior distribution and
use these in subsequent simulations to account for parameter uncertainty.
This approach IS computationally intensive because of having to re-fit the
model in NONMEM.

My feeling is that it is better to routinely do something to take into
account parameter uncertainty using the var-cov matrix and parametric
assumptions than to doing nothing to take into account parameter
uncertainty.  Hopefully, individuals working on utilities to automate this
parametric approach will be able to share them with us.

Ken

_______________________________________________________

From: Liping Zhang ZHANG_LIPING@lilly.com
Subject: RE: [NMusers] Simulating different populations
Date: Wed, 10 Aug 2005 12:20:37 -0500

Yes, Ken. I understand what you are saying, that is why I did not say
"computationally intensive", but says "computer-labor intensive",
indicating the long time it takes to get the final results and the custom
coding.

I did not follow up with the previous posting closely. I would be
interested to know who are the "he (Marc) and others" working on
developing utilities that will automate this process in NONMEM.

best regards, Liping

Liping Zhang, PhD

Global PK/PD Modeling and Trial Simulation
ELi Lilly & Company
Work: 317-277-8687       Fax: 317-433-6661
DC 0734
Email: zhang_liping@lilly.com
_______________________________________________________

Subject: RE : [NMusers] Simulating different populations
Date: Thu, 11 Aug 2005 12:43:13 +0200

Hi Marc and Max,

My 20 cents...

In order to do what Marc suggests, you make the assumptions:
(1) that NONMEM fixed effects posterior distribution is
multivariate normal (MVN)

(2) parameterized with mean equal to the THETA estimates and variance to
\$COV var-cov matrix from which you exclude random effects (co)variance
lines and columns, and

(3) that you have identified all covariates and/or multimodal distributions
using mixture models.

Those assumptions are probably not true for many parameters. For example clearance
is more likely to be log-normally distributed. F bioavailability is classically
bounded by 0 and 1. So in order to use your MVN assumption, you will need first
to re-parameterize your model with exponential, logistic functions, ... (for more
details, see http://www.page-meeting.org/page/page2002/PascalGirardPage2002.pdf )

Assuming that you have a function or subroutine to resample from MVN (R and Splus do
it very easily), second difficulty is to pass var-cov matrix from NONMEM output to your
MVN resampling function. Depending on the complexity of your model, this matrix can be
very large and keying all those numbers a real pain in ... and source of errors! So for
doing this, I have developed an INFN subroutine (to use with PREDPP or an equivalent for
\$PRED) that outputs THETAs and var-cov matrix in an ASCII file. Probably, "Marc and
others" have developed similar tools.

Anyway, I am still dreaming of an integrated tool doing all this without extra-coding...

Hope it helps,

Dr Pascal Girard
EA 3738, Ciblage Thérapeutique en Oncologie
Fac Médecine Lyon-Sud, BP12
69921 OULLINS Cedex
France
Tel  +33 (0)4 26 23 59 54 / Fax +33 (0)4 26 23 59 49
_______________________________________________________

From: "Gastonguay, Marc" marcg@metrumrg.com
Subject: Re: [NMusers] Simulating different populations
Date: Thu, 11 Aug 2005 08:59:42 -0400

Hi, Pascal.
Thanks for pointing out these limitations and assumptions. I think your point
about fixed effects parameters that have natural bounds is an excellent one; we do
we have to think carefully about the parameterization and choice of posterior
distributions. For this reason it may be simpler (but more computationally intensive)
to simulate from empirical distributions, such as those that might result from a
nonparametric bootstrap.

For the implementation part, I agree that an automated extraction of relevant components
is mandatory. It sure would be nice to have an all-in-one tool that did this for us,
but I think that this discussion also illustrates the wonderful flexibility of an ASCII
text-based, command-line program like NONMEM, when coupled with a useful programming
language (we use R and Perl).