Python Forum
The twelfth root of a matrix
Thread Rating:
  • 0 Vote(s) - 0 Average
  • 1
  • 2
  • 3
  • 4
  • 5
The twelfth root of a matrix
#1
[Image: 181125043338927257.png.html]

Hello everybody,

I need some help as I'm trying desperately to calculate the 12TH root of the 3 matrices in the image attached.
In fact, these are annual transition matrices and I would like to convert them in order to get monthly matrices. I know there is a way to calculate this using Python.
is someone familiar with this kind of calculation?

Thanks in advance for your help,

Best regards,
Reply
#2
Please, if this is not the right topic, tell me where I should post my thread to find my answers.

Best regards
Reply
#3
Typically, we expect people to post their Python code attempt (in code tags, not as an image) and we help them from there. (You might want to read the full rules: https://python-forum.io/misc.php?action=help)

Without code it's hard to say much. I would say that you can use the regular exponentiation operator (**) with 1/12.0 on numpy matrixes, I believe, and it would Just Work.
Reply
#4
Quote:I'm trying desperately to calculate the 12TH root of the 3 matrices
It depends on what you mean by that. If you mean the 12th root term by term, Mickseydel's answer should work. But if you mean in the sense of operator calculus of matrices (or functional calculus), the 12th root is a completely different mathematical operation.
Reply
#5
Hello Gribouillis and Micseydel,

Thanks a lot for your reply.
@Griboullis : you got the point, I don't want to calculate the 12TH root term by term (I can easily do that using excel).

I use excel MMULT function (12 times) to get matrix B as a result of matrix (A)^12.

But having matrix B, I'm not able to retrieve matrix A which would be matrix (B)^(1/12). Excel is not able to do that. Someone told me it could be possible using numpy matrices but I have no idea how to code this... I have never used Python before (@Micseydel: that's the reason why I did not attach a code attempt).

It would be wonderful if you could help me doing this!!!

Best Regards,
Reply
#6
Here is how I was able to do it for a particular matrix
import numpy as np
from numpy import linalg as la

A = np.array([
    [0.94328, 0.04779, 0.00892658],
    [0.367962, 0.52199275, 0.1100448852],
    [0, 0, 1],
    ])
# compute the eigenvalues and eigenvectors of A
w, v = la.eig(A)
# We are lucky, A is diagonalizable and its
# eigenvalues (the components of w) are non negative
# Aa = np.dot(v, np.dot(np.diag(w), la.inv(v)))
ww = w ** (1/12)
B = np.dot(v, np.dot(np.diag(ww), la.inv(v)))
Aaa = la.matrix_power(B, 12)
print(Aaa - A)
Output:
[[ -1.44328993e-15 -1.73472348e-16 -1.35308431e-16] [ -6.66133815e-16 -3.33066907e-16 1.80411242e-16] [ 0.00000000e+00 0.00000000e+00 0.00000000e+00]]
You need to figure out what to do if there are complex eigenvalues or if A is not diagonalizable. Usually, there won't be a unique matrix B solution of this problem, and usually the matrix B will be complex.

In the three matrices shown in the original post, the third line is 0, 0, 1. If this is always the case, it makes the problem simpler because 1 is always an eigenvalue, and we're basically working in dimension 2.
Reply
#7
Thanks a lot Gribouillis!!!

You're right regarding these 3 matrices : the third line is an absorbing state and will always be the case. So, as you mentioned, I hope this should be easier.
By the way, the output seems a bit odd to me, I would expect the result to be more similar to (I give you an example) :

For Matrix B as :

93,40% 6,00% 0,60%
18,00% 71,90% 10,10%
0,00% 0,00% 100,00%

I should get Matrix A :

99,376864% 0,601245% 0,0218908%
1,803736% 97,222401% 0,9738627%
0% 0% 100%

Have you got a solution to modify the code in order to get the same result?

Thank you so much for your help!

Best Regards,
Reply
#8
The format method can produce % output:
>>> x = 0.01803736
>>> print("{:%}".format(x))
1.803736%
Reply
#9
Thanks Gribouillis for the percentage format.
In fact, I was surprised by the output values, I'm not expecting negative values as a result.
Maybe it is the print (Aaa - A) which is misleading me?
I guess that Aa = np.dot(v, np.dot(np.diag(w), la.inv(v))) is materializing the relation A=PDP-1, is that right?
I think that if I understand your code, it is the value of B that I'm looking for.

Thanks a lot for your patience :-)
Reply
#10
Yes I should have used np.allclose(A, la.matrix_power(B, 12)) to check the difference. That said it remains the issue that the matrix B is not unique and it would also be a good idea to use the fact that the third row is 0 0 1. If you have code that tries to do this, you can post it here for further discussion.
Reply


Forum Jump:

User Panel Messages

Announcements
Announcement #1 8/1/2020
Announcement #2 8/2/2020
Announcement #3 8/6/2020