Rコードの最適化例：行列のクロス積
をテンプレートにして作成
[
トップ
] [
新規
|
一覧
|
検索
|
最終更新
|
ヘルプ
]
開始行:
SIZE(20){COLOR(magenta){Rコードの最適化例：行列のクロス積}} ([[Rコード最適化のコツと実例集]] に戻る)~
行列の積は順序が大事という例 (r-help 記事より、2004.10.07)~
失礼：順序の問題ではなく、クロス積を計算する二つの同値な関数 crossprod と %*% の違いというべきでした。実際、?crossprod によれば crossprod 関数の方が早いと書いてありました。どうも outer , %o% と違い単なる別名になっていないようです。~
二次形式の計算を様々に行ない比較する。結構計算速度が違うことがわかる。
> f1 <- function(a,X){ ignore <- t(a) %*% X %*% a }
> f2 <- function(a,X){ ignore <- crossprod(t(crossprod(a,X)),a) }
> f3 <- function(a,X){ ignore <- crossprod(a,X) %*% a }
> f4 <- function(a,X){ ignore <- (t(a) %*% X) %*% a }
> f5 <- function(a,X){ ignore <- t(a) %*% (X %*% a) }
> f6 <- function(a,X){ ignore <- crossprod(a,crossprod(X,a)) }
> a <- rnorm(100); X <- matrix(rnorm(10000),100,100)
> print(system.time( for(i in 1:10000){ a1<-f1(a,X)}))
[1] 3.98 0.02 4.09 0.00 0.00
> print(system.time( for(i in 1:10000){ a2<-f2(a,X)}))
[1] 2.40 0.00 2.45 0.00 0.00
> print(system.time( for(i in 1:10000){ a3<-f3(a,X)}))
[1] 1.31 0.00 1.34 0.00 0.00
> print(system.time( for(i in 1:10000){ a4<-f4(a,X)}))
[1] 4.01 0.00 4.01 0.00 0.00
> print(system.time( for(i in 1:10000){ a5<-f5(a,X)}))
[1] 4.06 0.00 4.09 0.00 0.00
> print(system.time( for(i in 1:10000){ a6<-f6(a,X)}))
[1] 1.19 0.00 1.19 0.00 0.00
> c(a1,a2,a3,a4,a5,a6) # 計算結果のチェック
[1] -56.73955 -56.73955 -56.73955 -56.73955 -56.73955 -56.73955
-crossprodも%*%も直接BLASを呼び出しますがtはごにょごにょしますんでその差かと. -- [[なかま]] &new{2004-10-07 (木) 11:06:32};
-中間さんのコメント参考にもう一つ実験。確かに COLOR(red){t(a)} を計算するオーバーヘッドがあるようですが、それだけでは説明できないような。
b <- t(a)
f0 <- function(a,b,X){ ignore <- b %*% X %*% a }
> print(system.time( for(i in 1:10000){ a1<-f0(a,b,X)}))
[1] 3.01 0.05 3.05 0.00 0.00
-よくよくみるとsrc/main/array.cには -- [[なかま]] &new{2004-10-07 (木) 22:08:41};
/* Don't trust the BLAS to handle NA/NaNs correctly: PR#4582
* The test is only O(n) here
*/
for (i = 0; i < nrx*ncx; i++)
if (ISNAN(x[i])) {have_na = TRUE; break;}
if (!have_na)
for (i = 0; i < nry*ncy; i++)
if (ISNAN(y[i])) {have_na = TRUE; break;}
などのチェック処理が存在するため,%*%の場合はdo_matprod->matprodと渡り、このチェック処理の為のオーバヘッドが入ります。
方やcrossprodの方はdo_matprod->crossprodと渡り、NANのチェック処理は行われません。
つまるところ安全に計算させるのと,チューニングされるであろう関数の存在と言うところでしょうか。
-いつも的確なコメントありがとうございます -> 中間さん。ソースを読まにゃあかんな。 -- &new{2004-10-08 (金) 07:39:44};
#comment
終了行:
SIZE(20){COLOR(magenta){Rコードの最適化例：行列のクロス積}} ([[Rコード最適化のコツと実例集]] に戻る)~
行列の積は順序が大事という例 (r-help 記事より、2004.10.07)~
失礼：順序の問題ではなく、クロス積を計算する二つの同値な関数 crossprod と %*% の違いというべきでした。実際、?crossprod によれば crossprod 関数の方が早いと書いてありました。どうも outer , %o% と違い単なる別名になっていないようです。~
二次形式の計算を様々に行ない比較する。結構計算速度が違うことがわかる。
> f1 <- function(a,X){ ignore <- t(a) %*% X %*% a }
> f2 <- function(a,X){ ignore <- crossprod(t(crossprod(a,X)),a) }
> f3 <- function(a,X){ ignore <- crossprod(a,X) %*% a }
> f4 <- function(a,X){ ignore <- (t(a) %*% X) %*% a }
> f5 <- function(a,X){ ignore <- t(a) %*% (X %*% a) }
> f6 <- function(a,X){ ignore <- crossprod(a,crossprod(X,a)) }
> a <- rnorm(100); X <- matrix(rnorm(10000),100,100)
> print(system.time( for(i in 1:10000){ a1<-f1(a,X)}))
[1] 3.98 0.02 4.09 0.00 0.00
> print(system.time( for(i in 1:10000){ a2<-f2(a,X)}))
[1] 2.40 0.00 2.45 0.00 0.00
> print(system.time( for(i in 1:10000){ a3<-f3(a,X)}))
[1] 1.31 0.00 1.34 0.00 0.00
> print(system.time( for(i in 1:10000){ a4<-f4(a,X)}))
[1] 4.01 0.00 4.01 0.00 0.00
> print(system.time( for(i in 1:10000){ a5<-f5(a,X)}))
[1] 4.06 0.00 4.09 0.00 0.00
> print(system.time( for(i in 1:10000){ a6<-f6(a,X)}))
[1] 1.19 0.00 1.19 0.00 0.00
> c(a1,a2,a3,a4,a5,a6) # 計算結果のチェック
[1] -56.73955 -56.73955 -56.73955 -56.73955 -56.73955 -56.73955
-crossprodも%*%も直接BLASを呼び出しますがtはごにょごにょしますんでその差かと. -- [[なかま]] &new{2004-10-07 (木) 11:06:32};
-中間さんのコメント参考にもう一つ実験。確かに COLOR(red){t(a)} を計算するオーバーヘッドがあるようですが、それだけでは説明できないような。
b <- t(a)
f0 <- function(a,b,X){ ignore <- b %*% X %*% a }
> print(system.time( for(i in 1:10000){ a1<-f0(a,b,X)}))
[1] 3.01 0.05 3.05 0.00 0.00
-よくよくみるとsrc/main/array.cには -- [[なかま]] &new{2004-10-07 (木) 22:08:41};
/* Don't trust the BLAS to handle NA/NaNs correctly: PR#4582
* The test is only O(n) here
*/
for (i = 0; i < nrx*ncx; i++)
if (ISNAN(x[i])) {have_na = TRUE; break;}
if (!have_na)
for (i = 0; i < nry*ncy; i++)
if (ISNAN(y[i])) {have_na = TRUE; break;}
などのチェック処理が存在するため,%*%の場合はdo_matprod->matprodと渡り、このチェック処理の為のオーバヘッドが入ります。
方やcrossprodの方はdo_matprod->crossprodと渡り、NANのチェック処理は行われません。
つまるところ安全に計算させるのと,チューニングされるであろう関数の存在と言うところでしょうか。
-いつも的確なコメントありがとうございます -> 中間さん。ソースを読まにゃあかんな。 -- &new{2004-10-08 (金) 07:39:44};
#comment
ページ名:
![[RjpWiki]](image/../Rlogo.png "[RjpWiki]")
