From: Stephen Duffull <*stephen.duffull*>

Date: Thu, 14 Jun 2012 18:51:00 +0000

It turns out that there are many symmetric matrices where off-diagonals are=

between -1 and 1 that have negative determinants. Just to add my 2c there =

are also non-positive definite symmetric matrices that have positive determ=

inants. Some are obvious and some are not. So if you were to parameterise=

using SD and corr you would need some care. I am pleased I don't have to c=

ode search algorithms that require positive definite matrices as a search s=

tep :-)

Sampling algorithms don't have this issue as there are good priors that gua=

rantee "legal" matrices.

Steve

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

*> From: owner-nmusers *

On

*> Behalf Of Bob Leary
*

*> Sent: Friday, 15 June 2012 5:11 a.m.
*

*> To: Bauer, Robert; Matt Hutmacher; 'Nick Holford'; nmusers *

*> Subject: RE: [NMusers] omega matrix - friendly suggestion
*

*>
*

*> Easy countereaxample to Nick's hypothesis that a covariance matrix with a=
*

n

*> underlying correlation matrix with all elements of magnitude <=1 is pos=
*

itive

*> definite
*

*>
*

*> Let all sd's be 1, and let the correlation matrix be c =
*

*> 1.0000 0.9900 -0.9900
*

*> 0.9900 1.0000 0.9900
*

*> -0.9900 0.9900 1.0000
*

*>
*

*> This is in fact an impossible correlation matrix, since it is not poitive
*

*> definite:
*

*>
*

*> > det(c)
*

*> ans =
*

*> -3.8809
*

*>
*

*> (THe determinant of a positive definite matrix must be positive)
*

*>
*

*>
*

*> The two most common ways to parameterize a symmetric pos def matrix that
*

*> guarantee positive definiteness are the cholesly factor method, and the
*

*> matrix exponential of arbitrary symmetric matrix . Neither one can be use=
*

d to

*> force an arbitrary element of a full block matrix to a given value while
*

*> preserving postive definiteness.
*

*> _______________________________________
*

*> From: owner-nmusers *

lf Of

*> Bauer, Robert [Robert.Bauer *

*> Sent: Wednesday, June 13, 2012 4:05 PM
*

*> To: Matt Hutmacher; 'Nick Holford'; nmusers *

*> Subject: RE: [NMusers] omega matrix - friendly suggestion
*

*>
*

*> Thank you all for your suggestions. The more flexible OMEGA element entr=
*

ies

*> will not be easy to implement in NONMEM's structure, but I will look into=
*

it.

*> A couple of things need to be noted.
*

*>
*

*> 1) Regarding Matt's suggestions, to enter initial Omega elements in Chole=
*

sky

*> format in NM 7.2, $OMEGA CHOLESKY may be used.
*

*> 2) The traditional methods (FO/FOCE/Laplace) actually vary the omega elem=
*

ents

*> in the Cholesky format for estimation, than transform back in the varianc=
*

e

*> format after each iteration. This assures that as NONMEM arbitrarily sea=
*

rches

*> throughout the Cholesky space (that is, as each cholesky element is
*

*> arbitrarily varied in search of the minimum OFV), the variance omega is
*

*> guaranteed to be positive definite. Thus, initially fixing an off-diagon=
*

al

*> variance element to a non-zero value within a block of variances cannot b=
*

e

*> honored by the search algorithm past the first iteration. Additionally, =
*

to my

*> best understanding, arbitrarily fixing an off-diagonal variance element t=
*

o 0

*> does not correspond to fixing a particular Cholesky off-diagonal element =
*

to 0,

*> except for banded matrix patterns, which is already allowed in NONMEM.
*

*> 3) EM algorithms such as Imp and SAEM update the Omega variance in the
*

*> variance domain, and also, allow a simple constraint filter to be super-
*

*> imposed on the resulting matrix after each iteration. In NONMEM, this ca=
*

n be

*> done using constraint.f90, wherein a person may code in fixing any omega
*

*> element to 0, or any other value, over-riding the EM algorithm's update. =
*

This

*> is okay in EM, as the omega elements are not primary search parameters, b=
*

ut

*> are constructed as a consequence of the expectation step, which provide
*

*> conditional means and variances from which the omega elements are evaluat=
*

ed.

*> User beware, should this constraint cause lack of positive definiteness. =
*

The

*> constraint.f90 routine is not used for the traditional update methods, as=
*

this

*> would upset their gradient search pattern by providing a dis-continuous
*

*> intrusion.
*

*>
*

*>
*

*>
*

*> 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: Robert.Bauer *

*> Web: www.iconplc.com
*

*>
*

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

*> From: owner-nmusers *

On

*> Behalf Of Matt Hutmacher
*

*> Sent: Wednesday, June 13, 2012 2:35 PM
*

*> To: 'Nick Holford'; nmusers *

*> Subject: RE: [NMusers] omega matrix - friendly suggestion
*

*>
*

*> Dear Paolo, Nick, all,
*

*>
*

*> A way to ensure a nonnegative definite (NND) variance-covariance matrix i=
*

s to

*> use a Cholesky decomposition of the matrix (Ref 1). I use this often and=
*

it

*> greatly stabilizes estimation such that I can often get better convergenc=
*

e and

*> improved $COV step success for larger OMEGA matrix. It also is great if =
*

you

*> sample from the $COV for simulation as an NND $OMEGA will be achieved wit=
*

h

*> each sampling. I currently implement this by hand such that all $OMEGA
*

*> becomes a Diagonal matrix with 1 FIXED as all of its entries.
*

*> Note that I have not evaluated whether NONMEM computes the appropriate WR=
*

ES or

*> CWRES for this parameterization (as I have calculated these myself - Bob =
*

may

*> be able to comment on whether this is a concern - ie putting the variance=
*

s

*> into linear combinations of thetas affecting the residuals). It is not en=
*

ough

*> to just model variances and correlations; such a parameterization will no=
*

t

*> guarantee NND OMEGA matrices.
*

*>
*

*> With respect to the original question which was to constrain certain
*

*> covarainces to 0, one might (I have not tried this) be able to do this by
*

*> writing the covariance matrix as wished, then solving for the Cholesky
*

*> components, and modeling with the derived Cholesky parameterization. Upo=
*

n

*> completion, one can then compute the $OMEGA matrix from the Cholesky
*

*> parameterization. I do not know how easy this would be to generalize in =
*

a

*> software application, but if possible, this could perhaps be done behind =
*

the

*> scenes as an option (may make standard errors a bit tricky) and would pro=
*

vide

*> a lot of flexibility for analysis (as opposed to hand coding discussed be=
*

low).

*>
*

*> There is also a modified Cholesky that can be used, which has the nice fe=
*

ature

*> that one can simply, conveniently remove all the variance components rela=
*

ted

*> to a certain random effect (Ref 2). This was recently used in the attach=
*

ed

*> reference to aid in selection of the $OMEGA structure with the fixed effe=
*

cts.

*>
*

*> Best, Matt
*

*>
*

*>
*

*> Ref.1: Pinhiero, J. and Bates, D. (1996). Unconstrained parameterizations=
*

for

*> variance-covariance matrices. Statistics and Computing 6, 289-286.
*

*>
*

*> Ref.2: Bondell HD, Krishna A, Ghosh (2010). Joint Variable Selection for =
*

Fixed

*> and Random Effects in Linear Mixed-Effects Model. Biometrics 66, 1069-10=
*

77.

*>
*

*>
*

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

*> From: owner-nmusers *

On

*> Behalf Of Nick Holford
*

*> Sent: Wednesday, June 13, 2012 1:20 PM
*

*> To: nmusers *

*> Subject: Re: [NMusers] omega matrix - friendly suggestion
*

*>
*

*> Paolo,
*

*>
*

*> You may remember this discussion about the history of the correlation fea=
*

ture

*> for OMEGA (and SIGMA) blocks.
*

*> http://www.cognigencorp.com/nonmem/current/2009-July/1804.html
*

*>
*

*> IMHO the only sensible way for a human to write a covariance matrix is us=
*

ing

*> the SD and correlation style introduced in NONMEM 7. If I understand thin=
*

gs

*> correctly (this kind of math is not my strong point) it should be impossi=
*

ble

*> to write non-positive definitive matrices provided the correlations are
*

*> -1<=0<=+1 (although NONMEM might complain at -1 or +1).
*

*>
*

*> Can anyone provide any good reason to prefer the variance-covariance
*

*> parameterization over the SD-correlation parametization?
*

*>
*

*> The other issue brought up by Pavel is a different one -- how to specify =
*

more

*> simply a zero correlation|covariance between random effects.
*

*>
*

*> Monolix has a 0|1 matrix that allows a user to fix any covariance matrix
*

*> element to zero so presumably the math exists somewhere to allow NONMEM h=
*

ow to

*> understand a matrix to be recognized in this format.
*

*>
*

*> Bob (Bauer) -- do you know how to do this? Perhaps you could talk to Ma=
*

rc

*> Lavielle and consider adding something similar to NM-TRAN. I'd prefer the
*

*> style shown by Pavel to having to write a separate matrix of
*

*> 0|1 elements as used in the Monolix GUI.
*

*>
*

*> Best wishes,
*

*>
*

*> Nick
*

*>
*

*>
*

*>
*

*> On 13/06/2012 5:16 p.m., Paolo Denti wrote:
*

*> > Dear Pavel and Ken,
*

*> > I also share the same pain, especially when coding correlations for
*

*> > between-occasion variability ETAs. This implies reorganizing them in
*

*> > blocks and renumbering everything. A real chore, and a very
*

*> > error-prone process.
*

*> >
*

*> > Maybe one can think of using the OMEGA in the correlation format,
*

*> > which should make it easier to write "legitimate" OMEGA matrices. Or
*

*> > NONMEM can check the positive definiteness of the initial estimate and
*

*> > complain if necessary (I believe it already does so).
*

*> >
*

*> > So, dear NONMEM developers, please count a +1 in the survey for this
*

*> > feature. :)
*

*> >
*

*> > Thank you and ciao,
*

*> > Paolo
*

*> >
*

*> > On 2012/06/12 17:08, Ken Kowalski wrote:
*

*> >>
*

*> >> Dear Pavel,
*

*> >>
*

*> >> I certainly feel your pain but you have to be careful how you fix
*

*> >> certain elements in Omega to ensure that you have a valid positive
*

*> >> definite covariance matrix. The starting values in your $OMEGA block
*

*> >> do not give rise to a valid covariance matrix. Note in particular
*

*> >> that the covariance between ETA3 and ETA4 is too large relative to
*

*> >> the variances for ETA3 and ETA4 such that the correlation is greater
*

*> >> than 1.0, i.e.,
*

*> >>
*

*> >> Corr(ETA3,ETA4) = 0.03/[SQRT(0.0166)*SQRT(0.0166)]=1.807 > 1
*

*> >>
*

*> >> You also have the same problem for Corr(ETA1,ETA3) > 1.
*

*> >>
*

*> >> Ken
*

*> >>
*

*> >> Kenneth G. Kowalski
*

*> >>
*

*> >> President & CEO
*

*> >>
*

*> >> A2PG - Ann Arbor Pharmacometrics Group, Inc.
*

*> >>
*

*> >> 110 Miller Ave., Garden Suite
*

*> >>
*

*> >> Ann Arbor, MI 48104
*

*> >>
*

*> >> Work: 734-274-8255
*

*> >>
*

*> >> Cell: 248-207-5082
*

*> >>
*

*> >> Fax: 734-913-0230
*

*> >>
*

*> >> ken.kowalski *

*> >>
*

*> >> www.a2pg.com <http://www.a2pg.com>
*

*> >>
*

*> >> *From:*owner-nmusers *

*> >> [mailto:owner-nmusers *

*> >> *nonmem *

*> >> *Sent:* Tuesday, June 12, 2012 10:24 AM
*

*> >> *To:* nmusers *

*> >> *Subject:* [NMusers] omega matrix - friendly suggestion
*

*> >>
*

*> >> Hello NONMEM Community,
*

*> >>
*

*> >> Sometimes it takes more time to choose the best omega matrix than to
*

*> >> develop a PD model. Selecting the omega matric is a tedious, time
*

*> >> consuming and less creative part of the model development. I hope
*

*> >> you feel my pain. Will it be helpful to rewrite the NONMEM software
*

*> >> so that any element of the omega matrix can be fixed to any value?
*

*> >> It may look like this:
*

*> >>
*

*> >> $OMEGA BLOCK (4) 3.60E-02 FIX
*

*> >>
*

*> >> 0.01 3.23E-02
*

*> >> 0.03 0 FIX 1.66E-02
*

*> >> 0.01 0 FIX 0.03 FIX 1.66E-02
*

*> >>
*

*> >> This change can make many NONMEM users happy.
*

*> >>
*

*> >> Thanks!
*

*> >>
*

*> >> Pavel
*

*> >>
*

*> >
*

*> > --
*

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

*> > Paolo Denti, PhD
*

*> > Junior Lecturer
*

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

*>
*

*> First World Conference on Pharmacometrics, 5-7 September 2012 Seoul, Kore=
*

a

*> http://www.go-wcop.org
*

*>
*

*> Dept Pharmacology& Clinical Pharmacology, Bldg 505 Room 202D 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
*

*>
*

*>
*

*>
*

*>
*

*> ICON plc made the following annotations.
*

*> -------------------------------------------------------------------------=
*

-----

*> This e-mail transmission may contain confidential or legally privileged
*

*> information that is intended only for the individual or entity named in t=
*

he e-

*> mail address. If you are not the intended recipient, you are hereby notif=
*

ied

*> that any disclosure, copying, distribution, or reliance upon the contents=
*

of

*> this e-mail is strictly prohibited. If you have received this e-mail
*

*> transmission in error, please reply to the sender, so that ICON plc can
*

*> arrange for proper delivery, and then please delete the message.
*

*> Thank You,
*

*> ICON plc
*

*> South County Business Park
*

*> Leopardstown
*

*> Dublin 18
*

*> Ireland
*

*> Registered number:
*

*> 145835_________________________________________________________________
*

*> NOTICE: The information contained in this electronic mail message is inte=
*

nded

*> only for the personal and confidential use of the designated recipient(s)
*

*> named above. This message may be an attorney-client communication, may be
*

*> protected by the work product doctrine, and may be subject to a protectiv=
*

e

*> order. As such, this message is privileged and confidential. If the reade=
*

r of

*> this message is not the intended recipient or an agent responsible for
*

*> delivering it to the intended recipient, you are hereby notified that you=
*

have

*> received this message in error and that any review, dissemination,
*

*> distribution, or copying of this message is strictly prohibited. If you h=
*

ave

*> received this communication in error, please notify us immediately by
*

*> telephone and e-mail and destroy any and all copies of this message in yo=
*

ur

*> possession (whether hard copies or electronically stored copies). Thank y=
*

ou.

*> *

Received on Thu Jun 14 2012 - 14:51:00 EDT

Date: Thu, 14 Jun 2012 18:51:00 +0000

It turns out that there are many symmetric matrices where off-diagonals are=

between -1 and 1 that have negative determinants. Just to add my 2c there =

are also non-positive definite symmetric matrices that have positive determ=

inants. Some are obvious and some are not. So if you were to parameterise=

using SD and corr you would need some care. I am pleased I don't have to c=

ode search algorithms that require positive definite matrices as a search s=

tep :-)

Sampling algorithms don't have this issue as there are good priors that gua=

rantee "legal" matrices.

Steve

On

n

itive

d to

lf Of

ies

it.

sky

ents

e

rches

al

e

to my

o 0

to 0,

n be

This

ut

ed.

The

this

On

s to

it

e and

you

h

ES or

may

s

ough

t

n

a

the

vide

low).

ature

ted

ed

cts.

for

Fixed

77.

On

ture

ing

gs

ble

more

ow to

rc

a

of

-----

he e-

ied

of

nded

e

r of

have

ave

ur

ou.

Received on Thu Jun 14 2012 - 14:51:00 EDT