On 25-04-2013, at 17:18, Ravi Varadhan <ravi.varad...@jhu.edu> wrote:
> Hi, > I am generating large Kac matrices (also known as Clement matrix). This a > tridiagonal matrix. I was wondering whether there is a vectorized solution > that avoids the `for' loops to the following code: > > n <- 1000 > > Kacmat <- matrix(0, n+1, n+1) > > for (i in 1:n) Kacmat[i, i+1] <- n - i + 1 > > for (i in 2:(n+1)) Kacmat[i, i-1] <- i-1 > > The above code is fast, but I am curious about vectorized ways to do this. You could vectorize like this Kacmat <- matrix(0, n+1, n+1) Kacmat[row(Kacmat)==col(Kacmat)-1] <- n -(1:n) + 1 Kacmat[row(Kacmat)==col(Kacmat)+1] <- 1:n But this show that your version is pretty quick f1 <- function(n) { Kacmat <- matrix(0, n+1, n+1) for (i in 1:n) Kacmat[i, i+1] <- n - i + 1 for (i in 2:(n+1)) Kacmat[i, i-1] <- i-1 Kacmat } f2 <- function(n) { Kacmat <- matrix(0, n+1, n+1) Kacmat[row(Kacmat)==col(Kacmat)-1] <- n -(1:n) + 1 Kacmat[row(Kacmat)==col(Kacmat)+1] <-1:n Kacmat } library(compiler) f1.c <- cmpfun(f1) f2.c <- cmpfun(f2) n <- 5000 system.time(K1 <- f1(n)) system.time(K2 <- f2(n)) system.time(K3 <- f1.c(n)) system.time(K4 <- f2.c(n)) identical(K2,K1) identical(K3,K1) identical(K4,K1) # > system.time(K1 <- f1(n)) # user system elapsed # 0.386 0.120 0.512 # > system.time(K2 <- f2(n)) # user system elapsed # 3.779 1.141 4.940 # > system.time(K3 <- f1.c(n)) # user system elapsed # 0.323 0.119 0.444 # > system.time(K4 <- f2.c(n)) # user system elapsed # 3.607 0.852 4.472 # > identical(K2,K1) # [1] TRUE # > identical(K3,K1) # [1] TRUE # > identical(K4,K1) # [1] TRUE Berend ______________________________________________ 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.