On specific data, write_fst is much slower than write_parquet (+0.9.8 is even slower than 0.9.0)

[arunsrinivasan]

For a specific dataset that we have, I found two issues:

1. To write compressed data (of about 120MB at compress=100 on disk) on approximately 670,000 rows and 359 cols with 3 char cols and rest all integer cols with quite low unique values), fst v0.9.0 takes about 1-2 mins at compress = 100 (1 min on Linux and 2 mins on Windows) and about 4-5x slower with v0.9.8.
2. write_parquet writes the same data in <5s with sizes smaller (by a small margin) than fst at compress=100.

I have tried to generate a synthetic data to reproduce. So the stats above won't perfectly match, but the trend is noticeable. Here are the timings/file sizes on the synthetic data. The code to generate the synthetic data is provided below. Note that fst performs better on this synthetic data of same dimensions which is still quite bad.

All of these were tested with 8 threads for both FST and PARQUET.
All packages were compiled from source with -O3 flag.
R version is 4.4.2 (although I think it should be reproducible on any R >= 4.0)
Note that for Parquet, I enable dictionary encoding only on character/factor cols, as can be seen in saveParquet function below.
For compress = 100 on this synthetic data, v0.9.8 is 3.3x slower than v0.9.0 (real case is 4-5x slower), but file compression is quite bad still. In real case the difference between fst and parquet is about 7MB with parquet being lesser, 120 vs 113MB).

Timings

     format                    version    compression user_time_sec system_time_sec elapsed_time_sec file_size_mb
 1:     fst                  fst 0.9.0  compress = 50          0.63            1.03             0.67       649.00
 2:     fst                  fst 0.9.0  compress = 90         55.98            0.29             8.36       307.42 
 3:     fst                  fst 0.9.0  compress = 95        106.24            0.49            15.92       264.52
 4:     fst                  fst 0.9.0  compress = 99        119.19            0.18            18.24       230.49
 5:     fst                  fst 0.9.0 compress = 100        145.39            0.21            22.99       220.02
 6:     fst fst 0.9.8 + fstcore 0.10.0  compress = 50          0.41            1.30             0.61       649.00
 7:     fst fst 0.9.8 + fstcore 0.10.0  compress = 90        184.62            0.31            27.66       304.97
 8:     fst fst 0.9.8 + fstcore 0.10.0  compress = 95        224.83            0.52            33.30       262.81
 9:     fst fst 0.9.8 + fstcore 0.10.0  compress = 99        293.43            0.30            46.25       228.54
10:     fst fst 0.9.8 + fstcore 0.10.0 compress = 100        447.81            0.17            76.03       218.01
11: parquet             arrow 19.0.1.1   zstd level 2          3.86            0.05             3.91       106.32

Data generation

require(data.table)

make_synthetic_data <- function(n = 669610L, seed = 1L) {
  set.seed(seed)

  DT <- data.table(.row_id = seq_len(n))

  add_na_runs <- function(x, p_na) {
    if (p_na <= 0) return(x)

    block <- sample(50:2000, 1L)
    starts <- seq(1L, length(x), by = block)

    na_blocks <- sample(
      starts,
      size = floor(length(starts) * p_na),
      replace = FALSE
    )

    for (s in na_blocks) {
      e <- min(length(x), s + block - 1L)
      x[s:e] <- NA_integer_
    }

    x
  }

  make_run_col <- function(n, k, min_val, max_val, avg_run, p_na = 0, p_zero = 0) {
    n_runs <- ceiling(n / avg_run)

    vals <- as.integer(round(runif(k, min_val, max_val)))
    run_vals <- sample(vals, n_runs, replace = TRUE)

    lens <- pmax(1L, rpois(n_runs, lambda = avg_run))
    x <- rep.int(run_vals, lens)
    x <- x[seq_len(n)]

    if (p_zero > 0) {
      idx <- sample.int(n, floor(n * p_zero))
      x[idx] <- 0L
    }

    add_na_runs(as.integer(x), p_na)
  }

  # 3 character columns, with locality
  dates <- as.character(seq.Date(as.Date("2020-01-01"),
                                 as.Date("2026-08-13"),
                                 by = "day"))

  DT[, symbol := rep(
    sample(sprintf("SYM%05d", 1:12000), ceiling(n / 50), replace = TRUE),
    each = 50
  )[1:n]]

  DT[, date := rep(
    sample(dates, ceiling(n / 300), replace = TRUE),
    each = 300
  )[1:n]]

  DT[, bucket := rep(
    sample(LETTERS[1:7], ceiling(n / 500), replace = TRUE),
    each = 500
  )[1:n]]

  col_id <- 0L

  # 7 constant integer columns
  for (j in seq_len(7L)) {
    col_id <- col_id + 1L
    DT[[sprintf("int_const_%03d", col_id)]] <-
      rep.int(sample(c(0L, 1000000L, -5000L, -87000L), 1L), n)
  }

  # 35 binary / tiny-cardinality integer columns
  for (j in seq_len(35L)) {
    col_id <- col_id + 1L
    DT[[sprintf("int_tiny_%03d", col_id)]] <-
      make_run_col(
        n       = n,
        k       = sample(2:8, 1L),
        min_val = -1000000L,
        max_val = 1000000L,
        avg_run = sample(200:3000, 1L),
        p_na    = runif(1, 0, 0.10),
        p_zero  = runif(1, 0.05, 0.50)
      )
  }

  # 90 low-cardinality integer columns
  for (j in seq_len(90L)) {
    col_id <- col_id + 1L
    DT[[sprintf("int_low_%03d", col_id)]] <-
      make_run_col(
        n       = n,
        k       = sample(100:6000, 1L),
        min_val = -5000000L,
        max_val = 5000000L,
        avg_run = sample(50:800, 1L),
        p_na    = runif(1, 0, 0.20),
        p_zero  = runif(1, 0, 0.20)
      )
  }

  # 154 medium-cardinality integer columns, tuned lighter than previous v3
  for (j in seq_len(154L)) {
    col_id <- col_id + 1L
    DT[[sprintf("int_med_%03d", col_id)]] <-
      make_run_col(
        n       = n,
        k       = sample(3000:30000, 1L),
        min_val = -50000000L,
        max_val = 50000000L,
        avg_run = sample(3:80, 1L),
        p_na    = runif(1, 0.02, 0.25),
        p_zero  = runif(1, 0, 0.06)
      )
  }

  # 70 high-cardinality integer columns, tuned lighter than previous v3
  for (j in seq_len(70L)) {
    col_id <- col_id + 1L
    DT[[sprintf("int_high_%03d", col_id)]] <-
      make_run_col(
        n       = n,
        k       = sample(150000:380000, 1L),
        min_val = -200000000L,
        max_val = 200000000L,
        avg_run = sample(1:8, 1L),
        p_na    = runif(1, 0, 0.18),
        p_zero  = runif(1, 0, 0.015)
      )
  }

  DT[, .row_id := NULL]

  stopifnot(nrow(DT) == n)
  stopifnot(ncol(DT) == 359L)
  stopifnot(sum(vapply(DT, is.integer, logical(1))) == 356L)
  stopifnot(sum(vapply(DT, is.character, logical(1))) == 3L)

  DT[]
}

x_syn <- make_synthetic_data()

Benchmarking FST

require(fst) # separately tested on v0.9.0 and v0.9.8+fstcore v0.10.0
threads_fst(8L)

for (cc in c(50L, 90L, 95L, 99L, 100L)) {
  out <- file.path(tempdir(), paste0("synthetic_fst_compress", cc, ".fst"))
  if (file.exists(out)) file.remove(out)

  gc()
  cat("\nwrite_fst synthetic compress =", cc, "\n")
  print(system.time(write_fst(x_syn, out, compress = cc)))
  cat("File MB:", round(file.info(out)$size / 1024^2, 2), "\n")
}

Benchmarking Parquet

require(arrow) # v19.0.1.1
saveParquet <- function(obj, f_name, threads = 8L) {
  arrow::set_cpu_count(threads)

  # Defaults for write_parquet arguments
  compression_algorithm <- "zstd"
  compression_level     <- 2L
  size_mb               <- as.numeric(object.size(obj) / 1024^2)
  target_groups         <- if (size_mb <= 32L)       1L
                            else if (size_mb <= 128L) 4L
                            else if (size_mb <= 512L) 8L
                            else                      12L
  chunk_size            <- as.integer(ceiling(nrow(obj) / target_groups))
  dictionary            <- sapply(obj, \(col) is.character(col) || is.factor(col))

  arrow::write_parquet(
    x                 = arrow::as_arrow_table(obj),
    sink              = f_name,
    compression       = compression_algorithm,
    compression_level = compression_level,
    chunk_size        = chunk_size,
    use_dictionary    = dictionary,
    write_statistics  = TRUE
  )
}

parquet_file <- file.path(tempdir(), "ppp.parquet")
if (file.exists(parquet_file)) file.remove(parquet_file)
system.time(saveParquet(x_syn, parquet_file))
cat("File MB:", round(file.info(parquet_file)$size / 1024^2, 2), "\n")