代码之家  ›  专栏  ›  技术社区  ›  Joseph Wood

R中的算术运算速度比整数快。发生什么事?

  •  9
  • Joseph Wood  · 技术社区  · 6 年前

    我正在将一些主要利用数字数据(即双倍)的代码转换为整数,并做了一个快速基准测试,看看我获得了多少效率。

    令我惊讶的是它慢了。。。大约20%。我以为我做错了什么,但最初的代码只是对中等大小的向量进行一些基本的算术运算,所以我知道不是这样的。也许我的环境一团糟?我重新开始,结果一样。。。整数效率较低。

    这开始了一系列的测试和跳进兔子洞。这是我的第一次考试。我们用基R求一百万个元素的和 sum . 注意,对于R版本 3.5.0 计时有点不同,对于v 3.5.1,计时也差不多(仍然不是人们所期望的那样):

    set.seed(123)
    int1e6 <- sample(1:10, 1e6, TRUE)
    dbl1e6 <- runif(1e6, 1, 10)
    
    head(int1e6)
    # [1] 5 3 6 8 6 2
    class(int1e6)
    # [1] "integer"
    
    head(dbl1e6)
    # [1] 5.060628 2.291397 2.992889 5.299649 5.217105 9.769613
    class(dbl1e6)
    #[1] "numeric"
    
    mean(dbl1e6)
    # [1] 5.502034
    mean(int1e6)
    # [1] 5.505185
    
    ## R 3.5.0
    library(microbenchmark)
    microbenchmark(intSum = sum(int1e6), dblSum = sum(dbl1e6), times = 1000)
    Unit: microseconds
      expr      min       lq      mean   median       uq      max neval
    intSum 1033.677 1043.991 1147.9711 1111.438 1200.725 2723.834  1000
    dblSum  817.719  835.486  945.6553  890.529  998.946 2736.024  1000
    
    ## R 3.5.1
    Unit: microseconds
      expr     min       lq      mean   median        uq      max neval
    intSum 836.243 877.7655  966.4443 950.1525  997.9025 2077.257  1000
    dblSum 866.939 904.7945 1015.3445 986.4770 1046.4120 2541.828  1000
    
    class(sum(int1e6))
    # [1] "integer"
    class(sum(dbl1e6))
    #[1] "numeric"
    

    从这里开始,3.5.0版和3.5.1版给出了几乎相同的结果。

    这是我们第一次潜水 rabbit hole . 以及 (见 ?sum ),我们看到了 总和 standardGeneric . 深入挖掘,我们看到它最终 R_execMethod here R_execClosure 被称为next,后面跟着许多不同的可能分支。我认为标准的做法是 eval 接下来,但我不确定。我的猜测是,最终,一个函数被调用 arithimetic.c 但我找不到任何能精确求和的向量。不管怎样,基于我对方法调度和 C 一般来说,我天真的假设是,调用如下所示的函数:

    template <typename T>
    T sum(vector<T> x) {
    
        T mySum = 0;
    
        for (std::size_t i = 0; i < x.size(); ++i)
            mySum += x[i];
    
        return mySum;
    }
    

    我知道在 C类 ,但你明白我的意思。我的信念是,最终,一堆相同类型的元素被添加到相同类型的元素中,并最终返回。在 Rcpp 我们要的是:

    template <typename typeReturn, typename typeRcpp>
    typeReturn sumRcpp(typeRcpp x) {
        typeReturn mySum = 0;
        unsigned long int mySize = x.size();
    
        for (std::size_t i = 0; i < mySize; ++i)
            mySum += x[i];
    
        return mySum;
    }
    
    // [[Rcpp::export]]
    SEXP mySumTest(SEXP Rx) {
        switch(TYPEOF(Rx)) {
            case INTSXP: {
                IntegerVector xInt = as<IntegerVector>(Rx);
                int resInt = sumRcpp<int>(xInt);
                return wrap(resInt);
            }
            case REALSXP: {
                NumericVector xNum = as<NumericVector>(Rx);
                double resDbl = sumRcpp<double>(xNum);
                return wrap(resDbl);
            }
            default: {
                Rcpp::stop("Only integers and numerics are supported");   
            }
        }
    }
    

    microbenchmark(mySumTest(int1e6), mySumTest(dbl1e6))
    Unit: microseconds
                 expr      min       lq      mean    median        uq      max neval
    mySumTest(int1e6)  103.455  160.776  185.2529  180.2505  200.3245  326.950   100
    mySumTest(dbl1e6) 1160.501 1166.032 1278.1622 1233.1575 1347.1660 1644.494   100
    

    二元运算符

    这让我进一步思考。也许只是复杂的事情 标准通用 奇怪的是 . 所以,让我们跳过所有的jazz,直接转到二进制运算符( +, -, *, /, %/% )

    set.seed(321)
    int1e6Two <- sample(1:10, 1e6, TRUE)
    dbl1e6Two <- runif(1e6, 1, 10)
    
    ## addition
    microbenchmark(intPlus = int1e6 + int1e6Two, 
                   dblPlus = dbl1e6 + dbl1e6Two, times = 1000)
    Unit: milliseconds
       expr      min       lq     mean   median       uq      max neval
    intPlus 2.531220 3.214673 3.970903 3.401631 3.668878 82.11871  1000
    dblPlus 1.299004 2.045720 3.074367 2.139489 2.275697 69.89538  1000
    
    ## subtraction
    microbenchmark(intSub = int1e6 - int1e6Two,
                   dblSub = dbl1e6 - dbl1e6Two, times = 1000)
    Unit: milliseconds
      expr      min       lq     mean   median       uq      max neval
    intSub 2.280881 2.985491 3.748759 3.166262 3.379755 79.03561  1000
    dblSub 1.302704 2.107817 3.252457 2.208293 2.382188 70.24451  1000
    
    ## multiplication
    microbenchmark(intMult = int1e6 * int1e6Two, 
                   dblMult = dbl1e6 * dbl1e6Two, times = 1000)
    Unit: milliseconds
       expr      min       lq     mean   median       uq      max neval
    intMult 2.913680 3.573557 4.380174 3.772987 4.077219 74.95485  1000
    dblMult 1.303688 2.020221 3.078500 2.119648 2.299145 10.86589  1000
    
    ## division
    microbenchmark(intDiv = int1e6 %/% int1e6Two,
                   dblDiv = dbl1e6 / dbl1e6Two, times = 1000)
    Unit: milliseconds
      expr      min       lq     mean   median       uq      max neval
    intDiv 2.892297 3.210666 3.720360 3.228242 3.373456 62.12020  1000
    dblDiv 1.228171 1.809902 2.558428 1.842272 1.990067 64.82425  1000
    

    unique(c(class(int1e6 + int1e6Two), class(int1e6 - int1e6Two),
             class(int1e6 * int1e6Two), class(int1e6 %/% int1e6Two)))
    # [1] "integer"
    
    unique(c(class(dbl1e6 + dbl1e6Two), class(dbl1e6 - dbl1e6Two),
             class(dbl1e6 * dbl1e6Two), class(dbl1e6 / dbl1e6Two)))
    # [1] "numeric"
    

    在每种情况下,我们都可以看到,对于数值数据类型,运算速度要快40%-70%。真正奇怪的是,当操作的两个向量相同时,我们会得到更大的差异:

    microbenchmark(intPlus = int1e6 + int1e6, 
                   dblPlus = dbl1e6 + dbl1e6, times = 1000)
    Unit: microseconds
       expr      min       lq     mean   median       uq      max neval
    intPlus 2522.774 3148.464 3894.723 3304.189 3531.310 73354.97  1000
    dblPlus  977.892 1703.865 2710.602 1767.801 1886.648 77738.47  1000
    
    microbenchmark(intSub = int1e6 - int1e6,
                   dblSub = dbl1e6 - dbl1e6, times = 1000)
    Unit: microseconds
      expr      min       lq     mean   median       uq      max neval
    intSub 2236.225 2854.068 3467.062 2994.091 3214.953 11202.06  1000
    dblSub  893.819 1658.032 2789.087 1730.981 1873.899 74034.62  1000
    
    microbenchmark(intMult = int1e6 * int1e6, 
                   dblMult = dbl1e6 * dbl1e6, times = 1000)
    Unit: microseconds
       expr      min       lq     mean   median       uq      max neval
    intMult 2852.285 3476.700 4222.726 3658.599 3926.264 78026.18  1000
    dblMult  973.640 1679.887 2638.551 1754.488 1875.058 10866.52  1000
    
    microbenchmark(intDiv = int1e6 %/% int1e6,
                   dblDiv = dbl1e6 / dbl1e6, times = 1000)
    Unit: microseconds
      expr      min       lq     mean   median       uq      max neval
    intDiv 2879.608 3355.015 4052.564 3531.762 3797.715 11781.39  1000
    dblDiv  945.519 1627.203 2706.435 1701.512 1829.869 72215.51  1000
    
    unique(c(class(int1e6 + int1e6), class(int1e6 - int1e6),
             class(int1e6 * int1e6), class(int1e6 %/% int1e6)))
    # [1] "integer"
    
    unique(c(class(dbl1e6 + dbl1e6), class(dbl1e6 - dbl1e6),
             class(dbl1e6 * dbl1e6), class(dbl1e6 / dbl1e6)))
    # [1] "numeric"
    

    这几乎是100%的增长与每种运营商类型!!!

    在基R中有一个正则for循环如何?

    funInt <- function(v) {
        mySumInt <- 0L
        for (element in v)
            mySumInt <- mySumInt + element
        mySumInt
    }
    
    funDbl <- function(v) {
        mySumDbl <- 0
        for (element in v)
            mySumDbl <- mySumDbl + element
        mySumDbl
    }
    
    microbenchmark(funInt(int1e6), funDbl(dbl1e6))
    Unit: milliseconds
              expr      min       lq     mean   median       uq      max neval
    funInt(int1e6) 25.44143 25.75075 26.81548 26.09486 27.60330 32.29436   100
    funDbl(dbl1e6) 24.48309 24.82219 25.68922 25.13742 26.49816 29.36190   100
    
    class(funInt(int1e6))
    # [1] "integer"
    class(funDbl(dbl1e6))
    # [1] "numeric"
    

    两者之间的差异并不惊人,但仍有人预计整数和的表现将优于双和。我真的不知道该怎么想。

    为什么数字数据类型在基R中的基本算术运算上的性能确切地优于整数数据类型?

    编辑。忘了提这个:

    sessionInfo()
    R version 3.5.1 (2018-07-02)
    Platform: x86_64-apple-darwin15.6.0 (64-bit)
    Running under: macOS High Sierra 10.13.6
    
    1 回复  |  直到 6 年前
        1
  •  9
  •   Ralf Stubner    6 年前

    F.Privé 在评论中的“随机猜测”真的很好!功能 do_arith arithmetic.c . 首先对于标量,我们看到 REALSXP is simple + 使用。为了 INTSXP 有一个 dispatch R_integer_plus ,它确实检查整数溢出:

    static R_INLINE int R_integer_plus(int x, int y, Rboolean *pnaflag)
    {
        if (x == NA_INTEGER || y == NA_INTEGER)
        return NA_INTEGER;
    
        if (((y > 0) && (x > (R_INT_MAX - y))) ||
        ((y < 0) && (x < (R_INT_MIN - y)))) {
        if (pnaflag != NULL)
            *pnaflag = TRUE;
        return NA_INTEGER;
        }
        return x + y;
    }
    

    类似于其他二进制操作。对于向量来说也是类似的。在 integer_binary 有一个 dispatch real_binary 使用标准操作时不进行任何检查。

    我们可以使用以下Rcpp代码看到这一点:

    #include <Rcpp.h>
    // [[Rcpp::plugins(cpp11)]]
    #include <cstdint>
    using namespace Rcpp;
    
    // [[Rcpp::export]]
    IntegerVector sumInt(IntegerVector a, IntegerVector b) {
      IntegerVector result(no_init(a.size()));
      std::transform(a.begin(), a.end(), b.begin(), result.begin(),
                     [] (int32_t x, int32_t y) {return x + y;});
      return result;
    }
    
    // [[Rcpp::export]]
    IntegerVector sumIntOverflow(IntegerVector a, IntegerVector b) {
      IntegerVector result(no_init(a.size()));
      std::transform(a.begin(), a.end(), b.begin(), result.begin(),
                     [] (int32_t x, int32_t y) {
        if (x == NA_INTEGER || y == NA_INTEGER)
          return NA_INTEGER;
        if (((y > 0) && (x > (INT32_MAX - y))) ||
            ((y < 0) && (x < (INT32_MIN - y))))
          return NA_INTEGER;
        return x + y;
      });
      return result;
    }
    
    // [[Rcpp::export]]
    NumericVector sumReal(NumericVector a, NumericVector b) {
      NumericVector result(no_init(a.size()));
      std::transform(a.begin(), a.end(), b.begin(), result.begin(),
                     [] (double x, double y) {return x + y;});
      return result;
    }
    
    /*** R
    set.seed(123)
    int1e6 <- sample(1:10, 1e6, TRUE)
    int1e6two <- sample(1:10, 1e6, TRUE)
    dbl1e6 <- runif(1e6, 1, 10)
    dbl1e6two <- runif(1e6, 1, 10)
    
    microbenchmark::microbenchmark(int1e6 + int1e6two,
                                   sumInt(int1e6, int1e6two),
                                   sumIntOverflow(int1e6, int1e6two),
                                   dbl1e6 + dbl1e6two,
                                   sumReal(dbl1e6, dbl1e6two),
                                   times = 1000)
    */
    

    结果:

    Unit: microseconds
                  expr      min        lq     mean    median       uq       max neval
    int1e6 + int1e6two 1999.698 2046.2025 2232.785 2061.7625 2126.970  5461.816  1000
                sumInt  812.560  846.1215 1128.826  861.9305  892.089 44723.313  1000
        sumIntOverflow 1664.351 1690.2455 1901.472 1702.6100 1760.218  4868.182  1000
    dbl1e6 + dbl1e6two 1444.172 1501.9100 1997.924 1526.0695 1641.103 47277.955  1000
               sumReal 1459.224 1505.2715 1887.869 1530.5995 1675.594  5124.468  1000
    

    将溢出检查引入C++代码中会显著降低性能。尽管没有标准差 . 因此,如果你知道你的整数是“好行为”,你可以通过跳过R的错误检查来直接获得C/C++的性能。这让我想起 another question 得出了类似的结论。由R执行的错误检查可能代价高昂。

    Unit: microseconds
               expr      min       lq     mean    median       uq       max neval
    int1e6 + int1e6 1761.285 2000.720 2191.541 2011.5710 2029.528 47397.029  1000
             sumInt  648.151  761.787 1002.662  767.9885  780.129 46673.632  1000
     sumIntOverflow 1408.109 1647.926 1835.325 1655.6705 1670.495 44958.840  1000
    dbl1e6 + dbl1e6 1081.079 1119.923 1443.582 1137.8360 1173.807 44469.509  1000
            sumReal 1076.791 1118.538 1456.917 1137.2025 1250.850  5141.558  1000
    

    双打(R和C++)都有显著的性能提升。对于整数,性能也有一些提高,但不如对双倍数的提高。