On page 22 of the R-introduction guide it's written: the quadratic form x^{'} A^{-1} x which is used in multivariate computations, should be computed by something like x%*%solve(A,x), rather than computing the inverse of A. Why isn't it good to compute t(x) %*% solve(A) %*% x? Thanks a lot for help!
Pedro: Solving the equation in this fashion can be computationally slow and unstable. There is R an article in R News that answers this question directly and is a pretty easy read with good examples. Check it out at the following link: http://cran.r-project.org/doc/Rnews/Rnews_2004-1.pdf -----Original Message----- From: r-help-bounces at stat.math.ethz.ch [mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Alvarez Pedro Sent: Thursday, November 03, 2005 8:02 AM To: r-help at stat.math.ethz.ch Subject: [R] quadratic form On page 22 of the R-introduction guide it's written: the quadratic form x^{'} A^{-1} x which is used in multivariate computations, should be computed by something like x%*%solve(A,x), rather than computing the inverse of A. Why isn't it good to compute t(x) %*% solve(A) %*% x? Thanks a lot for help! ______________________________________________ R-help at stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Alvarez Pedro wrote:> On page 22 of the R-introduction guide it's written: > > the quadratic form x^{'} A^{-1} x which is used in > multivariate computations, should be computed by > something like x%*%solve(A,x), rather than computing > the inverse of A. > > Why isn't it good to compute t(x) %*% solve(A) %*% x? > > Thanks a lot for help! > > ______________________________________________ > R-help at stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.htmlSee Bates, D. (2004): Least Squares Calculations in R. R News 4(1), 17-20, http://CRAN.R-project.org/doc/Rnews/ Uwe Ligges
On Thu, 3 Nov 2005, Alvarez Pedro wrote:> On page 22 of the R-introduction guide it's written: > > the quadratic form x^{'} A^{-1} x which is used in > multivariate computations, should be computed by > something like x%*%solve(A,x), rather than computing > the inverse of A. > > Why isn't it good to compute t(x) %*% solve(A) %*% x?The answer is only two lines above; Numerically, it is both inefficient and potentially unstable to compute @code{x <- solve(A) %*% b} instead of @code{solve(A,b)}. See also the footnote. -- Brian D. Ripley, ripley at stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UK Fax: +44 1865 272595
Alvarez Pedro <palvarez7777 at yahoo.es> writes:> On page 22 of the R-introduction guide it's written: > > the quadratic form x^{'} A^{-1} x which is used in > multivariate computations, should be computed by > something like x%*%solve(A,x), rather than computing > the inverse of A. > > Why isn't it good to compute t(x) %*% solve(A) %*% x?It's just a waste of CPU time. For k x k matrices, solution of Ax=b is of computational complexity O(k^2) whereas inversion of A is O(k^3). This is obviously more important for k=1000 than for k=5. -- O__ ---- Peter Dalgaard ??ster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk) FAX: (+45) 35327907
Hi Alvarez If you define quad.form.inv <- function (M, x) { drop(crossprod(x, solve(M, x))) } then you will avoid an expensive call to %*% as well. HTH Robin On 3 Nov 2005, at 13:01, Alvarez Pedro wrote:> On page 22 of the R-introduction guide it's written: > > the quadratic form x^{'} A^{-1} x which is used in > multivariate computations, should be computed by > something like x%*%solve(A,x), rather than computing > the inverse of A. > > Why isn't it good to compute t(x) %*% solve(A) %*% x? > > Thanks a lot for help! > > ______________________________________________ > R-help at stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting- > guide.html >-- Robin Hankin Uncertainty Analyst National Oceanography Centre, Southampton European Way, Southampton SO14 3ZH, UK tel 023-8059-7743
On 3 Nov 2005, at 17:25, Robin Hankin wrote:> Hi Alvarez > > > If you define > > quad.form.inv <- function (M, x) > { > drop(crossprod(x, solve(M, x))) > } > > then you will avoid an expensive call to %*% as well. >Is %*% really expensive in all platforms? I had a function that used QR decomposition instead of quadratic forms, but then I got a message from Canada suggesting that %*% would be faster. Indeed, it was in not-too-large data sets and in Mac (or powerpc). I run some tests with real applications, and found that my 800MHz iBook G4 run like a 2.5GHz Intel machine when %*% was used. This really was architecture dependent, since the performance boost was similar under OS X and Linux in the very same PowerPC. So it seems that %*% is very cheap if you have PowerPC, but it may be expensive in Intel. (I also run a test in Sun, and it was somewhere between Intel and PowerPC.) cheers, jari oksanen -- Jari Oksanen, Oulu, Finland
If you meant QR vs. inverting X'X for linear regression, the motivation for using QR is not speed, but numerical stability. There's no univerally good least squares algorithm that would be uniformly better than anything else for any kind of data. Andy> From: Jari Oksanen > > > On 3 Nov 2005, at 17:25, Robin Hankin wrote: > > > Hi Alvarez > > > > > > If you define > > > > quad.form.inv <- function (M, x) > > { > > drop(crossprod(x, solve(M, x))) > > } > > > > then you will avoid an expensive call to %*% as well. > > > Is %*% really expensive in all platforms? I had a function > that used QR > decomposition instead of quadratic forms, but then I got a > message from > Canada suggesting that %*% would be faster. Indeed, it was in > not-too-large data sets and in Mac (or powerpc). I run some > tests with > real applications, and found that my 800MHz iBook G4 run like > a 2.5GHz > Intel machine when %*% was used. This really was architecture > dependent, since the performance boost was similar under OS X > and Linux > in the very same PowerPC. So it seems that %*% is very cheap if you > have PowerPC, but it may be expensive in Intel. (I also run a test in > Sun, and it was somewhere between Intel and PowerPC.) > > cheers, jari oksanen > -- > Jari Oksanen, Oulu, Finland > > ______________________________________________ > R-help at stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! > http://www.R-project.org/posting-guide.html > >