14

vec長い (1E8 エントリから始まる)ベクトル があり、それを範囲 にバインドしたいとします[a,b]。確かにvec[vec < a] = aとをコーディングできますvec[vec > b] = bが、これには、データに対する 2 つのパスと、一時的な指標ベクトル (~800MB、2 回) のための大きな RAM 割り当てが必要です。メイン メモリからローカル キャッシュにデータを 1 回だけコピーすれば、より適切に実行できるため、2 つのパスの燃焼時間は短縮されます (メイン メモリへの呼び出しは、キャッシュ ミスと同様に良くありません)。そして、これが複数のスレッドでどれだけ改善できるかは誰にもわかりませんが、貪欲にならないようにしましょう。:)

ベースRまたは見落としているパッケージに優れた実装がありますか、それともこれはRcpp（または私の旧友data.table）の仕事ですか？

4

2 に答える 2

13

単純な C ソリューションは次のとおりです。

library(inline)

fun4 <-
    cfunction(c(x="numeric", a="numeric", b="numeric"), body4,
              language="C")
body4 <- "
    R_len_t len = Rf_length(x);
    SEXP result = Rf_allocVector(REALSXP, len);
    const double aa = REAL(a)[0], bb = REAL(b)[0], *xp = REAL(x);
    double *rp = REAL(result);

    for (int i = 0; i < len; ++i)
        if (xp[i] < aa)
            rp[i] = aa;
        else if (xp[i] > bb)
            rp[i] = bb;
        else
            rp[i] = xp[i];

    return result;
"
fun4 <-
    cfunction(c(x="numeric", a="numeric", b="numeric"), body4,
              language="C")

単純な並列バージョン (Dirk が指摘するように、これはCFLAGS = -fopenmp~/.R/Makevars 内、および openmp をサポートするプラットフォーム / コンパイラー上)

body5 <- "
    R_len_t len = Rf_length(x);
    const double aa = REAL(a)[0], bb = REAL(b)[0], *xp = REAL(x);
    SEXP result = Rf_allocVector(REALSXP, len);
    double *rp = REAL(result);

#pragma omp parallel for
    for (int i = 0; i < len; ++i)
        if (xp[i] < aa)
            rp[i] = aa;
        else if (xp[i] > bb)
            rp[i] = bb;
        else
            rp[i] = xp[i];

    return result;
"
fun5 <-
    cfunction(c(x="numeric", a="numeric", b="numeric"), body5,
              language="C")

そしてベンチマーク

> z <- runif(1e7)
> benchmark(fun1(z,0.25,0.75), fun4(z, .25, .75), fun5(z, .25, .75),
+           replications=10)
                 test replications elapsed  relative user.self sys.self
1 fun1(z, 0.25, 0.75)           10   9.087 14.609325     8.335    0.739
2 fun4(z, 0.25, 0.75)           10   1.505  2.419614     1.305    0.198
3 fun5(z, 0.25, 0.75)           10   0.622  1.000000     2.156    0.320
  user.child sys.child
1          0         0
2          0         0
3          0         0
> identical(res1 <- fun1(z,0.25,0.75), fun4(z,0.25,0.75))
[1] TRUE
> identical(res1, fun5(z, 0.25, 0.75))
[1] TRUE

私のクアッドコアラップトップで。数値入力、エラー チェックなし、NA 処理などを想定しています。

于 2012-05-07T03:56:22.800 に答える
3

pmin始めに：あなたのソリューションと/ソリューションの間に大きな違いはありませんpmax（私はせっかちなので、 n=1e8 ではなく n=1e7 で試してみてください）- pmin/pmaxは実際にはわずかに遅いです。

fun1 <- function(x,a,b) {x[x<a] <- a; x[x>b] <- b; x}
fun2 <- function(x,a,b) pmin(pmax(x,a),b)
library(rbenchmark)
z <- runif(1e7)

benchmark(fun1(z,0.25,0.75),fun2(z,0.25,0.75),rep=50)

                 test replications elapsed relative user.self sys.self
1 fun1(z, 0.25, 0.75)           10  21.607  1.00000     6.556   15.001
2 fun2(z, 0.25, 0.75)           10  23.336  1.08002     5.656   17.605
于 2012-05-06T23:07:48.260 に答える