1

次の関数を使用して、一部のデータのカーネル密度を推定しています

KDens = function(x,h,N) {
         fx = matrix(0,N,1)
         Kx = AK(x,h,N)
         for (i in 1:N) {
         fx[i] = sum(Kx[i,], na.rm=T)/N
         }
         return(fx) }

これがループの高速化に関する最初の質問ではないことはわかっています。サイトをチェックしたところ、一部の機能を使用した方が高速な場合があることがわかりましたがapply、ループを正しく設定できた場合は常にそうとは限りません。

上記のコードでは、「不要なもの」はすべてループから除外されています。これは、私の理解が正しければ、計算を高速化することが示唆されているためです。KDensしかし、上記の関数とdensityRにデフォルトで実装されている関数を比較してみました。まあ、私のマシンでは〜30秒必要ですdensityが、1〜2秒必要です。KDens

trywiththis <- rnorm(4800)
x = trywiththis
N = length(trywiththis)
h = 1.059*sd(trywiththis , na.rm=T)*(N^(-0.2))

編集:私が提供した情報は完全ではありませんでした

kerf = function(x){ return(dnorm(x)) }
ker = function(x,x0,h){
       temp = kerf((x-x0)/h)/h
       return(temp)
       }

AK = function(x,h,N) {
      K = array(0,c(N,N))                 
         for (i in 1:N) {
         for (j in 1:N) {
         K[i,j] = ker(x[i],x[j],h) 
       }}
       return(K) }

KDens 関数を高速化したい場合、どうすればよいでしょうか?

4

1 に答える 1

4

これを試してみてください... 元の長さ 4800 のデータセットの場合、2.5 秒かかります。

KDens2 = function(x,h,N) {
Kx <- outer( x , x , FUN = function(x,y) dnorm( ( x-y ) / h ) / h )
fx <- as.matrix( rowSums( Kx ) / N , ncol = 1 )
return( fx )
}

テスト

set.seed(1)
trywiththis <- rnorm(480)
x = trywiththis
N = length(trywiththis)
h = 1.059*sd(trywiththis , na.rm=T)*(N^(-0.2))

#Produce same result? (posibly not identical because of 'dnorm' function)
all.equal( KDens(x,h,N) , KDens2(x,h,N) )
[1] TRUE

#Rough timing on N=480 length data...
system.time( KDens2(x,h,N) )
#   user  system elapsed 
#   0.01    0.00    0.02 

system.time( KDens(x,h,N) )
#   user  system elapsed 
#    2.7     0.0     2.7 

そしていつN=4800...

system.time( KDens2(x,h,N) )
   user  system elapsed 
   2.33    0.19    2.51
于 2014-04-22T11:02:46.257 に答える