Re: [R] Exponentiate a matrix

From: Dimitrios Rizopoulos <Dimitris.Rizopoulos_at_med.kuleuven.be>
Date: Thu 21 Sep 2006 - 18:58:27 GMT

Quoting Douglas Bates <bates@stat.wisc.edu>:

> On 9/21/06, Dimitrios Rizopoulos <Dimitris.Rizopoulos@med.kuleuven.be> wrote:

>>
>> Quoting Duncan Murdoch <murdoch@stats.uwo.ca>:
>>
>> > 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:
>>
>> I think there was function mtx.exp() in the Malmig package, but it
>> seems that this package has been withdrawn from CRAN. An old version
>> appears to exist in:
>>
>> http://r.meteo.uni.wroc.pl/src/contrib/Descriptions/Malmig.html
>>
>> Best,
>> Dimitris
>
> Is that function for matrix powers or for the exponential of a matrix
> (which is what I initally thought that Harold wanted)?  There is a
> function expm in the Matrix package, patterned on the octave function
> of the same name, the calculates the matrix exponential for a square
> matrix.

this function calculates the n-th power of a matrix, and this is what I thought Harold wanted, i.e.,

P %*% P %*% P %*% P

should be equal to

mtx.exp(P, 4)

>>
>>
>> > "%^%" <- 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
>> >
>> > ______________________________________________
>> > R-help@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.
>> >
>> >
>>
>>
>>
>> Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm
>>
>> ______________________________________________
>> R-help@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.
>>

>
> ______________________________________________
> R-help@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.
>
>



Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm



R-help@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 Fri Sep 22 05:18:31 2006

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.1.8, at Fri 22 Sep 2006 - 08:30:06 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.