AI4S-YB · xuzhougeng · May 5, 2026 · May 4, 2026 · May 4, 2026 · May 4, 2026
diff --git a/.gitignore b/.gitignore
@@ -18,6 +18,9 @@ target
 # Test output directories (generated by cargo test / manual runs)
 tests/output*/
 
+# Git worktrees
+.worktrees/
+
 # Benchmark directory (results now live in https://github.com/seqeralabs/RustQC-benchmarks)
 benchmark/
 

diff --git a/CHANGELOG.md b/CHANGELOG.md
@@ -1,5 +1,29 @@
 # RustQC Changelog
 
+## [Unreleased]
+
+### Added
+
+- `rustqc atac` subcommand for single-pass ATAC-seq QC, replicating
+  ATACseqQC 1.36.0's `bamQC`, `fragSizeDist`, `TSSEscore`, `NFRscore`,
+  `PTscore`, and `estimateLibComplexity`. Numerical fidelity targets the
+  upstream R package (byte-identical / 1e-6 / 1e-3 per metric — see the
+  documentation site).
+- Tn5 +4/−5 shift and fixed-interval (NFR / mono / di / tri) length split
+  helpers — opt-in via `--emit-shifted-bam` / `--emit-split-bams` (file
+  writing reserved for a follow-up release).
+- Shared BAM/preseq infrastructure lifted out of `src/rna/` to crate root
+  (`src/bam_flags.rs`, `src/bam_io.rs`, `src/cpp_rng.rs`, `src/preseq.rs`).
+  No behavior change to existing `rustqc rna` outputs.
+
+### TODO
+
+- Factor-footprinting metrics (`factorFootprints`, `vPlot`, etc.) —
+  requires BSgenome / PWM motif DBs.
+- Random-forest split branch (only the fixed-interval split lands here).
+- `EmitWriters::open` actual file writing for `--emit-shifted-bam` /
+  `--emit-split-bams`.
+
 ## [Version 0.3.0](https://github.com/AI4S-YB/RustQC/releases/tag/v0.3.0) - 2026-04-22
 
 First release of the AI4S-YB fork. Focus: Windows support and a pure-Rust

diff --git a/Cargo.lock b/Cargo.lock
diff --git a/Cargo.toml b/Cargo.toml
@@ -77,6 +77,9 @@ plotters = { version = "0.3", features = ["fontconfig-dlopen"] }
 [target.'cfg(windows)'.dependencies]
 plotters = "0.3"
 
+[dev-dependencies]
+tempfile = "3"
+
 [build-dependencies]
 cc = "1"
 

diff --git a/README.md b/README.md
@@ -6,7 +6,7 @@
 </picture>
 </h1>
 
-<h4 align="center">Fast genomics quality control tools for sequencing data, written in Rust.</h4>
+<h4 align="center">Fast genomics quality control tools for RNA-seq and ATAC-seq data, written in Rust.</h4>
 
 <p align="center">
   <a href="https://github.com/seqeralabs/RustQC/actions/workflows/ci.yml"><img src="https://github.com/seqeralabs/RustQC/actions/workflows/ci.yml/badge.svg" alt="CI"></a>
@@ -36,7 +36,12 @@
 
 It currently includes:
 
-- `rustqc rna` is a single-command RNA-Seq QC tool that runs all QC analyses in one pass. Designed to slot into the [nf-core/rnaseq pipeline](https://nf-co.re/rnaseq/), but works anywhere:
+- `rustqc rna` — single-command RNA-Seq QC, designed to slot into the [nf-core/rnaseq pipeline](https://nf-co.re/rnaseq/), but works anywhere:
+- `rustqc atac` — single-command ATAC-seq QC, designed to slot into the [nf-core/atacseq pipeline](https://nf-co.re/atacseq/), but works anywhere.
+
+### `rustqc rna`
+
+`rustqc rna` is a single-command RNA-Seq QC tool that runs all QC analyses in one pass:
 
 | Tool                | Reimplements                                                                            | Description                                                           |
 | ------------------- | --------------------------------------------------------------------------------------- | --------------------------------------------------------------------- |
@@ -56,7 +61,22 @@ It currently includes:
 | idxstats            | [samtools](http://www.htslib.org/) `idxstats`                                           | Per-chromosome read counts                                            |
 | stats               | [samtools](http://www.htslib.org/) `stats`                                              | Full samtools stats output including all histogram sections           |
 
-All outputs are format- and numerically identical to the upstream tools, and compatible with [MultiQC](https://multiqc.info/) for reporting.
+All RNA outputs are format- and numerically identical to the upstream tools, and compatible with [MultiQC](https://multiqc.info/) for reporting.
+
+### `rustqc atac`
+
+`rustqc atac` is a single-command ATAC-seq QC tool that runs all QC analyses in one pass, replicating [ATACseqQC](https://bioconductor.org/packages/ATACseqQC/) 1.36.0:
+
+| Module | Reimplements | Description |
+| ------------------- | --------------------------------------------------------------------------------------- | --------------------------------------------------------------------- |
+| bamQC | [ATACseqQC](https://bioconductor.org/packages/ATACseqQC/) `bamQC` | Mapping rates, NRF, PBC1/2, MAPQ histogram |
+| fragSizeDist | [ATACseqQC](https://bioconductor.org/packages/ATACseqQC/) `fragSizeDist` | Fragment-length distribution + plot |
+| TSSEscore | [ATACseqQC](https://bioconductor.org/packages/ATACseqQC/) `TSSEscore` | TSS enrichment score (sliding window + loess smoothing) |
+| NFRscore | [ATACseqQC](https://bioconductor.org/packages/ATACseqQC/) `NFRscore` | Nucleosome-free region score per TSS |
+| PTscore | [ATACseqQC](https://bioconductor.org/packages/ATACseqQC/) `PTscore` | Promoter / transcript-body score per TSS |
+| estimateLibComplexity | [ATACseqQC](https://bioconductor.org/packages/ATACseqQC/) `estimateLibComplexity` | Library complexity extrapolation (preseq) |
+| Tn5 shift | [ATACseqQC](https://bioconductor.org/packages/ATACseqQC/) `shiftGAlignmentsList` | +4/−5 Tn5 shift |
+| Split (fixed) | [ATACseqQC](https://bioconductor.org/packages/ATACseqQC/) `splitGAlignmentsByCut` (fixed-interval) | NFR/mono/di/tri length split |
 
 ## Quick start
 

diff --git a/docs/astro.config.mjs b/docs/astro.config.mjs
@@ -71,6 +71,14 @@ export default defineConfig({
             { label: "Samtools", slug: "rna/samtools" },
           ],
         },
+        {
+          label: "ATAC",
+          items: [
+            { label: "Overview", slug: "atac/index" },
+            { label: "CLI Reference", slug: "atac/cli" },
+            { label: "Numerical Fidelity", slug: "atac/numerical-fidelity" },
+          ],
+        },
         {
           label: "About",
           items: [

diff --git a/docs/src/content/docs/atac/cli.mdx b/docs/src/content/docs/atac/cli.mdx
@@ -0,0 +1,82 @@
+---
+title: CLI Reference
+description: Command-line interface reference for rustqc atac.
+---
+
+## Synopsis
+
+```
+rustqc atac [OPTIONS] --gtf <GTF> <INPUT>...
+```
+
+`<INPUT>` is one or more paired-end BAM/SAM/CRAM files. Single-end inputs are
+rejected at startup.
+
+## Full help output
+
+```
+ATAC-Seq QC — single-pass analysis of paired-end BAM/SAM/CRAM files.
+
+Runs bamQC, fragSizeDist, TSSEscore, NFRscore, PTscore, and library complexity
+in one pass. Requires a GTF annotation and paired-end input. Optionally emits
+Tn5-shifted and length-split BAMs.
+```
+
+## Options
+
+### Input / Output
+
+| Flag | Env var | Default | Description |
+|------|---------|---------|-------------|
+| `-g, --gtf <GTF>` | `RUSTQC_GTF` | — | GTF gene annotation (plain or `.gz`); TSS coords source |
+| `-r, --reference <FASTA>` | `RUSTQC_REFERENCE` | — | Reference FASTA (required for CRAM) |
+| `-o, --outdir <OUTDIR>` | `RUSTQC_OUTDIR` | `.` | Output directory |
+| `--sample-name <NAME>` | `RUSTQC_SAMPLE_NAME` | derived from filename | Override sample name |
+| `--flat-output` | `RUSTQC_FLAT_OUTPUT` | false | Write outputs to a flat directory (no subdirs) |
+| `-c, --config <CONFIG>` | — | — | YAML configuration file |
+| `-j, --json-summary [<PATH>]` | `RUSTQC_JSON_SUMMARY` | — | JSON summary path (use `-` for stdout) |
+
+### ATAC-specific
+
+| Flag | Env var | Default | Description |
+|------|---------|---------|-------------|
+| `--mito-chrom <NAME>` | `RUSTQC_MITO_CHROM` | auto-detect `^chrM$\|^MT$\|^Mito$` | Mitochondrial chromosome name |
+| `--emit-shifted-bam` | `RUSTQC_EMIT_SHIFTED_BAM` | false | Emit +4/−5 Tn5-shifted BAM |
+| `--emit-split-bams` | `RUSTQC_EMIT_SPLIT_BAMS` | false | Emit NFR/mono/di/tri BAMs (fixed intervals) |
+| `--tsse-flank <N>` | `RUSTQC_TSSE_FLANK` | `1000` | TSSEscore flank (bp) |
+
+### General
+
+| Flag | Env var | Default | Description |
+|------|---------|---------|-------------|
+| `-t, --threads <THREADS>` | `RUSTQC_THREADS` | `1` | Number of threads |
+| `-Q, --mapq <MAPQ_CUT>` | `RUSTQC_MAPQ` | `30` | MAPQ cutoff |
+| `-q, --quiet` | `RUSTQC_QUIET` | false | Suppress output except warnings/errors |
+| `-v, --verbose` | — | false | Show additional detail |
+
+## Examples
+
+```bash
+# Basic run
+rustqc atac sample.markdup.bam --gtf genes.gtf --outdir results/
+
+# With Tn5 shift + length-split BAMs (Phase 14 reserved — not yet wired to disk)
+rustqc atac sample.markdup.bam --gtf genes.gtf --outdir results/ \
+    --emit-shifted-bam --emit-split-bams
+
+# Custom mito chromosome name
+rustqc atac sample.markdup.bam --gtf genes.gtf --mito-chrom Mito --outdir results/
+
+# With sample-name override and JSON summary to stdout
+rustqc atac sample.markdup.bam --gtf genes.gtf --outdir results/ \
+    --sample-name MyExperiment -j -
+```
+
+## Notes on `--emit-shifted-bam` / `--emit-split-bams`
+
+These flags are accepted by the CLI and the Tn5 shift / fixed-interval split
+logic runs in memory, but **file writing is reserved for a follow-up release**.
+The internal `EmitWriters::open` function is a stub in this version — the
+shifted/split BAMs are not written to disk. The flags are present so pipeline
+configurations can be written now and the output will appear automatically when
+file writing is wired up.
diff --git a/docs/src/content/docs/atac/index.mdx b/docs/src/content/docs/atac/index.mdx
@@ -0,0 +1,58 @@
+---
+title: ATAC-Seq QC
+description: Single-pass ATAC-seq quality control and Tn5 preprocessing in Rust.
+---
+
+`rustqc atac` is a single-command ATAC-seq QC tool that runs all QC analyses
+in one pass. Designed to slot into the [nf-core/atacseq pipeline][nfc] but works
+anywhere.
+
+[nfc]: https://nf-co.re/atacseq/
+
+## Modules
+
+| Module | Reimplements | Description |
+|--------|--------------|-------------|
+| bamQC | [ATACseqQC][a] `bamQC` | Mapping rates, NRF, PBC1/2, MAPQ histogram |
+| fragSizeDist | [ATACseqQC][a] `fragSizeDist` | Fragment-length distribution + plot |
+| TSSEscore | [ATACseqQC][a] `TSSEscore` | TSS enrichment score (sliding window + loess smoothing) |
+| NFRscore | [ATACseqQC][a] `NFRscore` | Nucleosome-free region score per TSS |
+| PTscore | [ATACseqQC][a] `PTscore` | Promoter / transcript-body score per TSS |
+| estimateLibComplexity | [ATACseqQC][a] `estimateLibComplexity` | Library complexity extrapolation (preseq) |
+| Tn5 shift | [ATACseqQC][a] `shiftGAlignmentsList` | +4/−5 Tn5 shift |
+| Split (fixed) | [ATACseqQC][a] `splitGAlignmentsByCut` (fixed-interval branch) | NFR/mono/di/tri split |
+
+[a]: https://bioconductor.org/packages/ATACseqQC/
+
+## Usage
+
+```bash
+rustqc atac sample.markdup.bam --gtf genes.gtf --outdir results/
+```
+
+ATAC-seq input must be paired-end. Single-end inputs are rejected at startup.
+
+## Outputs
+
+```
+outdir/
+  bamqc/<sample>.bamqc.tsv               # rates, NRF, PBC1/2
+  bamqc/<sample>.mapq.tsv                # MAPQ histogram
+  fragsize/<sample>.fragsize.{tsv,svg}   # fragment-length distribution
+  tsse/<sample>.tsse.{tsv,svg}           # TSS enrichment
+  nfr/<sample>.nfr.tsv                   # NFR score per TSS
+  pt/<sample>.pt.tsv                     # PT score per TSS
+  lib_complexity/<sample>.libcomplexity.{tsv,svg}
+  <sample>.atac.summary.json
+```
+
+## Scope and TODOs
+
+This release covers the core ATAC-seq QC metrics plus Tn5 +4/−5 shift and fixed-interval
+split (NFR < 100, mono 180–247, di 315–473, tri 558–615 bp). The following ATACseqQC
+features are **not yet implemented**:
+
+- Factor footprinting (`factorFootprints`, `footprintsScanner`, `vPlot`,
+  `enrichedFragments`, `distanceDyad`) — depends on BSgenome / PWM motif DBs.
+- Random-forest split branch — only the fixed-interval branch is implemented.
+- Multi-BAM batching — current driver requires exactly one BAM per invocation.
diff --git a/docs/src/content/docs/atac/numerical-fidelity.mdx b/docs/src/content/docs/atac/numerical-fidelity.mdx
@@ -0,0 +1,42 @@
+---
+title: Numerical Fidelity
+description: Validation strategy and per-metric tolerances against ATACseqQC 1.36.0.
+---
+
+`rustqc atac`'s outputs are validated against the upstream ATACseqQC R package
+(version 1.36.0, Bioconductor 3.23). The reference outputs are generated by the
+script at `tests/atac/golden/run_r_reference.R` and committed as fixtures
+(integration tests skip gracefully when goldens are absent — CI does not run R).
+
+## Per-metric acceptance bars
+
+| Metric | Bar |
+|--------|-----|
+| `bamQC` (rates / NRF / PBC1 / PBC2 / MAPQ histogram) | byte-identical (`abs diff <= 1e-12`) |
+| `fragSizeDist` | byte-identical histogram |
+| `NFRscore` (per TSS) | abs diff ≤ 1e-6 |
+| `PTscore` (per TSS) | abs diff ≤ 1e-6 |
+| `TSSEscore` (post-loess scalar) | abs diff ≤ 1e-3 |
+
+## Loess port
+
+ATACseqQC's `TSSEscore` calls R's `stats::loess.smooth(family = "gaussian")`. We
+port a minimal locally-weighted polynomial fit (tricube weights, degree 2,
+span 2/3, no robust reweighting) sufficient to land within `1e-3` of R on the
+GL1–GL4 fixtures. The pre-loess sliding-window vector (length 20) is byte-
+identical to R; only the loess smoothing introduces the 1e-3 tolerance.
+
+## Fixtures
+
+The GL1–GL4 BAMs at `tests/data/atac/` are the same files shipped inside
+ATACseqQC's `inst/extdata/` and align to chr1 of hg19. The `GL_tss.gtf` in the
+same directory is a synthetic 14-transcript subset on chr1 covering the read
+ranges of GL1 (genome-wide) and GL2-4 (~565kb–996kb dense cluster).
+
+## What this release does **not** validate
+
+- Factor footprinting (out of scope for v1).
+- Random-forest split branch (only fixed-interval split is implemented).
+- Tn5 shift / split BAM round-trip vs ATACseqQC outputs — `EmitWriters::open`
+  is a stub in this release; integration tests are reserved (`if true { return; }`)
+  until file writing lands.