R-alpha: matrix exponentiation

Jim Lindsey (jlindsey@luc.ac.be)
Wed, 20 Nov 1996 09:31:12 +0100 (MET)


From: Jim Lindsey <jlindsey@luc.ac.be>
Message-Id: <9611200831.AA00309@alpha.luc.ac.be>
Subject: R-alpha: matrix exponentiation
To: r-testers@stat.math.ethz.ch
Date: Wed, 20 Nov 1996 09:31:12 +0100 (MET)

1.With x a matrix,
  diag(x) <- 1
  gives Error: couldn't find function "diag<-"
2.For the function, sum(), the help says na.rm=T, but the function
  has F.
3.Here is a matrix exponentiation function for anyone working on
  stochastic processes. If anyone knows a better algorithm than brute
  force, please let me know. Jim
--------------------------------------------------------------------
mexp <- function(x, t=1, n=50){
	if(!is.matrix(x)){
		print("x must be a matrix")
		return(1)}
	if(length(dim(x))!=2){
		print("x must be a two dimensional matrix")
		return(2)}
	if(dim(x)[1]!=dim(x)[2]){
		print("x must be a square matrix")
		return(3)}
	p <- diag(dim(x)[1])
	tt <- p
	for(r in 1:n){
		tt <- x%*%tt*t/r
		p <- p+tt}
	p}
=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-
r-testers mailing list -- For info or help, send "info" or "help",
To [un]subscribe, send "[un]subscribe"
(in the "body", not the subject !)  To: r-testers-request@stat.math.ethz.ch
=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-