1

As I mentioned in the title, so far I'm using the following code to power a matrix in R. This method is using the expm package and its %^% operator.

Qmp <- expm::`%^%`( as.matrix(Qm), (as.numeric(L)-1) )

Last day, I'm trying to execute that for a kind big matrix which dimensions are

> dim(Qm)
[1] 17328 17328

and I realized that is taking too long time to complete. To be honest, I'm running that code for an hour and still didn't finish.

So I was wondering if there is any other way to approach such task and be quicker.

EDIT:

After the powering of the Qm I'm also performing and a multiplication of two matrices like

Pm <- Qmp %*% as.matrix(Rm)

where 

dim(Rm)
>17328    58
J. Doe
  • 619
  • 4
  • 16
  • And to what is the exponent in your case? Depending on its value, it might be that 1 hour is not that much. – SergGr Mar 23 '18 at 18:44
  • At this specific session, the power is 6 – J. Doe Mar 23 '18 at 18:46
  • Do you know any property beforehand? Rank, symmetry, spectral radius etc..? Where does the matrix come from? – Frostic Mar 23 '18 at 19:00
  • 1
    R can not even do `crossprod(Qm,Qm)` in a reasonable time for a matrix that big. You could use Rcpp matrix multiplication. I dont know it is worth using other technique for a low power like 6. https://stackoverflow.com/questions/35923787/fast-large-matrix-multiplication-in-r – Frostic Mar 23 '18 at 19:06
  • @MaxFt: Actually, Qm and Rm are specific parts of a transition matrix that describes a weighted graph. I think that it is worth to use another technique because the power could be anything from 2 to 50. It happened in that case to be 6. In the next loop, it might be 34. So It is important to find another way to approach it. So, you suggest me to look for Rcpp solution? – J. Doe Mar 23 '18 at 19:49
  • 1
    I am almost sure that you should use Rcpp for such size of matrix. I would be very careful since computers round your numbers up to a certain level of precision. This could break your matrix properties. The sums of the row could become dramatically different to one when you take a large power and do many computation. I will look for an answer until someone else comes up with something. Also, you could focus on computing low powers if your chain is ergodic and the powered matrix serie is supposed to converge – Frostic Mar 23 '18 at 20:15
  • @MaxFt: What if we try compressing the initial transition matrix because every row has multiple 0 values and then use the `expm()`? Or that approach is already taken into account from `expm`? – J. Doe Mar 26 '18 at 11:55

0 Answers0