yes there are. But Christofer doesn't need exp(A), but A^n.
But there is a matpow-function %^% in this package, which is a little
bit slower, I think:
library(expm)
states<-1000
tm <- matrix(runif(states^2),nrow=states) # random transition matrix
for illustration
tm <- t(apply(tm,1,function (x) x/sum(x))) # make its rows sum to 1
p0<-pm <- c(0.5,0.5,rep(0,states-2)) # half of cases start in state 1,
half in state 2
n<-10000
system.time({dd<-eigen(tm,symmetric=F)
as.real(p0%*% dd$vectors%*% diag(dd$values^n)%*%solve(dd$vectors))})
User System elapsed
15.20 0.09 15.57
system.time(p0%*%(tm%^%n))
User System elapsed
38.61 0.00 39.62
Am 01.02.2011 17:16, schrieb Ben Bolker:
> Eik Vettorazzi <E.Vettorazzi <at> uke.uni-hamburg.de> writes:
>
>>
>> if you have a homogeneous mc (= a constant transition matrix), your
>> state at time 10 is given by (chapman-kolmogorov)
>> p10=p0 %*% tm^(10)
>> so you need a matrix power function.
>
> There are matrix exponential functions in the Matrix and expm
> packages ... don't know about their speed
>
> ______________________________________________
> R-help@r-project.org 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.
--
Eik Vettorazzi
Institut für Medizinische Biometrie und Epidemiologie
Universitätsklinikum Hamburg-Eppendorf
Martinistr. 52
20246 Hamburg
T ++49/40/7410-58243
F ++49/40/7410-57790
______________________________________________
R-help@r-project.org 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.
