Re: [R] A function for raising a matrix to a power?

From: Ted Harding <ted.harding_at_nessie.mcc.ac.uk>
Date: Sun, 06 May 2007 11:35:41 +0100 (BST)


[Apologies: This will probably break the thread, but at the moment I cannot receive mail since my remote mail-server is down, and so am responding to the message as found on the R-help archives; hence this is not a "Reply"]

From:
Atte Tenkanen attenka at utu.fi
Sun May 6 11:07:07 CEST 2007

> Is there a function for raising a matrix to a power? For example if you
> like to compute A%*%A%*%A, is there an abbreviation similar to A^3?
>
> Atte Tenkanen
>
> > A=rbind(c(1,1),c(-1,-2))
> > A
> [,1] [,2]
> [1,] 1 1
> [2,] -1 -2
> > A^3
> [,1] [,2]
> [1,] 1 1
> [2,] -1 -8
>
> But:
>
> > A%*%A%*%A
> [,1] [,2]
> [1,] 1 2
> [2,] -2 -5

Of course, "^" alone acts element-wise on the matrix A, so the result of A^3 is the matrix B where B[i,j] = A[i,j]^3.

However, you can write your own, and call it for instance "%^%":

"%^%"<-function(A,n){
  if(n==1) A else {B<-A; for(i in (2:n)){A<-A%*%B}}; A   }

Then:

> A

     [,1] [,2]
[1,] 1 1
[2,] -1 -2

> A%^%1

     [,1] [,2]
[1,] 1 1
[2,] -1 -2

> A%^%2

     [,1] [,2]
[1,] 0 -1
[2,] 1 3

> A%^%3

     [,1] [,2]
[1,] 1 2
[2,] -2 -5

> A%^%100

              [,1] [,2]
[1,] -1.353019e+20 -3.542248e+20
[2,] 3.542248e+20 9.273727e+20

Hoping this helps!
Ted.



E-Mail: (Ted Harding) <ted.harding_at_nessie.mcc.ac.uk> Fax-to-email: +44 (0)870 094 0861
Date: 06-May-07                                       Time: 11:35:31
------------------------------ XFMail ------------------------------

______________________________________________
R-help_at_stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Sun 06 May 2007 - 10:40:15 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Sun 06 May 2007 - 12:31:27 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.