2016-10-24 2 views
1

AとBが次元10^8の5つの高さの高いスキニー行列であると考えてみましょう。 ;2つのスキニー・トール行列の間で行方向のドット積を計算する最速の方法R

r=10^8 
c=5 
A=matrix(runif(r*c,0,1),r,c) 
B=matrix(runif(r*c,0,1),r,c) 

Aの各行の内積を、Bの対応する行と比較したいと思います。

rowSums(A * B) 

しかし、それはかなり遅く、より速い方法があるかどうかを知りたいと思います。

答えて

2

これはあなたを失望させるかもしれませんが、Rレベルでは、これはあなた自身でCコードを書くことなく得ることができる最高のものです。問題はrowSums(A * B)をすることによって、あなたが効果的に

C <- A * B 
rowSums(C) 

を行っている最初の行は3つの大きな背の高い、薄い行列のフルスキャンを実行していることです。第2のラインは、1つの大きな背の高いマトリックスのフルスキャンを実行する。したがって、等価的に、背の高い薄いマトリックスを4回スキャンします(メモリーが豊富です)。

は、実際には、このような操作のために、最適なアルゴリズムが唯一行方向の外積行うことで、二回n * p背の高い薄型の行列をスキャンする必要があります。このように

rowsum <- numeric(n) 
for j = 1, 2, ... p 
    rowsum += A[,i] * B[,i] 

を、我々はまた、マトリックスCの発生を回避します。注意してください。上のコードは、有効なRコードまたはCコードではなく、偽のコードです。 C.


にアイデアが明確であり、我々はこれをプログラムしたいが、あなたの状況にアナロジーは、sum(x * y)crossprod(x, y)との速度差であるxを仮定し、yは、同じ長さの大きなベクトルで。

x <- runif(1e+7) 
y <- runif(1e+7) 

system.time(sum(x * y)) 
# user system elapsed 
# 0.124 0.032 0.158 

system.time(crossprod(x, y)) 
# user system elapsed 
# 0.036 0.000 0.056 

最初のケースでは、長いベクトルを4回スキャンし、2番目のケースでは2回スキャンします。統計計算で


関連

rowSums(A * B)

一般点ごとの予測分散に関連付けられた回帰計算に見られる、事実diag(tcrossprod(A, B))の効率的な評価です。例えば、モデル行列のQR分解からの薄い円を用いた通常の線形二乗回帰では、適合値の点ごとの分散はdiag(tcrossprod(Q))であり、rowSums(Q^2)によってより効率的に計算されます。しかし、すでに説明した理由から、これはまだ最速の評価ではありません。

関連する問題