cfDNA-immune-activation

Imported from vlachosg37/cfDNA-immune-activation (AGENTS.md). Install upstream with npx skills add vlachosg37/cfDNA-immune-activation. Copyright stays with the author.

Codex Agents: cfDNA-immune-activation

This document defines the project goal, the end-to-end implementation plan, and the engineering standards that all Codex agents must follow when writing or modifying code in this repository.

1) Project goal

Primary goal: Enable immune cell–resolved epigenetic profiling from cfDNA using Biomodal duet (6-base) sequencing data that includes 5mC and 5hmC calls.

We start from mixed cfDNA and aim to produce cell-type–specific readouts of immune activation and dysfunction by:

• estimating mixture fractions (composition),
• probabilistically labeling fragments/reads by immune cell-of-origin (cell type),
• quantifying promoter/marker 5mC and 5hmC within labeled immune fragments (state).

This repository should be understandable and compelling as a portfolio piece for bioinformatics engineering: clear interfaces, reproducibility, tests, and production-grade workflows.

2) End-to-end plan (target pipeline)

Inputs

• Biomodal duet modBAMs with 5mC/5hmC encoded via BAM tags (MM/ML).
• Reference FASTA (+ index).
• UXM marker resources + atlas (initially Loyfer et al. U markers/atlas).
• Genome annotation resources for promoter and regulatory-region summaries (see Section 4).

Outputs (target end-state)

• modcalls/{sample}_modcalls.tsv[.gz]
  • Per-cytosine per-read calls (5mC/5hmC/unmodified) with fragment and read metadata.
• uxm/mixture_fractions.tsv
  • Per-sample cell/tissue mixture fractions (NNLS).
• labels/{sample}_fragment_labels.tsv.gz
  • Per-fragment posterior over cell types and confidence metrics.
  • Optional: BAM tags (XC=top celltype, XP=max posterior, etc).
• promoters/{sample}.{celltype}.promoter_5mC_5hmC.tsv.gz
  • Cell-type–specific promoter-level summaries (5mC, 5hmC, coverage).
• ccres/{sample}.{celltype}.ccre_5mC_5hmC.tsv.gz
  • Cell-type–specific promoter-proximal cCRE summaries (5mC, 5hmC, coverage).
• scores/{sample}.{celltype}.immune_scores.tsv
  • Immune activation/dysfunction scores per immune cell type.
• HTML reports for QC and interpretability.

Stages

Stage A — Extraction (already implemented)

1. Parse duet BAM and write per-cytosine calls (modcalls.tsv).

Stage B — UXM deconvolution (already implemented) 2. Convert read-level calls to fragment-level UXM summaries at marker regions. 3. Solve atlas mixture fractions (NNLS; non-negative, sum-to-1).

Stage C — Per-fragment probabilistic cell-type labeling (to implement) 4. For each fragment that overlaps atlas markers:

• Compute evidence per marker (CpG counts; methylated counts; optional ML-weighted).
• Compute per-celltype likelihoods using atlas expected unmethylation probabilities.
• Convert to posteriors using priors (optionally blended with NNLS priors).
• Output top-1 label + posterior vector + LLR + evidence metrics.

5. Apply conservative confidence gating:
   • Min CpGs, min discriminative-marker requirement, posterior threshold, LLR threshold.
   • Fragments with insufficient evidence => label unknown.

Stage D — Cell-type–specific immune activation/dysfunction profiling (to implement) 6. Using labeled immune fragments (Blood-T, Blood-B, Blood-NK, Blood-Mono+Macro, Blood-Granul):

• Quantify promoter-level 5mC and 5hmC separately.
• Quantify promoter-proximal cCRE-level 5mC and 5hmC separately.
• Compute gene/panel scores (activation and dysfunction programs).
• Generate QC and summary outputs.

Notes on 5hmC

• Atlas-based labeling is initially based on WGBS-style 5mC expectations; for labeling we may treat 5hmC as methylated for compatibility.
• Downstream immune activation profiling should preserve 5mC vs 5hmC as separate signals.

3) Supported immune cell types (initial scope)

From deconvolution output populations, we focus on:

• Blood-T
• Blood-B
• Blood-NK
• Blood-Mono+Macro
• Blood-Granul

Other cell types may be present in mixture fractions but are not the initial focus of activation/dysfunction profiling.

4) Annotation resources and region definitions (hg38 / GRCh38)

This repo uses genome-wide annotation resources that are generated once and can later be filtered to gene panels via config.

4.1 Promoters (gene-level)

File

• resources/annotations/promoters.all.hg38.bed

Definition

• Gene-level promoters derived from a GENCODE GTF using feature == "gene" (one promoter per gene; avoids transcript explosion).
• Promoter interval: TSS ± 1000 bp (upstream=1000, downstream=1000).
• BED6 format with name = gene_name|gene_id.

Intended use

• Cell-type–specific promoter 5mC/5hmC quantification.
• Gene/panel scoring based on promoter methylation/hydroxymethylation.

4.2 Promoter-proximal ENCODE cCREs

Files

• resources/annotations/ccres_near_promoters.all.hg38.bed
• resources/annotations/ccres_near_promoters.all.hg38.links.tsv.gz

Definition

• ENCODE cCRE combined track (hg38), converted from UCSC bigBed to BED.
• “Near promoters” defined as cCREs within 2000 bp of any promoter interval (NEAR_W=2000).
• ccres_near_promoters.all.hg38.bed contains unique cCRE intervals (BED6).
• ccres_near_promoters.all.hg38.links.tsv.gz preserves gene-aware links.

Links file schema

• cCRE BED6 fields:
  1. ccre_chr
  2. ccre_start
  3. ccre_end
  4. ccre_id
  5. ccre_score
  6. ccre_strand (often .; treat as unstranded)
• Promoter mapping fields: 7) gene_name_gene_id (format: SYMBOL|ENSG...) 8) gene_strand 9) dist_bp (0 if overlapping promoter interval; otherwise gap in bp)

Intended use

• Summarize 5mC/5hmC in promoter-proximal regulatory regions.
• Filter cCREs to gene panels by using the links table (do not rely on cCRE intervals alone).

4.3 Panel filtering rule (future-proofing)

Agents should assume we will support selecting subsets later via config:

• Default: use *.all.hg38.*
• Panel mode: filter by a gene list (symbols) and regenerate:
  • resources/annotations/promoters.<panel>.hg38.bed
  • resources/annotations/ccres_near_promoters.<panel>.hg38.bed
  • resources/annotations/ccres_near_promoters.<panel>.hg38.links.tsv.gz

Panel filtering must be reproducible and should use explicit matching of gene_name|gene_id.

5) Core modeling requirements (fragment labeling)

Agents must implement labeling so it is explicit, probabilistic, and auditable.

Required properties

• Operates per fragment/read-pair using only overlapped marker regions.
• Produces a posterior distribution over candidate cell types.
• Provides explicit confidence metrics:
  • max posterior
  • LLR (best vs second-best)
  • number of markers used
  • number of CpGs used
• Allows an unknown label for fragments lacking evidence.

Recommended scoring model

• Prefer binomial or beta-binomial likelihood on CpG counts per marker (robust for small n).
• Avoid U/M/X hard-calling as the only model unless justified; if used, X must be modeled explicitly.
• Incorporate per-cytosine ML uncertainty if available (avoid losing information).
• Priors: allow optional use of NNLS mixture fractions, but avoid circularity by blending with uniform priors.

Marker discriminativeness

• Implement (or at least support) a marker “informativeness” filter (e.g., max Δp across immune types).
• Label only fragments overlapping sufficiently informative markers to improve precision.

6) Engineering standards for all agents (required)

6.1 Dignified Python (adapted)

Agents must follow these principles:

1. Look Before You Leap (LBYL)
   Prefer explicit checks over exception-driven control flow for:

   • dict keys, config keys, file existence, parsing branches, optional values.

2. Never swallow exceptions
   No except: pass. No broad except Exception unless:

   • adding context and re-raising, or
   • implementing a clearly documented fallback.

3. Defer import-time computation
   No I/O, filesystem scans, or environment parsing at import time.

4. Check existence before resolution
   Path.exists() before .resolve() or .is_relative_to().

5. Verify casts at runtime when cheap
   If using typing.cast, prefer an isinstance assertion when O(1).

6. Keyword-only args for complex functions
   For functions with many parameters (especially CLI-exposed), enforce keyword-only args.

7. Declare variables close to use
   Avoid long-lived intermediate variables unless reused or improves clarity.

6.2 Pipeline-specific rules (bioinformatics)

1. Determinism

   • Stable sorting, stable iteration. No implicit randomness.
   • If randomness is required, seed it and log it.

2. Streaming & memory discipline

   • Do not load entire BAM/TSV into memory without justification.
   • Prefer streaming iterators and chunked processing.
   • Support gzipped TSV where appropriate.

3. pysam correctness

   • Handle secondary/supplementary alignments explicitly.
   • Avoid assumptions about perfect pairing; validate flags.
   • Avoid repeated random access; prefer linear scans when possible.

4. Logging and observability

   • Use logging (not print) for scripts.
   • Emit counters: fragments processed, fragments labeled, unknown fraction, per-celltype counts.
   • Errors must be visible and actionable.

5. CLI hygiene

   • All scripts provide --help.
   • Validate inputs early and fail with informative messages.
   • Explicit exit codes.

6. No surprise side-effects

   • Do not write outside --outdir.
   • Do not auto-download resources.
   • Do not assume current working directory; accept absolute/relative paths robustly.

7. Reproducibility metadata

   • Record tool versions, parameters, and resource identifiers in outputs and/or logs.
   • Prefer writing a small run_metadata.json per run.

Agents must ensure that changes involving annotation resources also:

• document the definition (e.g., promoter window, NEAR_W),
• pin source versions/releases where possible,
• preserve schema stability.

7) Testing expectations

Agents should add tests when implementing new behavior.

Minimum expected tests:

• Unit tests for likelihood/posterior computation on synthetic markers.
• Edge cases:
  • no marker overlap ⇒ unknown
  • low CpG counts ⇒ unknown
  • two close cell types ⇒ low confidence
• Small integration test using tests/data/test1.bam (where feasible) validating:
  • output files exist
  • output schema columns and basic sanity counts

Performance sanity:

• Ensure labeling is vectorized where reasonable and does not scale quadratically with fragments.

8) Output formats (must be stable)

Fragment labels (TSV)

labels/{sample}_fragment_labels.tsv.gz must include at least:

• fragment_id
• chrom, start, end (or equivalent coordinate representation)
• celltype_top
• posterior_top
• llr_top_vs_second
• n_markers_used
• n_cpgs_used
• celltype_posteriors (JSON or delimited; documented clearly)

9) How agents should work in this repo

Agents must:

• Make changes in reviewable commits.
• Prefer incremental changes that keep the pipeline runnable.
• Update docs and parameter templates when adding new CLI/workflow options.
• Keep Nextflow and Snakemake behavior aligned where feasible.

A separate reviewer/verifier agent should be used to block:

• swallowed exceptions
• import-time side effects
• nondeterministic behavior
• large memory regressions
• undocumented output schema changes