Hello NONMEM Users,
As SAEM does not provide a useful objective function, the manuals recommend using IMP after SAEM. It works well in many cases when IMP works well. When IMP works well, SAEM is not always needed. SAEM is really needed when the other methods do not work well.
The issue is that there are hard cases when SAEM works very well and IMP does not work at all. SAEM provides meaningful and consistent PK/PD parameters across very different runs, while IMP provides objective function, which varies so greatly that it looks meaningless. Another potential issue with IMP is that even when it works well with a problem, it occasionally provides low values of objective functions after SAEM or as the first estimate (less frequently in the middle of IMP run) and then becomes unstable or jumps to much higher objective function and then converges to something between the low and the high values for a long time. It almost looks like IMP does not show some kind of integration/computation errors and keeps running providing a funny objective function.
It seems like we cannot estimate objective function when SAEM runs well and IMP does not. It reduces the value of SAEM. Is there a way around it?
Thanks,
Pavel
SAEM and IMP
Pavel:
Try using the latest nonmem 7.3 version if you can, and see if this still
occurs. Improvements are continually being made. Try also the following after
an SAEM step:
$EST METHOD=IMP EONLY=1 MAPITER=0 (...plus other options as neded)
Or
$EST METHOD=IMPMAP EONLY=1 (...plus other options as needed)
Robert J. Bauer, Ph.D.
Vice President, Pharmacometrics, R&D
ICON Development Solutions
7740 Milestone Parkway
Suite 150
Hanover, MD 21076
Tel: (215) 616-6428
Mob: (925) 286-0769
Email: [email protected]
Web: www.iconplc.com
Quoted reply history
-----Original Message-----
From: [email protected] [mailto:[email protected]] On
Behalf Of Pavel Belo
Sent: Thursday, May 15, 2014 11:02 AM
To: [email protected]
Subject: [NMusers] SAEM and IMP
Hello NONMEM Users,
As SAEM does not provide a useful objective function, the manuals recommend
using IMP after SAEM. It works well in many cases when IMP works well. When
IMP works well, SAEM is not always needed. SAEM is really needed when the
other methods do not work well.
The issue is that there are hard cases when SAEM works very well and IMP does
not work at all. SAEM provides meaningful and consistent PK/PD parameters
across very different runs, while IMP provides objective function, which varies
so greatly that it looks meaningless. Another potential issue with IMP is that
even when it works well with a problem, it occasionally provides low values of
objective functions after SAEM or as the first estimate (less frequently in the
middle of IMP run) and then becomes unstable or jumps to much higher objective
function and then converges to something between the low and the high values
for a long time. It almost looks like IMP does not show some kind of
integration/computation errors and keeps running providing a funny objective
function.
It seems like we cannot estimate objective function when SAEM runs well and IMP
does not. It reduces the value of SAEM. Is there a way around it?
Thanks,
Pavel
Pavel,
You might try IMP after SAEM with EONLY=1 (and a limited number of iterations)
so that the parameters remain stationary and only the expectation step is
executed.
Cheers... Brian
Quoted reply history
-----Original Message-----
From: [email protected] [mailto:[email protected]] On
Behalf Of Pavel Belo
Sent: Thursday, May 15, 2014 11:02 AM
To: [email protected]
Subject: [NMusers] SAEM and IMP
Hello NONMEM Users,
As SAEM does not provide a useful objective function, the manuals recommend
using IMP after SAEM. It works well in many cases when IMP works well. When
IMP works well, SAEM is not always needed. SAEM is really needed when the
other methods do not work well.
The issue is that there are hard cases when SAEM works very well and IMP does
not work at all. SAEM provides meaningful and consistent PK/PD parameters
across very different runs, while IMP provides objective function, which varies
so greatly that it looks meaningless. Another potential issue with IMP is that
even when it works well with a problem, it occasionally provides low values of
objective functions after SAEM or as the first estimate (less frequently in the
middle of IMP run) and then becomes unstable or jumps to much higher objective
function and then converges to something between the low and the high values
for a long time. It almost looks like IMP does not show some kind of
integration/computation errors and keeps running providing a funny objective
function.
It seems like we cannot estimate objective function when SAEM runs well and IMP
does not. It reduces the value of SAEM. Is there a way around it?
Thanks,
Pavel
Hi Pavel
I have experienced a similar problem. In my case, the
following code for IMP after SAEM (using NM7.3) greatly reduced the
Monte Carlo OFV noise from variations of about +/- 60 points to variations of
+/-
6 points (though still not good enough for covariate testing):
$EST METHOD=IMP LAPLACE INTER NITER=15 ISAMPLE=3000 EONLY=1 DF=7
IACCEPT=0.3
ISAMPEND=10000 STDOBJ=2 MAPITER=0 PRINT=1 SEED=123456 RANMETHOD=3S2
The settings are explained in the NM7.3 guide. If you are using
NM7.3, you can also try IACCEPT=0.0 whereupon "NONMEM will determine the
most appropriate IACCEPT level for each subject". Of course the settings
for DF and IACCEPT in the above code will depend on the type of data you have.
Which brings me to my own question. If I have both continous and categorical
DVs in the dataset (which would mean different optimal settings) and I
am using F_FLAG accordingly, what would the 'right' values of DF and
IACCEPT be? I have noticed that the DF automatically chosen by NONMEM for
individuals in the dataset can vary from 0-8 and this appears to be random.
Emmanuel
Quoted reply history
From: Pavel Belo <[email protected]>
>To: [email protected]
>Sent: Thursday, May 15, 2014 11:01 AM
>Subject: [NMusers] SAEM and IMP
>
>
>Hello NONMEM Users,
>
>As SAEM does not provide a useful objective function, the manuals
>recommend using IMP after SAEM. It works well in many cases when IMP
>works well. When IMP works well, SAEM is not always needed. SAEM is
>really needed when the other methods do not work well.
>
>The issue is that there are hard cases when SAEM works very well and IMP
>does not work at all. SAEM provides meaningful and consistent PK/PD
>parameters across very different runs, while IMP provides objective
>function, which varies so greatly that it looks meaningless. Another
>potential issue with IMP is that even when it works well with a problem,
>it occasionally provides low values of objective functions after SAEM or
>as the first estimate (less frequently in the middle of IMP run) and
>then becomes unstable or jumps to much higher objective function and
>then converges to something between the low and the high values for a
>long time. It almost looks like IMP does not show some kind of
>integration/computation errors and keeps running providing a funny
>objective function.
>
>It seems like we cannot estimate objective function when SAEM runs well
>and IMP does not. It reduces the value of SAEM. Is there a way around
>it?
>
>Thanks,
>Pavel
>
>
>
Hi Emmanuel,
While I am a strong advocate of using quasi-random rather than pseudo- random
sequences for importance sampling in EM methods like IMP, there is a
theoretical (and very real) problem with their use in the context you
suggested in your message, namely with a multivariate t distribution as the
importance sampling distribution. The 3S2 option implies you are using a Sobol
quasi-random sequence, while
the DF=7 implies the use of a multivariate T-distribution with 7 degrees of
freedom. The standard way of generating
a p-dimensional multivariate t -random variable with DF degrees of freedom is
to generate a p-dimensional multivariate normal and then divide by an
additional independent random variable which is basically the square root of a
1-d chi square random variable with DF degrees of freedom. Thus to generate a
p-dimensional importance sample, you actually need to use p+1 independent
random variables. If you simply use a p+1 dimensional Sobol vector as the
base quasi-random draw, the nonlinear mapping from p+1 dimensions to the final
p dimensional result destroys the low discrepancy property of the final
sequence in the p-dimensional space and in fact introduces a significant
amount of bias in the final result. The problem arises directly from the p+1
vs p dimensional mismatch.
There is no problem if the final p-dimensional result can be generated from a
p-dimensional quasi-random sequence, which is the case for multivariate normal
Importance samples. So quasi random sequences should really only be used for
the DF=0 multivariate normal importance sampling distribution case, not the
multivariate DF>0 multivariate t case.
I ran across this effect in testing the Sobol-based importance sampling EM
algorithm QRPEM in Phoenix NLME. It is very real and the net effect is to
introduce a significant bias. There is a partial fix that works but gives up
some of the benefit of using low-discrepancy sequences - namely use a
p-dimensional quasi-random vector to generate the p-dimensional multivariate
normal, but
then use a 1-d pseudo-random sequence to generate the chi-square random
variable.
Quoted reply history
From: [email protected] [mailto:[email protected]] On
Behalf Of Emmanuel Chigutsa
Sent: Thursday, May 15, 2014 1:03 PM
To: Pavel Belo; [email protected]
Subject: Re: [NMusers] SAEM and IMP
Hi Pavel
I have experienced a similar problem. In my case, the following code for IMP
after SAEM (using NM7.3) greatly reduced the Monte Carlo OFV noise from
variations of about +/- 60 points to variations of +/- 6 points (though still
not good enough for covariate testing):
$EST METHOD=IMP LAPLACE INTER NITER=15 ISAMPLE=3000 EONLY=1 DF=7 IACCEPT=0.3
ISAMPEND=10000 STDOBJ=2 MAPITER=0 PRINT=1 SEED=123456 RANMETHOD=3S2
The settings are explained in the NM7.3 guide. If you are using NM7.3, you can
also try IACCEPT=0.0 whereupon "NONMEM will determine the most appropriate
IACCEPT level for each subject". Of course the settings for DF and IACCEPT in
the above code will depend on the type of data you have. Which brings me to my
own question. If I have both continous and categorical DVs in the dataset
(which would mean different optimal settings) and I am using F_FLAG
accordingly, what would the 'right' values of DF and IACCEPT be? I have noticed
that the DF automatically chosen by NONMEM for individuals in the dataset can
vary from 0-8 and this appears to be random.
Dear Emmanuel and Pavel,
Further to Bob's answer, recall also that delta OFV in the likelihood ratio
test is only asymptotocaily chi squared distributed,and this is not the only
reason why you should not get too hung up on OFV to help choose your models.
For example, Lavielle 2010 in Biometrics showed nicely how the SAEM algorithm
can estimate parameters of complex differential equation models for joint HIV
viral load and CD4 counts. Using OFV-based metrics a latent model was chosen
whereby the majority of circulating T-cells were infected with virus - this
model also gave nice fits to the data. When immunologists look at CD4 cells of
HIV infected patients however, they find that much less than 1% (closer to
0.01%) of circulating T-cells contain virus (most of the virus making up the
latent reservoir is stuck to folicular cells in the periphery), so one would
have to question the meaning of the parameters identified as the best fit by
SAEM. By all means use SAEM to fit ODE models that don't run/converge with
other algorithms (I do), but choose models with parameters that make
mechanistic sense rather than relying too heavily on OFV-based metrics. A nice
VPC always goes down well too.
Joe
Joseph F Standing
MRC Fellow, UCL Institute of Child Health
Antimicrobial Pharmacist, Great Ormond Street Hospital
Tel: +44(0)207 905 2370
Mobile: +44(0)7970 572435
Quoted reply history
________________________________________
From: [email protected] [[email protected]] On Behalf Of
Bob Leary [[email protected]]
Sent: 15 May 2014 19:22
To: Emmanuel Chigutsa; Pavel Belo; [email protected]
Subject: RE: [NMusers] SAEM and IMP
Hi Emmanuel,
While I am a strong advocate of using quasi-random rather than pseudo- random
sequences for importance sampling in EM methods like IMP, there is a
theoretical (and very real) problem with their use in the context you
suggested in your message, namely with a multivariate t distribution as the
importance sampling distribution. The 3S2 option implies you are using a Sobol
quasi-random sequence, while
the DF=7 implies the use of a multivariate T-distribution with 7 degrees of
freedom. The standard way of generating
a p-dimensional multivariate t -random variable with DF degrees of freedom is
to generate a p-dimensional multivariate normal and then divide by an
additional independent random variable which is basically the square root of a
1-d chi square random variable with DF degrees of freedom. Thus to generate a
p-dimensional importance sample, you actually need to use p+1 independent
random variables. If you simply use a p+1 dimensional Sobol vector as the
base quasi-random draw, the nonlinear mapping from p+1 dimensions to the final
p dimensional result destroys the low discrepancy property of the final
sequence in the p-dimensional space and in fact introduces a significant
amount of bias in the final result. The problem arises directly from the p+1
vs p dimensional mismatch.
There is no problem if the final p-dimensional result can be generated from a
p-dimensional quasi-random sequence, which is the case for multivariate normal
Importance samples. So quasi random sequences should really only be used for
the DF=0 multivariate normal importance sampling distribution case, not the
multivariate DF>0 multivariate t case.
I ran across this effect in testing the Sobol-based importance sampling EM
algorithm QRPEM in Phoenix NLME. It is very real and the net effect is to
introduce a significant bias. There is a partial fix that works but gives up
some of the benefit of using low-discrepancy sequences – namely use a
p-dimensional quasi-random vector to generate the p-dimensional multivariate
normal, but
then use a 1-d pseudo-random sequence to generate the chi-square random
variable.
From: [email protected] [mailto:[email protected]] On
Behalf Of Emmanuel Chigutsa
Sent: Thursday, May 15, 2014 1:03 PM
To: Pavel Belo; [email protected]
Subject: Re: [NMusers] SAEM and IMP
Hi Pavel
I have experienced a similar problem. In my case, the following code for IMP
after SAEM (using NM7.3) greatly reduced the Monte Carlo OFV noise from
variations of about +/- 60 points to variations of +/- 6 points (though still
not good enough for covariate testing):
$EST METHOD=IMP LAPLACE INTER NITER=15 ISAMPLE=3000 EONLY=1 DF=7 IACCEPT=0.3
ISAMPEND=10000 STDOBJ=2 MAPITER=0 PRINT=1 SEED=123456 RANMETHOD=3S2
The settings are explained in the NM7.3 guide. If you are using NM7.3, you can
also try IACCEPT=0.0 whereupon "NONMEM will determine the most appropriate
IACCEPT level for each subject". Of course the settings for DF and IACCEPT in
the above code will depend on the type of data you have. Which brings me to my
own question. If I have both continous and categorical DVs in the dataset
(which would mean different optimal settings) and I am using F_FLAG
accordingly, what would the 'right' values of DF and IACCEPT be? I have noticed
that the DF automatically chosen by NONMEM for individuals in the dataset can
vary from 0-8 and this appears to be random.
Bob:
Yes, when I was developing the code for using the Sobol method in NONMEM, I
too found that using the t-distribution with Sobol/Quasi normal method resulted
in biased assessments of the objective function, when I used the unmodified
standard technique of generating t-samples. So, I experimented with various
alternative methods like you did, and the final algorithms used in NONMEM 7.2
and 7.3 appear to avoid the bias, at least when I tried some simple models,
such as the example1 problem in the NONMEM ..\examples directory. For such a
problem, I tried the following DF settings, with and without using the SOBOL
method (RANMETHOD=3 (default), or RANMETHOD=3S2), DF=0, 1, 2, 3, 4, 5, 6, 7, 8,
10, 11, and 12. I think beyond 12 DF is sufficiently normal like, and not much
use.
In all cases, the OBJF was within 1 unit of each other, and always within a STD
of the Monte Carlo noise to the OBJF (I used EONLY=1, 20 iterations to obtain
OBJFV replicates for each setting). By the way, the .cnv file publishes the
mean and standard deviation of the last CITER objective function values (see
nm730.pdf manual on root.cnv).
Certainly it is worth assuring this lack of bias under other kinds of, and more
complex, models, just to make sure, so it is worth while for modelers to use
caution when mixing Sobol with DF>0.
I should point out that odd-number DF’s such as 7 makes the corrective
algorithm for Sobol with t-distributions work harder, whereas even-numbered
DF’s are a little cleaner. So I recommend even-numbered DF’s, or 1, up to 10
to work with Sobol, at least as far as my algorithm is concerned (so
DF=1,2,4,6,8,10 seem most efficient to use with Sobol).
The algorithm I use for t-distributions is in the subroutine TDEV2, and is in
..source\general.f90. I publish it because I borrowed freely from the
literature, and I consider it therefore “open-source”.
Just a note for orienting you to the code, IRANM=5 refers to the Sobol method.
I am afraid my theoretical expertise in the Sobol/Quasi random methods is not
as good as yours, so please feel free to review the code, and perhaps let me
and the modeling community know if there may be any potential pitfalls in using
this t-distribution algorithm with Sobol. I am also happy to consider
improvements you may come up with on this matter.
Robert J. Bauer, Ph.D.
Vice President, Pharmacometrics, R&D
ICON Development Solutions
7740 Milestone Parkway
Suite 150
Hanover, MD 21076
Tel: (215) 616-6428
Mob: (925) 286-0769
Email: [email protected]<mailto:[email protected]>
Web: http://www.iconplc.com/
Bob, thanks for your reply and pointer to the source code for tdev2 . I have
looked at it but perhaps you can clarify what the source of tdev2 is
and how it is used (I cannot see any calls made to it in your source code, so I
am not really
sure how it is being used). I can see for DF =1, tdev2 returns a Cauchy
distributed random variable, which is the 1D t-distribution for DF=1, and for
t=2 it looks like it returns the 1D
t-distribution for DF=2. For the other DF values it is not so clear, but I am
guessing it returns a 1D t-distribution with the specified DF.
If you simply assemble a p-dimensional sample vector with components which are
independent 1D t random variables with a specified DF , indeed it will avoid
the dimensional mismatch
and work OK without bias in the QR case, but it will not be a true
multivariate t-distribution (the components of a multivariate t-distribution
are not independent).
The usual argument for using a multivariate t as an importance sampler is that
a) It has fatter tails than a multivariate normal
b) It preserves the radial symmetry of the multivariate normal in the case
where the covariance matrix is the unit matrix (i.e. the density is only a
function of the distance from the origin)
Assembling a p-dimensional sample vector whose individual components are 1D t
distributions with a specified DF satisfies a), but not b). It’s a perfectly
valid importance sampling distribution, and indeed has been suggested in the
literature as much simpler alternative to the multivariate t, particularly in
low dimensions, but it is not a multivariate t (which is the more usual meaning
of ‘t-distribution’ in the multivariate importance sampling context) . It is
not radially symmetric (which may or may not be important, depending on the
posterior multidimensional shape and
the dimensionality. There is a directionality bias which increases with
dimensionality)
The documentation is somewhat ambiguous –
DF=4
The proposal density is to be t distribution with 4 degrees of freedom.
Based on the code I have seen, I now guess this means that each component of
the proposal is an independent t random variable with DF=4.
But I am only guessing – perhaps you could clarify what distribution is really
being referred to by ‘t-distribution’ in the p-dimensional case – multivariate
p-dimensional t with specified DF or
Cartesian product of p independent 1D t random variables with specified DF?
Quoted reply history
From: [email protected] [mailto:[email protected]] On
Behalf Of Bauer, Robert
Sent: Thursday, May 15, 2014 9:57 PM
To: [email protected]
Subject: RE: [NMusers] SAEM and IMP
Bob:
Yes, when I was developing the code for using the Sobol method in NONMEM, I
too found that using the t-distribution with Sobol/Quasi normal method resulted
in biased assessments of the objective function, when I used the unmodified
standard technique of generating t-samples. So, I experimented with various
alternative methods like you did, and the final algorithms used in NONMEM 7.2
and 7.3 appear to avoid the bias, at least when I tried some simple models,
such as the example1 problem in the NONMEM ..\examples directory. For such a
problem, I tried the following DF settings, with and without using the SOBOL
method (RANMETHOD=3 (default), or RANMETHOD=3S2), DF=0, 1, 2, 3, 4, 5, 6, 7, 8,
10, 11, and 12. I think beyond 12 DF is sufficiently normal like, and not much
use.
In all cases, the OBJF was within 1 unit of each other, and always within a STD
of the Monte Carlo noise to the OBJF (I used EONLY=1, 20 iterations to obtain
OBJFV replicates for each setting). By the way, the .cnv file publishes the
mean and standard deviation of the last CITER objective function values (see
nm730.pdf manual on root.cnv).
Certainly it is worth assuring this lack of bias under other kinds of, and more
complex, models, just to make sure, so it is worth while for modelers to use
caution when mixing Sobol with DF>0.
I should point out that odd-number DF’s such as 7 makes the corrective
algorithm for Sobol with t-distributions work harder, whereas even-numbered
DF’s are a little cleaner. So I recommend even-numbered DF’s, or 1, up to 10
to work with Sobol, at least as far as my algorithm is concerned (so
DF=1,2,4,6,8,10 seem most efficient to use with Sobol).
The algorithm I use for t-distributions is in the subroutine TDEV2, and is in
..source\general.f90. I publish it because I borrowed freely from the
literature, and I consider it therefore “open-source”.
Just a note for orienting you to the code, IRANM=5 refers to the Sobol method.
I am afraid my theoretical expertise in the Sobol/Quasi random methods is not
as good as yours, so please feel free to review the code, and perhaps let me
and the modeling community know if there may be any potential pitfalls in using
this t-distribution algorithm with Sobol. I am also happy to consider
improvements you may come up with on this matter.
Robert J. Bauer, Ph.D.
Vice President, Pharmacometrics, R&D
ICON Development Solutions
7740 Milestone Parkway
Suite 150
Hanover, MD 21076
Tel: (215) 616-6428
Mob: (925) 286-0769
Email: [email protected]<mailto:[email protected]>
Web: http://www.iconplc.com/
Third attempt to send to nmusers. Ignore if you have already received this
message.
Bob:
First. My apologies in case you are getting repeated e-mails from me. I am
finding I need to attempt several times for nmusers list-serv system to publish
my e-mails. Something to do with >40000 characters. So I need to cut off some
of the trailing e-mail portions that have already been published.
In my original investigations, the normal and t-distribution were using
Box-Mueller type random sample generations, which of course did not work well,
so I switched to inverse CDF methods. For normal distribution (DF=0), the IMP
uses the modified GASDEV, where you see that when ISELECT=5 (Sobol), DINVNORM
is used, and the the Box-Mueller method is used when Sobol method is not used.
Robert J. Bauer, Ph.D.
Vice President, Pharmacometrics, R&D
ICON Development Solutions
7740 Milestone Parkway
Suite 150
Hanover, MD 21076
Tel: (215) 616-6428
Mob: (925) 286-0769
Email: [email protected]<mailto:[email protected]>
Web: http://www.iconplc.com/
Quoted reply history
From: Bob Leary [mailto:[email protected]]
Sent: Monday, May 19, 2014 2:36 PM
To: Bauer, Robert; [email protected]<mailto:[email protected]>
Subject: RE: [NMusers] SAEM and IMP
Thanks, Bob. This indeed has been an interesting and thought provoking
discussion.
I too took another look at this. I thought, without really doing any analysis,
that the directionality bias using independent univariate t-components
would get really bad as the dimensionality increased, because that’s the way it
works with using fat-tailed independent double exponentials (samples there tend
to fail disproportionately
near the coordinate axes, and this gets arbitrarily bad as the dimensionality
increases without bound). But when I actually did the analysis in the student
t-case, it’s not bad at all and simply goes to a limit as the dimensionality
increases. So in fact, there may be some real advantages, particularly in
the Sobol quasi-random case, to doing it the way you do with independent
univariate t’s rather than using the multivariate t.
So my original message to Emmanuel to discourage him from using Sobol in
combination with the t-distribution was based on the erroneous assumption that
IMP was using the multivariate-t rather than
Independent univariate t’s. So his suggestion is probably OK (except as you
noted, that an even value of DF might be cleaner than DF=7 given the way tdev2
works). But there is still the question of why you were seeing a bias with
your original method before going to tdev2 –was it simply a different method of
generating the same distribution?
Along somewhat similar lines, I did notice that your routine GASDEV for
generating normal N(0,1) components (presumably used in the DF=0 case)
actually offers two different ways of doing it –
Box-Muller and inverse cdf of a uniform 0-1. It is not clear which one you get
when you simply specify DF=0. The inverse cdf method of course is by
definition OK for the quasi-random case – the cdf of the random vectors
generated with the Sobol distribution using the inverse cdf method will
simply be that Sobol distribution. But for Box Muller, this certainly will not
be the case, and it is not clear what the discrepancy of the Box Muller
generated CDF will be like. Certainly Box Muller is obviously OK in the
pseudorandom case, and will have the same discrepancy as the inverse cdf method
using pseudorandom uniform 0-1 inputs,
but whether the CDF of a Box Muller – generated sequence of Normal random
vectors using Sobol vector inputs is actually a low discrepancy sequence (or
even a relatively low discrepancy sequence) is not at all clear. I know this
has been looked at in 2-dimensions and there Box Muller empirically looks OK,
but I don’t know if anything has actually been proved. The low discrepancy
property tends to be fragile with respect to nonlinear transformations. So it
is possible that tdev2 actually does preserve this property , or at least does
relatively little damage to it compared to whatever you were using before.
From: Bauer, Robert [mailto:[email protected]]
Sent: Monday, May 19, 2014 11:08 AM
To: Bob Leary; [email protected]<mailto:[email protected]>
Subject: RE: [NMusers] SAEM and IMP
Bob:
Yes, these are 1 Dimensional t-distribution random deviates generated by TDEV2,
followed by scaling p sets of them with an offset vector and cholesky matrix to
provide covariance and mean. It provides the tails as you say (property a),
and is conveniently used in the Sobol process.
Your discussion on the properties of various t-distribution sample creation
techniques and their possible impact in importance sampling is interesting.
With that in mind, I just completed testing example6, which has an 8
dimensional eta space, and found no bias in the objective function evaluation,
or in parameter assessment when setting DF=2 or 4 (I have not tried others),
with or without using Sobol (RANMETHOD=3S2). Also, using the 1 dimensional
t-distribution random deviates retained the Sobol method’s stochastic noise
reduction ability in this example. When I used the multivariate t-distribution
algorithm, Sobol’s stochastic noise reduction was not as good, confirming what
you related earlier.
As you point out as well, property b (radial symmetry) may or may not be needed
depending on the posterior density. I am inclined to think that since the
posterior density is not going to be perfectly t-distributed or normal
distributed anyway, then the sampling density matching the posterior density
regarding the radial symmetry property may be of less of relevance. The more
important aspect may be to choose a sampling density that has long tails in
cases where the posterior density also has long tails, to promote the general
efficiency of the sampling density for that posterior density.
It is possible that because the sampling density is also dynamically scaled to
best fit the posterior density, this reduces any inefficiency in fitting the
posterior density that might occur from the sampling density not having radial
symmetry.
Robert J. Bauer, Ph.D.
Vice President, Pharmacometrics, R&D
ICON Development Solutions
7740 Milestone Parkway
Suite 150
Hanover, MD 21076
Tel: (215) 616-6428
Mob: (925) 286-0769
Email: [email protected]<mailto:[email protected]>
Web: http://www.iconplc.com/