[R] Coding matrix equation
Petr Savicky
savicky at praha1.ff.cuni.cz
Mon Apr 11 14:17:09 CEST 2011
On Mon, Apr 11, 2011 at 08:43:03AM +0100, matthew.r.robinson at sheffield.ac.uk wrote:
> Hi all,
>
> I have two matrices:
>
> G<-matrix(c(2.0, 0.5, 0.5, 0.5, 2.0, 0.5, 0.5, 0.5,2.0),3,3)
> P<-matrix(c(1.0, 0.5, 0.5, 0.5, 1.0, 0.5, 0.5, 0.5,1.0),3,3)
>
> and I want to run this equation to get a new matrix F:
>
> F = [P+2G]^-1/2 P [P+2G]^-1/2
>
> Could someone please tell me how to code this in R?
Hi.
Try the following.
sqrtSymMat <- function(A)
{
out <- eigen(A)
D <- diag(out$values)
U <- out$vectors
U %*% sqrt(D) %*% t(U)
}
A <- sqrtSymMat(solve(P + 2*G))
F <- A %*% P %*% A
See also the function svd() and its help ?svd.
The operator A^(1/2) works component-wise. There may be a function
computing the square root of a positive semidefinite matrix in some
of the CRAN extension packages
http://cran.at.r-project.org/web/packages/index.html
but i am not sure. The package mvtnorm
http://cran.at.r-project.org/web/packages/mvtnorm/index.html
computes the square root of a matrix internally. See the help
of the function ?rmvnorm.
Hope this helps.
Petr Savicky.
More information about the R-help
mailing list