NEW FASTQ QC statistics and qcstats TSV (milestone 4)

Collect per-file FASTQ statistics (base composition, GC / non-ATCG fractions,
sequence-length range, encoding) as fastq/paired/preprocess run, in registration
order. qcstats({fastq}) serialises them to the transposed TSV produced by
writeOutputTSV and returns a counts file; write copies it out.

This required a faithful port of Haskell's `show :: Double -> String` (fixed-point
for exponents 0..7, scientific otherwise, e.g. 3.896103896103896e-2) in
values::show_double, reusing Rust's shortest-round-tripping digits.

Also register qcstats as an always-loaded builtin (it comes from the stats
module), add the file-backed NGOCounts value, and label preprocess output stats
with the statement line number (preproc.lnoN.*).

With this, tests/preprocess3 passes against the Rust binary: paired preprocess
with gzip output (compared by decompressed content via check.sh) and a
byte-identical output.fqstats.tsv.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>

Author: luispedro

rust-migration.md

@@ -29,17 +29,21 @@
 >   files (encoding-checked, empty singles dropped), `preprocess` processes mates in lockstep
 >   (both survive → pair, one survives → singleton via `keep_singles`), and `write` derives
 >   `pair.1`/`pair.2`/`singles` names (`_formatFQOname`) and concatenates per-slot files.
+>   FASTQ QC statistics are collected as `fastq`/`paired`/`preprocess` run (per-file base
+>   composition, GC/non-ATCG fractions, sequence-length range, encoding), and `qcstats({fastq})`
+>   serialises them to the transposed TSV (mirroring `writeOutputTSV`); `write` of the resulting
+>   counts file copies it out. This needed a faithful port of Haskell's `show :: Double -> String`
+>   (fixed vs. scientific notation, e.g. `3.896103896103896e-2`) in `values::show_double`.
 >   **First functional tests now pass against the Rust binary** (identical output to the
 >   committed `expected.*`): `tests/write_fq`, `tests/write_fq_inline`, `tests/preprocess`
->   (all six outputs: copy, substrim, endstrim, smoothtrim, `avg_quality` filter, and
->   `n_to_zero_quality`), `tests/regression-fqgz` (gz input → uncompressed output), and
->   `tests/preprocess3_empty_singles` (paired preprocess → `pair.1`/`pair.2`). Their
+>   (copy, substrim, endstrim, smoothtrim, `avg_quality` filter, `n_to_zero_quality`),
+>   `tests/regression-fqgz` (gz input), `tests/preprocess3_empty_singles` (paired preprocess),
+>   and `tests/preprocess3` (paired preprocess + gz output + `qcstats` TSV). Their
 >   `ngless "1.1"` headers were bumped to `"1.5"` (these features have only minimum-version
 >   checks, no version-conditional behavior, so Haskell output is unchanged).
 >   Simplifications to lift next: files are read whole rather than streamed (no
->   `FileOrStream`/bounded queues), bzip2/zstd compression, and FASTQ QC statistics
->   (`qcstats`, needed by `tests/preprocess3`). Still not started: module loading and the
->   `map`/`count`/SAM subsystems.
+>   `FileOrStream`/bounded queues), bzip2/zstd compression, and per-position quality
+>   percentiles. Still not started: module loading and the `map`/`count`/SAM subsystems.
 
 ## Context
 

src/fastq.rs

@@ -52,6 +52,72 @@ impl FastQEncoding {
             FastQEncoding::Solexa => 64,
         }
     }
+
+    /// Human-readable name, matching `encodingName` (used in the qcstats TSV).
+    pub fn name(self) -> &'static str {
+        match self {
+            FastQEncoding::Sanger => "Sanger (33 offset)",
+            FastQEncoding::Solexa => "Solexa (64 offset)",
+        }
+    }
+}
+
+/// Per-file FASTQ statistics, mirroring the relevant fields of `FQStatistics` plus the derived
+/// `gcFraction`/`nonATCGFrac`/`nBasepairs`.
+#[derive(Clone, Debug, PartialEq)]
+pub struct FastQStats {
+    pub n_seq: i64,
+    pub min_len: i64,
+    pub max_len: i64,
+    /// Base counts: (A, C, G, T, other), each case-insensitive.
+    pub bp: (i64, i64, i64, i64, i64),
+}
+
+impl FastQStats {
+    pub fn num_basepairs(&self) -> i64 {
+        let (a, c, g, t, o) = self.bp;
+        a + c + g + t + o
+    }
+
+    /// GC fraction over A/C/G/T only (mirrors `gcFraction`).
+    pub fn gc_fraction(&self) -> f64 {
+        let (a, c, g, t, _) = self.bp;
+        (c + g) as f64 / (a + c + g + t) as f64
+    }
+
+    /// Fraction of non-ATCG bases over all bases (mirrors `nonATCGFrac`).
+    pub fn non_atcg_fraction(&self) -> f64 {
+        let (a, c, g, t, o) = self.bp;
+        o as f64 / (a + c + t + g + o) as f64
+    }
+}
+
+/// Compute FASTQ statistics from decoded reads (mirrors `fqStatsC`). For an empty input the
+/// sequence-length range is `(maxBound, 0)` as in the Haskell code.
+pub fn stats_from_reads(reads: &[ShortRead]) -> FastQStats {
+    let (mut a, mut c, mut g, mut t, mut o) = (0i64, 0i64, 0i64, 0i64, 0i64);
+    let mut min_len = i64::MAX;
+    let mut max_len = 0i64;
+    for r in reads {
+        let len = r.sequence.len() as i64;
+        min_len = min_len.min(len);
+        max_len = max_len.max(len);
+        for &b in r.sequence.as_bytes() {
+            match b.to_ascii_lowercase() {
+                b'a' => a += 1,
+                b'c' => c += 1,
+                b'g' => g += 1,
+                b't' => t += 1,
+                _ => o += 1,
+            }
+        }
+    }
+    FastQStats {
+        n_seq: reads.len() as i64,
+        min_len,
+        max_len,
+        bp: (a, c, g, t, o),
+    }
 }
 
 impl ShortRead {
@@ -435,6 +501,25 @@ mod tests {
         );
     }
 
+    #[test]
+    fn stats_from_reads_cases() {
+        // Two reads, case-insensitive base counting, one N (non-ATCG).
+        let reads = vec![
+            ShortRead::new("a", "ACGT", vec![30, 30, 30, 30]),
+            ShortRead::new("b", "acgN", vec![20, 20, 20, 20]),
+        ];
+        let st = stats_from_reads(&reads);
+        assert_eq!(st.n_seq, 2);
+        assert_eq!(st.min_len, 4);
+        assert_eq!(st.max_len, 4);
+        assert_eq!(st.bp, (2, 2, 2, 1, 1)); // A, C, G, T, other
+        assert_eq!(st.num_basepairs(), 8);
+        // GC over ATCG only: (C+G)/(A+C+G+T) = 4/7.
+        assert_eq!(st.gc_fraction(), 4.0 / 7.0);
+        // non-ATCG over all: 1/8.
+        assert_eq!(st.non_atcg_fraction(), 1.0 / 8.0);
+    }
+
     #[test]
     fn read_quality_methods() {
         // avg_quality: mean of the quality values.