From: Jurgen Bulitta <*jbulitta*>

Date: Mon, 10 Aug 2009 19:54:42 -0400

Dear Nick,

Dear Paolo,

Dear Leonid,

Dear All,

Depending on the application, I think there are the following reasons for c=

onstraining ka and ke:

1) Standard-two-stage (STS) approach: As discussed previously, calculating =

average and SD from STS estimates is likely to yield biased averages and SD=

s and this bias may have a great impact when these descriptive statistics a=

re used as means and variances in a Monte Carlo simulation. Thus, manually =

switching of individual ka & ke (if needed) makes perfect sense for STS.

2) Nonparametric population PK with nonparametric simulation: As the interi=

ndividual distribution of PK parameters can take any shape, it will make no=

difference for a non-parametric simulation based on the list of support po=

ints whether ka and ke are constrained or not. Even if a support point (wit=

h e.g. 10% probability and ka>ke) splits into two equivalent support points=

with 5% probability each (5% for ka>ke and 5% for ke>ka) this will not aff=

ect the results of a nonparametric simulation. Thus for a nonparametric sim=

ulation, constraining is not needed.

3) Nonparametric population PK with parametric simulation: Once one compute=

s averages and SDs (etc.) from the list of support points, the same issue o=

f biased averages and SDs arises as for STS, if ka and ke are not constrain=

ed. As averages and SDs are commonly calculated for nonparametric pop PK mo=

dels, constraining at the individual parameter estimate level seems most re=

asonable, as proposed by Vladimir Piotrovsky (option 2 in Paolo's email). I=

think this was the approach taken at the 2005 PAGE software comparison.

I am not aware of a method of constraining a nonparametric model at the pop=

ulation level.

4) Parametric EM-algorithms: Assume true mean ka and true mean ke are withi=

n 10% of each other and BSV has CV of 30% each. Clearly, there must be subs=

tantial overlap in the true individual estimates of ka and ke.

4.1) As one EM algorithm, the MC-PEM algorithm assumes that the interindivi=

dual distribution of individual model parameters is normal (or lognormal, o=

r can be transformed to be normal). Otherwise, the population mean and popu=

lation variance-covariance matrix obtained from the conditional means and c=

ondition variance matrices of each subject (see eq. 23 and 27 in [1] and re=

f. [2]) does not necessarily minimize the overall log-likelihood. Truncatio=

n of the interindividual distribution of PK parameters, for example by use =

of the EXIT command (option 3 in Paolo's email), seems to violate the norma=

l distribution assumption.

4.2) The conditional means calculated in the expectation-step can be biased=

due to e.g. a bi-modal distribution caused by the flip-flop (see also Serg=

e's comment at the end of the discussion in 2003). This might be a more sev=

ere concern compared to 4.1).

If such bias in the conditional means and in the conditional variance matri=

ces occurs, the computation of the population means and population variance=

-covariance matrix via eq. (21) and (22) (see Bauer et al. [1]) does no lon=

ger guarantee to minimize the objective function, I think. Constraining ka =

and ke at the population level (option 1 in Paolo's email) is unlikely to s=

olve this problem, if mean ka and mean ke are close.

If one does not constrain ka and ke at the individual parameter level, bias=

of conditional means is likely to occur, if:

a) the mean ka and mean ke are close given their interindividual variabilit=

y

b) the data are sparse and one therefore uses a wide envelope function to a=

pproximate an individual's conditional density (this problem is more likely=

to occur if one uses low values of gefficiency [please see S-ADAPT manual]=

or if one uses e.g. a t-distribution with small degrees of freedom as enve=

lope function)

Therefore, constraining at the individual parameter level via estimation of=

CL, V, and (ka-kel) (i.e. the Vladimir Piotrovsky method) and use of a ful=

l var-cov matrix (as proposed by Leonid) seems the best choice for the MC-P=

EM algorithm. I think this also applies to the SAEM and MCMC algorithm. The=

se algorithms usually have no problem with estimating a full var-cov matrix=

and do not require excessively longer runtimes to estimate a full var-cov =

matrix, as opposed to FOCE.

SUMMARY:

1) There is no need to constrain ka and ke for a nonparametric estimation w=

ith nonparametric simulation. However, as constraining does not hurt, I wou=

ld still constrain ka and ke, since interpretation of mean parameters is ea=

sier.

2) For nonparametric estimation, ka and ke should be constrained at the ind=

ividual level, if means and SDs (etc.) are to be computed from the list of =

support points (as commonly done).

3) Constraining ka and ke at the individual level and estimation of a full =

var-cov matrix seems most appropriate for parametric EM algorithms and for =

nonparametric estimation with subsequent parametric Monte Carlo simulations=

.

Hope these comments are helpful despite their length.

Nick, I would agree that in many cases, fortunately the issues mentioned ab=

ove only have a minor practical impact.

Best wishes

Juergen

References:

[1] Bauer RJ et al. AAPS Journal, 2007, 9(1): E60-E83. See bottom left and =

right column on page E64.

[2] Schumitzky A . EM algorithms and two stage methods in pharmacokinetics =

population analysis. In: D'Argenio DZ , ed. Advanced Methods of Pharmacokin=

etic and Pharmacodynamic Systems Analysis. vol. 2. Boston, MA : Kluwer Acad=

emic Publishers ; 1995 :145- 160.

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

Jurgen Bulitta, Ph.D., Senior Scientist

Ordway Research Institute, Albany, NY, USA

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

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

From: owner-nmusers

Behalf Of Nick Holford

Sent: Friday, August 07, 2009 5:03 PM

To: nmusers

Subject: Re: [NMusers] Code to avoid flip-flop kinetics

Hi,

I agree with Leonid -- especially this:

*> If population K and KA are close, then it is unclear why individual K
*

*> and KA should be consistently ordered as K < KA and not vice versa for
*

*> some subjects. If so, why not to allow the model to decide whether to
*

*> have flip flop or not?
*

I do not understand what the concern is with flip-flop kinetics and why

people want to force K<KA. If you have a good reason then various

methods are available (as reviewed again by Paulo).

I wonder if someone would like to describe why it is of interest to

force K<KA?

Nick

Leonid Gibiansky wrote:

*> Paolo,
*

*> Option 1 is very helpful. Option 2 is not attractive for the reason
*

*> that you stated, especially, extra correlation introduced by this
*

*> trick. If used, it should be used with the full OMEGA block (CL-V-KA).
*

*> I used it once but only as the last resort when nothing else worked. I
*

*> have not tried option 3 but this option artificially restricts
*

*> distributions, so I am not sure whether it is good or not even when it
*

*> is working.
*

*>
*

*> On the other hand, if population K and KA are sufficiently far apart,
*

*> you are unlikely to get individual K and KA flip-flopped.
*

*>
*

*> If population K and KA are close, then it is unclear why individual K
*

*> and KA should be consistently ordered as K < KA and not vice versa for
*

*> some subjects. If so, why not to allow the model to decide whether to
*

*> have flip flop or not?
*

*>
*

*> I would try to use option (1), and if you like the model, diagnostic
*

*> plots, etc, I would not worry about individual K and KA relation. One
*

*> of the diagnostics could be the fraction of patients with the
*

*> flip-flop. If it is small, this would justify the approach. Another
*

*> possible diagnostic is ETA_KA - ETA_K differences, If it is close to
*

*> normal you are OK, if it has two mirror-symmetric (relative to the
*

*> y-axis) peaks, then flip flop is interfering with the model.
*

*>
*

*> Thanks
*

*> Leonid
*

*>
*

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

*> Leonid Gibiansky, Ph.D.
*

*> President, QuantPharm LLC
*

*> web: www.quantpharm.com
*

*> e-mail: LGibiansky at quantpharm.com
*

*> tel: (301) 767 5566
*

*>
*

*>
*

*>
*

*>
*

*> Paolo Denti wrote:
*

*>> Dear NMUsers,
*

*>> having to deal with the flip-flop kinetics phenomenon, I had a look
*

*>> at previous posts on the NMusers list.
*

*>> I found this post particularly enlightening:
*

*>> http://www.cognigencorp.com/nonmem/nm/99aug072003.html
*

*>>
*

*>> Some code was proposed to avoid the flip-flop at population and
*

*>> individual level. Here's a not-so-brief summary.
*

*>>
*

*>> This parameterization solves the issue at population level:
*

*>> CL=THETA(1)*EXP(ETA(1))
*

*>> V=THETA(2)*EXP(ETA(2))
*

*>> TVKE=THETA(1)/THETA(2)
*

*>> TVKA=TVKE+THETA(3)
*

*>> KA=TVKA*EXP(ETA(3))
*

*>>
*

*>> However, it does not prevent the phenomenon occurring at individual
*

*>> level. Vlamidir Piotrovsky proposed the code below, which does solve
*

*>> the problem at individual level, but makes the interpretation of the
*

*>> results a bit awkward and introduces correlation among the model
*

*>> parameters. In particular, variance of ETA3 was greatly increased.
*

*>> CL=THETA(1)*EXP(ETA(1))
*

*>> V=THETA(2)*EXP(ETA(2))
*

*>> KE=CL/V
*

*>> KA=KE+THETA(3)*EXP(ETA(3))
*

*>>
*

*>> Another approach, suggested by Nick Holford, implements error
*

*>> recovery using EXIT 1. The code is reported below:
*

*>> CL=THETA(cl)*EXP(ETA(cl))
*

*>> V=THETA(v)*EXP(ETA(v))
*

*>> KA=THETA(ka)*EXP(ETA(ka))
*

*>> K=CL/V
*

*>> IF (KA.LE.K) EXIT 1 101 ; try again (PREDERR message error code 101)
*

*>>
*

*>> As far as I understand, this interrupts the computation whenever the
*

*>> flip-flop occurs and lets NONMEM restart. However, if such an error
*

*>> arises at initialization, NONMEM does not recover and the run goes no
*

*>> further. Nick probably experienced something similar, but apparently
*

*>> received no answer
*

*>> http://www.cognigencorp.com/nonmem/nm/99oct072004.html
*

*>>
*

*>> Does anyone know of a way around this drawback? Or have other code to
*

*>> deal with flip-flop kinetics?
*

*>>
*

*>> Thank you in advance,
*

*>> Paolo
*

*>>
*

*>>
*

*>> --
*

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

*>> Paolo Denti, Post-Doctoral Fellow
*

*>> Division of Clinical Pharmacology
*

*>> Department of Medicine
*

*>> University of Cape Town
*

*>>
*

*>> K45 Old Main Building
*

*>> Groote Schuur Hospital
*

*>> Observatory, Cape Town
*

*>> 7925 South Africa
*

*>> phone: +27 21 404 7719
*

*>> fax: +27 21 448 1989
*

*>> email: paolo.denti *

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

--

Nick Holford, Professor Clinical Pharmacology

Dept Pharmacology & Clinical Pharmacology

University of Auckland, 85 Park Rd, Private Bag 92019, Auckland, New Zealan=

d

n.holford

mobile: +64 21 46 23 53

http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford

Received on Mon Aug 10 2009 - 19:54:42 EDT

Date: Mon, 10 Aug 2009 19:54:42 -0400

Dear Nick,

Dear Paolo,

Dear Leonid,

Dear All,

Depending on the application, I think there are the following reasons for c=

onstraining ka and ke:

1) Standard-two-stage (STS) approach: As discussed previously, calculating =

average and SD from STS estimates is likely to yield biased averages and SD=

s and this bias may have a great impact when these descriptive statistics a=

re used as means and variances in a Monte Carlo simulation. Thus, manually =

switching of individual ka & ke (if needed) makes perfect sense for STS.

2) Nonparametric population PK with nonparametric simulation: As the interi=

ndividual distribution of PK parameters can take any shape, it will make no=

difference for a non-parametric simulation based on the list of support po=

ints whether ka and ke are constrained or not. Even if a support point (wit=

h e.g. 10% probability and ka>ke) splits into two equivalent support points=

with 5% probability each (5% for ka>ke and 5% for ke>ka) this will not aff=

ect the results of a nonparametric simulation. Thus for a nonparametric sim=

ulation, constraining is not needed.

3) Nonparametric population PK with parametric simulation: Once one compute=

s averages and SDs (etc.) from the list of support points, the same issue o=

f biased averages and SDs arises as for STS, if ka and ke are not constrain=

ed. As averages and SDs are commonly calculated for nonparametric pop PK mo=

dels, constraining at the individual parameter estimate level seems most re=

asonable, as proposed by Vladimir Piotrovsky (option 2 in Paolo's email). I=

think this was the approach taken at the 2005 PAGE software comparison.

I am not aware of a method of constraining a nonparametric model at the pop=

ulation level.

4) Parametric EM-algorithms: Assume true mean ka and true mean ke are withi=

n 10% of each other and BSV has CV of 30% each. Clearly, there must be subs=

tantial overlap in the true individual estimates of ka and ke.

4.1) As one EM algorithm, the MC-PEM algorithm assumes that the interindivi=

dual distribution of individual model parameters is normal (or lognormal, o=

r can be transformed to be normal). Otherwise, the population mean and popu=

lation variance-covariance matrix obtained from the conditional means and c=

ondition variance matrices of each subject (see eq. 23 and 27 in [1] and re=

f. [2]) does not necessarily minimize the overall log-likelihood. Truncatio=

n of the interindividual distribution of PK parameters, for example by use =

of the EXIT command (option 3 in Paolo's email), seems to violate the norma=

l distribution assumption.

4.2) The conditional means calculated in the expectation-step can be biased=

due to e.g. a bi-modal distribution caused by the flip-flop (see also Serg=

e's comment at the end of the discussion in 2003). This might be a more sev=

ere concern compared to 4.1).

If such bias in the conditional means and in the conditional variance matri=

ces occurs, the computation of the population means and population variance=

-covariance matrix via eq. (21) and (22) (see Bauer et al. [1]) does no lon=

ger guarantee to minimize the objective function, I think. Constraining ka =

and ke at the population level (option 1 in Paolo's email) is unlikely to s=

olve this problem, if mean ka and mean ke are close.

If one does not constrain ka and ke at the individual parameter level, bias=

of conditional means is likely to occur, if:

a) the mean ka and mean ke are close given their interindividual variabilit=

y

b) the data are sparse and one therefore uses a wide envelope function to a=

pproximate an individual's conditional density (this problem is more likely=

to occur if one uses low values of gefficiency [please see S-ADAPT manual]=

or if one uses e.g. a t-distribution with small degrees of freedom as enve=

lope function)

Therefore, constraining at the individual parameter level via estimation of=

CL, V, and (ka-kel) (i.e. the Vladimir Piotrovsky method) and use of a ful=

l var-cov matrix (as proposed by Leonid) seems the best choice for the MC-P=

EM algorithm. I think this also applies to the SAEM and MCMC algorithm. The=

se algorithms usually have no problem with estimating a full var-cov matrix=

and do not require excessively longer runtimes to estimate a full var-cov =

matrix, as opposed to FOCE.

SUMMARY:

1) There is no need to constrain ka and ke for a nonparametric estimation w=

ith nonparametric simulation. However, as constraining does not hurt, I wou=

ld still constrain ka and ke, since interpretation of mean parameters is ea=

sier.

2) For nonparametric estimation, ka and ke should be constrained at the ind=

ividual level, if means and SDs (etc.) are to be computed from the list of =

support points (as commonly done).

3) Constraining ka and ke at the individual level and estimation of a full =

var-cov matrix seems most appropriate for parametric EM algorithms and for =

nonparametric estimation with subsequent parametric Monte Carlo simulations=

.

Hope these comments are helpful despite their length.

Nick, I would agree that in many cases, fortunately the issues mentioned ab=

ove only have a minor practical impact.

Best wishes

Juergen

References:

[1] Bauer RJ et al. AAPS Journal, 2007, 9(1): E60-E83. See bottom left and =

right column on page E64.

[2] Schumitzky A . EM algorithms and two stage methods in pharmacokinetics =

population analysis. In: D'Argenio DZ , ed. Advanced Methods of Pharmacokin=

etic and Pharmacodynamic Systems Analysis. vol. 2. Boston, MA : Kluwer Acad=

emic Publishers ; 1995 :145- 160.

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

Jurgen Bulitta, Ph.D., Senior Scientist

Ordway Research Institute, Albany, NY, USA

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

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

From: owner-nmusers

Behalf Of Nick Holford

Sent: Friday, August 07, 2009 5:03 PM

To: nmusers

Subject: Re: [NMusers] Code to avoid flip-flop kinetics

Hi,

I agree with Leonid -- especially this:

I do not understand what the concern is with flip-flop kinetics and why

people want to force K<KA. If you have a good reason then various

methods are available (as reviewed again by Paulo).

I wonder if someone would like to describe why it is of interest to

force K<KA?

Nick

Leonid Gibiansky wrote:

--

Nick Holford, Professor Clinical Pharmacology

Dept Pharmacology & Clinical Pharmacology

University of Auckland, 85 Park Rd, Private Bag 92019, Auckland, New Zealan=

d

n.holford

mobile: +64 21 46 23 53

http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford

Received on Mon Aug 10 2009 - 19:54:42 EDT