diag_Omegas function

Simultaneously diagonalize two covariance matrices

diag_Omegas Simultaneously diagonalizes two covariance matrices using eigenvalue decomposition.

diag_Omegas(Omega1, Omega2)

Arguments

• Omega1: a positive definite $(dxd)$ covariance matrix $(d>1)$
• Omega2: another positive definite $(dxd)$ covariance matrix

Returns

Returns a length $d^2 + d$ vector where the first $d^2$ elements are $vec(W)$ with the columns of $W$ being (specific) eigenvectors of the matrix $\Omega_2\Omega_1^{-1}$ and the rest $d$ elements are the corresponding eigenvalues "lambdas". The result satisfies $WW' = Omega1$ and $Wdiag(lambdas)W' = Omega2$.

If Omega2 is not supplied, returns a vectorized symmetric (and pos. def.) square root matrix of Omega1.

Details

See the return value and Muirhead (1982), Theorem A9.9 for details.

Warning

No argument checks! Does not work with dimension $d=1$!

Examples

# Create two (2x2) coviance matrices using the parameters W and lambdas:
d <- 2 # The dimension
W0 <- matrix(1:(d^2), nrow=2) # W
lambdas0 <- 1:d # The eigenvalues
(Omg1 <- W0%*%t(W0)) # The first covariance matrix
(Omg2 <- W0%*%diag(lambdas0)%*%t(W0)) # The second covariance matrix

# Then simultaneously diagonalize the covariance matrices:
res <- diag_Omegas(Omg1, Omg2)

# Recover W:
W <- matrix(res[1:(d^2)], nrow=d, byrow=FALSE)
tcrossprod(W) # == Omg1, the first covariance matrix

# Recover lambdas:
lambdas <- res[(d^2 + 1):(d^2 + d)]
W%*%diag(lambdas)%*%t(W) # == Omg2, the second covariance matrix

References

• Muirhead R.J. 1982. Aspects of Multivariate Statistical Theory, Wiley.