Reading fathom csv benchmarks

[mhpob]

File from actel issue #164

download.file(
    'https://github.com/user-attachments/files/18086120/NexTrak-R1.801032.2024-10-25.093307.csv',
    destfile = 'NexTrak-R1.801032.2024-10-25.093307.csv'
)

fathom_csv <- 'NexTrak-R1.801032.2024-10-25.093307.csv'

data.table version:

read_fathom_dt <- function(fathom_csv, data_type = "DET"){
  
  # Update: this is needed as grep is not in default Windows
  search_fun <- ifelse(.Platform$OS.type == "windows", 'FINDSTR', "grep")
  
  data.table::fread(
    cmd = paste(
      search_fun, data_type, fathom_csv
    ),
    header = F,
    col.names = {
      as.character(
        data.table::fread(cmd = paste(
          search_fun, paste0(data_type, "_DESC"), fathom_csv
        ),
        header = F)
      )
    }
  )
}

readr version:

read_fathom_readr <- function(fathom_csv, data_type = "DET") {
  
  search_fun <- ifelse(.Platform$OS.type == "windows", 'FINDSTR', "grep")
  
  # cheating here. Could/should be switched from system
  header <- system2(
    search_fun,
    c(paste0(data_type, "_DESC"),
      fathom_csv
    ),
    stdout = T) |> 
    strsplit(',') |> 
    unlist()
  
  select_detections <- function(x, pos) x[x[,1] == data_type,]
  
  readr::read_csv_chunked(
    fathom_csv,
    readr::DataFrameCallback$new(select_detections),
    col_names = header,
    col_types = readr::cols()
  ) |> 
    readr::type_convert() 
}

Drag race:

microbenchmark::microbenchmark(
  readr = read_fathom_readr(fathom_csv),
  dt = read_fathom_dt(fathom_csv),
  times = 20
)

Unit: milliseconds
  expr       min        lq      mean    median        uq       max neval cld
 readr 2311.1762 2474.3908 2670.6241 2616.3300 2717.2186 4093.1096    20  a 
    dt  299.9092  309.8798  320.0476  313.5691  328.8836  361.5658    20   b

Note that data.table can read from HTTP if not on Windows, if desired. readr, as coded cannot.

read_fathom_readr('https://github.com/user-attachments/files/18086120/NexTrak-R1.801032.2024-10-25.093307.csv')
read_fathom_dt('https://github.com/user-attachments/files/18086120/NexTrak-R1.801032.2024-10-25.093307.csv')

[jdpye]

Determining which of these methods is quickest, most efficient against large datasets, and potentially which ones do and do not preserve data integrity, will be the 'work' of this ticket. Winners get to be the way that surimi reads fathom output!

[mhpob]

I'm still trying to see if R ships with some version of grep in Windows.

It looks like data.table::fread(cmd) writes stdout to tempdir(), then reads that in. cmd uses R's system() in UNIX-alike OSs and shell() in Windows, selecting first what is in Sys.getenv("R_SHELL"), then what is in Sys.getenv("COMSPEC"), which usually is cmd.exe.

When I run shell("which grep"), it returns /usr/bin/grep which makes me think it's hitting my WSL2 instance. This is also what is returned when I run which grep directly in cmd.exe. The GH Actions Windows server also has WSL, so I'm not sure if that will help. Going to post in r-support to see what happens.

[mhpob]

Very much thinking that we do need a non-grep workaround for most Windows users:

• https://bsky.app/profile/obrien.page/post/3ltn4ulusnk2j
• https://otn.slack.com/archives/CRH7PGQ77/p1752177308673329

[chrisholbrook]

Very much thinking that we do need a non-grep workaround for most Windows users:

• https://bsky.app/profile/obrien.page/post/3ltn4ulusnk2j
• https://otn.slack.com/archives/CRH7PGQ77/p1752177308673329

Early version of glatos::read_vdat_csv() used similar grep implementation and was abandoned for same reason.

[mhpob]

Looks like findstr was introduced in Windows 2000, so the OS switch to findstr from grep should be adequate to address this on Windows and UNIX-alike. As alluded to above, grep seems to be able to handle files over HTTP and findstr cannot, but @jdpye notes this is an edge case that be ignored at the moment.

[mhpob]

There is some suggestion that the grep/findstr method may be slower than bringing the file into R and messing with it there, but I'm not finding this on initial tests.

read_fathom_dt_all <- function(fathom_csv, data_type = "DET"){
  
  dat_all <- data.table::fread(fathom_csv, header = F, fill = T)
  the_names <- as.character(dat_all[V1 == 'DET_DESC'])
  the_names <- the_names[the_names != '']
  
  dat <- dat_all[V1 == 'DET', 1:length(the_names)]
  data.table::setnames(dat, the_names)
  dat[,lapply(.SD, type.convert, as.is = T)] # note that this does not convert date-times from character
}

> microbenchmark::microbenchmark(
   dt = read_fathom_dt(fathom_csv),
   dat_all = read_fathom_dt_all(fathom_csv),
   times = 20
 )

Unit: milliseconds
    expr      min       lq     mean   median       uq      max neval cld
      dt 282.7064 287.6445 304.7159 301.1494 310.4046 355.3872    20  a 
 dat_all 378.5805 391.8400 431.2868 424.2040 430.9142 680.1297    20   b

[chrisholbrook]

I am really glad to see this out in the open in an issue here. I confronted many of these same questions for read_vdat_csv behind the curtain but am finding it not easy to recall how we landed where we did.

[mhpob]

Note that we tried pipe('grep...') |> readr::read_csv() as an fread(cmd) analog, but it would get ahead of itself and corrupt the imported data (not the original) in a manner that cannot be reproduced. The location of the corrupted data changed after each run.

[chrisholbrook]

There is some suggestion that the grep/findstr method may be slower than bringing the file into R and messing with it there, but I'm not finding this on initial tests.

I'm also finding this not to be the case. I ran your benchmarks with a second detection file (2.6M detections vs 176k) and also added glatos::read_vdat_csv, which is similar to your read_fathom_dt_all(). (I should have done this on your dt_all version, but it came in just after.

  library(glatos)
  
  # second file with 2.6M detection records (2.7M total rows)
  #  vs 157k detection records (178k total rows) for fathom_csv defined above
  fathom_csv2 <- "HR2-180_461742_2022-10-04_140923.csv"

  bm <- microbenchmark::microbenchmark(
    readr = read_fathom_readr(fathom_csv),
    readr2 = read_fathom_readr(fathom_csv2),
    dt = read_fathom_dt(fathom_csv),
    dt2 = read_fathom_dt(fathom_csv2),
    glat = glatos::read_vdat_csv(fathom_csv, record_types = "DET")$DET,
    glat2 = glatos::read_vdat_csv(fathom_csv2, record_types = "DET")$DET,
    times = 3
  )

bm

> Unit: milliseconds
>    expr        min         lq       mean     median        uq       max neval
>   readr  1583.0863  1937.4624  2149.0847  2291.8386  2432.084  2572.329     3
>  readr2 37494.3226 46105.5722 60324.6046 54716.8218 71739.746 88762.669     3
>      dt   311.4566   319.7693   323.6569   328.0820   329.757   331.432     3
>     dt2  4066.5124  4264.0864  7598.6681  4461.6605  9364.746 14267.832     3
>    glat   630.1305   671.6214  1898.7126   713.1123  2533.004  4352.895     3
>   glat2 15876.5372 17383.8898 22126.3439 18891.2425 25251.247 31611.252     3

Your original data.table (read_fathom_dt) outperforms on both files, with the glatos function coming in between readr and dt.

Note that the dt function also scales better... the larger file (~15x larger) took 13x longer to parse than the smaller; vs. 23-26x longer with the readr and glatos methods.

[chrisholbrook]

Have you considered if you want to enforce data types on the result?

For glatos, we decided to strictly enforce data types, which I believe is one reason why read_vdat_csv pulls the entire object into R and then separates by record type. This may be worth thinking about now because it may exclude some candidate methods.

We chose to enforce data types because if you don't, you can get inconsistent types for the same column across different objects (input files). E.g., using your functions above, if you don't have any sensor detections in a given file, you'll end up with a sensor value column of type logical. Another file might end up with sensor type numeric. Such differences can present challenges when combining data across files. In our view, enforcing up front to ensure consistency in the result is worth the cost (longer processing time is not prohibitive). So you might want to think ahead to that part of the workflow before you go too far into the benchmarking functions.

[mhpob]

All really fair points, which I'll leave to @jackVanish to chew on and answer as this project is his baby. Part of the reason we didn't put glatos' function in the comparison is because of the the reasons you outlined: it's safer and goes much farther than a simple reader so it's not really fair to compare it on the grounds of speed.. These are much closer to read.csv for Fathom files and so have all the same pitfalls, like columns that we know are one type being labeled as "logical" when empty

[chrisholbrook]

I'm 100% supportive of a simple-as-possible, read-it-as-it-comes, function for this. Maybe at this stage that's the focus and the check/clean/standardize stuff can be done in another step. E.g., a workflow like:
dt <- read_fathom(fathom_csv)
which might be a data.frame, tibble, or data.table. and then:
dt2 <- as_fathom_detections(dt)
returning an object of class fathom_detections after doing the check/clean/standardize.

And fwiw, my intent is not to endorse glatos::read_vdat_csv. I'm glad we have it but I've always felt there must be a better way (I could never have predicted we would end up doing it the way we did) and an independent effort is best. But I suppose I can't help throw it in the mix for curiosity sake.

[franksmithxyz]

I think that avoiding use of an external tool like grep is the way to go, that way you don't need to deal with these messy OS dependencies.

One approach could be something like

• read the fathom csv file into memory in chunks of N lines at a time using the readLines n= parameter, e.g. do it with 10,000 lines at a time.
  • this helps limit memory usage and will be advantageous when dealing with large HR files
• then for each chunk
  • use an R grep-type function to filter the lines that start with the record type you're looking for
  • pass just those lines to data.table::fread and produce a data table
• after processing all the chunks, merge the chunk-based data tables together.

[mhpob]

@franksmithxyz, readr::read_csv_chunked and data.table::fread are "chunked" by default; very fair point about system dependencies. I had cheated before by using grep/findstr to find the header row; the readr code no longer needs the system dependencies after refactoring to remove the cheat. I also refactored the data.table-only implementation using fread(skip = "DET"), which does a fast scan and only starts reading when it finds the given character. Seems to speed up a lot.

Benchmarking: refactored readr code, glatos (cc @chrisholbrook), various data.table implementations: 1) dt_sys using an OS switch between FINDSTR (Windows) and grep (otherwise); 2) dt_sysflag the same, but using a literal matching flag rather than regex which should theoretically be faster; 3/4) explicitly using grep and grep -F (dt_grep, dt_grepflag); 5/6) explicitly using FINDSTR and FINDSTR /l (dt_findstr, dt_findstrflag); and, finally, 7) all R dt_nodep.

I'm on Windows 11, 32 GB RAM, with Rtools (I think this is where grep is coming from?). readr is using 8 threads (4 logical), data.table is using 4 threads.

1. Pure-R implementation in data.table is fastest, but not significantly different from FINDSTR (dt_sys and dt_findstr) on my Windows machine. While not significant due to some long-running iterations, it's usually ~2x faster.
2. All flavors of pre-processing are not significantly different on my machine, though FINDSTR tended to be faster than grep
3. glatos::rad_vdat_csv data.table-based implementation came in after this (caveats noted above re: type checks, etc., so not really apples:apples)
4. readr implementation came in last

microbenchmark::microbenchmark(
  readr_nodep = readr_no_sysdep(fathom_csv),
  dt_sys = dt_sys(fathom_csv),
  dt_grep = dt_grep(fathom_csv),
  dt_findstr = dt_findstr(fathom_csv),
  dt_sysflag = dt_sysflag(fathom_csv),
  dt_grepflag = dt_grepflag(fathom_csv),
  dt_findstrflag = dt_findstrflag(fathom_csv),
  dt_nodep = dt_nosysdep(fathom_csv),
  glatos = glatos::read_vdat_csv(fathom_csv, record_types = "DET")$DET,
  times = 10
)

Unit: milliseconds
           expr       min        lq      mean    median        uq       max neval  cld
    readr_nodep 1938.3969 1982.5617 2063.3968 2016.1806 2101.9407 2320.5900    10 a   
         dt_sys  284.0313  301.0708  340.3796  319.7423  373.2678  436.4023    10  bc 
        dt_grep  395.6941  401.3058  461.4487  457.6349  508.9658  552.9412    10   c 
     dt_findstr  287.9290  299.8394  331.9545  318.0616  377.0791  404.5935    10  bc 
     dt_sysflag  299.1494  303.1353  353.6239  323.9110  406.8065  489.3555    10   c 
    dt_grepflag  310.4210  401.1186  422.5574  415.5152  433.5766  584.5359    10   c 
 dt_findstrflag  288.0306  298.8895  348.2724  357.3393  382.5329  411.8606    10   c 
       dt_nodep  111.9785  142.2629  217.0439  162.7199  334.9004  465.0917    10  b  
         glatos  651.8871  772.3790  879.3177  875.3200 1046.3406 1064.0245    10    d

data.table

System pre-processing

regex

dt_sys <- function(fathom_csv, data_type = "DET"){
  search_fun <- ifelse(.Platform$OS.type == "windows", 'FINDSTR', "grep")
  
  data.table::fread(
    cmd = paste(
      search_fun, data_type, fathom_csv
    ),
    header = F,
    col.names = {
      as.character(
        data.table::fread(cmd = paste(
          search_fun, paste0(data_type, "_DESC"), fathom_csv
        ),
        header = F)
      )
    }
  )
}
dt_grep <- function(fathom_csv, data_type = "DET"){
  search_fun <- "grep"
  
  data.table::fread(
    cmd = paste(
      search_fun, data_type, fathom_csv
    ),
    header = F,
    col.names = {
      as.character(
        data.table::fread(cmd = paste(
          search_fun, paste0(data_type, "_DESC"), fathom_csv
        ),
        header = F)
      )
    }
  )
}

dt_findstr <- function(fathom_csv, data_type = "DET"){
  search_fun <- "FINDSTR"
  
  data.table::fread(
    cmd = paste(
      search_fun, data_type, fathom_csv
    ),
    header = F,
    col.names = {
      as.character(
        data.table::fread(cmd = paste(
          search_fun, paste0(data_type, "_DESC"), fathom_csv
        ),
        header = F)
      )
    }
  )
}

literal (no regex)

dt_sysflag <- function(fathom_csv, data_type = "DET"){
  search_fun <- ifelse(.Platform$OS.type == "windows", 'FINDSTR /l', "grep -F")
  
  data.table::fread(
    cmd = paste(
      search_fun, data_type, fathom_csv
    ),
    header = F,
    col.names = {
      as.character(
        data.table::fread(cmd = paste(
          search_fun, paste0(data_type, "_DESC"), fathom_csv
        ),
        header = F)
      )
    }
  )
}
dt_grepflag <- function(fathom_csv, data_type = "DET"){
  search_fun <- "grep -F"
  
  data.table::fread(
    cmd = paste(
      search_fun, data_type, fathom_csv
    ),
    header = F,
    col.names = {
      as.character(
        data.table::fread(cmd = paste(
          search_fun, paste0(data_type, "_DESC"), fathom_csv
        ),
        header = F)
      )
    }
  )
}
dt_findstrflag <- function(fathom_csv, data_type = "DET"){
  search_fun <- "FINDSTR /l"
  
  data.table::fread(
    cmd = paste(
      search_fun, data_type, fathom_csv
    ),
    header = F,
    col.names = {
      as.character(
        data.table::fread(cmd = paste(
          search_fun, paste0(data_type, "_DESC"), fathom_csv
        ),
        header = F)
      )
    }
  )
}

no pre-processing

dt_nosysdep <- function(fathom_csv, data_type = "DET"){
  the_names <- as.character(
    data.table::fread(
      fathom_csv,
      skip = paste0(data_type, "_DESC"),
      nrows = 1, header = F
    )
  )
  
  
  dat_all <- data.table::fread(
    fathom_csv,
    skip = paste0(data_type, ","),
    header = F,
    fill = T
  )[V1 == data_type, 1:length(the_names)]
  
  
  data.table::setnames(dat_all, the_names)
  dat_all
}

readr

readr_no_sysdep <- function(fathom_csv, data_type = "DET") {
  header_id <- paste0(data_type, "_DESC")
  
  header <- readr::read_lines(fathom_csv, n_max = 50)
  header <- unlist(strsplit(grep(header_id, header, value = TRUE), ","))
  
  select_data <- function(x, pos){
    x[grepl(
      paste0("^", data_type, "$"),
      eval(parse(text = paste0("x$", header_id)))
    ),]
  }
  
  dat <- readr::read_csv_chunked(
    fathom_csv,
    readr::DataFrameCallback$new(select_data),
    col_names = header,
    show_col_types = FALSE
  )
  suppressMessages(readr::type_convert(dat))
}

[mhpob]

If we stick with the pure-R data.table function, class/type checks/conversions happen before we filter for the desired data type. This means it's likely many numeric/integer data will be coerced into character if one of the other data types has a character in that column. Because of this, we'll have to do a type check.

Luckily, the first three columns will always have the same type: character (data type), POSIXct (device time), and POSIXct (time). This is good, as conversion to/from datetime is the thing that takes the longest in this process. My thought is to skip the first three columns and convert the rest via type.convert. We could go out to the 7th column as offset, correction, model, and serial number will all be the same type, BUT the 4th column is "Time Offset (h)" which will usually be empty (NAs, so logical class). This should be integer if it eve has data, so starting at the 4th column will allow it to be caught.

dt_nosysdep <- function(fathom_csv, data_type = "DET"){
  the_names <- as.character(
    data.table::fread(
      fathom_csv,
      skip = paste0(data_type, "_DESC"),
      nrows = 1, header = F
    )
  )
  
  all_data <- data.table::fread(
    fathom_csv,
    skip = paste0(data_type, ","),
    header = F,
    fill = T
  )[V1 == data_type, 1:length(the_names)]
  
  data.table::setnames(all_data, the_names)
  all_data[, names(.SD) := lapply(.SD, type.convert, as.is=T),
          .SDcols = 4:length(the_names)]
  all_data[]
}

# "Memory Remaining (%)" column would be character without type conversion
dt_nosysdep(fathom_csv, "HEALTH_NEXTRAK")
      HEALTH_NEXTRAK_DESC   Device Time (UTC)                Time Time Offset (h)
                   <char>              <POSc>              <POSc>          <lgcl>
   1:      HEALTH_NEXTRAK 2025-05-14 06:00:00 2025-05-14 06:00:00              NA
   2:      HEALTH_NEXTRAK 2025-05-14 07:00:00 2025-05-14 07:00:00              NA
   3:      HEALTH_NEXTRAK 2025-05-14 08:00:00 2025-05-14 08:00:00              NA
   4:      HEALTH_NEXTRAK 2025-05-14 09:00:00 2025-05-14 09:00:00              NA
   5:      HEALTH_NEXTRAK 2025-05-14 10:00:00 2025-05-14 10:00:00              NA
  ---                                                                            
1180:      HEALTH_NEXTRAK 2025-07-02 09:00:00 2025-07-02 09:02:26              NA
1181:      HEALTH_NEXTRAK 2025-07-02 10:00:00 2025-07-02 10:02:26              NA
1182:      HEALTH_NEXTRAK 2025-07-02 11:00:00 2025-07-02 11:02:26              NA
1183:      HEALTH_NEXTRAK 2025-07-02 12:00:00 2025-07-02 12:02:26              NA
1184:      HEALTH_NEXTRAK 2025-07-02 13:00:00 2025-07-02 13:02:26              NA
      Time Correction (s)      Model Serial Number Memory Remaining (%)
                    <num>     <char>         <int>                <num>
   1:           0.1071032 NexTrak-R1        802463              99.9999
   2:           0.2312409 NexTrak-R1        802463              99.9997
   3:           0.3553786 NexTrak-R1        802463              99.9995
   4:           0.4795162 NexTrak-R1        802463              99.9994
   5:           0.6036539 NexTrak-R1        802463              99.9992
  ---                                                                  
1180:         146.4654149 NexTrak-R1        802463              99.7746
1181:         146.5895526 NexTrak-R1        802463              99.7744
1182:         146.7136902 NexTrak-R1        802463              99.7742
1183:         146.8378279 NexTrak-R1        802463              99.7741
1184:         146.9619656 NexTrak-R1        802463              99.7739

dt_nosysdep(fathom_csv, "DET")
       DET_DESC   Device Time (UTC)                Time Time Offset (h) Time Correction (s)
         <char>              <POSc>              <POSc>          <lgcl>               <num>
    1:      DET 2025-05-14 05:09:55 2025-05-14 05:09:55              NA          0.00348413
    2:      DET 2025-05-14 05:19:57 2025-05-14 05:19:57              NA          0.02424409
    3:      DET 2025-05-14 05:25:43 2025-05-14 05:25:43              NA          0.03620294
    4:      DET 2025-05-14 05:30:57 2025-05-14 05:30:57              NA          0.04700404
    5:      DET 2025-05-14 05:36:39 2025-05-14 05:36:40              NA          0.05882634
   ---                                                                                     
26453:      DET 2025-07-02 12:32:56 2025-07-02 12:35:23              NA        146.90598363
26454:      DET 2025-07-02 12:43:56 2025-07-02 12:46:23              NA        146.92874358
26455:      DET 2025-07-02 12:54:55 2025-07-02 12:57:22              NA        146.95146905
26456:      DET 2025-07-02 13:05:50 2025-07-02 13:08:17              NA        146.97405659
26457:      DET 2025-07-02 13:16:37 2025-07-02 13:19:04              NA        146.99636826
            Model Serial Number Channel Detection Type        Full ID    ID Raw Data
           <char>         <int>   <int>         <char>         <char> <int>   <lgcl>
    1: NexTrak-R1        802463       1       PPM_SELF A69-1601-16952 16952       NA
    2: NexTrak-R1        802463       1       PPM_SELF A69-1601-16952 16952       NA
    3: NexTrak-R1        802463       1            PPM A69-1601-16950 16950       NA
    4: NexTrak-R1        802463       1       PPM_SELF A69-1601-16952 16952       NA
    5: NexTrak-R1        802463       1            PPM A69-1601-16950 16950       NA
   ---                                                                              
26453: NexTrak-R1        802463       1       PPM_SELF A69-1601-16952 16952       NA
26454: NexTrak-R1        802463       1       PPM_SELF A69-1601-16952 16952       NA
26455: NexTrak-R1        802463       1       PPM_SELF A69-1601-16952 16952       NA
26456: NexTrak-R1        802463       1       PPM_SELF A69-1601-16952 16952       NA
26457: NexTrak-R1        802463       1       PPM_SELF A69-1601-16952 16952       NA
       Transmitter Serial Signal Strength (dB) Noise (dB) Gain (dB) Quality Score Station Name
                   <lgcl>                <num>      <num>    <lgcl>        <lgcl>       <lgcl>
    1:                 NA                   NA         NA        NA            NA           NA
    2:                 NA                   NA         NA        NA            NA           NA
    3:                 NA                 83.5       37.0        NA            NA           NA
    4:                 NA                   NA         NA        NA            NA           NA
    5:                 NA                 83.5       36.5        NA            NA           NA
   ---                                                                                        
26453:                 NA                   NA         NA        NA            NA           NA
26454:                 NA                   NA         NA        NA            NA           NA
26455:                 NA                   NA         NA        NA            NA           NA
26456:                 NA                   NA         NA        NA            NA           NA
26457:                 NA                   NA         NA        NA            NA           NA
       Latitude Longitude GPS HDOP
         <lgcl>    <lgcl>   <lgcl>
    1:       NA        NA       NA
    2:       NA        NA       NA
    3:       NA        NA       NA
    4:       NA        NA       NA
    5:       NA        NA       NA
   ---                            
26453:       NA        NA       NA
26454:       NA        NA       NA
26455:       NA        NA       NA
26456:       NA        NA       NA
26457:       NA        NA       NA