Venus (启明星)

Venus is a tumor multi-omics pipeline built on oxo-flow, covering the full path from raw FASTQ to a cohort clinical report: somatic SNV/indel calling, copy-number analysis, bulk RNA-seq, single-cell RNA-seq, and an auto-generated methods section.

What it does

Module	Tools	Outputs
QC	fastp, FastQC, MultiQC	trimmed FASTQs, qc/multiqc/multiqc_report.html
QC (WGS align)	samtools flagstat/stats, mosdepth	per-sample qc/alignqc/ metrics + wgs_alignqc_summary.tsv (incl. sex check)
DNA alignment	BWA-MEM2, Picard MarkDuplicates, GATK BQSR	sorted BAM + BAI per sample
Somatic calling	GATK Mutect2 (per-chromosome scatter) + FilterMutectCalls with orientation-bias priors	per-pair VCF (PASS), stats, MAF
Annotation	vep-rs (JSON cache) → VCF + MAF	deliver/{pair}.somatic.pass.vep.vcf.gz, deliver/{pair}.maf
CNV	CNVkit (cnv_method, default wgs; --reference for tumor-only)	per-pair CNV profiles
ASCAT	ASCAT 3.2 + alleleCounter (paired only)	purity/ploidy fit per pair, ascat/cohort_ascat.tsv
ecDNA	AmpliconSuite-pipeline (CNVkit seed intervals → AmpliconArchitect → AmpliconClassifier; tumor BAM only)	ecdna/<pair>/ classification TSVs, report/cohort_ecdna.tsv
Clinical interpretation	bcftools joins (gnomAD AF, ClinVar) + CIViC monthly dump tiering (AMP/ASCO/CAP via VICC mapping) + curated fusion/CNV panels + NMPA/CSCO therapy context; optional OncoKB API (token)	report/clinical/<pair>.{snv_evidence,cnv_evidence,sv_evidence,biomarkers,therapies}.tsv, report/clinical_cohort.tsv
SV	Manta (paired somatic + tumor-only candidates; --exome for WES)	deliver/{pair}.sv.vcf.gz (PASS, paired) or deliver/{pair}.sv.candidate.vcf.gz (unscored, tumor-only)
MSI	MSIsensor-pro (paired only)	deliver/{pair}.msi.tsv
Signatures	in-house SBS96 counter + NNLS vs COSMIC v3.1	report/signatures/{pair}.sbs96.counts.tsv, report/signatures/{pair}.exposures.tsv
Bulk RNA	2-pass STAR, featureCounts	BAM, gene-count matrix, alignment metrics, fusion screening candidates
Single-cell	dnbc4tools (--chemistry {config.scrna_chemistry}); QC summaries in dependency-light Python	count matrix, QC/clustering metrics
Report	venus helper scripts	cohort TSV, TMB, clinical report with methods generated from what actually ran

Analysis modes:

• Paired tumor–normal (pairs.tsv with a control sample) — Mutect2 tumor/normal, CNVkit --normal.
• Tumor-only (empty control column) — Mutect2 without -normal, CNVkit with a prebuilt normal .cnn reference. Gating is per-pair via when = "wildcard.control != ''", so a single cohort may mix both modes. Tumor-only calls still consume the germline resource AND the configured PoN (config.pon, default the bundled 1000G panel — override with a cohort PoN via --arg pon=... once ≥10 normals are available); declare any PoN used in the methods report.
• Tumor-only cohort (every pair has an empty control) — set --arg tumor_only=true. The paired-only cohort rules then redirect to header-only tumor-only companions sharing their output paths: ascat_cohort → ascat_cohort_tumoronly, msi_cohort → msi_cohort_tumoronly (ASCAT and MSIsensor-pro are paired-only analyses; the report renders their sections with explicit "not run" statements rather than empty tables). Without the flag the cohort rules fail on unproducible per-pair inputs (ascat_paired / msi_paired are never instantiated). The clinical module self-gates the same way per pair: clinical_sv_fusions_{paired,tumoronly} and clinical_biomarkers_{paired,tumoronly} consume the scored or candidate SV delivery respectively, and the interpretation labels every tumor-only call "putative somatic" with MSI marked not_run.
• Mixed cohort caveat — some pairs with a control, some without: per-pair rules self-gate via wildcard.control, but the cohort aggregations (ascat_cohort, msi_cohort) fan out over config.pairs_list and still reference the unproducible per-pair inputs of the tumor-only pairs. There is no injected list of paired-only pairs to expand over (engine limitation), so a mixed cohort currently requires running the paired and tumor-only chains as separate cohorts or waiting on engine support for paired-only pair lists.
• WGS vs WES/panel — set target_bed to the capture intervals for WES/panel, plus mutect_target_args, cnv_method, cnv_targets_args, sv_exome_args, and tmb_target_args per the checklist below; leave empty for WGS.

Groups and assay gating

config/groups.tsv maps sample ids to assay modules via three columns: name (group name), samples (comma-separated ids), and assay (a group metadata column that fans out per sample as the wildcard.assay value). Rules are gated per sample instance with when:

• when = "wildcard.assay == 'wgs'" — DNA chain: fastp/fastqc run for DNA and bulk RNA (any non-scRNA assay), but BWA-MEM2/MarkDuplicates/BQSR only for DNA.
• when = "wildcard.assay == 'rnaseq'" — STAR + featureCounts.
• when = "wildcard.assay == 'scrna'" — dnbc4tools counting (raw FASTQs; no trimming).
• Gates may combine config flags: when = "config.run_bqsr == 'true' && wildcard.assay == 'wgs'".

The assay name is arbitrary — any metadata column in groups.tsv becomes a wildcard.<key> available to rules. A rule that mentions {sample} fans out over every group; the when gate then prunes non-matching instances at plan time, so star_align never tries to align a DNA sample. Rules that aggregate with expand_inputs should expand over per-group lists (config.samples_DNA, config.samples_bulkRNA, config.samples_scRNA — auto-injected from groups.tsv) rather than the global config.samples_list, which now mixes assays.

Note: a sample id declared both in groups.tsv and in pairs.tsv (tumor/control) is allowed — the engine warns about duplicate owners but the rule paths here are group/pair-disjoint, so outputs never collide.

Platforms: BGI (DNBSEQ/MGISEQ) and Illumina FASTQs are interchangeable downstream of fastp. Single-cell libraries are counted by dnbc4tools 3.0; the chemistry is pinned per run via the scrna_chemistry config key (default scRNAv2HT; accepted values auto, scRNAv1HT, scRNAv2HT, scRNAv3HT, scRNAv5P — MGI stomics/DNBSEQ kits). 10x Genomics libraries are not covered by the current scRNA module and would need an additional counting tool.

Library strandedness

featureCounts strand mode is configurable via featurecounts_strandedness in venus.oxoflow (default "0" = unstranded; "1" = sense/forward, "2" = antisense/reverse; override per run with --arg featurecounts_strandedness=2). The unstranded default was verified empirically for this cohort rather than assumed — on a 200k read-pair subset of one sample aligned with STAR (93.4% unique, 147k spliced reads), three independent read-geometry checks against GENCODE exon structure all landed at ~50/50:

Check	Same-strand as R1	Opposite	Verdict
R1 orientation vs GTF exon strand	50.4%	49.6%	unstranded
R1 vs STAR intronMotif XS tag	50.3%	49.7%	unstranded
R1 vs splice-donor GT/CT motif from FASTA	49.4%	50.6%	unstranded

Unstranded counting is strand-agnostic and remains correct regardless of the true library type. If your own library is known stranded, set the config key accordingly — and verify first: a stranded library counted with the wrong -s silently halves or zeroes counts. Changing the value invalidates the featurecounts instances and every downstream consumer of the count matrix (merge_counts, neoantigen expression filter) on the next run.

WGS alignment QC and sex check (per sample)

Every DNA sample gets samtools flagstat (mapping rate / duplicates), samtools stats (insert-size distribution / GC) and a mosdepth depth profile (--by 1000000 --thresholds 10,20: 1Mb windows + ≥10x/20x base fractions — no whole-genome per-base file on disk) on the analysis-ready BQSR BAM. scripts/collect_alignqc.py aggregates them into qc/alignqc/wgs_alignqc_summary.tsv (one row per sample), which feeds MultiQC, report/cohort_summary.tsv and the clinical report's WGS QC section.

The summary also carries the sex check (predicted_sex): mosdepth per-chromosome means vs autosome mean — M if X/A<0.75 and Y/A≥0.05, F if X/A≥0.85 and Y/A≤0.02, else ambiguous (constants disclosed in the collector; report renders them too). The cohort has no recorded sex metadata, so treat every call as a coverage-based hypothesis to review, not a sample-sheet correction. Cross-sample contamination is NOT assessed — VerifyBamID2/freemix requires a 1000G genotype reference panel (SVDP bundle) and meaningful cohort density; deferred to P2. The contamination_freemix column carries an explicit not_assessed marker so downstream tables never present a blank as a measurement.

RNA fusion screening (screening-level)

star_align runs with --chimSegmentMin 12 --chimOutType WithinBAM SoftClip Junctions --chimJunctionOverhangMin 12 (STAR-Fusion sensitivity settings, unique-mapping chimeras only) and keeps Chimeric.out.junction per sample. scripts/fusion_screen.py aggregates split reads per junction, recomputes anchor lengths from the p-extended CIGARs, annotates both breakpoints against the GENCODE gene spans and flags hits against the curated actionable_fusion_pairs.tsv. Output: rna/fusion/{sample}.fusion_candidates.tsv → report/fusion_screen_cohort.tsv → clinical report "RNA fusion screening" section.

This is screening-level, not a fusion caller — no spanning-fragment support, no transcript reconstruction, no RT-artifact filtering. Every row carries the screening disclosure; positives require orthogonal RNA confirmation (FISH/RT-PCR/panel). It COMPLEMENTS the Manta BND screen (RNA junctions vs gDNA breakpoints see different evidence classes), it does not replace it. Full reconstruction (STAR-Fusion) needs its own CTAT genome lib index (~tens of GB) — P2.

scRNA ambient/doublet: not assessed

dnbc4tools 3.0 does not report ambient-RNA or doublet metrics, so report/scrna_metrics.tsv carries explicit not_assessed(tool_lacks) markers in the ambient_pct / doublet_pct columns (never blanks or zeros — report consumers can tell "never measured" from "measured 0"). The protected per-sample raw_matrix/ trios are kept precisely so that a P2 pass with scDblFinder/DecontX (scanpy ecosystem) can compute both without re-running the dnbc4tools count.

Repository layout

venus.oxoflow        # main workflow: config + include list + chromosome scatter
rules/*.oxoflow      # rule fragments (qc, align, varcall, varcall_merge, annotation,
                     #   cnv, ascat, sv, msi, signature, rna, scrna, ecdna, clinical,
                     #   neoantigen, report)
scripts/             # python/R helpers (MAF conversion, TMB, SBS96 fit, cohort tables,
                     #   clinical report, …) + clean_intermediates.sh (disk cleanup)
docs/                # method rationale, reference data, cleanup guide +
                     #   artifacts.md (full output inventory & cleanup classes)
envs/<mod>/pixi.toml # pixi environments, one dir per module (env specs in rules point here)
resources/           # vendored data tables (COSMIC v3.1 SBS96 GRCh38 signature matrix)
config/pairs.tsv     # tumor/control pairing (control may be empty for tumor-only)
config/groups.tsv    # sample groups + assay metadata (DNA, bulkRNA, scRNA)

Environments: materialize with bootstrap_envs.sh (issue #26)

The repo ships only envs/<mod>/pixi.toml + pixi.lock (~1MB); the materialized environments (envs/<mod>/.pixi/envs/default, several GB per module) are not part of the repo and must not be committed, copied, or rsync'd around with it. Lifecycle:

# fresh clone / fresh vendored workdir — materialize all module envs:
scripts/bootstrap_envs.sh
# after any pixi.toml / pixi.lock change — idempotent re-sync:
scripts/bootstrap_envs.sh            # or one module: scripts/bootstrap_envs.sh amplicon

Rules:

• pixi.lock is the single source of truth. Never pip install into a module env at rule runtime — that drifts the env from its lock and the engine's env fingerprint then re-invalidates the rule on every run (observed with mosek in the amplicon env, issue #23; mosek now lives in the lock as a [pypi-dependencies] entry).
• Run workdirs need materialized envs to execute; checkouts (the portal copy, mirrors, backups) do not — keep checkouts lock-only and bootstrap only where you run. A checkout that accumulated a materialized .pixi tree should have it deleted (re-bootstrap on demand); the same class of env-drift risk motivated the deletion-by-default policy.
• .gitignore covers envs/*/.pixi/ and the repo-root .pixi/ cache so git status stays clean even where envs were materialized.

Extending venus with a new module

A module is one rule fragment plus one environment; nothing else changes:

1. envs/<mod>/pixi.toml — the tool's conda-forge/bioconda deps. Keep it minimal; one env per module avoids the version conflicts that a shared env accumulates (e.g. Manta needs python 2.7 while delivery filtering needs modern bcftools — hence the separate envs/sv and envs/varcall).
2. rules/<mod>.oxoflow — a [workflow] header + [[rules]]. Conventions the rest of the pipeline relies on:
   • inputs reference upstream deliverables by path (bqsr/{experiment}/…, deliver/{pair_id}.maf) — the engine resolves producers automatically;
   • per-pair deliverables land in deliver/{pair_id}.<mod>.<ext>, per-sample analyses in report/<mod>/{pair_id}.<ext>; run logs go to logs/;
   • gate paired-only rules with when = "wildcard.control != ''" so a tumor-only pair simply skips the instance instead of failing;
   • gate COHORT rules that aggregate paired-only products with when = "config.tumor_only != 'true'", and add a tumor-only companion rule sharing the same output path with the complementary gate when = "config.tumor_only == 'true'" (see ascat_cohort / ascat_cohort_tumoronly in rules/ascat.oxoflow) so downstream consumers keep a hard input to resolve;
   • tunables live as [config] keys in venus.oxoflow (overridable at run time with --arg KEY=VALUE), not as literals inside rules.
3. Wire it up: [[include]] path = "rules/<mod>.oxoflow" in venus.oxoflow plus any new [config] keys. One-off data tables (signature matrices, panels) go in resources/ and are declared as rule inputs so checksums cover them.

WES/panel compatibility checklist

WGS is the default; a WES/panel run needs these config keys set:

Key	Purpose
target_bed	capture intervals — the documented anchor for the interval set; tools receive it via the *_args fragments below
mutect_target_args	e.g. -L /path/regions.bed — restricts Mutect2 to the capture targets
cnv_method	hybrid for WES/panel capture (default wgs)
cnv_targets_args	e.g. --targets /path/regions.bed for CNVkit
cnv_annotate_args	optional, paired mode (reference-building) only: e.g. --annotate /path/gene.bed so the built reference .cnn and downstream .cnr/.cns carry gene names (default "" = gene column "-"). In tumor-only mode batch rejects it (cnvkit ≥0.9 flag conflict) — annotate the reference .cnn itself with scripts/annotate_cnv_reference.py instead
sv_exome_args	set --exome so Manta's off-target depth filters keep real WES events
tmb_target_args	e.g. --target-bed /path/regions.bed so all_somatic_mut_per_mb uses merged capture Mb (the clinical tmb_mut_per_mb uses coding Mb for both WGS and WES)
msi_coverage	20 for WES depth floors — capture-restricted MSI is NOT supported (msisensor-pro has no BED flag)

Operating on a run workdir (vendored workflow copy)

A run workdir (processed-wsx/) contains its own copy of venus.oxoflow, rules/, scripts/, envs/, resources/. oxo-flow binds environment manifests and checkpoint fingerprints to workdir-resolved paths, so:

• run and dry-run from the workdir against its own copy (oxo-flow dry-run venus.oxoflow), never against the checkout elsewhere — otherwise every env manifest resolves to the checkout and everything looks stale;
• after editing the repo, sync changed files into the workdir first (rsync -a rules/ scripts/ envs/ venus.oxoflow <workdir>/);
• a config edit only invalidates rules that reference the changed keys; oxo-flow dry-run venus.oxoflow reports the exact stale/skip/completed split before anything runs (the Summary: N rules line, by contrast, counts total graph instances, not stale ones).

Cleaning intermediates

Rule definitions carry two cleanup markers (see docs/cleanup.md):

• temporary = true on pure-scratch rules (trim FASTQs, per-chunk VCFs, …) — once the engine supports wildcard cleanup, a fully successful run reclaims them automatically;
• protected_output = [...] on delivery/report/ASCAT/CNV/RNA-count/scRNA rules — these artifacts survive oxo-flow clean.

scripts/clean_intermediates.sh                    # preview what a full clean deletes
scripts/clean_intermediates.sh --mode all --apply # backup deliver/report/ascat/msi, then clean

Note the trade-off: clean deletes files but not checkpoint bookkeeping, so a later oxo-flow run lazily regenerates only what downstream rules need. All detection results and deliverables (MSI scores included — rules/msi.oxoflow protects deliver/{pair_id}.msi.tsv and msi/reference.list) are protected_output; the script additionally tars the undeclared evidence-grade side products (per-locus MSI detail, Manta germline VCFs, scRNA per-cell results) before deleting anything.

Every artifact the pipeline produces, its cleanup class, and how each cleanup strategy treats it is catalogued in docs/artifacts.md.

Requirements

• Linux x86-64, oxo-flow CLI ≥ 0.24 (the [engine] disk-pressure guard additionally needs engine commit a456c796+ — PR #845)
• pixi on PATH (workflow environments are resolved per rule into .oxo-flow/pixi/)
• Reference data (see venus.oxoflow [config] for the exact keys): GRCh38 FASTA + BWA-MEM2 index, STAR index + GENCODE GTF, GATK resource bundle (dbsnp/mills/gnomad/PoN), vep-rs binary + JSON cache, single-cell reference for dnbc4tools, CNVkit normal reference for tumor-only CNV.

Quick start

git clone https://github.com/WangLabCSU/oxo-flow-venus.git
cd oxo-flow-venus

# 1. Stage inputs (relative to your run workdir)
#    raw/{sample}_R1.fastq.gz / raw/{sample}_R2.fastq.gz  (WGS + bulk RNA)
#    single-cell FASTQs keep their library naming — see rules/scrna.oxoflow
# 2. Edit config/pairs.tsv and config/groups.tsv for your cohort
#    (groups.tsv assigns each sample id an assay: wgs / rnaseq / scrna)
# 3. Point the absolute paths in venus.oxoflow [config] at your references

# Validate, then run
oxo-flow validate
oxo-flow run --background

# Reports land in report/ ; delivery files in deliver/

Every config key can be overridden per run without editing the file:

oxo-flow run --arg target_bed=/data/exome.bed \
  --arg mutect_target_args="-L /data/exome.bed" --arg cnv_method=hybrid \
  --arg cnv_targets_args="--targets /data/exome.bed" --arg sv_exome_args="--exome" \
  --arg tmb_target_args="--target-bed /data/exome.bed" --arg msi_coverage=20

Note: run_bqsr=false is NOT a runnable end-state yet — every downstream caller consumes bqsr/{sample}/{sample}.bqsr.bam and no companion rule publishes the un-recalibrated BAM under that path, so the flag would leave the whole somatic chain without a producer (see venus.oxoflow).

Outputs of note

• deliver/{pair_id}.somatic.pass.vcf.gz — PASS somatic variants (unannotated)
• deliver/{pair_id}.somatic.pass.vep.vcf.gz — VEP-annotated VCF
• deliver/{pair_id}.maf — MAF (2.4.1 subset) for downstream TMB/mutational-signature tools
• deliver/{pair_id}.sv.vcf.gz — Manta somatic SVs filtered to PASS (paired); deliver/{pair_id}.sv.candidate.vcf.gz — unscored candidates (tumor-only)
• deliver/{pair_id}.msi.tsv — MSIsensor-pro MSI score (paired runs only; tumor-only MSI needs a panel-of-normals baseline)
• report/signatures/{pair_id}.sbs96.counts.tsv / .exposures.tsv — SBS96 channel counts and NNLS-fitted COSMIC v3.1 exposures (with cosine similarity in the run log)
• report/tmb/{pair_id}.tmb.tsv — TMB under two conventions: tmb_mut_per_mb (coding non-silent per coding Mb — clinical, ≥10 mut/Mb cutoff) and all_somatic_mut_per_mb (all PASS per target Mb — WGS literature convention)
• report/venus_clinical_report.html — cohort report; its methods section is generated from rule_runs in .oxo-flow/checkpoint.json, so it always describes the commands that actually executed
• qc/multiqc/multiqc_report.html, rna/qc/rna_qc_summary.tsv, report/scrna_metrics.tsv
• qc/alignqc/wgs_alignqc_summary.tsv — per-sample WGS mapping rate, duplicate rate, depth, ≥10x/20x fractions, insert-size median, sex check (screening-level; contamination explicitly not_assessed)
• report/fusion_screen_cohort.tsv — RNA fusion screening candidates (panel hits first; screening-level, needs orthogonal confirmation)

License

MIT — see LICENSE.