4

R で特定の種類のスパース マトリックスを設定しようとしています。次のコードでは、必要な結果が得られますが、非常に遅いです。

library(Matrix)

f <- function(x){

  out <- rbind(head(x, -1), tail(x, -1))
  out <- bdiag(split(out, col(out)))
  return(out)

}#END f 

x <- outer(1:250, (1:5)/10, '+')
do.call(rBind, apply(x, 1, f))

私が取り組んでいるシミュレーション研究では、これを何千回も行う必要があるため、これはかなり深刻なボトルネックです。この場合の Rprof() の出力はかなり紛らわしいものでした。高速化する方法について、ご提案をいただければ幸いです。

お時間をいただきありがとうございます。

4

1 に答える 1

9

このコードは、非常に遅いことをすべて回避するため、はるかに高速に実行されます (私のボックスでは 3.36 秒に対して <0.01 秒) rBind。重要なのは、最初に行インデックス、列インデックス、およびゼロ以外のセルの値を準備することです。を 1 回呼び出すsparseMatrix(i,j,x)だけで、 を 1 回も呼び出す必要なく疎行列が作成されrBind()ます。

library(Matrix)
A <- 1:250
B <- (1:5)/10
x <- outer(A, B, '+')

f2 <- function(x){
    n <- length(x)
    rep(x, each=2)[-c(1, 2*n)]
}

system.time({
  val <- as.vector(apply(x,1,f2))
  n <- length(val)
  i <- seq_len(n)
  j <- rep(rep(seq_len(length(B)-1), each=2), length.out=n)
  outVectorized <- sparseMatrix(i = i, j = j, x = val)
})
#    user  system elapsed 
#       0       0       0 

結果が同じであることを示すためだけに:

## Your approach
f <- function(x){
    out <- rbind(head(x, -1), tail(x, -1))
    out <- bdiag(split(out, col(out)))
    return(out)
}

system.time(outRBinded <- do.call(rBind, apply(x, 1, f)))
#    user  system elapsed 
#    3.36    0.00    3.36 

identical(outVectorized, outRBinded)
# [1] TRUE
于 2012-05-30T19:36:46.503 に答える