Более быстрый способ найти первое ИСТИННОЕ значение в векторе

В одной функции мне очень часто нужно использовать код, например:

which(x==1)[1]
which(x>1)[1]
x[x>10][1]

где x - числовой вектор. summaryRprof() показывает, что я трачу > 80% времени на реляционные операторы. Интересно, есть ли функция, которая выполняет сравнение только до тех пор, пока не будет достигнуто первое значение TRUE, чтобы ускорить мой код. Для цикла медленнее, чем приведенные выше параметры.

Ответы

Ответ 1

База R предоставляет Position и Find для поиска первого индекса и значения, соответственно, для которых предикат возвращает истинное значение. Эти функции более высокого порядка возвращаются сразу после первого удара.

f<-function(x) {
  r<-vector("list",3)
  r[[1]]<-which(x==1)[1]
  r[[2]]<-which(x>1)[1]
  r[[3]]<-x[x>10][1]
  return(r)
}

p<-function(f,b) function(a) f(a,b)
g<-function(x) {
  r<-vector("list",3)
  r[[1]]<-Position(p(`==`,1),x)
  r[[2]]<-Position(p(`>`,1),x)
  r[[3]]<-Find(p(`>`,10),x)
  return(r)
}

Относительная производительность сильно зависит от вероятности нахождения раннего раннего значения относительно стоимости предиката и накладных расходов Position/Find.

library(microbenchmark)
set.seed(1)
x<-sample(1:100,1e5,replace=TRUE)
microbenchmark(f(x),g(x))

Unit: microseconds
 expr      min        lq     mean    median        uq      max neval cld
 f(x) 5034.283 5410.1205 6313.861 5798.4780 6948.5675 26735.52   100   b
 g(x)  587.463  650.4795 1013.183  734.6375  950.9845 20285.33   100  a

y<-rep(0,1e5)
microbenchmark(f(y),g(y))

Unit: milliseconds
 expr        min         lq       mean     median         uq        max neval cld
 f(y)   3.470179   3.604831   3.791592   3.718752   3.866952   4.831073   100  a 
 g(y) 131.250981 133.687454 137.199230 134.846369 136.193307 177.082128   100   b

Ответ 2

Я не знаю о чистом R-способе сделать это, поэтому я написал функцию C > , чтобы сделать это для пакета Quantstrat. Эта функция была написана с определенной целью, поэтому она не такая общая, как хотелось бы. Например, вы можете заметить, что он работает только с реальными/двойными/числовыми данными, поэтому обязательно вызовите Data до вызова функции .firstCross.

#include <R.h>
#include <Rinternals.h>

SEXP firstCross(SEXP x, SEXP th, SEXP rel, SEXP start)
{
    int i, int_rel, int_start;
    double *real_x=NULL, real_th;

    if(ncols(x) > 1)
        error("only univariate data allowed");

    /* this currently only works for real x and th arguments
     * support for other types may be added later */
    real_th = asReal(th);
    int_rel = asInteger(rel);
    int_start = asInteger(start)-1;

    switch(int_rel) {
        case 1:  /* >  */
            real_x = REAL(x);
            for(i=int_start; i<nrows(x); i++)
                if(real_x[i] >  real_th)
                    return(ScalarInteger(i+1));
            break;
        case 2:  /* <  */
            real_x = REAL(x);
            for(i=int_start; i<nrows(x); i++)
                if(real_x[i] <  real_th)
                    return(ScalarInteger(i+1));
            break;
        case 3:  /* == */
            real_x = REAL(x);
            for(i=int_start; i<nrows(x); i++)
                if(real_x[i] == real_th)
                    return(ScalarInteger(i+1));
            break;
        case 4:  /* >= */
            real_x = REAL(x);
            for(i=int_start; i<nrows(x); i++)
                if(real_x[i] >= real_th)
                    return(ScalarInteger(i+1));
            break;
        case 5:  /* <= */
            real_x = REAL(x);
            for(i=int_start; i<nrows(x); i++)
                if(real_x[i] <= real_th)
                    return(ScalarInteger(i+1));
            break;
        default:
            error("unsupported relationship operator");
  }
  /* return number of observations if relationship is never TRUE */
  return(ScalarInteger(nrows(x)));
}

И здесь функция R, которая ее вызывает:

.firstCross <- function(Data, threshold=0, relationship, start=1) {
    rel <- switch(relationship[1],
            '>'    =  ,
            'gt'   = 1,
            '<'    =  ,
            'lt'   = 2,
            '=='   =  ,
            'eq'   = 3,
            '>='   =  ,
            'gte'  =  ,
            'gteq' =  ,
            'ge'   = 4,
            '<='   =  ,
            'lte'  =  ,
            'lteq' =  ,
            'le'   = 5)
    .Call('firstCross', Data, threshold, rel, start)
}

Некоторые тесты, просто для удовольствия.

> library(quantstrat)
> library(microbenchmark)
> firstCross <- quantstrat:::.firstCross
> set.seed(21)
> x <- rnorm(1e6)
> microbenchmark(which(x > 3)[1], firstCross(x,3,">"), times=10)
Unit: microseconds
                  expr      min       lq    median       uq      max neval
       which(x > 3)[1] 9482.081 9578.072 9597.3870 9690.448 9820.176    10
 firstCross(x, 3, ">")   11.370   11.675   31.9135   34.443   38.614    10
> which(x>3)[1]
[1] 919
> firstCross(x,3,">")
[1] 919

Обратите внимание, что firstCross даст большее относительное ускорение, чем больше Data (потому что R-реляционные операторы должны закончить сравнение всего вектора).

> x <- rnorm(1e7)
> microbenchmark(which(x > 3)[1], firstCross(x,3,">"), times=10)
Unit: microseconds
                  expr      min        lq    median        uq        max neval
       which(x > 3)[1] 94536.21 94851.944 95799.857 96154.756 113962.794    10
 firstCross(x, 3, ">")     5.08     5.507    25.845    32.164     34.183    10
> which(x>3)[1]
[1] 97
> firstCross(x,3,">")
[1] 97

... и не будет заметно быстрее, если первое значение TRUE находится рядом с концом вектора.

> microbenchmark(which(x==last(x))[1], firstCross(x,last(x),"eq"),times=10)
Unit: milliseconds
                         expr      min       lq   median       uq       max neval
       which(x == last(x))[1] 92.56311 93.85415 94.38338 98.18422 106.35253    10
 firstCross(x, last(x), "eq") 86.55415 86.70980 86.98269 88.32168  92.97403    10
> which(x==last(x))[1]
[1] 10000000
> firstCross(x,last(x),"eq")
[1] 10000000

Ответ 3

Хороший вопрос и ответ... просто добавить any() не быстрее, чем which() или match(), но оба быстрее, чем [], которые, как я думаю, могут создать большой вектор, бесполезный T, F. Поэтому я догадываюсь, что No..short ответа выше.

    v=rep('A', 10e6)
    v[5e6]='B'
    v[10e6]='B'

    microbenchmark(which(v=='B')[1])
    Unit: milliseconds
                   expr      min       lq   median       uq      max neval
     which(v == "B")[1] 332.3788 337.6718 344.4076 347.1194 503.4022   100

    microbenchmark(any(v=='B'))
    Unit: milliseconds
              expr      min      lq   median       uq      max neval
     any(v == "B") 334.4466 335.114 335.6714 347.5474 356.0261   100

    microbenchmark(v[v=='B'][1])
    Unit: milliseconds
               expr      min       lq  median       uq      max neval
     v[v == "B"][1] 601.5923 605.3331 609.191 612.0689 707.1409   100

    microbenchmark(match("B", v))
    Unit: milliseconds
               expr      min       lq   median       uq      max neval
    match("B", v) 339.2872 344.7648 350.5444 359.6746 915.6446   100

Любые другие идеи?