3

重複の可能性:
Rでループ操作を高速化する

ループについていくつか質問があります。R はベクトル化された計算でより高速に動作することを知っています。これを利用するために以下のコードを変更したいと思います。フォーラムで他の回答を調べると、sapply 関数は for ループの内側を置き換えることができるようですが、ゼロのベクトルを生成しているため、エラーが発生します。Tao は 1000 のままで、これが問題を引き起こしていると思います。

私の主な関心事は速度です。アルゴリズム全体にループを作成し、さらに分析するために異なる V および n サイズでプロットする必要があるためです。

ご協力いただきありがとうございます

代替ループ

tao = 1000
L = (tao - 1)   
n = 10      
V = 5               
I = 10000                       
V_s = matrix(rnorm(I), I, 1)
V_b = matrix(rnorm(I), I, 1)

signal <- matrix(0, L, 1)  

for( j in (n:L)){

    sapply(((j-n+1):j),function (tao) signal[j] = signal[j] + abs(V_s[tao] - V_b[tao]))

    signal[j] = (signal[j] / (n * V) )

} 

オリジナルループ

tao = 1000
L = (tao - 1)   
n = 10      
V = 5               
I = 10000                       
V_s = matrix(rnorm(I), I, 1)
V_b = matrix(rnorm(I), I, 1)

signal <- matrix(0, L, 1)  

for( j in (n:L)){

    for( tao in ((j-n+1):j))    {

        signal[j] = (signal[j] + abs(V_s[tao] - V_b[tao]))

    }
        signal[j] = (signal[j] / (n * V) )

}
4

3 に答える 3

12

フィルターを使用すると、ループがなくても計算を実行できます (それsapplyは隠れたループにすぎません)。

absdif <- abs(V_s - V_b)
signal <- filter(absdif[1:L], rep(1/(n*V), n), sides=1)
signal[is.na(signal)] <- 0

ただし、フィルター処理に慣れていない場合、2 行目で何が起こっているかを理解するのは簡単ではありません。詳しく見てみましょう。

まず、ループで頻繁に使用するV_sとの絶対差を計算します。V_b続いてフィルターです。あなたの計算はn、各時間 value で過去の値を合計するだけですj。したがって、次のようなものがあります

signal[j] <- sum(absdif[j-n+1:j])

これは、畳み込みフィルターが行うこととまったく同じです-いくつかの値を合計します-一般的な形式で、重みを掛けて乗算します。重みとして1/(n*V)、すべての値に対して を選択します。これは、外側のループで行う正規化に対応します。最後の引数はsides=1、過去の値のみを取得するようにフィルターに指示するだけです (つまり、sides=2を意味しますsum(absdif[(j-n/2):(j+n/2)]))。

最後の行は、最初の値を埋めているだけNAです (フィルターには合計を計算するのに十分なデータがありません - これは最初のn値をスキップすることと同じです)。

最後に、いくつかのタイミング:

フルループ ソリューション:

   User      System       total 
  0.037       0.000       0.037 

ジュバの解決策:

   User      System       total 
  0.007       0.000       0.008 

フィルターを使用したソリューション:

   User      System       total 
  0.000       0.000       0.001 

フィルターの概念は非常によく研究されており、信じられないほど高速に実行できることに注意してください。

編集: に記載されているように、R は標準コマンド?filterで高速フーリエ変換を使用しません。filter通常、FFT は畳み込みを実装する最も効率的な方法です。ただし、これでも filter コマンドを次のように置き換えることで実行できます。

signal <- convolve(absdif[1:L], rep(1/(n*V), n), type='filter')

nに設定する代わりに、最初のエントリが削除されることに注意してくださいNA。ただし、結果は同じです。この時間のタイミングは役に立ちません - 合計時間はsystem.time... の 3 桁の出力を下回っています。ただし、 の R ヘルプにある次の注意事項に注意してくださいfilter

convolve(, type="filter") は計算に FFT を使用するため、単変量系列の長いフィルターの場合は高速になる可能性がありますが、時系列を返さず (したがって、時間の配置が不明確になります)、欠損値を処理しません。 . たとえば、長さ 1000 の一連のフィルタに対して長さ 100 のフィルタを使用すると、filter の方が高速になります。

于 2013-01-25T10:55:41.963 に答える
3

計算をベクトル化することは、常に *apply 関数を使用することを意味するわけではありません。

たとえば、2 番目の for ループを vector indexing に置き換えることで、物事を単純化して高速化できます。

for(j in (n:L)){
  sel <- (j-n+1):j
  signal[j] <- sum(abs(V_s[sel] - V_b[sel])) / (n*V)
}

このソリューションの場合、私のシステムでの実行時間は次のとおりです。

utilisateur     système      écoulé 
      0.008       0.004       0.009 

あなたのforループの場合は次のとおりです。

utilisateur     système      écoulé 
       0.06        0.00        0.06 

ところで、tao名前を 2 つの異なるものに使用しないでください。

于 2013-01-25T10:48:18.077 に答える
0

明示的なループが正しい計算であると仮定して、これを試してください:

 signal[j]<- signal[j] + 
              sapply((j-n+1):j, 
                   FUN = function(iter){ 
                           abs(V_s[iter] - V_b[iter])
                   }, V_s = V_s, V_b = V_b)

sapply は、V_s と V_b の間の iter 番目のインデックスの絶対差を返していることに注意してください。次に、これが signal[j] に追加されます

于 2013-01-25T10:27:52.417 に答える