R help - Sep 2006 - Exponentiate a matrix

Suppose I have a square matrix P

P <- matrix(c(.3,.7, .7, .3), ncol=2)

I know that 
> P * P 
Returns the element by element product, whereas
> P%*%P
Returns the matrix product.

Now, P^2 also returns the element by element product. But, is there a
slick way to write

P %*% P %*% P

Obviously, P^3 does not return the result I expect.

Thanks,
Harold






	[[alternative HTML version deleted]]

On 9/21/2006 10:40 AM, Doran, Harold wrote:
> Suppose I have a square matrix P
> 
> P <- matrix(c(.3,.7, .7, .3), ncol=2)
> 
> I know that 
> 
>> P * P 
> 
> Returns the element by element product, whereas
> 
>> P%*%P
> 
> Returns the matrix product.
> 
> Now, P^2 also returns the element by element product. But, is there a
> slick way to write
> 
> P %*% P %*% P
> 
> Obviously, P^3 does not return the result I expect.

I don't think there's anything built in, but it's easy to write your
own:

"%^%" <- function(mat, pow) {
   stopifnot(length(pow) == 1, all.equal(pow, round(pow)), nrow(mat) == 
ncol(mat))
   pow <- round(pow)
   if (pow < 0) {
     mat <- solve(mat)
     pow <- abs(pow)
   }
   result <- diag(nrow(mat))
   while (pow > 0) {
     result <- result %*% mat
     pow <- pow - 1
   }
   result
}

Now P %^% 3 will give you the matrix cube.

Duncan Murdoch

I am getting a bit rusty on some of these things, but I seem to recall 
that there is a numerical advantage (speed and/or accuracy?) to 
diagonalizing:

 > expM <- function(X,e) { v <- La.svd(X); v$u %*% diag(v$d^e) %*% v$vt
}

 > P <- matrix(c(.3,.7, .7, .3), ncol=2)
 > P %*% P %*% P
      [,1]  [,2]
[1,] 0.468 0.532
[2,] 0.532 0.468
 > expM(P,3)
      [,1]  [,2]
[1,] 0.468 0.532
[2,] 0.532 0.468

I think this also works for non-integer, negative, large, and complex 
exponents:

 > expM(P, 1.5)
          [,1]      [,2]
[1,] 0.3735089 0.6264911
[2,] 0.6264911 0.3735089
 > expM(P, 1000)
     [,1] [,2]
[1,]  0.5  0.5
[2,]  0.5  0.5
 > expM(P, -3)
        [,1]    [,2]
[1,] -7.3125  8.3125
[2,]  8.3125 -7.3125
 > expM(P, 3+.5i)
                  [,1]              [,2]
[1,] 0.4713+0.0141531i 0.5287-0.0141531i
[2,] 0.5287-0.0141531i 0.4713+0.0141531i
 >

Paul Gilbert

Doran, Harold wrote:
>Suppose I have a square matrix P
>
>P <- matrix(c(.3,.7, .7, .3), ncol=2)
>
>I know that 
>
>  
>
>>P * P 
>>    
>>
>
>Returns the element by element product, whereas
>
>  
>
>>P%*%P
>>    
>>
>
>Returns the matrix product.
>
>Now, P^2 also returns the element by element product. But, is there a
>slick way to write
>
>P %*% P %*% P
>
>Obviously, P^3 does not return the result I expect.
>
>Thanks,
>Harold
>
>
>
>
>
>
>	[[alternative HTML version deleted]]
>
>______________________________________________
>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.
>  
>
===================================================================================
La version fran?aise suit le texte anglais.

------------------------------------------------------------------------------------

This email may contain privileged and/or confidential inform...{{dropped}}