17

If I am given a matrix, for example,

$$A = \begin{bmatrix} 0.7 & 0.2 & 0.1 \\[0.3em] 0.2 & 0.5 & 0.3 \\[0.3em] 0 & 0 & 1 \end{bmatrix}$$

how do I calculate the fractional powers of $A$, say, $A^{\frac{1}{2}}$ or $A^{\frac{3}{2}}$?

Rodrigo de Azevedo
  • 23,223
user121692
  • 457
  • 1
    Why do you even think it is defined? – Git Gud Mar 30 '14 at 13:59
  • 3
    @GitGud Why wouldn't it be? – fgp Mar 30 '14 at 14:01
  • @fgp Why would it be? Matrices necessarily having square roots? Square roots over what field? Furthermore the square root need not be unique, so using the symbol $\sqrt A$ leads to problems. – Git Gud Mar 30 '14 at 14:03
  • And for what's it's worth, this is a stochastic matrix, so non-real solutions are probably not acceptable. – Git Gud Mar 30 '14 at 14:04
  • I tried typing this into my graphing calculator and I am getting an error, so does this mean the answer is not unique? – user121692 Mar 30 '14 at 14:05
  • 1
    @user121692 No, it means it doesn't make sense. You need to define the symbol $A^{1/2}$. If you're just looking for a matrix $X$ such that $X^2=A$, then you should specify so, but this isn't always possible. – Git Gud Mar 30 '14 at 14:06
  • @GitGud Well, you don't need a field to define square roots - for every ring $R$ there's the ring of polynomials $R[x]$ and you can thus ask for zeros of $s_a(x) = x^2 - a$ for any $a \in R$. Thus, $A^\frac{p}{q}$ is naturally interpreted as "find a matrix $B$ such that $A^p = B^q$". Whether or not there's a unique solution is a different question... – fgp Mar 30 '14 at 14:08
  • @fgp The field thing was just to set a background. And no, you cannot define a symbol as having more than one meaning in a fixed context. – Git Gud Mar 30 '14 at 14:10
  • 1
    @Git Gud Well for my specific purposes, the matrix I am working with are Markov chains transition matrix and the matrix A I wrote above is the annual transition probability. I need to calculate semi-annual, I suppose what I am looking for then is as you said $X$ such that $X^2 = A?$ – user121692 Mar 30 '14 at 14:11
  • @user121692 OK, that has been answered. It's possible here. – Git Gud Mar 30 '14 at 14:11

2 Answers2

21

If a matrix is diagonalizable, then diagonalize it, $A=PDP^{-1}$ and apply the power to the diagonal

$$A^n=PD^n P^{-1}$$

The diagonal values are acted on individually.

octave gives:

$$P=\begin{bmatrix} 0.85065 & -0.52573 & 0.57735\\ 0.52573 & 0.85065 & 0.57735\\ 0.00000 & 0.00000 & 0.57735\end{bmatrix}$$

$$D=\operatorname{diag}(0.82361, 0.37639,1)$$ I realize this is a numerical uglyness but I don't have a symbolic manipulation software at hand from this computer. However, the eigenvalues are different so this is a diagonalization.

The square root is $$\sqrt{A}= \begin{bmatrix}0.82626 & 0.13149 & 0.04225\\ 0.13149 & 0.69477 & 0.17374\\ 0.00000 & 0.00000 & 1.00000\end{bmatrix}$$

This definition satisfies the requirement for roots that $(A^{1/p})^p=A$ for positive definite matrices (just like with $\sqrt{x}$ for scalars).

In a similar way, you can define functions on matrices through their power series. For instance, $e^A=P \exp(D)P^{-1}$ is perfectly well defined for diagonalizable matrices.

The convergence criteria and domain of these functions gets generalized and usually involves conditions for eigenvalues, positive-definiteness, symmetry, ortogonality and so on.

Note that the term square root of a matrix is sometimes used to represent a Cholesky decomposition, which instead works as $A=LL^T$ where $L$ is a lower triangular matrix. This is not the square root in the strictest sense, but it works like one for some numerical procedures.

orion
  • 16,154
  • 1
  • 34
  • 45
  • Did you read the question? – Git Gud Mar 30 '14 at 14:01
  • 2
    Sure, the matrix is diagonalizable, use $n=\frac{1}{2}$ and calculate. – orion Mar 30 '14 at 14:02
  • 'The'${{{{}}}}$? – Git Gud Mar 30 '14 at 14:10
  • If there are 8 "square roots", is only one of them a valid transition matrix? – Neil W Mar 30 '14 at 14:25
  • 1
    @Neil the same applies as for reals: you define it to be the one with positive eigenvalues. That's why you start having problems selecting the branch when the original matrix is not positive definite. – orion Mar 30 '14 at 14:30
  • Thanks for the answer but $D^\frac{1}{2}$ still gives me an error? – user121692 Mar 30 '14 at 14:37
  • What kind of an error? How are you calculating this? – orion Mar 30 '14 at 14:38
  • So you are just square rooting each individual entry? – user121692 Mar 30 '14 at 14:40
  • 1
    Exactly. A calculator probably doesn't allow this operation on a matrix, but serious software for matrix manipulation swallows this without a problem. For instance, in octave/matlab A^(1/2) works directly (and gives the same result). And yes, just take the square root of the eigenvalues. – orion Mar 30 '14 at 14:42
  • Thanks, let me try this again. – user121692 Mar 30 '14 at 14:44
  • The definition of exponential of a matrix is actually defined even for non-diagonalizable matrices, at least over the field of complex numbers, and can be calculated using Jordan Canonical Form for matrices non-diagonalizable matrices, since the Maclaurin Series of exp(M) always terminates if M is in Jordan Canonical Form. The rules that you can use any basis to take powers of, add, and multiply scalars by matrices explains why the Taylor series definitions of functions allow you to work on similar matrices like the JCF, RCF, or diagonalization to evaluate such functions. – Mr. Nichan Feb 08 '25 at 20:24
0

I notice there's a scipy module that supposedly does this (scipy.linalg.fractional_matrix_power):

https://docs.scipy.org/doc/scipy/reference/generated/scipy.linalg.fractional_matrix_power.html

The citation there is:

Nicholas J. Higham and Lijing lin (2011) “A Schur-Pade Algorithm for Fractional Powers of a Matrix.” SIAM Journal on Matrix Analysis and Applications, 32 (3). pp. 1056-1078. ISSN 0895-4798

Someone on GitHub commented that those authors have improved their algorithm since then:

Nicholas J. Higham and Lijing lin (2013) "An Improved Schur--Padé Algorithm for Fractional Powers of a Matrix and Their Fréchet Derivatives" SIAM Journal on Matrix Analysis and Applications (Impact Factor: 1.59). 07/2013; 34(3).

I found a version of the latter on something called "Amanote Research".

I should note that, if you define fractional powers of a matrix X like so:

$X^{p/q} = \{R^p |R^q=X\}$

then sets like $X^{1/2}$ can definitely have multiple elements. E.g., if Z is a zero matrix in a dimension greater than one, then $Z^{1/2}$ has uncountably infinitely many distinct elements.

However, similar to how square roots of real numbers are often defined to only exist for positive numbers and to only be the positive root (even though there is also a negative one), I have seen, in the description of Singular Value Decomposition, a convention where square roots are defined only for matrices with an eigenbasis and all positive real eigenvalues, and are defined in exactly the way orion says (take the positive square root of each element in the diagonalization, then undo the diagonalization, so you get a matrix with same eigenvectors, but all the eigenvalues square-rooted). For such matrices, this is one element of the $X^{1/2}$ set I described above. I haven't checked what that scipy module does, but maybe it only handles cases where the definition orion gave works.

Mr. Nichan
  • 263