Finding a length n needle in a haystack

Problem: Find index of a short vector (needle) inside a larger vector (haystack)

What’s a good way to find the location of a needle in a haystack?

haystack <- sample(0:12, size = 2000, replace = TRUE)
needle   <- c(2L, 10L, 8L) 

• The haystack may be large.
• The needle is definitely a vector
• Only the index for the first occurrence of needle within haystack is needed
• Assuming that there is always at least 1 occurrence

Solution 1: for loop

• Naively test every position in a for loop
• Advantage of this method over a vectorised approach, is that it can terminate early i.e. as soon as a value is found

forloop_find <- function(needle, haystack) {
  n <- length(needle) - 1L
  for (i in seq(haystack)) {
    if (identical(haystack[i:(i+n)], needle)) {
      return(i)
    }
  }
}

forloop_find(needle, haystack)

## [1] 931

Solution 2: vapply

vapply_find <- function(needle, haystack) {
  n <- length(needle) - 1L
  which(vapply(seq(haystack), function(i) {identical(haystack[i:(i+n)], needle)}, FUN.VALUE = logical(1L)))[1L]
}

vapply_find(needle, haystack)

## [1] 931

Solution 3: zoo::rollapply

rollapply_find <- function(needle, haystack) {
  which(zoo::rollapply(haystack, length(needle), identical, needle))[1L]
}

rollapply_find(needle, haystack)

## [1] 931

Solution 4: stacked lead() finding

• for each item in the needle
  • shift the haystack and check for the needle element
  • stack the shifted haystacks and sum them
  • check where sum == length(needle)

lead_find <- function(needle, haystack) {
  v  <- haystack == needle[1]
  for (i in seq(2, length(needle))) {
    v <- v + (lead(haystack, i-1L) == needle[i])
  }
  
  which(v == length(needle))[1L]
}

lead_find(needle, haystack)

## [1] 931

Solution 5: shift

• data.table::shift is a vectorised version of dplyr::lead()
• Adapted from user “Jaap” on stackoverflow

shift_find <- function(needle, haystack) {
  shifted_haystack <- data.table::shift(haystack, type='lead', 0:(length(needle)-1))
  
  v <- Map('==', shifted_haystack, needle)
  v <- Reduce('+', v)
  
  which(v == length(needle))[1]
}

shift_find(needle, haystack)

## [1] 931

Solution 5: Rcpp for loop

• update 2018-04-04
  • Shaved a few % points off the runtime by using the value of j instead of an explicit match count.

Rcpp::cppFunction('int rcpp_find(NumericVector needle, NumericVector haystack) {
    int needle_length   = needle.size();
    int haystack_length = haystack.size();

    for (int i = 0; i < (haystack_length - needle_length); i++ ) {
        int j;
        for (j = 0; j < needle_length; j++ ) {
                if (needle[j] != haystack[i + j]) {
                    break;
                }
        }

        if (j == needle_length) {
            return(i+1);
        }
    }

    return(0);

    }')


rcpp_find(needle, haystack)

## [1] 931

Solution 6: pile find (by dmytro)

• dmytro on twitter offered this one.

pile_find <- function(needle, haystack) { 
  index <- seq_along(haystack) 
  for (i in seq_along(needle)) { 
    pile <- haystack[index] 
    index <- index[pile==needle[i]] + 1L
  } 
  (index-length(needle))[1]
} 

pile_find(needle, haystack)

## [1] 931

Solution 7: Sieved find (by carroll_jono)

• carroll_jono on twitter offered this one.
• This is the fastest of the pure-R methods.

sieved_find <- function(needle, haystack) { 
  sieved <- which(haystack == needle[1L]) 
  for (i in seq.int(1L, length(needle)-1L)) {
    sieved <- sieved[haystack[sieved + i] == needle[i + 1L]]
  }
  sieved[1L] 
}

sieved_find(needle, haystack)

## [1] 931

Solution 8: Coerce to string and use grep

• Convert each integer into a string e.g. c(1, 2, 3) becomes “001.002.003”
• Current tests are with integers less than 1000, so I’m safe to use %03i as a format string.

grep_find <- function(needle, haystack) {
  expanded_width <- 4  # 3 digits plus a "."
  res <- stringr::str_locate_all(paste(sprintf("%03i", haystack), collapse="."), paste(sprintf("%03i", needle), collapse="."))[[1]][,1][[1]]
  
  (res-1)/expanded_width + 1
}

grep_find(needle, haystack)

## [1] 931

Solution 9: Idomatic C++ solution hrbrmstr

Rcpp::cppFunction('int idiomaticrcpp_find(NumericVector needle, NumericVector haystack) {
  NumericVector::iterator it;
  it = std::search(haystack.begin(), haystack.end(), needle.begin(), needle.end());
  int pos = it - haystack.begin() + 1;
  if (pos > haystack.size()) pos = -1;
  return(pos);
}')

idiomaticrcpp_find(needle, haystack)

## [1] 931

Solution 10: Rcpp boyer-moore jimhester

• Boyer-Moore implementation in Rcpp by jimhester

# https://gist.github.com/jimhester/2d02d63d0459d9fdadcf15e08c3afc0a
Rcpp::sourceCpp("cache/2018-04-03-boyer-moore.cpp")
boyermoorecpp_find(needle, haystack)

## [1] 931

needle is in the middle of the haystack

needle is very close to beginning of haystack

The for loop and rcpp are the only methods to benefit from a needle that is thought to be at the beginning of the haystack i.e. they can both return early as soon as they find match.

When the needle is close to the beginning of the haystack, then the searching time is very short, and the for loop is quite competitive with custom Rcpp code.

Conclusion

Summary

• Rcpp is the fastest way. No real difference between my C-style of C++ and idiomatic C++
• Boyer-Moore in Rcpp no faster than a vanilla solution for current benchmark sizes
• For a pure R solution, the sieve_find() method is fastest and very competitive with Rcpp