Clever Geek Handbook

Potentiation matrix without looping

I don't see a function for my problem. I have matrix A.

 A<-matrix(c(1,0,0,1,1,0,0,1,1),3,3)

I want to calculate A ^ 13, not element with element, but row with column. Example for A ^ 2

A%*%A

I can do with a loop, but my question is, is there an easy way or function to calculate this?

My result for A ^ 13 should be:

[1 13 78]
[0  1 13]
[0  0  1]
+3


source to share


2 answers


you can use expm 

 library(expm)
 A%^%13
 #   [,1] [,2] [,3]
 #[1,]    1   13   78
 #[2,]    0    1   13
 #[3,]    0    0    1

Benchmarks



On a smaller matrix

 set.seed(42)
 A <- matrix(sample(0:1, 100*100, replace=TRUE), 100, 100)
 f1 <- function() Reduce("%*%", replicate(200, A, FALSE))
 f2 <- function() A%^%200

 library(microbenchmark)
 microbenchmark(f1(), f2(), unit="relative", times=25L)
 # Unit: relative
 #expr      min       lq   median      uq      max neval
 #f1() 237.8381 237.6059 235.6788 239.053 225.3443    25
 #f2()   1.0000   1.0000   1.0000   1.000   1.0000    25
+2


source


You may be looking for Reduce

:



Reduce("%*%", replicate(13, A, FALSE))
#      [,1] [,2] [,3]
# [1,]    1   13   78
# [2,]    0    1   13
# [3,]    0    0    1
+3


source






More articles:


All Articles