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
omega matrix - friendly suggestion
10 messages
9 people
Latest: Jun 15, 2012
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
<mailto:[email protected]> [email protected]
http://www.a2pg.com www.a2pg.com
Quoted reply history
From: [email protected] [mailto:[email protected]] On
Behalf Of [email protected]
Sent: Tuesday, June 12, 2012 10:24 AM
To: [email protected]
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
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
Quoted reply history
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
>
> [email protected] <mailto:[email protected]>
>
> www.a2pg.com http://www.a2pg.com
>
> *From:* [email protected] [ mailto: [email protected] ] *On Behalf Of * [email protected]
>
> *Sent:* Tuesday, June 12, 2012 10:24 AM
> *To:* [email protected]
> *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: [email protected]
Paolo,
You may remember this discussion about the history of the correlation feature 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 using the SD and correlation style introduced in NONMEM 7. If I understand things correctly (this kind of math is not my strong point) it should be impossible 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 how to understand a matrix to be recognized in this format.
Bob (Bauer) -- do you know how to do this? Perhaps you could talk to Marc 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
Quoted reply history
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
> >
> > [email protected] <mailto:[email protected]>
> >
> > www.a2pg.com http://www.a2pg.com
> >
> > *From:* [email protected] [ mailto: [email protected] ] *On Behalf Of * [email protected]
> >
> > *Sent:* Tuesday, June 12, 2012 10:24 AM
> > *To:* [email protected]
> > *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:[email protected]
> ------------------------------------------------
--
Nick Holford, Professor Clinical Pharmacology
First World Conference on Pharmacometrics, 5-7 September 2012
Seoul, Korea 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: [email protected]
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
Dear Paolo, Nick, all,
A way to ensure a nonnegative definite (NND) variance-covariance matrix is
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
convergence 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 with 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 WRES
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
variances into linear combinations of thetas affecting the residuals). It is
not enough to just model variances and correlations; such a parameterization
will not 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. Upon
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
provide a lot of flexibility for analysis (as opposed to hand coding
discussed below).
There is also a modified Cholesky that can be used, which has the nice
feature that one can simply, conveniently remove all the variance components
related to a certain random effect (Ref 2). This was recently used in the
attached reference to aid in selection of the $OMEGA structure with the
fixed effects.
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-1077.
Quoted reply history
-----Original Message-----
From: [email protected] [mailto:[email protected]] On
Behalf Of Nick Holford
Sent: Wednesday, June 13, 2012 1:20 PM
To: [email protected]
Subject: Re: [NMusers] omega matrix - friendly suggestion
Paolo,
You may remember this discussion about the history of the correlation
feature 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 using
the SD and correlation style introduced in NONMEM 7. If I understand things
correctly (this kind of math is not my strong point) it should be impossible
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 how
to understand a matrix to be recognized in this format.
Bob (Bauer) -- do you know how to do this? Perhaps you could talk to Marc
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
>>
>> [email protected] <mailto:[email protected]>
>>
>> www.a2pg.com http://www.a2pg.com
>>
>> *From:*[email protected]
>> [mailto:[email protected]] *On Behalf Of
>> *[email protected]
>> *Sent:* Tuesday, June 12, 2012 10:24 AM
>> *To:* [email protected]
>> *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:[email protected]
> ------------------------------------------------
--
Nick Holford, Professor Clinical Pharmacology
First World Conference on Pharmacometrics, 5-7 September 2012 Seoul, Korea
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: [email protected]
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
Thank you all for your suggestions. The more flexible OMEGA element entries
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 Cholesky
format in NM 7.2,
$OMEGA CHOLESKY
may be used.
2) The traditional methods (FO/FOCE/Laplace) actually vary the omega elements
in the Cholesky format for estimation, than transform back in the variance
format after each iteration. This assures that as NONMEM arbitrarily searches
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-diagonal variance element to
a non-zero value within a block of variances cannot be honored by the search
algorithm past the first iteration. Additionally, to my best understanding,
arbitrarily fixing an off-diagonal variance element to 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 can 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, but are constructed as a
consequence of the expectation step, which provide conditional means and
variances from which the omega elements are evaluated. 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: [email protected]
Web: www.iconplc.com
Quoted reply history
-----Original Message-----
From: [email protected] [mailto:[email protected]] On
Behalf Of Matt Hutmacher
Sent: Wednesday, June 13, 2012 2:35 PM
To: 'Nick Holford'; [email protected]
Subject: RE: [NMusers] omega matrix - friendly suggestion
Dear Paolo, Nick, all,
A way to ensure a nonnegative definite (NND) variance-covariance matrix is
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
convergence 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 with 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 WRES
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
variances into linear combinations of thetas affecting the residuals). It is
not enough to just model variances and correlations; such a parameterization
will not 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. Upon
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
provide a lot of flexibility for analysis (as opposed to hand coding
discussed below).
There is also a modified Cholesky that can be used, which has the nice
feature that one can simply, conveniently remove all the variance components
related to a certain random effect (Ref 2). This was recently used in the
attached reference to aid in selection of the $OMEGA structure with the
fixed effects.
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-1077.
-----Original Message-----
From: [email protected] [mailto:[email protected]] On
Behalf Of Nick Holford
Sent: Wednesday, June 13, 2012 1:20 PM
To: [email protected]
Subject: Re: [NMusers] omega matrix - friendly suggestion
Paolo,
You may remember this discussion about the history of the correlation
feature 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 using
the SD and correlation style introduced in NONMEM 7. If I understand things
correctly (this kind of math is not my strong point) it should be impossible
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 how
to understand a matrix to be recognized in this format.
Bob (Bauer) -- do you know how to do this? Perhaps you could talk to Marc
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
>>
>> [email protected] <mailto:[email protected]>
>>
>> www.a2pg.com http://www.a2pg.com
>>
>> *From:*[email protected]
>> [mailto:[email protected]] *On Behalf Of
>> *[email protected]
>> *Sent:* Tuesday, June 12, 2012 10:24 AM
>> *To:* [email protected]
>> *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:[email protected]
> ------------------------------------------------
--
Nick Holford, Professor Clinical Pharmacology
First World Conference on Pharmacometrics, 5-7 September 2012 Seoul, Korea
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: [email protected]
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
Nick,
There is also no guarantee that you can get a positive-definite matrix by
specifying variances>0
and "correlations" between -1 and 1.
Say A is a 3x3 matrix with all 1's in diagonal and all -0.6 in off-diagonal.
It is not positive definite. det(A)=-0.512.
Siva
Quoted reply history
________________________________________
From: [email protected] [[email protected]] On Behalf Of
Matt Hutmacher [[email protected]]
Sent: Wednesday, June 13, 2012 2:35 PM
To: 'Nick Holford'; [email protected]
Subject: RE: [NMusers] omega matrix - friendly suggestion
Dear Paolo, Nick, all,
A way to ensure a nonnegative definite (NND) variance-covariance matrix is
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
convergence 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 with 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 WRES
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
variances into linear combinations of thetas affecting the residuals). It is
not enough to just model variances and correlations; such a parameterization
will not 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. Upon
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
provide a lot of flexibility for analysis (as opposed to hand coding
discussed below).
There is also a modified Cholesky that can be used, which has the nice
feature that one can simply, conveniently remove all the variance components
related to a certain random effect (Ref 2). This was recently used in the
attached reference to aid in selection of the $OMEGA structure with the
fixed effects.
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-1077.
-----Original Message-----
From: [email protected] [mailto:[email protected]] On
Behalf Of Nick Holford
Sent: Wednesday, June 13, 2012 1:20 PM
To: [email protected]
Subject: Re: [NMusers] omega matrix - friendly suggestion
Paolo,
You may remember this discussion about the history of the correlation
feature 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 using
the SD and correlation style introduced in NONMEM 7. If I understand things
correctly (this kind of math is not my strong point) it should be impossible
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 how
to understand a matrix to be recognized in this format.
Bob (Bauer) -- do you know how to do this? Perhaps you could talk to Marc
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
>>
>> [email protected] <mailto:[email protected]>
>>
>> www.a2pg.com http://www.a2pg.com
>>
>> *From:*[email protected]
>> [mailto:[email protected]] *On Behalf Of
>> *[email protected]
>> *Sent:* Tuesday, June 12, 2012 10:24 AM
>> *To:* [email protected]
>> *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:[email protected]
> ------------------------------------------------
--
Nick Holford, Professor Clinical Pharmacology
First World Conference on Pharmacometrics, 5-7 September 2012 Seoul, Korea
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: [email protected]
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
Easy countereaxample to Nick's hypothesis that a covariance matrix with an
underlying correlation matrix with all elements of magnitude <=1
is positive 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 used to force an
arbitrary element of a full block matrix to a given value while preserving
postive definiteness.
_______________________________________
Quoted reply history
From: [email protected] [[email protected]] On Behalf Of
Bauer, Robert [[email protected]]
Sent: Wednesday, June 13, 2012 4:05 PM
To: Matt Hutmacher; 'Nick Holford'; [email protected]
Subject: RE: [NMusers] omega matrix - friendly suggestion
Thank you all for your suggestions. The more flexible OMEGA element entries
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 Cholesky
format in NM 7.2,
$OMEGA CHOLESKY
may be used.
2) The traditional methods (FO/FOCE/Laplace) actually vary the omega elements
in the Cholesky format for estimation, than transform back in the variance
format after each iteration. This assures that as NONMEM arbitrarily searches
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-diagonal variance element to
a non-zero value within a block of variances cannot be honored by the search
algorithm past the first iteration. Additionally, to my best understanding,
arbitrarily fixing an off-diagonal variance element to 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 can 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, but are constructed as a
consequence of the expectation step, which provide conditional means and
variances from which the omega elements are evaluated. 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: [email protected]
Web: www.iconplc.com
-----Original Message-----
From: [email protected] [mailto:[email protected]] On
Behalf Of Matt Hutmacher
Sent: Wednesday, June 13, 2012 2:35 PM
To: 'Nick Holford'; [email protected]
Subject: RE: [NMusers] omega matrix - friendly suggestion
Dear Paolo, Nick, all,
A way to ensure a nonnegative definite (NND) variance-covariance matrix is
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
convergence 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 with 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 WRES
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
variances into linear combinations of thetas affecting the residuals). It is
not enough to just model variances and correlations; such a parameterization
will not 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. Upon
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
provide a lot of flexibility for analysis (as opposed to hand coding
discussed below).
There is also a modified Cholesky that can be used, which has the nice
feature that one can simply, conveniently remove all the variance components
related to a certain random effect (Ref 2). This was recently used in the
attached reference to aid in selection of the $OMEGA structure with the
fixed effects.
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-1077.
-----Original Message-----
From: [email protected] [mailto:[email protected]] On
Behalf Of Nick Holford
Sent: Wednesday, June 13, 2012 1:20 PM
To: [email protected]
Subject: Re: [NMusers] omega matrix - friendly suggestion
Paolo,
You may remember this discussion about the history of the correlation
feature 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 using
the SD and correlation style introduced in NONMEM 7. If I understand things
correctly (this kind of math is not my strong point) it should be impossible
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 how
to understand a matrix to be recognized in this format.
Bob (Bauer) -- do you know how to do this? Perhaps you could talk to Marc
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
>>
>> [email protected] <mailto:[email protected]>
>>
>> www.a2pg.com http://www.a2pg.com
>>
>> *From:*[email protected]
>> [mailto:[email protected]] *On Behalf Of
>> *[email protected]
>> *Sent:* Tuesday, June 12, 2012 10:24 AM
>> *To:* [email protected]
>> *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:[email protected]
> ------------------------------------------------
--
Nick Holford, Professor Clinical Pharmacology
First World Conference on Pharmacometrics, 5-7 September 2012 Seoul, Korea
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: [email protected]
http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
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 determinants.
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 code search
algorithms that require positive definite matrices as a search step :-)
Sampling algorithms don't have this issue as there are good priors that
guarantee "legal" matrices.
Steve
Quoted reply history
> -----Original Message-----
> From: [email protected] [mailto:[email protected]] On
> Behalf Of Bob Leary
> Sent: Friday, 15 June 2012 5:11 a.m.
> To: Bauer, Robert; Matt Hutmacher; 'Nick Holford'; [email protected]
> Subject: RE: [NMusers] omega matrix - friendly suggestion
>
> Easy countereaxample to Nick's hypothesis that a covariance matrix with an
> underlying correlation matrix with all elements of magnitude <=1 is positive
> 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 used to
> force an arbitrary element of a full block matrix to a given value while
> preserving postive definiteness.
> _______________________________________
> From: [email protected] [[email protected]] On Behalf Of
> Bauer, Robert [[email protected]]
> Sent: Wednesday, June 13, 2012 4:05 PM
> To: Matt Hutmacher; 'Nick Holford'; [email protected]
> Subject: RE: [NMusers] omega matrix - friendly suggestion
>
> Thank you all for your suggestions. The more flexible OMEGA element entries
> 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 Cholesky
> format in NM 7.2, $OMEGA CHOLESKY may be used.
> 2) The traditional methods (FO/FOCE/Laplace) actually vary the omega elements
> in the Cholesky format for estimation, than transform back in the variance
> format after each iteration. This assures that as NONMEM arbitrarily searches
> 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-diagonal
> variance element to a non-zero value within a block of variances cannot be
> honored by the search algorithm past the first iteration. Additionally, to my
> best understanding, arbitrarily fixing an off-diagonal variance element to 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 can 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, but
> are constructed as a consequence of the expectation step, which provide
> conditional means and variances from which the omega elements are evaluated.
> 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: [email protected]
> Web: www.iconplc.com
>
> -----Original Message-----
> From: [email protected] [mailto:[email protected]] On
> Behalf Of Matt Hutmacher
> Sent: Wednesday, June 13, 2012 2:35 PM
> To: 'Nick Holford'; [email protected]
> Subject: RE: [NMusers] omega matrix - friendly suggestion
>
> Dear Paolo, Nick, all,
>
> A way to ensure a nonnegative definite (NND) variance-covariance matrix is 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 convergence 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 with
> 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 WRES 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 variances
> into linear combinations of thetas affecting the residuals). It is not enough
> to just model variances and correlations; such a parameterization will not
> 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. Upon
> 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 provide
> a lot of flexibility for analysis (as opposed to hand coding discussed below).
>
> There is also a modified Cholesky that can be used, which has the nice feature
> that one can simply, conveniently remove all the variance components related
> to a certain random effect (Ref 2). This was recently used in the attached
> reference to aid in selection of the $OMEGA structure with the fixed effects.
>
> 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-1077.
>
>
> -----Original Message-----
> From: [email protected] [mailto:[email protected]] On
> Behalf Of Nick Holford
> Sent: Wednesday, June 13, 2012 1:20 PM
> To: [email protected]
> Subject: Re: [NMusers] omega matrix - friendly suggestion
>
> Paolo,
>
> You may remember this discussion about the history of the correlation feature
> 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 using
> the SD and correlation style introduced in NONMEM 7. If I understand things
> correctly (this kind of math is not my strong point) it should be impossible
> 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 how to
> understand a matrix to be recognized in this format.
>
> Bob (Bauer) -- do you know how to do this? Perhaps you could talk to Marc
> 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
> >>
> >> [email protected] <mailto:[email protected]>
> >>
> >> www.a2pg.com http://www.a2pg.com
> >>
> >> *From:*[email protected]
> >> [mailto:[email protected]] *On Behalf Of
> >> *[email protected]
> >> *Sent:* Tuesday, June 12, 2012 10:24 AM
> >> *To:* [email protected]
> >> *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:[email protected]
> > ------------------------------------------------
>
> --
> Nick Holford, Professor Clinical Pharmacology
>
> First World Conference on Pharmacometrics, 5-7 September 2012 Seoul, Korea
> 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: [email protected]
> http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
>
>
>
>
>
In terms of determinants, the necessary and sufficient condtion for positive
definiteness is that all
determinants of the leading principal submatrices (which consist of the
elements of the leading 1 by1 , 2 by 2, 3 by 3,etc square submatrices starting
with the (1,1) element.) Another equivalent condition is that all eigenvalues
of the matrix are positive.
Quoted reply history
________________________________________
From: Stephen Duffull [[email protected]]
Sent: Thursday, June 14, 2012 1:51 PM
To: Bob Leary; Bauer, Robert; Matt Hutmacher; 'Nick Holford';
[email protected]
Subject: RE: [NMusers] omega matrix - friendly suggestion
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 determinants.
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 code search
algorithms that require positive definite matrices as a search step :-)
Sampling algorithms don't have this issue as there are good priors that
guarantee "legal" matrices.
Steve
> -----Original Message-----
> From: [email protected] [mailto:[email protected]] On
> Behalf Of Bob Leary
> Sent: Friday, 15 June 2012 5:11 a.m.
> To: Bauer, Robert; Matt Hutmacher; 'Nick Holford'; [email protected]
> Subject: RE: [NMusers] omega matrix - friendly suggestion
>
> Easy countereaxample to Nick's hypothesis that a covariance matrix with an
> underlying correlation matrix with all elements of magnitude <=1 is positive
> 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 used to
> force an arbitrary element of a full block matrix to a given value while
> preserving postive definiteness.
> _______________________________________
> From: [email protected] [[email protected]] On Behalf Of
> Bauer, Robert [[email protected]]
> Sent: Wednesday, June 13, 2012 4:05 PM
> To: Matt Hutmacher; 'Nick Holford'; [email protected]
> Subject: RE: [NMusers] omega matrix - friendly suggestion
>
> Thank you all for your suggestions. The more flexible OMEGA element entries
> 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 Cholesky
> format in NM 7.2, $OMEGA CHOLESKY may be used.
> 2) The traditional methods (FO/FOCE/Laplace) actually vary the omega elements
> in the Cholesky format for estimation, than transform back in the variance
> format after each iteration. This assures that as NONMEM arbitrarily searches
> 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-diagonal
> variance element to a non-zero value within a block of variances cannot be
> honored by the search algorithm past the first iteration. Additionally, to my
> best understanding, arbitrarily fixing an off-diagonal variance element to 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 can 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, but
> are constructed as a consequence of the expectation step, which provide
> conditional means and variances from which the omega elements are evaluated.
> 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: [email protected]
> Web: www.iconplc.com
>
> -----Original Message-----
> From: [email protected] [mailto:[email protected]] On
> Behalf Of Matt Hutmacher
> Sent: Wednesday, June 13, 2012 2:35 PM
> To: 'Nick Holford'; [email protected]
> Subject: RE: [NMusers] omega matrix - friendly suggestion
>
> Dear Paolo, Nick, all,
>
> A way to ensure a nonnegative definite (NND) variance-covariance matrix is 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 convergence 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 with
> 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 WRES 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 variances
> into linear combinations of thetas affecting the residuals). It is not enough
> to just model variances and correlations; such a parameterization will not
> 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. Upon
> 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 provide
> a lot of flexibility for analysis (as opposed to hand coding discussed below).
>
> There is also a modified Cholesky that can be used, which has the nice feature
> that one can simply, conveniently remove all the variance components related
> to a certain random effect (Ref 2). This was recently used in the attached
> reference to aid in selection of the $OMEGA structure with the fixed effects.
>
> 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-1077.
>
>
> -----Original Message-----
> From: [email protected] [mailto:[email protected]] On
> Behalf Of Nick Holford
> Sent: Wednesday, June 13, 2012 1:20 PM
> To: [email protected]
> Subject: Re: [NMusers] omega matrix - friendly suggestion
>
> Paolo,
>
> You may remember this discussion about the history of the correlation feature
> 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 using
> the SD and correlation style introduced in NONMEM 7. If I understand things
> correctly (this kind of math is not my strong point) it should be impossible
> 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 how to
> understand a matrix to be recognized in this format.
>
> Bob (Bauer) -- do you know how to do this? Perhaps you could talk to Marc
> 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
> >>
> >> [email protected] <mailto:[email protected]>
> >>
> >> www.a2pg.com http://www.a2pg.com
> >>
> >> *From:*[email protected]
> >> [mailto:[email protected]] *On Behalf Of
> >> *[email protected]
> >> *Sent:* Tuesday, June 12, 2012 10:24 AM
> >> *To:* [email protected]
> >> *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:[email protected]
> > ------------------------------------------------
>
> --
> Nick Holford, Professor Clinical Pharmacology
>
> First World Conference on Pharmacometrics, 5-7 September 2012 Seoul, Korea
> 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: [email protected]
> http://www.fmhs.auckland.ac.nz/sms/pharmacology/holford
>
>
>
>
>