rhistogramcut

Getting same output as cut() using speedier hist() or findInterval()?


I read this article http://www.r-bloggers.com/comparing-hist-and-cut-r-functions/ and tested hist() to be faster than cut() by ~4 times on my PC. My script loops through cut() many times so the time-saving would be significant. I thus tried to switch to the speedier function but am having difficulties getting the exact output as per cut().

From the sample code below:

data <- rnorm(10, mean=0, sd=1)  #generate data
my_breaks <- seq(-6, 6, by=1)  #create a vector that specifies my break points
cut(data, breaks=my_breaks)

I wish to get a vector comprising levels that each element of data is assigned to using my breakpoints, i.e. the exact output of cut:

 [1] (1,2]   (-1,0]  (0,1]   (1,2]   (0,1]   (-1,0]  (-1,0]  (0,1]   (-2,-1] (0,1]  
Levels: (-6,-5] (-5,-4] (-4,-3] (-3,-2] (-2,-1] (-1,0] (0,1] (1,2] (2,3] (3,4] (4,5] (5,6]
> 

My question: How do I use elements of the hist() output (i.e. breaks, counts, density, mids, etc) or findInterval to reach my objective?

Separately, I found an example from https://stackoverflow.com/questions/12379128/r-switch-statement-on-comparisons using findInterval, but this requires me to create the interval labels beforehand, which is not what I want.

Any help would be appreciated. Thanks in advance!


Solution

  • Here is an implementation based on your findInterval suggestion which is 5-6 times faster than classical cut:

    cut2 <- function(x, breaks) {
      labels <- paste0("(",  breaks[-length(breaks)], ",", breaks[-1L], "]")
      return(factor(labels[findInterval(x, breaks)], levels=labels))
    }
    
    library(microbenchmark)
    
    set.seed(1)
    data <- rnorm(1e4, mean=0, sd=1)
    
    microbenchmark(cut.default(data, my_breaks), cut2(data, my_breaks))
    
    # Unit: microseconds
    #                         expr      min        lq    median        uq      max neval
    # cut.default(data, my_breaks) 3011.932 3031.1705 3046.5245 3075.3085 4119.147   100
    #        cut2(data, my_breaks)  453.761  459.8045  464.0755  469.4605 1462.020   100
    
    identical(cut(data, my_breaks), cut2(data, my_breaks))
    # TRUE