Given a symmetric, positive definite, covariance matrix $S$, how can I calculate a matrix $C$ such that
$CSC^T$ is the Identity matrix $I$?
$\begingroup$
$\endgroup$
3
-
2$\begingroup$ Any invertible matrix for which $AA^T = S$ with $C=A^{-1}$ will work, so for example, Cholesky decomposition is one way - it's a simple algorithm for finding a lower triangular $A$ (or equivalently, in many implementations, an upper triangular $A^T$). That is, it finds $S = LL^T$, whence $L^{-1}S(L^{-1})^T = I$. But a number of other choices are available. $\endgroup$– Glen_bCommented Aug 16, 2013 at 0:57
-
2$\begingroup$ To appreciate @Glen_b's last point, observe that there are many orthogonal matrices $\mathbb{O}$: by definition, these preserve the identity under conjugation, $\mathbb{OIO}'=\mathbb{I}$. If $\mathbb{CSC}'=\mathbb{I}$, then clearly $\mathbb{I=OIO'=OCSC'O'=(OC)S(OC)'}$, whence $\mathbb{OC}$ works just as well. Conversely, $\mathbb{CSC'=I=BSB'}$ implies $\mathbb{B}$ is invertible and $\mathbb{B^{-1}C}$ is orthogonal. (The set of $n$ by $n$ orthogonal matrices is naturally a submanifold of $\mathbb{R}^{n^2}$ of dimension $n(n-1)/2$: that's the sense in which there are quite a lot of them.) $\endgroup$– whuber ♦Commented Aug 16, 2013 at 1:57
-
$\begingroup$ I expected this question would already have an answer for me to point to but I'm not seeing anything exactly on topic here. If nobody points a duplicate out or posts an answer, I may expand on my comment a little and turn it into one so that the question is answered at least once (and that expansion would include an explanation about how it's usually the case that if $S = X'X$ we do something else). $\endgroup$– Glen_bCommented Aug 16, 2013 at 2:03
Add a comment
|
1 Answer
$\begingroup$
$\endgroup$
Since I haven't found a duplicate and nobody else has pointed one out, I have expanded my comment into an answer and added an example to show you how it works.
Any invertible matrix $A$ for which $AA^T=S$ will work by taking $C=A^{−1}$ .
So for example, Cholesky decomposition is one common way to do it - it's a simple algorithm for finding a lower triangular $A$ (or equivalently, in many implementations, an upper triangular $A^T$). That is, it finds $S=LL^T$, whence $L^{−1}S(L^{−1})^T=I$.
Example:
Consider the matrix $$S=\left(\begin{matrix} 4 & 1 \\ 1 & 2.5 \end{matrix}\right) $$
We write $S = LL^T$, which implies that
$$s_{11} = l_{11}^2$$
That is, $l_{11}=2$ (we'll take positive square roots each time to make our diagonals on $L$ positive).
It also implies that $s_{21} = l_{21}l_{11}$, so $l_{21} = s_{21}/l_{11}$, i.e. $l_{21} = 1/2$. Then we see that $s_{22} = l_{21}^2+l_{22}^2$, so $l_{22}^2 = s_{22} - l_{21}^2 = 2.5 - 0.5^2 = 2.25$, in which case $l_{22} = 1.5$, so
$$L=\left(\begin{matrix} 2 & 0 \\ 0.5 & 1.5 \end{matrix}\right) $$
You can confirm that
$$C = L^{-1}=\left(\begin{matrix} \frac{1}{2} & 0 \\ -\frac{1}{6} & \frac{2}{3} \end{matrix}\right) $$
and that $CSC^T = I$.
In practice we don't invert $L$, however, but solve a linear system, it's more stable.
Another possibility is to use Singular Value Decomposition to find a matrix to do the diagonalizing, and then scale by the singular values.
But a number of other choices are available. (As whuber points out, a very large number.)
In practice, frequently we know $S$ is of the form $X^TX$ for some $X$. In that case, it's usually better to work on $X$ than $S$, via the QR decomposition.