6

R には大きな data.table があります。すべての行について、同様の値の x1 (+/- 許容範囲、tol) を持つ行をカウントしたいと考えています。adply を使用してこれを機能させることはできますが、遅すぎます。data.table が適しているように思えます-実際、計算の一部にすでに data.table を使用しています。

これを完全に data.table で行う方法はありますか? 次に例を示します。

library(data.table)
library(plyr)
my.df = data.table(x1 = 1:1000,
                   x2 = 4:1003)
tol = 3
adply(my.df, 1, function(df) my.df[x1 > (df$x1 - tol) & x1 < (df$x1 + tol), .N])

結果：

        x1   x2 V1
   1:    1    4  3
   2:    2    5  4
   3:    3    6  5
   4:    4    7  5
   5:    5    8  5
  ---             
 996:  996  999  5
 997:  997 1000  5
 998:  998 1001  5
 999:  999 1002  4
1000: 1000 1003  3

アップデート：

これは、実際のデータに少し近いサンプル データセットです。

set.seed(10)
x = seq(1,100000000,100000)
x = x + sample(1:50000, length(x), replace=T)
x2 = x + sample(1:50000, length(x), replace=T)
my.df = data.table(x1 = x,
                   x2 = x2)
setkey(my.df,x1)
tol = 100000

og = function(my.df) {
  adply(my.df, 1, function(df) my.df[x1 > (df$x1 - tol) & x1 < (df$x1 + tol), .N])
}

microbenchmark(r_ed <- ed(copy(my.df)),
               r_ar <- ar(copy(my.df)),
               r_og <- og(copy(my.df)),
               times = 1)

Unit: milliseconds
                    expr         min          lq      median          uq         max neval
 r_ed <- ed(copy(my.df))    8.553137    8.553137    8.553137    8.553137    8.553137     1
 r_ar <- ar(copy(my.df))   10.229438   10.229438   10.229438   10.229438   10.229438     1
 r_og <- og(copy(my.df)) 1424.472844 1424.472844 1424.472844 1424.472844 1424.472844     1

明らかに、@eddi と @Arun の両方からのソリューションは、私のソリューションよりもはるかに高速です。あとは、ロールを理解しようとするだけです。

4

4 に答える 4

9

（この特定の問題に対する）より高速な解決策については、@eddiの回答を参照してください。x1が整数でない場合にも機能します。

探しているアルゴリズムはInterval Treeです。そして、このタスクを実行するIRangesと呼ばれるバイオコンダクター パッケージがあります。それを打ち負かすのは難しいです。

require(IRanges)
require(data.table)
my.df[, res := countOverlaps(IRanges(my.df$x1, width=1), 
           IRanges(my.df$x1-tol+1, my.df$x1+tol-1))]

いくつかの説明：

コードを分解すると、次の 3 行で記述できます。

ir1 <- IRanges(my.df$x1, width=1)
ir2 <- IRanges(my.df$x1-tol+1, my.df$x1+tol-1)
cnt <- countOverlaps(ir1, ir2)

基本的に行うことは、2 つの「範囲」を作成することです (タイプir1しir2て、それらがどのようになっているかを確認するだけです)。次に、エントリごとに、ir1それらがいくつ重複しているかを尋ねますir2(これが「間隔ツリー」の部分です)。そして、これは非常に効率的です。デフォルトでは、暗黙的にtypeへの引数は「type = any」です。countOverlaps必要に応じて、他のタイプを調べることができます。とても便利です。また、関連するのはfindOverlaps関数です。

注: ir1 の幅 = 1 であるこの特定のケースでは、より高速なソリューションが存在する可能性があります (実際にはあります。@eddi を参照してください)。しかし、幅が可変または > 1 である問題については、これが最速のはずです。

ベンチマーク:

ag <- function(my.df) my.df[, res := sum(abs(my.df$x1-x1) < tol), by=x1]
ro <- function(my.df) {
            my.df[,res:= { y = my.df$x1
            sum(y > (x1 - tol) & y < (x1 + tol))
            }, by=x1]
      }
ar <- function(my.df) {
           my.df[, res := countOverlaps(IRanges(my.df$x1, width=1), 
            IRanges(my.df$x1-tol+1, my.df$x1+tol-1))]
      }


require(microbenchmark)
microbenchmark(r1 <- ag(copy(my.df)), r2 <- ro(copy(my.df)), 
               r3 <- ar(copy(my.df)), times=100)

Unit: milliseconds
                  expr      min       lq   median       uq       max neval
 r1 <- ag(copy(my.df)) 33.15940 39.63531 41.61555 44.56616 208.99067   100
 r2 <- ro(copy(my.df)) 69.35311 76.66642 80.23917 84.67419 344.82031   100
 r3 <- ar(copy(my.df)) 11.22027 12.14113 13.21196 14.72830  48.61417   100 <~~~

identical(r1, r2) # TRUE
identical(r1, r3) # TRUE
于 2013-08-08T12:04:48.430 に答える
4

これがより速いdata.table解決策です。のローリング マージ機能を使用するという考え方ですがdata.table、その前に、データを少し変更して、列x1を整数ではなく数値にする必要があります。これは、OP が厳密な不等式を使用しており、ローリング ジョインを使用するには、許容範囲を少し減らして浮動小数点数にする必要があるためです。

my.df[, x1 := as.numeric(x1)]

# set the key to x1 for the merges and to sort
# (note, if data already sorted can make this step instantaneous using setattr)
setkey(my.df, x1)

# and now we're going to do two rolling merges, one with the upper bound
# and one with lower, then get the index of the match and subtract the ends
# (+1, to get the count)
my.df[, res := my.df[J(x1 + tol - 1e-6), list(ind = .I), roll = Inf]$ind -
               my.df[J(x1 - tol + 1e-6), list(ind = .I), roll = -Inf]$ind + 1]


# and here's the bench vs @Arun's solution
ed = function(my.df) {
  my.df[, x1 := as.numeric(x1)]
  setkey(my.df, x1)
  my.df[, res := my.df[J(x1 + tol - 1e-6), list(ind = .I), roll = Inf]$ind -
                 my.df[J(x1 - tol + 1e-6), list(ind = .I), roll = -Inf]$ind + 1]
}

microbenchmark(ed(copy(my.df)), ar(copy(my.df)))
#Unit: milliseconds
#            expr       min       lq   median       uq      max neval
# ed(copy(my.df))  7.297928 10.09947 10.87561 11.80083 23.05907   100
# ar(copy(my.df)) 10.825521 15.38151 16.36115 18.15350 21.98761   100

注: Arun と Matthew の両方が指摘したように、x1が整数の場合、数値に変換してから少量を減算する必要はなく、上記の代わりにtol使用できます。tol - 1Ltol - 1e-6

于 2013-08-08T15:01:53.637 に答える
2

という事実を利用して

 abs(x-y) < tol ~    y-tol <= x <= y+ tol

パフォーマンスを 2 倍向上させることができます。

## wrap codes in 2 function for benchmarking
library(data.table)
library(plyr)
my.df = data.table(x1 = 1:1000,
                   x2 = 4:1003)
tol = 3
ag <- function()
my.df[, res := sum(abs(my.df$x1-x1) < tol), by=x1]
ro <- function()
  my.df[,res:= { y = my.df$x1
          sum(y > (x1 - tol) & y < (x1 + tol))
          }, by=x1]
## check equal results
identical(ag(),ro())
TRUE
library(microbenchmark)
## benchmarks 
microbenchmark(ag(),
               ro(),times=1)

Unit: milliseconds
 expr      min       lq   median       uq      max neval
 ag() 32.75638 32.75638 32.75638 32.75638 32.75638     1
 ro() 63.50043 63.50043 63.50043 63.50043 63.50043     1
于 2013-08-08T08:05:31.303 に答える
2

純粋な data.table ソリューションは次のとおりです。

my.df[, res:=sum(my.df$x1 > (x1 - tol) & my.df$x1 < (x1 + tol)), by=x1]

my.df <- adply(my.df, 1, 
           function(df) my.df[x1 > (df$x1 - tol) & x1 < (df$x1 + tol), .N])

identical(my.df[,res],my.df[,V1])
#[1] TRUE

ただし、一意の が多数ある場合、これはまだ比較的遅くなりますx1。結局のところ、膨大な数の比較を行う必要があり、それを回避する方法は今のところ思いつきません。

于 2013-08-08T07:36:42.453 に答える