3

約 300 個のファイルがあり、それぞれに 1000 個の時系列の実現が含まれています (各ファイルに約 76 MB)。

300000 の実現の完全なセットから各時間ステップで分位数 (0.05、0.50、0.95) を計算したいと思います。

ファイルが大きくなりすぎるため、実現を 1 つのファイルにまとめることはできません。

これを行う最も効率的な方法は何ですか?

各マトリックスはモデルを実行することによって生成されますが、乱数を含むサンプルを次に示します。

x <- matrix(rexp(10000000, rate=.1), nrow=1000)
4

1 に答える 1

4

少なくとも 3 つのオプションがあります。

  1. それは完全なセットからのものでなければなりませんか？ここでは、10% のサンプルが非常に適切な近似値となるはずです。
  2. 300k 要素はそれほど大きなベクトルではありませんが、300k x 100+ 列の行列は大きいです。行列全体ではなく、必要な列だけをメモリにプルします (必要に応じて、すべての列で繰り返すことができます)。
  3. 適切な範囲で開始できるように、場合によってはより小さなサンプルと組み合わせて、順番に実行してください。5 パーセンタイルについては、現在の推測を上回っている項目と下回っている項目の数を知る必要があります。次のようなものです：
    1. 1% のサンプルを取り、その 5 パーセンタイルを見つけます。正確な 5 パーセンタイルがその範囲内にあることを確認できるように、許容範囲を上下に移動します。
    2. チャンクで行列を読み取ります。チャンクごとに、範囲を超える観測値と範囲を下回る観測値の数を数えます。次に、範囲内にあるすべての観測を保持します。
    3. 最後のチャンクを読み取ると、3 つの情報 (上に数え、下に数え、内部の観測ベクトル) が得られます。変位値を取得する 1 つの方法は、ベクトル全体を並べ替えて n 番目の観測を見つけることです。上記の情報を使用して、範囲内の観測を並べ替え、(n-count_below) 番目を見つけます。

編集：（3）の例。

私はチャンピオンのアルゴリズム設計者ではなく、誰かがこのためのより優れたアルゴリズムをほぼ確実に設計したことに注意してください。また、この実装は特に効率的ではありません。速度が重要な場合は、Rcpp を検討するか、さらに最適化された R を検討してください。大量のリストを作成してから値を抽出するのはあまりスマートではありませんが、この方法でプロトタイプを作成するのは簡単だったので、それを採用しました。

library(plyr)

set.seed(1)

# -- Configuration -- #
desiredQuantile <- .25

# -- Generate sample data -- #

# Use some algorithm (sampling, iteration, or something else to come up with a range you're sure the true value lies within)
guessedrange <- c( .2, .3 )
# Group the observations to correspond to the OP's files
dat <- data.frame( group = rep( seq(100), each=100 ), value = runif(10000) )

# -- Apply the algorithm -- #

# Count the number above/below and return the values within the range, by group
res <- dlply( dat, .( group ), function( x, guessedrange ) {
  above <- x$value > guessedrange[2]
  below <- x$value < guessedrange[1]
  list(
    aboveCount  = sum( above ),
    belowCount = sum( below ),
    withinValues = x$value[ !above & !below ]
  )
}, guessedrange = guessedrange )
# Exract the count of values below and the values within the range
belowCount <- sum( sapply( res, function(x) x$belowCount ) )
belowCount
withinValues <- do.call( c, sapply( res, function(x) x$withinValues ) )
str(withinValues)
# Count up until we find the within value we want
desiredQuantileCount <- floor( desiredQuantile * nrow(dat) ) #! Should fix this so it averages when there's a tie
sort(withinValues)[ desiredQuantileCount - belowCount + 1 ]
# Compare to exact value
quantile( dat$value, desiredQuantile )

最終的に、値は正確なバージョンから少しずれています。1 つまたはいくつかの同じくらいばかげた説明によって、私はずれているのではないかと思いますが、根本的な何かが欠けているのかもしれません。

于 2014-02-24T11:01:10.153 に答える