Compression algorithm analysis on real genomic data (HG00096, chr1, 4.7GB SAM)

[swetank18]

Motivation

Following the suggestion to benchmark compression on large real-world genomic data,
I ran a full comparison of all supported compression algorithms on chromosome 1 of
HG00096 from the 1000 Genomes Project (4.7GB SAM, 11.2M reads).

Environment

• OS: Ubuntu 24.04
• ROOT: 6.36.04
• GCC: 13.3
• Dataset: HG00096.bam (1000 Genomes Project, chr1 extracted)
• SAM size: 4,776 MB | BAM size: 15 GB | Reads: 11,277,558

Conversion Benchmark

Algorithm	Output Size	Compression Ratio vs SAM	Compression Ratio vs BAM	Time
BAM (baseline)	15,000 MB	0.32x	1x	—
SAM (baseline)	4,776 MB	1x	0.32x	—
ZLIB (101)	1,021 MB	4.68x	14.7x	144s
LZ4 (404)	1,284 MB	3.72x	11.7x	213s
LZMA (505)	909 MB	5.26x	16.9x	192s
ZSTD (606)	923 MB	5.18x	16.6x	171s

Region Query Benchmark (1:1000000-2000000, 45,296 records)

Format	Time	Records
samtools (BAM)	0.07s	45,491
RAMTools ZLIB	0.86s	45,296
RAMTools LZ4	0.90s	45,296
RAMTools LZMA	0.92s	45,296
RAMTools ZSTD	0.92s	45,296

Key Findings

1. ZSTD is the best overall algorithm

ZSTD achieves near-LZMA compression (5.18x vs 5.26x) while being 12% faster.
It should be considered as the default instead of LZMA.

2. LZ4 underperforms expectations

LZ4 is typically the fastest algorithm but here it is both the slowest AND has
the worst compression ratio. This is likely because LZ4 is optimized for
decompression speed, not compression speed. For write-once genomic archives,
ZSTD or LZMA are better choices.

3. Region query speed gap vs samtools

RAMTools takes ~0.9s vs samtools 0.07s for the same region. This is a 12x gap.
Possible causes:

• RNTuple page size is not optimized for genomic access patterns
• The index granularity (one entry per ~1000 reads) may be too coarse
• samtools uses a highly optimized BAI index with bin-based lookup

4. Record count discrepancy

RAMTools returns 45,296 records vs samtools 45,491 for the same region.
This 195-read difference needs investigation — likely secondary/supplementary
reads are being filtered differently.

5. Quality score compression opportunity

Quality scores are the largest field in a SAM file (~30% of data).
Current implementation stores them as encoded strings. Opportunities:

• Illumina binning — reduce quality alphabet from 40 values to 8 bins
• Reference-based quality compression — similar to CRAM
• Quality score dropping — for applications that don't need per-base quality

6. Field-specific compression

Different fields have very different compressibility:

• seq — 2-bit encoding already applied, good
• qual — high entropy, hardest to compress, biggest opportunity
• cigar — highly repetitive, compresses well
• qname — structured patterns, could use delta encoding
• pos — sorted integers, ideal for delta encoding

Proposed Optimizations

1. Switch default to ZSTD — better speed/ratio tradeoff than LZMA
2. Per-field compression — apply different algorithms per column
   (RNTuple supports this natively via field-level compression settings)
3. Delta encoding for pos — sorted positions compress much better with delta
4. Investigate index granularity — finer index could close the query speed gap
5. Quality score binning — implement Illumina 8-bin compression as default option
6. Investigate record count discrepancy — align filtering logic with samtools

References

• CRAM compression paper: Hsi-Yang Fritz et al. (2011)
• ZSTD vs LZMA benchmarks: Facebook engineering blog
• RNTuple per-field compression: ROOT documentation

[swetank18]

@vgvassilev compression benchmark on HG00096 chr1 (4.7GB SAM, 11.2M reads). ZSTD appears to be a better default than LZMA. Also found a 12x query speed gap vs samtools and a record count discrepancy worth investigating.

[swetank18]

Updated Benchmarks (with latest develop branch + CRAM comparison)

Dataset: HG00096 chr1 | SAM: 4,776 MB | Reads: 11,277,558

Storage Comparison

Format	Size	Ratio vs SAM	Notes
BAM	15,000 MB	0.32x	Standard binary format
SAM	4,776 MB	1x	Uncompressed baseline
CRAM	1,800 MB	2.65x	Reference-based compression
RAM LZ4	1,284 MB	3.72x
RAM ZLIB	1,021 MB	4.68x
RAM ZSTD	923 MB	5.18x	Best speed/ratio tradeoff
RAM LZMA	909 MB	5.26x	Best compression

RAM (LZMA) achieves 2x better compression than CRAM and 16.5x better than BAM.

Region Query Speed (1:1000000-2000000, ~45k reads)

Format	Time	Records returned
BAM (samtools)	0.07s	45,491
CRAM (samtools)	0.21s	45,491
RAM ZLIB	0.86s	45,296
RAM LZ4	0.90s	45,296
RAM LZMA	0.92s	45,296
RAM ZSTD	0.92s	45,296

Key Observations

1. RAM compression beats CRAM by 2x — significant achievement
2. Query speed gap — RAM is 4x slower than CRAM, 13x slower than BAM
3. Record discrepancy — RAM consistently returns 195 fewer reads than BAM/CRAM for the same region. This needs investigation — likely secondary/supplementary reads being filtered differently.
4. ZSTD vs LZMA — ZSTD is only 1.5% worse compression but 12% faster conversion. Worth considering as default.
5. LZ4 underperforms — worse compression AND slower than ZLIB for this workload.

[swetank18]

Record Count Discrepancy — Resolved

The 195-read difference between RAMTools (45,296) and samtools default (45,491)
is explained by RAMTools applying FLAG_FILTER = UNMAPPED | SECONDARY | SUPPLEMENTARY
during region queries.

Verification:

• Total reads in region: 45,491
• Reads filtered by FLAG_FILTER: 195
• Expected RAMTools output: 45,296 (matches exactly)

[AdityaPandeyCN]

Are you really sure RNTuple ROOT provided more compression than CRAM for this particular sample.

[swetank18]

Are you really sure RNTuple ROOT provided more compression than CRAM for this particular sample.

You raise a valid point. Looking back at my benchmark, the CRAM was generated from the BAM file (which is already compressed) while RAM was generated from the uncompressed SAM. This may not be a fair comparison.

Also HG00096 is a low-coverage sample (~4-6x)
CRAM's reference-based compression typically performs much better on high coverage data where read sequences cluster heavily around reference positions.

I'll redo the comparison more carefully:

1. Generate CRAM directly from SAM (same source as RAM)
2. Verify CRAM is using full reference-based compression
3. Test on a higher-coverage sample if available

I'll update issue #34 with corrected results.

[swetank18]

Update: Corrected CRAM Comparison

Following @AdityaPandeyCN's question, I redid the CRAM comparison generating it directly from SAM (same source as RAM) using samtools 1.19.2 with the hs37d5 reference.

Result: CRAM from SAM = 1.8GB identical to CRAM from BAM

This confirms the original comparison was valid. CRAM size is independent of whether the source is BAM or SAM since it uses reference-based compression against hs37d5.fa either way.

Corrected Comparison Table (all generated from same SAM source)

Format	Size	Notes
BAM	15 GB	Original download
CRAM (from SAM, ref-based)	1.8 GB	samtools 1.19.2, hs37d5 reference
RAM LZMA	909 MB	2x better than CRAM
RAM ZSTD	923 MB	2x better than CRAM

Sample Context

HG00096 is a low-coverage (~4-6x) whole genome sample with 45M reads, 99.64% mapped. CRAM is known to perform better on high-coverage data where reference-based compression has more to exploit — the result may not generalize.

The 2x compression advantage of RAM over CRAM appears to hold for this low-coverage sample. Whether this generalizes to high-coverage data remains an open question and would be worth investigating with a high-coverage sample (e.g. GIAB HG001/NA12878 at 30x+).