10
$\begingroup$

Can we decompose real matrices in a way that is analogous to polar decomposition in the complex case?

What I mean is: given an invertible real matrix $M$, can we always write: $$ M = OP, $$

maybe uniquely, where $O$ is orthogonal, and $P$ is symmetric positive-definite?

Thanks.

$\endgroup$
4
  • $\begingroup$ en.wikipedia.org/wiki/Polar_decomposition $\endgroup$
    – Bill Cook
    Sep 29, 2015 at 11:56
  • $\begingroup$ @BillCook More or less... oh wait, I forgot symmetry! $\endgroup$
    – geodude
    Sep 29, 2015 at 11:57
  • $\begingroup$ The answer is yes. This is a consequence of the spectral theorem, whereby $M^TM$ can be orthogonally diagonalized. $\endgroup$ Sep 29, 2015 at 14:37
  • $\begingroup$ @Omnomnomnom Wow. Can you explain more? (Or make it into an answer?) $\endgroup$
    – geodude
    Sep 29, 2015 at 14:39

3 Answers 3

19
$\begingroup$

The answer is yes. In fact, we don't need the spectral theorem to prove it.

Suppose that $M$ is a real invertible matrix. Then $M^TM$ is positive definite and has the unique positive semidefinite square root $P = \sqrt{M^TM}$. Now, note that $P$ has the property $\|Px\| = \|Mx\|$ for all vectors $x$.

If $M$ is invertible, then $P$ is invertible as well, and we have $M = (MP^{-1})P$. We note that $MP^{-1}$ is orthogonal since for all $y = Px$, we have $$ \|MP^{-1}y\| = \|y\| $$ Thus, we have a polar decomposition with $O =( MP^{-1})$ .

We note that this decomposition is unique. In particular, suppose $M = OP$ for orthogonal $O$ and positive $P$. then $$ M^TM = PO^TOP = P^2 $$ And so, by the uniqueness of positive definite square roots, $P$ is uniquely determined. Then we can rearrange $M = OP$ to find $O = MP^{-1}$ is also unique determined.


Things get a bit trickier when $M$ is not invertible, but we can still guarantee a (non-unique) polar decomposition.

$\endgroup$
6
  • $\begingroup$ This can also be derived easily from the singular value decomposition. $\endgroup$ Sep 29, 2015 at 15:01
  • $\begingroup$ Wait...how do you prove that $PM^{-1}P=M$? $\endgroup$
    – geodude
    Sep 30, 2015 at 12:36
  • $\begingroup$ I think I meant to write $MP^{-1}$... the proof is pretty much the same, though. The only other thing we need to show is that $M$ invertible $\iff$ $M^TM$ invertible $\iff$ $P$ invertible, which isn't too hard. $\endgroup$ Sep 30, 2015 at 13:10
  • $\begingroup$ See my latest edit. $\endgroup$ Sep 30, 2015 at 13:12
  • $\begingroup$ Oh yeah, like this it works! Thanks. $\endgroup$
    – geodude
    Sep 30, 2015 at 13:29
0
$\begingroup$

A derivation of the polar decomposition for a 2x2 matrix can be found on Polar Decomposition

I also elaborate the relation with SVD. In short, shuffling around with the SVD decompostion:

$A= U \Sigma V^T$ = $(UV^T) (V \Sigma V^T) $= $R_{\theta} S_{AT}$

$S_{AT}$ scales along the principal axes of the ellipse defined by $A^T(unitcircle)$

$R_{\theta}$ rotates over the angle between the principal axes of $A$ and $A^T$

$A^{-1}=\left({\ \ U\ \mathrm{\Sigma}}^{-1}U^{-1}\right)^{-1}R_\theta$ = $S_A R_\theta$

$R_{\theta}$ rotates over the angle between the principal axes of $A$ and $A^T$

$S_{A}$ scales along the principal axes of the ellipse defined by $A(unitcircle)$

$\endgroup$
1
  • $\begingroup$ Is this voted down because it is wrong? not helpful? $\endgroup$ Nov 30, 2019 at 12:47
0
$\begingroup$

Theory of Polar Decomposition is described in Wikipedia For the sake of simplicity we shall restrict our attention to real-valued non-singular square (invertible) matrices.
Theorem. Every such a matrix B can be decomposed as follows and the decomposition is unique: $$ B = U Q $$ Where $U =$ orthogonal matrix and $Q=$ symmetric positive definite matrix.
Instead of providing still another proof of the above theorem, I have decided to be practical and describe how to numerically obtain the decomposition.
First define the transpose of the matrix multiplied with the original: $$ P = B^T B = \left(B^TB\right)^T = P^T \quad \mbox{: symmetric} $$ It follows that $P$ is symmetric and positive definite. Then what we need is the square root of the matrix $P$. In order to understand how to obtain it, let's take a look at Newton's method for obtaining the square root of a real positive number $p$: $$ f(x) = x^2 - p\quad \Longrightarrow \\ x_{n+1} = x_n - \frac{x_n^2-p}{2x_n} = \left(x_n + p\,x_n^{-1}\right)/2 $$ Let's do the same for our matrix $P$, iterations starting with the unit matrix: $$ X_0 = I \quad ; \quad X_{n+1} = \left(X_n + P X_n^{-1}\right)/2 \quad ; \quad \sqrt{P} = \lim_{n\to\infty} X_n $$ According to a MSE reference $X_n^{-1}$ is positive definite and symmetric, provided that $X_n$ is positive definite and symmetric. Products and sums of positive definite and symmetric matrices are positive definite and symmetric too. Therefore, by induction to $n$, it can be proved that $\sqrt{P}$ is positive definite and symmetric as well. Now define $Q=\sqrt{P}=Q^T$. At last define $U = B Q^{-1}$ and prove that $U$ is orthogonal: $$ U^{T}U = \left(B Q^{-1}\right)^T\left(B Q^{-1}\right) = \left(Q^{-1}\right)^T\left(B^TB\right)Q^{-1} = \\ \left(Q^T\right)^{-1} Q^2 Q^{-1} = \left(Q^{-1}Q\right)\left(QQ^{-1}\right) = II = I $$ The conclusion looks like trivial: $$ B = \left[B\left(B^TB\right)^{-1/2}\right]\left[\left(B^TB\right)^{1/2}\right] $$ So far so good about theory. I want to talk about practice now, which in modern times is computer programming:

Accompanying software

My favorite programming language is (Delphi) Pascal. First we need some standard matrix manipulation routines and a Newton-Raphson routine for calculating the matrix square root. These are implemented in a Unit called 'Wiskunde'. The software is completed by testing if theory works in practice as expected. The number of iterations in procedure 'Newton' has been determined as follows. A suitable stopping criterion is: $$ \left|\left| X^2 - B^TB\,\right|\right| < \epsilon $$ with $\epsilon = 10^{-10}$ and the norm of a $n \times n$ matrix $A$ defined as: $$ \left|\left|A\right|\right| = \sqrt{\sum_{i=1}^n \sum_{j=1}^n A_{ij}^2} $$ Convergence is very fast. When cast in MathJax format the output looks like this: $$ \small \begin{bmatrix} 0.361048 & -0.297419 & -0.227079 \\ 0.171654 & -0.181309 & -0.338205 \\ -0.127762 & -0.074326 & -0.417988 \end{bmatrix} = \\ \small \begin{bmatrix} 0.491571 & -0.868969 & -0.057015 \\ 0.599413 & 0.385125 & -0.701701 \\ -0.631714 & -0.310760 & -0.710187 \end{bmatrix} \begin{bmatrix} 0.361082 & -0.207928 & -0.050301 \\ -0.207928 & 0.211719 & 0.196967 \\ -0.050301 & 0.196967 & 0.547115 \end{bmatrix} $$
$\endgroup$

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .