Extract R analysis into modular backend abstraction - #153
Open
sanghoonio wants to merge 16 commits into
Open
Conversation
Create bedboss/bedstat/backends/ package with StatBackend ABC and RStatBackend implementation. Move R logic from bedstat.py into r_backend.py (no rewrite, just extraction). bedstat() becomes a thin dispatcher that delegates to the configured backend via factory. bedboss.py reads backend from bbagent.config.config.analysis.backend and passes it through. All R scripts and r_service.py untouched. Identical behavior when backend="r" (the default). Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Move magic strings "r" and "gtars" to BACKEND_R / BACKEND_GTARS constants in const.py. Also fixes black formatting in __init__.py. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
New files: - gtars_backend.py: GtarsStatBackend wrapping gtars genomicdist CLI - compress_distributions.py: fixed-size compression (histograms, KDE, dense arrays) - ref_utils.py: reference file resolution via refgenie + seqcol fallback Changes: - Register GtarsStatBackend in factory (replaces NotImplementedError) - Only create RServiceManager when backend is "r" - Skip R bedset plots when backend is "gtars" - Add --backend CLI override for run_all and run_stats - Add DEFAULT_PRECISION constant Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Reorder chrom.sizes resolution: local refgenie cache → seqcol API (~50KB) → refgenie FASTA pull (~3GB). Avoids downloading a full genome FASTA just to get chromosome sizes on cold start. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
5 tasks
R stores partition percentages as fractions (0.0615), not percent (6.082). Remove the * 100 multiplier so gtars output matches the existing database convention. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Derives median absolute neighbor distance from the gtars CLI's raw neighbor_distances list, stores it in the data dict so it gets written to the new BedStats.median_neighbor_distance column (added in bbconf PR #113). Computed BEFORE compress_distributions replaces the flat list with a KDE. The per-file KDE is still stored in the distributions JSONB blob for single-file views; only the scalar is used for bedset aggregation (mean ± sd across files). Replaces what used to be an aggregated neighbor_distances KDE at the bedset level — per earlier discussion, per-file KDE variance is low within bedsets (assay type dominates), and a scalar median gives enough signal at collection level. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
New 'gtars-py' backend that uses gtars Python bindings directly, no
subprocess. Exists alongside the existing 'gtars' CLI-subprocess
backend so we can benchmark both on real data and keep only the
better performer before merge.
Architecture:
- GtarsPyStatBackend in backends/gtars_py_backend.py
- GenomeRefs dataclass caches reference data per genome on the backend
instance (FASTA, chrom_sizes, GeneModel, PartitionList, TssIndex,
SignalMatrix). Loaded lazily on first compute() call for a given
genome, reused for all subsequent files.
- Same compute() signature as GtarsStatBackend
- Same output dict shape (compress_distributions + DB insertion unchanged)
- Same pypiper integration for status tracking (pm.timestamp markers
between Python steps, no subprocess wrapping)
- median_neighbor_distance scalar derived same way
Factory dispatch: create_backend("gtars-py") returns the new backend.
Extended bbconf's analysis.backend Literal to accept "gtars-py" (bbconf
commit on modular-backend-logic branch).
Output shape normalization: the gtars Python binding returns
calc_partitions as {partition: [names], count: [counts], total: n} but
the CLI JSON schema has {counts: [[name, count], ...], total: n}.
_normalize_partitions() converts between them so downstream code sees
the same shape regardless of backend.
Parity verified on a 126K-region hg38 ENCODE BED file:
- Scalars identical (number_of_regions, mean_region_width,
median_tss_dist, median_neighbor_distance, gc_content)
- All 14 partition fields (frequency+percentage for 7 categories) match
- widths, tss, neighbor_distances, region_distribution bins byte-match
after compress_distributions
Performance (5 ENCODE files, hg38):
- gtars (CLI): total 11.15s, median per-file 2.55s
- gtars-py: total 5.71s, median per-file 1.15s
- Speedup: ~2x batch, ~2.2x per-file median
- Biggest wins: FASTA loaded once (not per file), no subprocess startup,
no JSON round-trip to disk
Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Architecture change: batch orchestrators build a stats backend once and
pass it through the call chain as a StatBackend instance, instead of
creating a new backend per file inside bedstat().
Before: bedstat(backend="r", r_service=..., ...) reconstructed the
backend on every call. Fine for R because r_service was passed IN and
the RStatBackend constructor was cheap, but broken for future backends
that need per-instance caches (GtarsPyStatBackend's reference cache
would be wiped between files).
After: callers build one backend for their batch, pass it to bedstat()
for each file, and cleanup() at the end. Backend-internal resources
(R service, cached references) are held on the backend instance and
amortized across the whole batch.
Changes:
- bedstat(): signature takes `backend: StatBackend` instead of
`backend: str + r_service`; no longer calls create_backend internally
- Add build_backend(name) helper to backends/__init__.py that hides
backend-specific prerequisites (e.g. starts RServiceManager for "r")
- StatBackend base: add __enter__/__exit__ so backends work as context
managers (with build_backend("r") as backend: ...)
- run_all / insert_pep / reprocess_all / upload_all / _upload_gse:
build backend at batch-orchestrator level, pass through, cleanup
when done
- cli.py `bedstat` standalone command: build+use via `with` block
- run_all accepts optional backend=None for single-file standalone use;
builds + cleans up locally if caller didn't provide one
This is a breaking change — no backward compat shim. The backend
abstraction is pre-merge (PR #153).
Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Before: GtarsStatBackend was "mixed" — called gtars CLI for most stats but used gtars Python bindings for two things inside the backend: 1. RegionSet(bedfile) to compute bed_digest fallback 2. calc_gc_content(rs, assembly, ...) for GC content This blurred the line between the "CLI" and "Python bindings" backends. After: GtarsStatBackend is a pure CLI subprocess backend. - Pass --fasta to gtars CLI (uses PR #249's --fasta flag that outputs gc_content in the JSON). Read gc_mean + per_region GC from the CLI output instead of calling gtars Python bindings. - Require bed_digest from the caller. Backends no longer parse BED files. Resolution moves up to bedstat() (the orchestration layer), which computes the digest via gtars Python bindings once when absent. - Remove imports of gtars.models.RegionSet and calculate_gc_content. - Add get_fasta_path() helper to ref_utils.py (pure refgenie, no gtars dependency at import time) so the backend can resolve FASTA paths for the --fasta flag. The three backends now have cleanly separated execution domains: - RStatBackend: R subprocess (+ Python bindings helpers for file loading and GC content — kept as-is, R's native GC calc is slow) - GtarsStatBackend: gtars CLI subprocess ONLY - GtarsPyStatBackend: gtars Python bindings ONLY (no subprocess) Parity verified: pure-CLI gtars output matches gtars-py output exactly on all 19 scalar + partition fields for a 126K-region hg38 file. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
- Remove GtarsPyStatBackend (gtars CLI with .fab is faster and simpler) - Remove BACKEND_GTARS_PY constant and all references - Add get_fab_path() to ref_utils: auto-compiles .fab from FASTA on first use via gtars prep --fasta, cached forever - GtarsStatBackend prefers .fab over plain .fa for GC content - Fix reprocess_bedset: pass backend from config to run_bedbuncher - Update CLI help text and bedbuncher docstrings Two backends remain: "r" (RStatBackend) and "gtars" (GtarsStatBackend). Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
compress_region_distribution now takes chrom_sizes to compute correct per-chromosome array lengths (longest chr = n_bins, shorter chrs proportionally fewer). Without this, arrays were sized by max observed rid which varied by region coverage. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Add gtars genomicdist backend
Nothing consumes it. BedStatsModel has never declared the field, so extra="ignore" discarded the value before it could reach bed_stats, and the bbconf aggregation that was meant to read it has been dropped along with the collection-level histogram it fed. Removing the derivation also drops a median over the full neighbor_distances list on every file. The list itself is untouched and still reaches the distributions blob, which is what the per-file density plot renders. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01PRGhWjGiUXFoygd35QsKQ7
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
2 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Summary
Refactor bedboss bedstat into a modular backend system. Zero behavior change — all existing R-based analysis works identically.
What changed
New package:
bedboss/bedstat/backends/base.py:StatBackendABC withcompute()andcleanup()methods, context manager protocolr_backend.py:RStatBackend— extracts current R logic frombedstat.py(move, not rewrite)__init__.py:create_backend(name)low-level factory,build_backend(name)high-level constructor that handles backend-specific prerequisites (e.g. startsRServiceManagerfor R)Refactored orchestration:
bedstat.py: signature changed frombedstat(backend="r", r_service=..., ...)tobedstat(backend: StatBackend, ...)— takes an instance, not a stringbedboss.py,insert_pep,reprocess_all,upload_all: build one backend at batch-orchestrator level viabuild_backend(), pass through tobedstat()for each file, cleanup when donecli.py: standalonebedstatcommand useswithblock for automatic cleanupWhy this matters:
Backend instances hold resources across a batch (R keeps its subprocess alive, future backends can cache reference data). The old pattern reconstructed backends per file, which prevented cross-file amortization.
Constants
BACKEND_R = "r"andBACKEND_GTARS = "gtars"inconst.py(gtars backend comes in PR Add gtars genomicdist backend #155)Test plan
RServiceManagerlifecycle unchanged (created inbuild_backend, cleaned up by caller)🤖 Generated with Claude Code