From da57ff202d67ca254ad257c7d6c47f5fc7917562 Mon Sep 17 00:00:00 2001 From: gbouras13 Date: Thu, 13 Aug 2026 23:02:38 +0930 Subject: [PATCH 01/12] Finish gracefully when no genomes are resolved (coverm zero-alignment panic) A sample where phables resolves no genomes and has no unresolved phage-like edges produces an empty genomes_and_unresolved_edges.fasta, so coverm_map_genomes emits a valid but empty BAM. CoverM 0.7.0 panics on such a BAM ('index out of bounds: the len is 0 but the index is 0', coverage_printer.rs:467) instead of printing an empty table, killing the whole run at the very last stage after all the expensive work had already succeeded. Guard the rule: if the BAM has no alignments, write coverm's header with no data rows rather than invoking coverm. Everything downstream already handles an empty table (verified end-to-end: coverm_combine_ genomes writes header-only output, format_koverage_results.py's readlines()[1:] yields no rows and pandas writes header-only report TSVs), so the run now completes normally with empty report tables. The fallback header reproduces CoverM's own exactly -- verified against CoverM 0.7.0 source rather than guessed: coverage_printer.rs joins ' ' with a space, mosdepth_genome_coverage_estimators.rs ::column_headers defines the metric strings, and bin/coverm.rs builds the estimator list by iterating the -m flags in order. That also disproved an existing comment in format_koverage_results.py claiming the covered_fraction column reads 'Covered_fraction' -- it is actually 'Covered Fraction', two space-separated words, and count is 'Read Count'. Corrected; the IDX_* indices it documents were already right. --- phables/workflow/rules/postprocess.smk | 43 +++++++++++++++++-- .../scripts/format_koverage_results.py | 16 +++++-- 2 files changed, 51 insertions(+), 8 deletions(-) diff --git a/phables/workflow/rules/postprocess.smk b/phables/workflow/rules/postprocess.smk index d481e5c..1d33b02 100644 --- a/phables/workflow/rules/postprocess.smk +++ b/phables/workflow/rules/postprocess.smk @@ -90,7 +90,34 @@ rule coverm_map_genomes: rule coverm_bam2counts_genomes: - """Per-sample coverage stats over the resolved genomes.""" + """Per-sample coverage stats over the resolved genomes. + + Guarded against the zero-alignment case, which is a NORMAL outcome, not an + error: a sample where phables resolved no genomes AND had no unresolved + phage-like edges produces an empty genomes_and_unresolved_edges.fasta, so + coverm_map_genomes maps against an empty reference and emits a valid but + empty BAM. CoverM 0.7.0 panics outright on such a BAM rather than printing + an empty table -- + + [WARN coverm::contig] No primary alignments were observed for sample X + thread 'main' panicked at src/coverage_printer.rs:467:61: + index out of bounds: the len is 0 but the index is 0 + + -- which killed the whole run at the very last stage, after all the + expensive work had already succeeded. Everything downstream of here already + handles an empty table correctly (coverm_combine_genomes writes header-only + output; format_koverage_results.py's `readlines()[1:]` yields no rows and + pandas writes header-only report TSVs), so emitting the header ourselves is + all that's needed for the run to finish normally with empty report tables. + + The header below reproduces CoverM's own exactly -- " ", + space-joined, in the order the -m flags are given (confirmed against CoverM + 0.7.0 source: coverage_printer.rs writes "\\t{stoit_name} {estimator_header}", + mosdepth_genome_coverage_estimators.rs::column_headers defines the metric + strings, and bin/coverm.rs builds the estimator list by iterating the -m + flags in order). Note "Covered Fraction" is two space-separated words. + Stoit name is the BAM's basename without .bam, i.e. exactly {sample}. + """ input: os.path.join(OUTDIR, "postprocess", "temp", "{sample}.bam") output: @@ -103,9 +130,17 @@ rule coverm_bam2counts_genomes: os.path.join(LOGSDIR, "coverm_bam2counts_genomes.{sample}.log") shell: """ - coverm contig -b {input} \ - -m count -m rpkm -m tpm -m mean -m covered_fraction -m variance \ - > {output} 2> {log} + n_aln=$(samtools view -c {input}) + if [ "$n_aln" -eq 0 ]; then + echo "No alignments in {input} -- no genomes were resolved for this sample (and no unresolved phage-like edges), so there is nothing to compute coverage over. Writing a header-only coverage table instead of running coverm, which panics on a zero-alignment BAM. The run continues and finishes normally; the per-genome report tables will be empty." > {log} + S={wildcards.sample} + printf 'Contig\\t%s Read Count\\t%s RPKM\\t%s TPM\\t%s Mean\\t%s Covered Fraction\\t%s Variance\\n' \ + "$S" "$S" "$S" "$S" "$S" "$S" > {output} + else + coverm contig -b {input} \ + -m count -m rpkm -m tpm -m mean -m covered_fraction -m variance \ + > {output} 2> {log} + fi """ diff --git a/phables/workflow/scripts/format_koverage_results.py b/phables/workflow/scripts/format_koverage_results.py index e0f9c5e..9367c00 100644 --- a/phables/workflow/scripts/format_koverage_results.py +++ b/phables/workflow/scripts/format_koverage_results.py @@ -41,11 +41,19 @@ # / coverm_bam2counts_genomes / coverm_combine_genomes), completing PLAN.md # §4.8. The indices below therefore describe the CoverM-mode header: # -# Sample Contig Count RPKM TPM Mean Covered_fraction Variance -# 0 1 2 3 4 5 6 7 +# Sample Contig Read Count RPKM TPM Mean Covered Fraction Variance +# 0 1 2 3 4 5 6 7 # -# ("Covered_fraction" is coverm's `covered_fraction` method as it appears after -# coverm_combine_genomes strips the per-column " " prefix.) +# (Those are CoverM's own literal header strings as they appear after +# coverm_combine_genomes strips the per-column " " prefix -- verified +# against CoverM 0.7.0 source, mosdepth_genome_coverage_estimators.rs:: +# column_headers. The `count` method's header is "Read Count" and +# `covered_fraction`'s is "Covered Fraction" -- two space-separated words, NOT +# "Covered_fraction" as this comment previously claimed. Column ORDER follows +# the order the -m flags are passed in postprocess.smk, not any canonical +# order: bin/coverm.rs builds the estimator list by iterating the flags with +# .enumerate(). Only the indices below actually matter to this script; the +# names are documentation.) # # HISTORICAL, for anyone reading old output or an older checkout: this script # previously parsed Koverage's *native* "map" mode file, whose header is From fe8b5137c5eef0f21cacc1d4ce5a3865949953c2 Mon Sep 17 00:00:00 2001 From: gbouras13 Date: Thu, 13 Aug 2026 23:04:08 +0930 Subject: [PATCH 02/12] prostt5-rocm: bump torch pin to 2.9.1 torch==2.9.1+rocm6.3 confirmed to exist on the pinned ROCm 6.3 index (download.pytorch.org/whl/rocm6.3/torch/), with cp310-cp314 wheels -- so every Python the env's own python>=3.10 pin can resolve is covered (2.9.1 drops the cp39 wheel 2.7.1 shipped, which that pin already excludes anyway). The --extra-index-url is unchanged; rocm6.3 still carries this version. This intentionally diverges the conda-env ROCm path from the container path, which still builds on quay.io/pawsey/pytorch:2.7.1-rocm6.3.3 and asserts torch 2.7.1 -- noted inline in the yaml so the mismatch reads as deliberate rather than an oversight. --- phables/workflow/envs/prostt5-rocm.yaml | 39 ++++++++++++++++--------- 1 file changed, 25 insertions(+), 14 deletions(-) diff --git a/phables/workflow/envs/prostt5-rocm.yaml b/phables/workflow/envs/prostt5-rocm.yaml index 6b33dc9..02df6d5 100644 --- a/phables/workflow/envs/prostt5-rocm.yaml +++ b/phables/workflow/envs/prostt5-rocm.yaml @@ -8,21 +8,32 @@ dependencies: # Unlike CUDA, ROCm is never the default pip torch build -- this index # override is required, not just for reproducibility. - --extra-index-url https://download.pytorch.org/whl/rocm6.3 - # torch==2.7.1+rocm6.3 pinned exactly, not just floored -- matches Setonix's - # own verified-working Pawsey container (quay.io/pawsey/pytorch:2.7.1-rocm6.3.3), - # confirmed real via https://download.pytorch.org/whl/rocm6.3/torch/ (the - # exact wheel exists for cp39-cp313). Switched to this conda-env path instead - # of the --prostt5-container route after repeatedly hitting container-specific - # problems on Setonix that had nothing to do with phables/pholdlib itself: - # a broken torchaudio needing libomp.so, --home overriding the container's - # apparent CWD so relative output paths silently resolved to the wrong place, - # and a missing MIOpen AI-heuristic kernel-selection file for gfx90a. All real, - # all diagnosed, none of them things a plain conda env run directly on the host - # (no container namespace, no bind-mount surface) is exposed to at all. + # torch==2.9.1+rocm6.3 pinned exactly, not just floored. Confirmed real on + # https://download.pytorch.org/whl/rocm6.3/torch/ (wheels for cp310-cp314, + # so every Python the `python>=3.10` pin above can resolve is covered -- + # 2.9.1 drops the cp39 wheel 2.7.1 had, which that pin already excludes). + # + # NOTE this deliberately no longer matches the Pawsey base image the + # CONTAINER path uses (quay.io/pawsey/pytorch:2.7.1-rocm6.3.3, see + # container/Dockerfile, which bakes in and asserts torch 2.7.1). The two + # ROCm paths are now on different torch versions on purpose: the container + # takes whatever Pawsey verified, this conda env takes a newer torch. If + # you need them identical again, either rebuild the container from a newer + # Pawsey image or drop this back to 2.7.1. + # + # Switched to this conda-env path instead of the --prostt5-container route + # after repeatedly hitting container-specific problems on Setonix that had + # nothing to do with phables/pholdlib itself: a broken torchaudio needing + # libomp.so, --home overriding the container's apparent CWD so relative + # output paths silently resolved to the wrong place, and a missing MIOpen + # AI-heuristic kernel-selection file for gfx90a. All real, all diagnosed, + # none of them things a plain conda env run directly on the host (no + # container namespace, no bind-mount surface) is exposed to at all. # >=2.6 would already avoid transformers' checkpoint loader refusing # torch.load() on older torch outright (CVE-2025-32434 hardening in # check_torch_load_is_safe, which broke loading Rostlab/ProstT5_fp16 -- no - # safetensors file); pinning to exactly 2.7.1 goes further, matching a - # version already confirmed to work on this exact hardware (MI250X, gfx90a). - - torch==2.7.1 + # safetensors file); pinning exactly goes further, keeping the version + # deterministic across every array task rather than drifting as new + # releases land mid-run. + - torch==2.9.1 - pholdlib>=0.1.2 From 7ea315977804343ffca5e02e9dfdccd9eee33092 Mon Sep 17 00:00:00 2001 From: gbouras13 Date: Fri, 14 Aug 2026 08:19:56 +0930 Subject: [PATCH 03/12] Replace per-rule container flags with one monolithic Setonix image Removes --container and --prostt5-container entirely: both CLI options, both config keys, the CONTAINER_IMAGE/PROSTT5_CONTAINER preflight vars, all 15 per-rule 'container:' directives, and genes.smk's PROSTT5_CONTAINER if/else split around predict_3di. Every rule is back to a plain 'conda:' directive. In their place, container/Dockerfile is now one monolithic image built the way hybracter's is: Pawsey's ROCm base for the gfx90a userspace, miniforge installed inside, and EVERY per-rule conda env pre-built into the image by container/prebuild_envs.sh. Convert it to a single .sif and 'phables run' inside it finds its envs already there -- which is what actually kills the conda-env race condition (many array tasks each being their own Snakemake process, all creating the same shared env at once), since at runtime there is nothing left to create. Envs are built at phables' DEFAULT --conda-prefix on purpose: Snakemake hashes an env file's content together with the prefix path, so building and running at the same prefix in the same image makes the hashes match and nothing is rebuilt -- important because a .sif is read-only, so a rebuild attempt is a hard failure rather than a slow path. Databases stay OUT of the image. prebuild_envs.sh creates zero-byte placeholder DB files purely so the DAG can resolve, then runs --conda-create-envs-only once per flag combination (gene caller, detection mode, GPU backend, tree, install). Verified for real: the DAG resolves against zero-byte placeholder databases, and all six flag combinations dry-run clean after the container removal. container/test_image.sh smoke-tests the built image and fails the build if any pre-built env is missing a binary the rules actually invoke, or if torch didn't come from the prostt5-rocm env rather than the base image's own unused 2.7.1. Not yet built: no docker build or Setonix Apptainer run has happened. The pinned conda packages were confirmed present for linux-64/noarch (mmseqs2=13.45111, foldseek, coverm, fraggenescan, cogent3<2026.7 on bioconda/noarch; piqtree is pip-only, which phylotree.yaml already reflects). --- .github/workflows/build_container.yaml | 14 ++ container/Dockerfile | 173 ++++++----------- container/prebuild_envs.sh | 78 ++++++++ container/test_image.sh | 94 +++++++++ docs/container.md | 181 ++++++++++-------- docs/usage.md | 33 +--- phables/__main__.py | 26 --- phables/config/config.yaml | 27 --- phables/workflow/envs/prostt5-rocm.yaml | 31 +-- .../workflow/rules/02_phables_preflight.smk | 15 -- phables/workflow/rules/coverage.smk | 8 +- phables/workflow/rules/genes.smk | 121 ++++-------- phables/workflow/rules/gfa2fasta.smk | 4 +- phables/workflow/rules/phables.smk | 4 +- phables/workflow/rules/phylotree.smk | 8 +- phables/workflow/rules/postprocess.smk | 12 +- 16 files changed, 415 insertions(+), 414 deletions(-) create mode 100755 container/prebuild_envs.sh create mode 100755 container/test_image.sh diff --git a/.github/workflows/build_container.yaml b/.github/workflows/build_container.yaml index 23100bd..a5ef161 100644 --- a/.github/workflows/build_container.yaml +++ b/.github/workflows/build_container.yaml @@ -24,6 +24,15 @@ jobs: # pre-pulled docker images), freeing ~25-30GB. ``tool-cache: false`` # keeps the runner's node/python toolcache so other actions (checkout, # docker login) still work fast. + # + # NOTE this image is now MONOLITHIC -- on top of that base it installs + # miniforge and pre-builds every per-rule conda env, including two torch + # builds (rocm + cpu). That is several GB more than the base alone, and + # is close to what a hosted runner can do even after this cleanup. If a + # build starts failing on "no space left on device", the honest fix is + # building on a machine with real disk and pushing manually, not shaving + # more off here. The ``df -h`` steps below exist to make that diagnosis + # obvious rather than a guess. - name: Free up disk space uses: jlumbroso/free-disk-space@main with: @@ -46,6 +55,11 @@ jobs: IMAGE=quay.io/${{ secrets.QUAY_USERNAME }}/phables:${{ github.sha }} docker build -f container/Dockerfile -t $IMAGE . + - name: Show image size and free disk after build + run: | + docker image ls + df -h + - name: Push Docker image run: | IMAGE=quay.io/${{ secrets.QUAY_USERNAME }}/phables:${{ github.sha }} diff --git a/container/Dockerfile b/container/Dockerfile index c5bcb5c..5040106 100644 --- a/container/Dockerfile +++ b/container/Dockerfile @@ -1,137 +1,72 @@ # -# phables -- massive all-in-one image for Setonix (ROCm), replacing per-rule -# --use-conda env creation entirely. See ../docs/container.md for how to run -# it (--container --use-singularity --no-use-conda) and why this -# exists (a real, hit-in-production race condition: every sample's `phables -# run` is its own separate process, so many array tasks needing the same -# shared --use-conda env for the first time at once can corrupt each other's -# `mamba env create`). +# phables -- ONE monolithic image for Setonix: ROCm base + conda + every +# per-rule conda env pre-built INSIDE the image. Build it once, convert it to a +# single .sif, and every `phables run` inside it finds its conda envs already +# there. Same strategy as hybracter's own container (miniforge installed into +# the image, envs materialised at build time), NOT the per-rule `container:` +# directive approach -- phables has no --container/--prostt5-container flags +# any more, and its Snakemake rules carry plain `conda:` directives only. +# +# Why this exists: a real, hit-in-production race condition. Every sample's +# `phables run` is its own separate process (one SLURM array task each), and +# Snakemake only serialises conda-env creation WITHIN a single process's DAG. +# Many array tasks all needing the same not-yet-built shared env at once +# corrupt each other's `mamba env create` ("Fatal Python error: +# init_fs_encoding ... no codec search functions registered"). With every env +# already baked into the image there is nothing left to create at runtime, and +# so nothing left to race. +# +# Base is Pawsey's own verified-working ROCm image, as requested: it supplies +# the ROCm userspace matching Setonix's MI250X (gfx90a) kernel driver. NOTE +# that the base's own system-python torch (2.7.1) is NOT what the workflow +# uses -- predict_3di runs inside the prostt5-rocm conda env built below, which +# brings its own torch (2.9.1+rocm6.3, see workflow/envs/prostt5-rocm.yaml). +# The base image's torch is left untouched and simply goes unused. +# +# See ../docs/container.md for how to build the .sif and run it on Setonix. # -# Base is Pawsey's own verified-working ROCm+PyTorch image -- the SAME one -# phold's own container (../../phold/container/hpci/Dockerfile) builds from, -# deliberately: torch==2.7.1+rocm6.3 baked in here is exactly what -# envs/prostt5-rocm.yaml pins, confirmed working on real Setonix MI250X -# hardware this session. Nothing below reinstalls or touches torch -- pip -# installs pholdlib and every other pure-Python dependency directly into this -# same system Python, so they see (and don't conflict with) the base image's -# already-correct torch/ROCm stack rather than resolving a second, possibly -# different one. FROM quay.io/pawsey/pytorch:2.7.1-rocm6.3.3 ARG DEBIAN_FRONTEND="noninteractive" RUN apt-get update && apt-get install -y --no-install-recommends \ - wget \ - tar \ - build-essential \ + wget \ + bzip2 \ + ca-certificates \ + git \ + procps \ && rm -rf /var/lib/apt/lists/* -# the base image only provides python3, not python -- some tooling (and -# Snakemake's own script: directive) expects `python` on PATH -RUN ln -sf "$(command -v python3)" /usr/local/bin/python - # -# foldseek -- same install as phold's own Dockerfile (direct binary download, -# not conda/pip), since that's the reference build this image otherwise -# mirrors for the ROCm/torch stack +# Miniforge inside the image -- hybracter's containers/HPC/Dockerfile does +# exactly this (wget the installer, `bash ... -b -p `) rather than +# relying on a conda that may or may not exist on the host. conda-forge is +# miniforge's default channel, so no extra channel config is needed. # -RUN wget https://mmseqs.com/foldseek/foldseek-linux-avx2.tar.gz && \ - tar -xzf foldseek-linux-avx2.tar.gz && \ - rm foldseek-linux-avx2.tar.gz && \ - mv foldseek /opt/foldseek -ENV PATH="/opt/foldseek/bin:${PATH}" +ENV CONDA_DIR=/opt/miniforge3 +RUN wget -q https://github.com/conda-forge/miniforge/releases/latest/download/Miniforge3-Linux-x86_64.sh -O /tmp/miniforge.sh && \ + bash /tmp/miniforge.sh -b -p "${CONDA_DIR}" && \ + rm /tmp/miniforge.sh +ENV PATH="${CONDA_DIR}/bin:${PATH}" -# -# Every other per-rule tool that isn't a pure-Python package: bootstrap a -# standalone micromamba (no full miniforge needed) and install them all into -# one prefix, /opt/conda -- kept entirely separate from the base image's own -# system Python (added to PATH, but conda never pulls in a python of its own -# here since none of these packages are Python libraries, only standalone -# binaries -- confirmed by their real recipes: minimap2/samtools/mmseqs2/ -# fraggenescan/hmmer/mafft have no python dependency; coverm's own recipe -# pulls minimap2+samtools transitively too, listed here explicitly anyway for -# a single deterministic solve rather than relying on that side effect). -# mmseqs2 pinned to the exact version envs/mmseqs.yaml already uses. -RUN curl -Ls https://micro.mamba.pm/api/micromamba/linux-64/latest | tar -xvj bin/micromamba && \ - mv bin/micromamba /usr/local/bin/micromamba && \ - micromamba create -y -p /opt/conda -c conda-forge -c bioconda \ - minimap2 \ - samtools \ - "coverm>=0.6" \ - "mmseqs2=13.45111" \ - fraggenescan \ - hmmer \ - mafft \ - && micromamba clean -a -y -ENV PATH="/opt/conda/bin:${PATH}" - -# -# Everything else -- phables' own pure-Python deps, ProstT5/pholdlib, and -# phables itself -- pip installed directly into the base image's system -# Python. Deliberately NOT a separate venv/conda env: the whole point is for -# all of this to see the base image's already-installed, already-verified -# torch/ROCm build rather than resolving (and risking conflicting with) a -# second one. -# -# pholdlib's own declared torch requirement is unconstrained (confirmed -# earlier this session against envs/prostt5-rocm.yaml/-cpu.yaml/-cuda.yaml) -- -# `pip install pholdlib` alone pulls in its transitive deps (transformers, -# sentencepiece, h5py, loguru, etc.) via pip's own resolver without touching -# the already-installed torch. -RUN pip install --ignore-installed PyYAML && \ - python -c "import torch; print('torch (pre-existing, untouched):', torch.__version__)" +# Snakemake lives in the base conda env alongside phables itself -- it is the +# thing that READS the pre-built per-rule envs, so it must not be one of them. +RUN conda install -y -c conda-forge -c bioconda \ + "snakemake>=8" \ + && conda clean -a -y # -# NOT included: cogent3/piqtree (phylotree.yaml's own deps, for the optional -# --phylotree/build_tree rule -- not part of this workflow's usage). -# Deliberately left out rather than fought into pip (every cogent3 release on -# PyPI is an unreleased alpha, which breaks a plain version-range pip -# install). Note: --phylotree is therefore NOT supported under --container -- -# build_tree still declares container: CONTAINER_IMAGE (unchanged), so it -# will fail if actually invoked against this image. -# -RUN pip install \ - "pholdlib>=0.1.2" \ - "pyrodigal-gv>=0.3.1" \ - biopython \ - python-igraph \ - pysam \ - "networkx>=2.8.6" \ - scipy \ - numpy \ - pandas \ - tqdm \ - click \ - "metasnek>=0.0.3" \ - "flowpaths>=0.2.20" - -RUN python -c "import torch; print('torch (post pip installs, should be unchanged):', torch.__version__)" - -# -# phables itself -- installed from the exact source in this build context -# (this Dockerfile is meant to be built by CI on every commit, tagged by -# commit SHA, so it's the checked-out commit's own tree, not a fetched -# git+https ref pinned separately). +# phables from the exact source in this build context. CI builds this image per +# commit and tags it by commit SHA, so this is the checked-out tree -- not a +# separately-pinned git ref that could drift from the rules and env files the +# image bakes in. # COPY . /opt/phables_src RUN pip install /opt/phables_src -# -# tests -- fail the build loudly if any of the above didn't actually work, -# rather than shipping a broken image to quay.io -# -RUN foldseek --help -RUN minimap2 --version -RUN samtools --version -RUN mmseqs version -RUN coverm --version -RUN hmmsearch -h -# the bioconda package's actual binary is run_FragGeneScan.pl, not -# `fraggenescan` -- confirmed against genes.smk's own real invocation -RUN run_FragGeneScan.pl 2>&1 | head -5 || true -RUN mafft --version -RUN python --version -RUN python -c "import pholdlib; print('pholdlib OK')" -RUN python -c "import torch; assert torch.__version__.startswith('2.7.1'), torch.__version__" -RUN phables --version -RUN phables run -h +# Pre-build every per-rule conda env. See the script for the full rationale +# (default --conda-prefix, placeholder databases, one pass per flag combo). +RUN bash /opt/phables_src/container/prebuild_envs.sh /opt/phables_src + +# Fail the build loudly rather than pushing a broken image to quay.io. +RUN bash /opt/phables_src/container/test_image.sh diff --git a/container/prebuild_envs.sh b/container/prebuild_envs.sh new file mode 100755 index 0000000..d0a2f35 --- /dev/null +++ b/container/prebuild_envs.sh @@ -0,0 +1,78 @@ +#!/bin/bash +# Pre-builds EVERY per-rule conda env into the image, at Docker build time. +# Run from container/Dockerfile; not meant to be run on a host. +# +# Envs go into phables' own DEFAULT --conda-prefix (snake_base("workflow/conda"), +# i.e. inside the installed package) deliberately: a user inside the resulting +# container then runs plain `phables run ...` with no --conda-prefix flag and +# Snakemake resolves these exact envs. Snakemake names an env by hashing its +# file content together with the conda prefix path -- both identical at build +# time and run time here, since they're the same paths in the same image -- so +# the hashes match and nothing is rebuilt at runtime. That matters more than +# it sounds: a .sif is READ-ONLY when running, so an attempted rebuild is a +# hard failure, not just a slow path. +# +# --conda-create-envs-only builds a DAG's envs without running any of it. The +# DAG still has to RESOLVE, which requires the database files to EXIST -- but +# only to exist, since no job runs. Empty placeholder files are therefore +# enough. That was verified for real (dry-running every flag combination below +# against zero-byte placeholder DB files) before this script was written, and +# it's what keeps the multi-GB PHROGs/hallmark databases OUT of the image: +# mount the real ones at runtime via --databases, exactly as outside a +# container. +# +# One invocation per flag combination, because which envs a DAG needs depends +# on the flags -- gene caller, phage-detection mode and GPU backend each select +# different rules and env files. Together these cover every env under +# workflow/envs/ that any `phables run` (or `phables install`) can reach. + +set -euxo pipefail + +SRC="${1:-/opt/phables_src}" +DB=/tmp/placeholder_db + +mkdir -p "$DB/phrogs_mmseqs_db" "$DB/hallmark_db" +touch "$DB/marker.hmm" \ + "$DB/phrog_annot_v4.tsv" \ + "$DB/phrogs_mmseqs_db/phrogs_profile_db" \ + "$DB/hallmark_db/hallmark_db" \ + "$DB/hallmark_db/hallmark_categories.tsv" + +GFA="${SRC}/tests/data/ERR1301161/assembly_graph_after_simplification.gfa" +READS="${SRC}/tests/data/ERR1301161" +COMMON=(--input "$GFA" --reads "$READS" --databases "$DB" --threads 1) + +# 1. default path -> coverm, genecall (FragGeneScan), smg (HMMER), mmseqs, phables +phables run "${COMMON[@]}" --output /tmp/envbuild1 --conda-create-envs-only + +# 2. the other gene caller +phables run "${COMMON[@]}" --output /tmp/envbuild2 \ + --genecaller pyrodigal-gv --conda-create-envs-only + +# 3. ProstT5 + foldseek detection on ROCm -- the reason this image exists +phables run "${COMMON[@]}" --output /tmp/envbuild3 \ + --phagedetection prostt5-foldseek --gpu-backend rocm --conda-create-envs-only + +# 4. Same, CPU backend. Built because gpu_backend's own config default is `cpu`: +# without this, forgetting `--gpu-backend rocm` inside the container would hit +# a missing env on a read-only filesystem. (cuda is deliberately NOT built -- +# this is a ROCm image for Setonix, and the CUDA wheels would add several GB +# that could never be used on this hardware.) +phables run "${COMMON[@]}" --output /tmp/envbuild4 \ + --phagedetection prostt5-foldseek --gpu-backend cpu --conda-create-envs-only + +# 5. Optional phylogenetic tree -> phylotree env (MAFFT + cogent3/piqtree). +# Note conda resolves cogent3 fine, unlike pip, where every published +# release is a prerelease and a plain version range matches nothing. +phables run "${COMMON[@]}" --output /tmp/envbuild5 \ + --build-tree --conda-create-envs-only + +# 6. install.smk's own env (curl), so `phables install` works inside here too +phables install --output /tmp/envbuild6 --databases "$DB" --conda-create-envs-only + +rm -rf /tmp/envbuild1 /tmp/envbuild2 /tmp/envbuild3 /tmp/envbuild4 /tmp/envbuild5 /tmp/envbuild6 "$DB" +conda clean -a -y + +echo "=== pre-built conda envs ===" +CONDA_PREFIX_DIR="$(python -c 'import phables, os; print(os.path.join(os.path.dirname(phables.__file__), "workflow", "conda"))')" +ls -1 "$CONDA_PREFIX_DIR" diff --git a/container/test_image.sh b/container/test_image.sh new file mode 100755 index 0000000..6737a30 --- /dev/null +++ b/container/test_image.sh @@ -0,0 +1,94 @@ +#!/bin/bash +# Build-time smoke tests for the monolithic image -- fail the build loudly +# rather than pushing something broken to quay.io. Run from container/Dockerfile. +# +# Every per-rule tool is checked INSIDE whichever pre-built conda env provides +# it, not on the base PATH: none of them are on the base PATH, and that's the +# point of the per-rule env layout. Envs are located by searching the prefix +# rather than by hardcoded directory names, because Snakemake names each env +# by a content hash that this script has no business predicting. + +set -euo pipefail + +PREFIX="$(python -c 'import phables, os; print(os.path.join(os.path.dirname(phables.__file__), "workflow", "conda"))')" +echo "conda prefix: $PREFIX" +test -d "$PREFIX" + +echo "=== phables CLI ===" +phables --version +phables run -h > /dev/null +phables install -h > /dev/null +echo "CLI OK" + +echo "=== the --container / --prostt5-container flags must be GONE ===" +if phables run -h 2>&1 | grep -qE '\-\-(prostt5-)?container'; then + echo "ERROR: a container flag is still present in the CLI" >&2 + exit 1 +fi +echo "confirmed absent" + +echo "=== pre-built env count ===" +n=$(find "$PREFIX" -maxdepth 1 -mindepth 1 -type d | wc -l) +echo "found $n env directories" +find "$PREFIX" -maxdepth 1 -mindepth 1 -type d -exec basename {} \; +# 6 distinct env files are reachable at minimum (coverm, genecall, smg, mmseqs, +# phables, curl) before counting foldseek/prostt5-rocm/prostt5-cpu/phylotree. +test "$n" -ge 6 + +# Finds an executable in any pre-built env; fails if no env provides it. +check_bin() { + local want="$1" e + for e in "$PREFIX"/*/; do + if [ -x "${e}bin/${want}" ]; then + echo "OK: $want -> $e" + return 0 + fi + done + echo "MISSING from every pre-built env: $want" >&2 + return 1 +} + +echo "=== per-rule binaries ===" +check_bin minimap2 +check_bin samtools +check_bin coverm +check_bin mmseqs +check_bin foldseek +check_bin hmmsearch +# the bioconda package's real binary name, confirmed against genes.smk's own +# invocation -- NOT `fraggenescan` +check_bin run_FragGeneScan.pl +check_bin mafft +check_bin curl + +echo "=== torch must come from the prostt5-rocm env, not the base image ===" +# The base image ships its own system-python torch 2.7.1, which the workflow +# does NOT use. predict_3di runs inside the prostt5-rocm conda env, whose torch +# is 2.9.1+rocm6.3 (workflow/envs/prostt5-rocm.yaml). If nothing here reports +# 2.9.1+rocm, that env didn't build properly and predict_3di would be running +# on the wrong stack. +found=0 +for e in "$PREFIX"/*/; do + if [ -x "${e}bin/python" ]; then + if "${e}bin/python" - <<'PY' 2>/dev/null +import sys +try: + import torch +except Exception: + sys.exit(1) +sys.exit(0 if torch.__version__.startswith("2.9.1") and "rocm" in torch.__version__ else 1) +PY + then + echo "OK: torch 2.9.1+rocm in $e" + "${e}bin/python" -c "import pholdlib; print('pholdlib OK')" + found=1 + break + fi + fi +done +if [ "$found" -ne 1 ]; then + echo "ERROR: no pre-built env provides torch 2.9.1+rocm" >&2 + exit 1 +fi + +echo "=== all image tests passed ===" diff --git a/docs/container.md b/docs/container.md index 4adc96f..83cfabd 100644 --- a/docs/container.md +++ b/docs/container.md @@ -1,75 +1,93 @@ # Running phables from a single container -`--use-conda`'s per-rule environment creation is convenient but has a real -failure mode at HPC scale: every rule's conda env is named by a hash of its -env file, so every sample's `phables run` needs the *same* env directory — -but each sample is its own separate Snakemake process (a SLURM array task, -say), and Snakemake only serializes env creation *within* one process's DAG. -Two array tasks both needing the same not-yet-built env at the same moment -can race `mamba env create` into the same directory, corrupting it -(`Fatal Python error: init_fs_encoding ... no codec search functions -registered` is what that looks like when it happens). - -`container/Dockerfile` builds a single image with every per-rule tool this -workflow needs already installed — `--container` then points the whole -workflow at it, replacing conda entirely so there's nothing left to race. +`container/Dockerfile` builds **one monolithic image**: a ROCm base, conda +installed inside it, and every per-rule conda env this workflow can reach +already built into the image. Convert it to a single `.sif` and each +`phables run` inside it finds its environments already there. + +There are no container-related CLI flags. Earlier versions had `--container` +and `--prostt5-container`, which pointed individual Snakemake rules at images +via `container:` directives; both are **removed**. The rules carry plain +`conda:` directives only, and the container satisfies them by having the envs +pre-built rather than by bypassing conda. + +## Why + +`--use-conda`'s per-rule env creation has a real failure mode at HPC scale. +Every rule's conda env is named by a hash of its env file, so every sample's +`phables run` needs the *same* env directory — but each sample is its own +separate Snakemake process (a SLURM array task), and Snakemake only serialises +env creation *within* one process's DAG. Two array tasks both needing the same +not-yet-built env at the same moment can race `mamba env create` into the same +directory and corrupt it (`Fatal Python error: init_fs_encoding ... no codec +search functions registered` is what that looks like in the wild). + +Baking the envs into the image removes the failure mode entirely: at runtime +there is nothing left to create, so nothing left to race. ## Getting an image -CI (`.github/workflows/build_container.yaml`) builds and pushes an image to -`quay.io//phables:` on every push and pull request -— tag by the exact commit you want, there's no floating `latest`. +CI (`.github/workflows/build_container.yaml`) builds and pushes to +`quay.io//phables:` on every push and pull request — +tag by the exact commit you want; there is no floating `latest`. -## Running with it +On Setonix, pull it once as a `.sif`: ```bash -phables run --input assembly_graph.gfa --reads fastq \ - --container quay.io/gbouras13/phables: \ - --no-use-conda \ - --use-singularity +module load singularity/4.1.0-slurm +singularity pull phables.sif docker://quay.io/gbouras13/phables: ``` -`--use-singularity` is Snakemake's own flag for running rules inside -Apptainer/Singularity (the standard container runtime on HPC, incl. Setonix) -— pass it through phables' existing `snake_args` passthrough, same as -`--dry-run`/`--keep-going`/etc. `--use-conda` alone does **not** activate -container execution; explicitly turn it off (`--no-use-conda`) when using -`--container`, since combining the two makes Snakemake build a *separate* -conda env *inside* the container instead of using what's already -installed there — exactly the redundant-torch-reinstall problem this image -exists to avoid. - -## What's actually in the image - -Base: `quay.io/pawsey/pytorch:2.7.1-rocm6.3.3` — the same verified-working -ROCm+PyTorch build `envs/prostt5-rocm.yaml` pins, and the same base phold's -own container (`../../phold/container/hpci/Dockerfile`) builds from. Nothing -in phables' own Dockerfile reinstalls or touches torch — `pip install -pholdlib` and everything else pure-Python go straight into that same system -Python, so they use the base image's already-correct torch/ROCm stack rather -than resolving a second, possibly conflicting one. This is the direct answer -to "don't rebuild the ProstT5 conda env for this container" — there's no -separate ProstT5 env in the container at all; predict_3di runs against the -same Python everything else does. - -Everything else — foldseek (direct binary, matching phold's own install), -minimap2/samtools/coverm/mmseqs2/FragGeneScan/HMMER/MAFFT (via a standalone -micromamba install into `/opt/conda`, kept deliberately separate from the -system Python so it can't touch it), and phables' own Python dependencies -(pyrodigal-gv, biopython, python-igraph, pysam, flowpaths, cogent3/piqtree, -...) — covers every `envs/*.yaml` this workflow has, so the whole DAG can run -from this one image with `--no-use-conda`. - -## `--container` vs `--prostt5-container` - -`--prostt5-container` already existed for pointing *just* `predict_3di` at a -container (e.g. phold's own image, which also has pholdlib+torch). It still -works exactly as before, and if set, it wins for that one rule. `--container` -is new and broader: it's the default container for **every** rule, including -predict_3di if `--prostt5-container` isn't also given. In practice, setting -`--container` alone is enough — there's no reason to use `--prostt5-container` -separately unless you specifically want predict_3di on a *different* image -than the rest of the workflow. +## Running it + +Run phables *inside* the container. Nothing special is passed to phables +itself — `--use-conda` is its own default and the envs are already present: + +```bash +singularity exec --rocm \ + -B /scratch/pawsey1018:/scratch/pawsey1018 \ + phables.sif \ + phables run --input assembly_graph.gfa --reads fastq \ + --output phables_out \ + --databases /scratch/.../all_databases/databases \ + --phagedetection prostt5-foldseek --gpu-backend rocm \ + --prostt5-checkpoint /scratch/.../model.pt \ + --threads 8 +``` + +Notes that matter on Setonix: + +- **`--rocm`** exposes the host's GPU devices to the container. Without it, + ProstT5 silently falls back to CPU. +- **Bind-mount your scratch** so databases, inputs and outputs are visible. + Databases are deliberately *not* in the image (see below). +- **Don't pass `--conda-prefix`.** The envs were built at the default prefix + (inside the installed package), and Snakemake resolves an env by hashing its + file content *together with the prefix path* — changing the prefix changes + the hash, and Snakemake would try to rebuild into a read-only filesystem. +- **`--gpu-backend rocm`** selects the prostt5-rocm env. The `cpu` env is also + prebuilt (it's the config default, so forgetting the flag still works); + `cuda` deliberately is not — this is a ROCm image. + +## What's in the image + +- **Base**: `quay.io/pawsey/pytorch:2.7.1-rocm6.3.3`, Pawsey's own verified + ROCm build, supplying the ROCm userspace matching Setonix's MI250X (gfx90a). + Its system-python torch (2.7.1) is *not* what the workflow uses and goes + untouched. +- **Miniforge** at `/opt/miniforge3`, plus Snakemake — the thing that reads the + prebuilt envs, so deliberately not one of them. +- **phables**, pip-installed from the build context (the commit CI tagged). +- **Every per-rule conda env**, prebuilt by `container/prebuild_envs.sh`: + coverm (minimap2/samtools/CoverM), genecall (FragGeneScan), pyrodigal-gv, + smg (HMMER), mmseqs, foldseek, phylotree (MAFFT + cogent3/piqtree), + prostt5-rocm (torch 2.9.1+rocm6.3 + pholdlib), prostt5-cpu, and curl for + `phables install`. + +**Databases are not included** — PHROGs and the hallmark DB are multi-GB and +separately versioned. Mount them and point `--databases` at them, exactly as +outside a container. `phables install` also works inside the image if you'd +rather fetch them from there. ## Building it yourself @@ -77,16 +95,29 @@ than the rest of the workflow. docker build -f container/Dockerfile -t phables:local . ``` -The base image is large (~14GB compressed) — building on a standard GitHub -Actions runner needs the disk-cleanup step already in -`build_container.yaml` (`jlumbroso/free-disk-space`), or the very first -`FROM` line fails with `no space left on device`. Building locally needs -proportionate free disk space too. - -**Not yet built or run for real** — the Dockerfile and the workflow-wide -`--container` plumbing were written and the Snakemake DAG wiring was verified -via real `--dry-run`s (with and without `--container` set, confirming both -resolve identically at the structural level), but no actual `docker build` -or real Apptainer/Singularity execution against Setonix hardware has -happened yet. Build it, push a first tag, and run one real sample through it -before trusting this for production batches. +`container/prebuild_envs.sh` does the env pre-building; it creates zero-byte +placeholder database files purely so the Snakemake DAG can *resolve*, then runs +`--conda-create-envs-only` once per flag combination (gene caller, detection +mode, GPU backend, tree) so that every reachable env gets built. The +placeholder trick is what keeps the databases out of the image; it was verified +by dry-running each of those combinations against zero-byte DB files. +`container/test_image.sh` then smoke-tests the result and fails the build if any +env is missing a binary the rules actually invoke, or if torch didn't come from +the prostt5-rocm env. + +**Disk**: this image is large — a ~14GB compressed ROCm base plus conda envs +including two torch builds. The CI workflow runs `jlumbroso/free-disk-space` +first because a stock GitHub Actions runner has only ~14GB free and the base +alone won't fit. Even with that, this is close to the limit of what a hosted +runner can build; if CI starts failing on `no space left on device`, building +on a machine with real disk and pushing manually is the fallback. + +## Status + +Not yet built or run for real. The Snakemake wiring was verified with real +dry-runs of every flag combination (and of the DAG resolving against +placeholder databases), the pinned conda packages were confirmed to exist for +`linux-64`/`noarch`, and both build scripts are syntax-checked — but no +`docker build` and no Setonix Apptainer run has happened yet. Build it, push a +tag, and put one real sample through it before trusting it for production +batches. diff --git a/docs/usage.md b/docs/usage.md index ee99c90..a5894de 100644 --- a/docs/usage.md +++ b/docs/usage.md @@ -18,14 +18,6 @@ Options: --use-conda / --no-use-conda Use conda for Snakemake rules [default: use-conda] --conda-prefix PATH Custom conda env directory - --container PATH container image with every per-rule tool - already installed (container/Dockerfile), - replacing --use-conda's per-rule env - creation for the WHOLE workflow -- not just - predict_3di (--prostt5-container). Needs - --use-singularity passed as a trailing - snakemake arg; don't also pass --use-conda - alongside this. --profile TEXT Snakemake profile --snake-default TEXT Customise Snakemake runtime args [default: --rerun-incomplete, --printshellcmds, @@ -112,13 +104,6 @@ Options: batch immediately [default: 4000] --prostt5-max-batch INTEGER max sequences per ProstT5 batch -- device- specific, tune per GPU [default: 20] - --prostt5-container TEXT container image with pholdlib + torch - already installed (e.g. phold's own image), - used instead of a conda env for predict_3di. - Needs --use-singularity passed as a trailing - snakemake arg -- --use-conda alone won't - honour it. Overrides --gpu-backend for this - rule. --evalue FLOAT maximum e-value for phrog annotations [default: 1e-10] --seqidentity FLOAT minimum sequence identity for phrog @@ -132,26 +117,26 @@ Options: --prefix TEXT prefix for genome identifier -h, --help Show this message and exit. - +  If you use Phables in your work, please cite Phables as, - +  Vijini Mallawaarachchi, Michael J Roach, Przemyslaw Decewicz, Bhavya Papudeshi, Sarah K Giles, Susanna R Grigson, George Bouras, Ryan D Hesse, Laura K Inglis, Abbey L K Hutton, Elizabeth A Dinsdale, Robert A Edwards, Phables: from fragmented assemblies to high-quality bacteriophage genomes, Bioinformatics, Volume 39, Issue 10, October 2023, btad586, https://doi.org/10.1093/bioinformatics/btad586 - - +  +  For more information on Phables please visit: https://phables.readthedocs.io/ - - +  +  CLUSTER EXECUTION: phables run ... --profile [profile] For information on Snakemake profiles see: https://snakemake.readthedocs.io/en/stable/executing/cli.html#profiles - +  RUN EXAMPLES: Required: phables run --input [assembly graph file] Specify threads: phables run ... --threads [threads] @@ -190,7 +175,6 @@ Options: * `--databases` - path to the databases directory [default: wherever `phables install` put them] * `--use-conda` / `--no-use-conda` - use conda for Snakemake rules [default: `use-conda`] * `--conda-prefix` - custom conda env directory -* `--container` - run the whole workflow from a single container image instead of per-rule conda envs (needs `--use-singularity`, and `--no-use-conda`) -- see [Running from a single container](container.md) * `--snake-default` - customise Snakemake runtime args [default: `--rerun-incomplete, --printshellcmds, --nolock, --show-failed-logs`] ### Phage-gene detection: `--phagedetection` @@ -213,7 +197,8 @@ Two independent methods for finding phage-like genes on unitigs, feeding the sam * `--prostt5-max-residues`, `--prostt5-max-seq-len`, `--prostt5-max-batch` - ProstT5 batching knobs, device-specific -- tune per GPU rather than trusting the defaults on unfamiliar hardware [defaults: 4000, 4000, 20] * `--gpu-backend` - which PyTorch build the ProstT5 conda env solves against: `cpu`, `cuda`, or `rocm` (e.g. Setonix's MI250X). Independent of `--prostt5-cpu`, which forces the CPU device at runtime even inside a GPU-capable env -- this controls which env gets *built* [default: `cpu`] * `--foldseek-gpu` - use Foldseek's CUDA GPU search mode for the hallmark scan. Requires `--gpu-backend cuda`, a CUDA-capable Foldseek build on `PATH` (not the plain bioconda package), and a `*_gpu`-suffixed, `makepaddedseqdb`-prepared hallmark DB -* `--prostt5-container` - use a container with `pholdlib`+torch already installed (e.g. phold's own image) instead of a conda env for the ProstT5 step. Needs `--use-singularity` passed as a trailing Snakemake arg -- `--use-conda` alone won't honour it. Overrides `--gpu-backend` for this rule. + +To run all of this from one prebuilt image on HPC (Setonix), see [Running from a single container](container.md) — the image ships every conda env already built, so there are no per-rule container flags to pass. ## Example usage diff --git a/phables/__main__.py b/phables/__main__.py index 4b575ae..3922a7a 100644 --- a/phables/__main__.py +++ b/phables/__main__.py @@ -93,20 +93,6 @@ def common_options(func): type=click.Path(), show_default=False, ), - click.option( - "--container", - default=None, - required=False, - help=( - "container image with every per-rule tool already installed " - "(container/Dockerfile), replacing --use-conda's per-rule env " - "creation for the WHOLE workflow -- not just predict_3di " - "(--prostt5-container). Needs --use-singularity passed as a " - "trailing snakemake arg; don't also pass --use-conda alongside " - "this." - ), - type=click.Path(), - ), click.option( "--profile", help="Snakemake profile", default=None, show_default=False ), @@ -355,18 +341,6 @@ def run_options(func): type=int, show_default=True, ), - click.option( - "--prostt5-container", - default=None, - required=False, - help=( - "container image with pholdlib + torch already installed (e.g. " - "phold's own image), used instead of a conda env for predict_3di. " - "Needs --use-singularity passed as a trailing snakemake arg -- " - "--use-conda alone won't honour it. Overrides --gpu-backend for " - "this rule." - ), - ), click.option( "--evalue", default=1e-10, diff --git a/phables/config/config.yaml b/phables/config/config.yaml index 5d6e637..658a157 100644 --- a/phables/config/config.yaml +++ b/phables/config/config.yaml @@ -59,33 +59,6 @@ prostt5_cpu: False # regardless of gpu_backend unless explicitly requested and gpu_backend=cuda. gpu_backend: cpu foldseek_gpu: False -# Optional: a single container image with EVERY per-rule tool this workflow -# needs already installed (e.g. "docker://quay.io/gbouras13/phables:" -- -# see container/Dockerfile). When set, every rule in the workflow runs with -# ONLY `container:` -- no `conda:` -- eliminating Snakemake's own per-rule -# --use-conda env creation entirely (and the real race condition that comes -# with it when many samples run as separate concurrent processes on an HPC -# array: multiple jobs needing the same not-yet-built shared env at once can -# corrupt each other's env-creation). Deliberately not combined with a -# `conda:` directive for the same reason as prostt5_container below: -# Snakemake's documented "conda inside container" feature builds a *separate* -# env inside the container rather than using what's already installed, which -# defeats the entire point here. When empty, every rule falls back to its own -# conda: env exactly as before. -# Needs `--use-singularity` (Apptainer/Singularity on most HPC, incl. Setonix) -# passed through phables' existing snake_args passthrough, since --use-conda -# alone won't honour the container: directive -- and --use-conda should -# generally NOT also be passed alongside this, for the reason above. -container: -# Optional: a container image (e.g. "docker://quay.io/gbouras13/phold:") -# that already has pholdlib + torch/ROCm baked in, used for JUST the -# predict_3di rule -- e.g. to point predict_3di at a different image than the -# rest of the workflow. Falls back to `container` above when unset, so -# setting `container` alone is enough to route predict_3di through it too; -# only set this separately if predict_3di specifically needs a different -# image. When both are empty, falls back to the conda envs selected by -# gpu_backend above. Same --use-singularity requirement as `container` above. -prostt5_container: # Batch limits are device-specific -- the defaults are conservative CPU/laptop-MPS # values, not tuned for any particular GPU. Retune per device (phold's autotune.py # is the reference approach) rather than trusting these at scale. diff --git a/phables/workflow/envs/prostt5-rocm.yaml b/phables/workflow/envs/prostt5-rocm.yaml index 02df6d5..0361d7d 100644 --- a/phables/workflow/envs/prostt5-rocm.yaml +++ b/phables/workflow/envs/prostt5-rocm.yaml @@ -13,22 +13,23 @@ dependencies: # so every Python the `python>=3.10` pin above can resolve is covered -- # 2.9.1 drops the cp39 wheel 2.7.1 had, which that pin already excludes). # - # NOTE this deliberately no longer matches the Pawsey base image the - # CONTAINER path uses (quay.io/pawsey/pytorch:2.7.1-rocm6.3.3, see - # container/Dockerfile, which bakes in and asserts torch 2.7.1). The two - # ROCm paths are now on different torch versions on purpose: the container - # takes whatever Pawsey verified, this conda env takes a newer torch. If - # you need them identical again, either rebuild the container from a newer - # Pawsey image or drop this back to 2.7.1. + # This env is the ONLY torch the workflow uses on ROCm, in both ways of + # running it. container/Dockerfile builds on Pawsey's ROCm base + # (quay.io/pawsey/pytorch:2.7.1-rocm6.3.3) for its ROCm userspace, but + # pre-builds THIS env inside the image and leaves the base's own + # system-python torch 2.7.1 untouched and unused -- so predict_3di runs on + # 2.9.1+rocm6.3 whether or not you're in the container. Keep it that way: + # the base image's torch version is an implementation detail of the base, + # not a version this workflow is pinned to. # - # Switched to this conda-env path instead of the --prostt5-container route - # after repeatedly hitting container-specific problems on Setonix that had - # nothing to do with phables/pholdlib itself: a broken torchaudio needing - # libomp.so, --home overriding the container's apparent CWD so relative - # output paths silently resolved to the wrong place, and a missing MIOpen - # AI-heuristic kernel-selection file for gfx90a. All real, all diagnosed, - # none of them things a plain conda env run directly on the host (no - # container namespace, no bind-mount surface) is exposed to at all. + # Phables drives ProstT5 through a plain conda env rather than pointing + # rules at a prebuilt torch container, after repeatedly hitting + # container-specific problems on Setonix that had nothing to do with + # phables/pholdlib itself: a broken torchaudio needing libomp.so, --home + # overriding the container's apparent CWD so relative output paths + # silently resolved to the wrong place, and a missing MIOpen AI-heuristic + # kernel-selection file for gfx90a. All real, all diagnosed, none of them + # things a conda env is exposed to. # >=2.6 would already avoid transformers' checkpoint loader refusing # torch.load() on older torch outright (CVE-2025-32434 hardening in # check_torch_load_is_safe, which broke loading Rostlab/ProstT5_fp16 -- no diff --git a/phables/workflow/rules/02_phables_preflight.smk b/phables/workflow/rules/02_phables_preflight.smk index 4d59bc1..dc6ef6d 100644 --- a/phables/workflow/rules/02_phables_preflight.smk +++ b/phables/workflow/rules/02_phables_preflight.smk @@ -36,21 +36,6 @@ GC = config['genecaller'] PD = config['phagedetection'] GPU_BACKEND = config['gpu_backend'] FOLDSEEK_GPU = config['foldseek_gpu'] -# CONTAINER_IMAGE: workflow-wide container (e.g. a phables release image with -# every per-rule tool baked in), replacing --use-conda's per-rule env creation -# entirely when set -- every rule below does -# `conda: None if CONTAINER_IMAGE else os.path.join(...)` + -# `container: CONTAINER_IMAGE`, so a real value here disables conda for ALL of -# them at once (not just predict_3di). Run with --use-singularity, no -# --use-conda, once this is set. -# -# PROSTT5_CONTAINER falls back to CONTAINER_IMAGE when not set explicitly -- -# --prostt5-container remains available to point JUST predict_3di at a -# different image (e.g. phold's own) than the rest of the workflow, but -# --container alone is now sufficient to route everything, predict_3di -# included, through one image. -PROSTT5_CONTAINER = config['prostt5_container'] or config['container'] -CONTAINER_IMAGE = config['container'] EV = config['evalue'] SI = config['seqidentity'] CT = config['covtol'] diff --git a/phables/workflow/rules/coverage.smk b/phables/workflow/rules/coverage.smk index 9d3908a..30be549 100644 --- a/phables/workflow/rules/coverage.smk +++ b/phables/workflow/rules/coverage.smk @@ -96,9 +96,7 @@ rule coverm_map: resources: mem_mb = config["resources"]["jobMem"] conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "coverm.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "coverm.yaml") log: os.path.join(LOGSDIR, "coverm_map.{sample}.log") shell: @@ -119,9 +117,7 @@ rule coverm_bam2counts: output: os.path.join(OUTDIR, "preprocess", "temp", "{sample}.cov") conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "coverm.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "coverm.yaml") log: os.path.join(LOGSDIR, "coverm_bam2counts.{sample}.log") shell: diff --git a/phables/workflow/rules/genes.smk b/phables/workflow/rules/genes.smk index 2dce368..1a2afad 100644 --- a/phables/workflow/rules/genes.smk +++ b/phables/workflow/rules/genes.smk @@ -25,9 +25,7 @@ if GC == "pyrodigal-gv": log: os.path.join(LOGSDIR, "gene_call_pyrodigal_gv.log") conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "genecall.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "genecall.yaml") script: os.path.join("..", "scripts", "gene_caller.py") @@ -48,9 +46,7 @@ else: out = os.path.join(LOGSDIR, "gene_call_fraggenescan_out.log"), err = os.path.join(LOGSDIR, "gene_call_fraggenescan_err.log"), conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "smg.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "smg.yaml") shell: """ run_FragGeneScan.pl -genome={input.genome} -out={params.frag} -complete=0 -train=complete -thread={threads} 1>{log.out} 2>{log.err} @@ -71,9 +67,7 @@ rule scan_smg: hmm_out=os.path.join(LOGSDIR, "smg_scan_hmm_out.log"), hmm_err=os.path.join(LOGSDIR, "smg_scan_hmm_err.log") conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "smg.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "smg.yaml") shell: """ hmmsearch --domtblout {output.hmmout} --cut_tc --cpu {threads} {input.hmm} {input.faa} 1>{log.hmm_out} 2> {log.hmm_err} @@ -98,9 +92,7 @@ rule scan_phrogs: log: os.path.join(LOGSDIR, "phrogs_scan.log") conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "mmseqs.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "mmseqs.yaml") shell: """ mkdir -p {params.out_path} @@ -135,74 +127,35 @@ if PD == "prostt5-foldseek": QUERY_3DI = os.path.join(OUTDIR, "preprocess", "hallmark", "proteins_3di.fasta") - if PROSTT5_CONTAINER: - - # container:-only, deliberately with no conda: alongside it. Snakemake's - # documented way to combine the two ("Ad-hoc combination of Conda package - # management with containers") builds a fresh, isolated conda env *inside* - # the container rather than exposing what the image already has installed - # -- which would just reinstall a second, redundant torch and ignore the - # image's own verified-working one (e.g. phold's own container, which - # already bundles pholdlib + a working torch/ROCm stack for Setonix, per - # the same logic phold's own Snakemake rules use). Needs - # `--use-singularity` passed through phables' snake_args passthrough -- - # --use-conda alone won't honour this directive. - rule predict_3di: - input: - faa = PROTEINS_FILE, - threads: - config["resources"]["jobCPU"] - resources: - mem_mb = config["resources"]["jobMem"] - output: - threedi = QUERY_3DI - params: - checkpoint = config["prostt5_checkpoint"], - model_name = config["prostt5_model"], - model_dir = config["prostt5_model_dir"], - half_precision = config["prostt5_half_precision"], - cpu = config["prostt5_cpu"], - max_residues = config["prostt5_max_residues"], - max_seq_len = config["prostt5_max_seq_len"], - max_batch = config["prostt5_max_batch"], - log: - os.path.join(LOGSDIR, "predict_3di.log") - container: - PROSTT5_CONTAINER - script: - os.path.join("..", "scripts", "predict_3di.py") - - else: - - rule predict_3di: - input: - faa = PROTEINS_FILE, - threads: - config["resources"]["jobCPU"] - resources: - mem_mb = config["resources"]["jobMem"] - output: - threedi = QUERY_3DI - params: - checkpoint = config["prostt5_checkpoint"], - model_name = config["prostt5_model"], - model_dir = config["prostt5_model_dir"], - half_precision = config["prostt5_half_precision"], - cpu = config["prostt5_cpu"], - max_residues = config["prostt5_max_residues"], - max_seq_len = config["prostt5_max_seq_len"], - max_batch = config["prostt5_max_batch"], - log: - os.path.join(LOGSDIR, "predict_3di.log") - conda: - # gpu_backend selects which torch build this env solves against -- - # cpu/cuda/rocm need different PyTorch wheels (conda envs are - # solved once from a static file, so this has to be three files, - # not one file with a runtime switch). See the individual env - # files for what each backend actually needs and why. - os.path.join("..", "envs", f"prostt5-{GPU_BACKEND}.yaml") - script: - os.path.join("..", "scripts", "predict_3di.py") + rule predict_3di: + input: + faa = PROTEINS_FILE, + threads: + config["resources"]["jobCPU"] + resources: + mem_mb = config["resources"]["jobMem"] + output: + threedi = QUERY_3DI + params: + checkpoint = config["prostt5_checkpoint"], + model_name = config["prostt5_model"], + model_dir = config["prostt5_model_dir"], + half_precision = config["prostt5_half_precision"], + cpu = config["prostt5_cpu"], + max_residues = config["prostt5_max_residues"], + max_seq_len = config["prostt5_max_seq_len"], + max_batch = config["prostt5_max_batch"], + log: + os.path.join(LOGSDIR, "predict_3di.log") + conda: + # gpu_backend selects which torch build this env solves against -- + # cpu/cuda/rocm need different PyTorch wheels (conda envs are + # solved once from a static file, so this has to be three files, + # not one file with a runtime switch). See the individual env + # files for what each backend actually needs and why. + os.path.join("..", "envs", f"prostt5-{GPU_BACKEND}.yaml") + script: + os.path.join("..", "scripts", "predict_3di.py") rule build_hallmark_query_db: @@ -216,9 +169,7 @@ if PD == "prostt5-foldseek": log: os.path.join(LOGSDIR, "build_hallmark_query_db.log") conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "foldseek.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "foldseek.yaml") script: os.path.join("..", "scripts", "build_foldseek_query_db.py") @@ -271,9 +222,7 @@ if PD == "prostt5-foldseek": log: os.path.join(LOGSDIR, "scan_hallmark.log") conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "foldseek.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "foldseek.yaml") shell: """ foldseek search {params.query_prefix} {input.hallmark_db} {params.result} {params.tmp} \ diff --git a/phables/workflow/rules/gfa2fasta.smk b/phables/workflow/rules/gfa2fasta.smk index 531aa2c..2e43c07 100644 --- a/phables/workflow/rules/gfa2fasta.smk +++ b/phables/workflow/rules/gfa2fasta.smk @@ -15,8 +15,6 @@ rule run_gfa2fasta: log: os.path.join(LOGSDIR, "gfa2fasta.log") conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "phables.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "phables.yaml") script: os.path.join('..', 'scripts', 'gfa2fasta.py') \ No newline at end of file diff --git a/phables/workflow/rules/phables.smk b/phables/workflow/rules/phables.smk index 9a1ca07..549d454 100644 --- a/phables/workflow/rules/phables.smk +++ b/phables/workflow/rules/phables.smk @@ -42,8 +42,6 @@ rule run_phables: log: os.path.join(LOGSDIR, "phables_output.log") conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "phables.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "phables.yaml") script: os.path.join("..", "scripts", "phables.py") diff --git a/phables/workflow/rules/phylotree.smk b/phables/workflow/rules/phylotree.smk index 0f61961..caae230 100644 --- a/phables/workflow/rules/phylotree.smk +++ b/phables/workflow/rules/phylotree.smk @@ -11,9 +11,7 @@ rule build_msa: log: os.path.join(LOGSDIR, "mafft_output.log") conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "phylotree.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "phylotree.yaml") shell: """ mafft --auto --thread {threads} {input} > {output} @@ -36,8 +34,6 @@ rule build_tree: log: os.path.join(LOGSDIR, "piqtree_output.log") conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "phylotree.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "phylotree.yaml") script: os.path.join("..", "scripts", "phylotree.py") \ No newline at end of file diff --git a/phables/workflow/rules/postprocess.smk b/phables/workflow/rules/postprocess.smk index 1d33b02..438f55a 100644 --- a/phables/workflow/rules/postprocess.smk +++ b/phables/workflow/rules/postprocess.smk @@ -75,9 +75,7 @@ rule coverm_map_genomes: resources: mem_mb = config["resources"]["jobMem"] conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "coverm.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "coverm.yaml") log: os.path.join(LOGSDIR, "coverm_map_genomes.{sample}.log") shell: @@ -123,9 +121,7 @@ rule coverm_bam2counts_genomes: output: temp(os.path.join(OUTDIR, "postprocess", "temp", "{sample}.cov")) conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "coverm.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "coverm.yaml") log: os.path.join(LOGSDIR, "coverm_bam2counts_genomes.{sample}.log") shell: @@ -186,8 +182,6 @@ rule format_genome_coverage: log: os.path.join(LOGSDIR, "format_koverage_results_output.log") conda: - None if CONTAINER_IMAGE else os.path.join("..", "envs", "phables.yaml") - container: - CONTAINER_IMAGE + os.path.join("..", "envs", "phables.yaml") script: os.path.join("..", "scripts", "format_koverage_results.py") \ No newline at end of file From b2781ba45142159428a9ee169382ebc63e48c225 Mon Sep 17 00:00:00 2001 From: gbouras13 Date: Fri, 14 Aug 2026 08:34:52 +0930 Subject: [PATCH 04/12] container: fix two real build failures in env pre-build 1. tests/data/ is gitignored, so it doesn't exist in a fresh clone or CI checkout -- the build died at the pre-build step with 'Invalid value for --reads: Path ... does not exist'. prebuild_envs.sh now generates its own throwaway inputs (a two-segment GFA and a pair of tiny gzipped FASTQs) instead of depending on untracked data. Verified these produce an identical DAG to the real test data. 2. 'phables install' was pointed at the placeholder databases dir, whose files satisfy install.smk's own download targets -- so Snakemake said 'Nothing to be done', the DAG was empty, and the curl env would have been silently missing from the image (with test_image.sh's own 'check_bin curl' the only thing that would have caught it). It now uses an empty databases dir so the four *_download rules are in the DAG and their env actually gets built. prebuild_envs.sh has now been run end-to-end with its --conda-create-envs-only calls swapped for dry-runs: all six invocations succeed, the synthetic inputs resolve, and cleanup fires. The conda solves themselves and test_image.sh against a real image remain unverified. --- container/Dockerfile | 2 +- container/prebuild_envs.sh | 38 +++++++++++++++++++++++++++++------- docs/container.md | 40 +++++++++++++++++++++++++------------- 3 files changed, 59 insertions(+), 21 deletions(-) diff --git a/container/Dockerfile b/container/Dockerfile index 5040106..56f5c55 100644 --- a/container/Dockerfile +++ b/container/Dockerfile @@ -66,7 +66,7 @@ RUN pip install /opt/phables_src # Pre-build every per-rule conda env. See the script for the full rationale # (default --conda-prefix, placeholder databases, one pass per flag combo). -RUN bash /opt/phables_src/container/prebuild_envs.sh /opt/phables_src +RUN bash /opt/phables_src/container/prebuild_envs.sh # Fail the build loudly rather than pushing a broken image to quay.io. RUN bash /opt/phables_src/container/test_image.sh diff --git a/container/prebuild_envs.sh b/container/prebuild_envs.sh index d0a2f35..9411498 100755 --- a/container/prebuild_envs.sh +++ b/container/prebuild_envs.sh @@ -1,6 +1,8 @@ #!/bin/bash # Pre-builds EVERY per-rule conda env into the image, at Docker build time. -# Run from container/Dockerfile; not meant to be run on a host. +# Run from container/Dockerfile; not meant to be run on a host. Takes no +# arguments and reads nothing from the repo -- it generates its own throwaway +# inputs (see below), so it works in a bare CI checkout. # # Envs go into phables' own DEFAULT --conda-prefix (snake_base("workflow/conda"), # i.e. inside the installed package) deliberately: a user inside the resulting @@ -28,8 +30,8 @@ set -euxo pipefail -SRC="${1:-/opt/phables_src}" DB=/tmp/placeholder_db +WORK=/tmp/envbuild_inputs mkdir -p "$DB/phrogs_mmseqs_db" "$DB/hallmark_db" touch "$DB/marker.hmm" \ @@ -38,8 +40,22 @@ touch "$DB/marker.hmm" \ "$DB/hallmark_db/hallmark_db" \ "$DB/hallmark_db/hallmark_categories.tsv" -GFA="${SRC}/tests/data/ERR1301161/assembly_graph_after_simplification.gfa" -READS="${SRC}/tests/data/ERR1301161" +# Synthetic minimal inputs, generated here rather than taken from +# tests/data/: that directory is in .gitignore, so it does NOT exist in a +# fresh clone or in a CI checkout -- depending on it made the image build +# fail with "Invalid value for '--reads': Path ... does not exist". Nothing +# below is ever actually processed (no job runs under +# --conda-create-envs-only); these files exist purely so click's exists=True +# checks pass and the DAG can resolve. Verified to produce the same DAG as +# the real test data. +mkdir -p "$WORK/reads" +printf 'H\tVN:Z:1.0\nS\tedge_1\tACGTACGTACGTACGTACGTACGTACGTACGT\tLN:i:32\nS\tedge_2\tTTTTGGGGCCCCAAAATTTTGGGGCCCCAAAA\tLN:i:32\nL\tedge_1\t+\tedge_2\t+\t0M\n' \ + > "$WORK/assembly_graph.gfa" +printf '@r1\nACGT\n+\nIIII\n' | gzip > "$WORK/reads/sample1_R1.fastq.gz" +printf '@r1\nACGT\n+\nIIII\n' | gzip > "$WORK/reads/sample1_R2.fastq.gz" + +GFA="$WORK/assembly_graph.gfa" +READS="$WORK/reads" COMMON=(--input "$GFA" --reads "$READS" --databases "$DB" --threads 1) # 1. default path -> coverm, genecall (FragGeneScan), smg (HMMER), mmseqs, phables @@ -67,10 +83,18 @@ phables run "${COMMON[@]}" --output /tmp/envbuild4 \ phables run "${COMMON[@]}" --output /tmp/envbuild5 \ --build-tree --conda-create-envs-only -# 6. install.smk's own env (curl), so `phables install` works inside here too -phables install --output /tmp/envbuild6 --databases "$DB" --conda-create-envs-only +# 6. install.smk's own env (curl), so `phables install` works inside here too. +# Deliberately pointed at an EMPTY databases dir, not "$DB": the placeholder +# files in $DB satisfy install.smk's own download targets, so Snakemake says +# "Nothing to be done", the DAG is empty, and NO env gets created -- the +# curl env would silently be missing from the image. An empty dir puts the +# four *_download rules in the DAG so their conda env actually gets built. +# (--conda-create-envs-only still downloads nothing.) +mkdir -p /tmp/empty_db +phables install --output /tmp/envbuild6 --databases /tmp/empty_db --conda-create-envs-only -rm -rf /tmp/envbuild1 /tmp/envbuild2 /tmp/envbuild3 /tmp/envbuild4 /tmp/envbuild5 /tmp/envbuild6 "$DB" +rm -rf /tmp/envbuild1 /tmp/envbuild2 /tmp/envbuild3 /tmp/envbuild4 /tmp/envbuild5 /tmp/envbuild6 \ + "$DB" "$WORK" /tmp/empty_db conda clean -a -y echo "=== pre-built conda envs ===" diff --git a/docs/container.md b/docs/container.md index 83cfabd..b30fa77 100644 --- a/docs/container.md +++ b/docs/container.md @@ -95,12 +95,22 @@ rather fetch them from there. docker build -f container/Dockerfile -t phables:local . ``` -`container/prebuild_envs.sh` does the env pre-building; it creates zero-byte -placeholder database files purely so the Snakemake DAG can *resolve*, then runs -`--conda-create-envs-only` once per flag combination (gene caller, detection -mode, GPU backend, tree) so that every reachable env gets built. The -placeholder trick is what keeps the databases out of the image; it was verified -by dry-running each of those combinations against zero-byte DB files. +`container/prebuild_envs.sh` does the env pre-building. It generates its own +throwaway inputs — a two-segment GFA and a pair of tiny gzipped FASTQs, plus +zero-byte placeholder database files — purely so the Snakemake DAG can +*resolve*, then runs `--conda-create-envs-only` once per flag combination (gene +caller, detection mode, GPU backend, tree) so every reachable env gets built. +Nothing is ever processed, since no job runs. + +Two things it deliberately does **not** do, both of which broke a real build: + +- It doesn't use `tests/data/`. That directory is in `.gitignore`, so it + doesn't exist in a fresh clone or a CI checkout — depending on it failed with + `Invalid value for '--reads': Path ... does not exist`. +- It points `phables install` at an *empty* databases directory, not the + placeholder one. The placeholders satisfy install.smk's own download targets, + so Snakemake reports "Nothing to be done", the DAG is empty and the `curl` + env is silently never built. `container/test_image.sh` then smoke-tests the result and fails the build if any env is missing a binary the rules actually invoke, or if torch didn't come from the prostt5-rocm env. @@ -114,10 +124,14 @@ on a machine with real disk and pushing manually is the fallback. ## Status -Not yet built or run for real. The Snakemake wiring was verified with real -dry-runs of every flag combination (and of the DAG resolving against -placeholder databases), the pinned conda packages were confirmed to exist for -`linux-64`/`noarch`, and both build scripts are syntax-checked — but no -`docker build` and no Setonix Apptainer run has happened yet. Build it, push a -tag, and put one real sample through it before trusting it for production -batches. +A real `docker build` got as far as the env pre-build step before failing on +the two issues listed above; both are fixed, and `prebuild_envs.sh` has since +been run end-to-end (with its `--conda-create-envs-only` calls swapped for +dry-runs) so that all six invocations, the synthetic input generation and the +cleanup are known to work as written. The pinned conda packages were confirmed +to exist for `linux-64`/`noarch`. + +Still unverified: a complete `docker build` (the conda solves themselves, and +`test_image.sh` against a real image), and any Setonix Apptainer run. Build it, +push a tag, and put one real sample through it before trusting it for +production batches. From f3a717e2ed1439ec8c76f133f81389a1a5eeb01e Mon Sep 17 00:00:00 2001 From: gbouras13 Date: Fri, 14 Aug 2026 08:51:50 +0930 Subject: [PATCH 05/12] Reuse the base image's ROCm torch instead of installing a second one Adds --gpu-backend system: predict_3di declares NO conda env and runs against the torch + pholdlib already in the ambient python. Conda envs are isolated, so a rule declaring conda: can never see a torch installed outside its env -- omitting the directive is the only way to reuse an existing, known-good GPU torch. Verified via --list-conda-envs: system declares 5 envs, rocm declares 6 (the extra being prostt5-rocm). The container now uses it, which required inverting how the image is built: phables, Snakemake and pholdlib install into the BASE IMAGE'S python (the one holding Pawsey's verified ROCm torch), not miniforge's. Snakemake runs a script: rule with no conda env using its own interpreter, so that interpreter has to be the one with torch. Miniforge is now APPENDED to PATH, never prepended, so it cannot shadow the base python -- prepending it is exactly why the previous build ran phables under a torch-less python and needed its own torch env. Consequences: no prostt5-* env is built at all (the prostt5-cpu build is gone too), the image drops several GB and one whole conda solve, and it runs the torch already verified on gfx90a rather than an unverified 2.9.1. The image sets its own config default to gpu_backend: system so a plain works with no extra flag -- without that, the config default (cpu) would look for an env the image deliberately lacks, which is fatal on a read-only .sif. The build asserts torch's version is unchanged across the pip install, and test_image.sh now fails if any pre-built env contains a torch of its own -- i.e. if a second copy got installed after all. prostt5-rocm.yaml is untouched and still pins 2.9.1 for non-container use. --- container/Dockerfile | 86 ++++++++++++++++++++++---------- container/prebuild_envs.sh | 27 +++++----- container/test_image.sh | 55 ++++++++++---------- docs/container.md | 70 +++++++++++++++----------- docs/usage.md | 15 ++++-- phables/__main__.py | 9 +++- phables/config/config.yaml | 9 ++++ phables/workflow/rules/genes.smk | 19 ++++--- 8 files changed, 180 insertions(+), 110 deletions(-) diff --git a/container/Dockerfile b/container/Dockerfile index 56f5c55..8d29e96 100644 --- a/container/Dockerfile +++ b/container/Dockerfile @@ -1,9 +1,9 @@ # # phables -- ONE monolithic image for Setonix: ROCm base + conda + every # per-rule conda env pre-built INSIDE the image. Build it once, convert it to a -# single .sif, and every `phables run` inside it finds its conda envs already -# there. Same strategy as hybracter's own container (miniforge installed into -# the image, envs materialised at build time), NOT the per-rule `container:` +# single .sif, and every `phables run` inside it finds its envs already there. +# Same strategy as hybracter's own container (miniforge installed into the +# image, envs materialised at build time), NOT the per-rule `container:` # directive approach -- phables has no --container/--prostt5-container flags # any more, and its Snakemake rules carry plain `conda:` directives only. # @@ -13,15 +13,24 @@ # Many array tasks all needing the same not-yet-built shared env at once # corrupt each other's `mamba env create` ("Fatal Python error: # init_fs_encoding ... no codec search functions registered"). With every env -# already baked into the image there is nothing left to create at runtime, and -# so nothing left to race. +# baked into the image there is nothing left to create at runtime, and so +# nothing left to race. # -# Base is Pawsey's own verified-working ROCm image, as requested: it supplies -# the ROCm userspace matching Setonix's MI250X (gfx90a) kernel driver. NOTE -# that the base's own system-python torch (2.7.1) is NOT what the workflow -# uses -- predict_3di runs inside the prostt5-rocm conda env built below, which -# brings its own torch (2.9.1+rocm6.3, see workflow/envs/prostt5-rocm.yaml). -# The base image's torch is left untouched and simply goes unused. +# TORCH: this image does NOT install a second PyTorch. The Pawsey ROCm base +# already ships a torch that is known to work on Setonix's MI250X (gfx90a), so +# predict_3di reuses it via `--gpu-backend system` (see below) instead of a +# prostt5-* conda env. That is the ONLY way to reuse it: conda envs are +# isolated, so a rule declaring `conda:` can never see a torch installed +# outside its env -- which is exactly why an earlier version of this image +# ended up downloading multiple GB of a second, unverified torch. +# +# The consequence for how the image is built: phables and Snakemake are +# installed into the BASE IMAGE'S python, not into miniforge's. Snakemake runs +# `script:` rules that declare no conda env using its own interpreter, so that +# interpreter has to be the one holding torch. Miniforge is installed only to +# provide the `conda` binary that BUILDS the other per-rule envs, and is +# appended to PATH rather than prepended so it can never shadow the base +# python. # # See ../docs/container.md for how to build the .sif and run it on Setonix. # @@ -37,32 +46,55 @@ RUN apt-get update && apt-get install -y --no-install-recommends \ procps \ && rm -rf /var/lib/apt/lists/* +# Pin down the base image's python BEFORE miniforge exists, so later steps can +# target it unambiguously no matter where the base keeps it. Everything phables +# needs goes into this interpreter. +RUN command -v python3 > /opt/base_python_path && \ + ln -sf "$(cat /opt/base_python_path)" /usr/local/bin/python && \ + echo "base python: $(cat /opt/base_python_path)" && \ + "$(cat /opt/base_python_path)" -c "import torch; print('base torch:', torch.__version__)" + # -# Miniforge inside the image -- hybracter's containers/HPC/Dockerfile does -# exactly this (wget the installer, `bash ... -b -p `) rather than -# relying on a conda that may or may not exist on the host. conda-forge is -# miniforge's default channel, so no extra channel config is needed. +# Miniforge -- installed ONLY for its `conda`/`mamba` binaries, which Snakemake +# shells out to when building the per-rule envs. Deliberately APPENDED to PATH: +# prepending it would make miniforge's python the default, and phables would +# then run under an interpreter with no torch in it, silently defeating the +# whole point of building on a ROCm base. # ENV CONDA_DIR=/opt/miniforge3 RUN wget -q https://github.com/conda-forge/miniforge/releases/latest/download/Miniforge3-Linux-x86_64.sh -O /tmp/miniforge.sh && \ bash /tmp/miniforge.sh -b -p "${CONDA_DIR}" && \ rm /tmp/miniforge.sh -ENV PATH="${CONDA_DIR}/bin:${PATH}" - -# Snakemake lives in the base conda env alongside phables itself -- it is the -# thing that READS the pre-built per-rule envs, so it must not be one of them. -RUN conda install -y -c conda-forge -c bioconda \ - "snakemake>=8" \ - && conda clean -a -y +ENV PATH="${PATH}:${CONDA_DIR}/bin" # -# phables from the exact source in this build context. CI builds this image per -# commit and tags it by commit SHA, so this is the checked-out tree -- not a -# separately-pinned git ref that could drift from the rules and env files the -# image bakes in. +# phables + Snakemake + pholdlib into the BASE python. pholdlib declares no +# torch version constraint, so pip pulls its own dependencies (transformers, +# sentencepiece, h5py, ...) without touching the already-installed torch -- the +# assertion below fails the build if that ever stops being true, since a +# silently-replaced torch is precisely the failure this image exists to avoid. # COPY . /opt/phables_src -RUN pip install /opt/phables_src +RUN BASE_PY="$(cat /opt/base_python_path)" && \ + TORCH_BEFORE="$("$BASE_PY" -c 'import torch; print(torch.__version__)')" && \ + "$BASE_PY" -m pip install --no-cache-dir "snakemake>=8" pholdlib /opt/phables_src && \ + TORCH_AFTER="$("$BASE_PY" -c 'import torch; print(torch.__version__)')" && \ + echo "torch before=${TORCH_BEFORE} after=${TORCH_AFTER}" && \ + test "$TORCH_BEFORE" = "$TORCH_AFTER" + +# +# Make `system` the image's own default GPU backend, so a plain +# `phables run --phagedetection prostt5-foldseek` inside the container reuses +# the base torch with no extra flag. Without this the config default (`cpu`) +# would send predict_3di looking for a prostt5-cpu env that this image +# deliberately does not contain -- and a .sif is read-only, so that is a hard +# failure rather than a slow rebuild. +# +RUN BASE_PY="$(cat /opt/base_python_path)" && \ + CONFIG="$("$BASE_PY" -c 'import phables, os; print(os.path.join(os.path.dirname(phables.__file__), "config", "config.yaml"))')" && \ + sed -i 's/^gpu_backend: .*/gpu_backend: system/' "$CONFIG" && \ + grep -q '^gpu_backend: system$' "$CONFIG" && \ + echo "image default gpu_backend set to system in $CONFIG" # Pre-build every per-rule conda env. See the script for the full rationale # (default --conda-prefix, placeholder databases, one pass per flag combo). diff --git a/container/prebuild_envs.sh b/container/prebuild_envs.sh index 9411498..abbd310 100755 --- a/container/prebuild_envs.sh +++ b/container/prebuild_envs.sh @@ -65,25 +65,22 @@ phables run "${COMMON[@]}" --output /tmp/envbuild1 --conda-create-envs-only phables run "${COMMON[@]}" --output /tmp/envbuild2 \ --genecaller pyrodigal-gv --conda-create-envs-only -# 3. ProstT5 + foldseek detection on ROCm -- the reason this image exists +# 3. ProstT5 + foldseek detection -> the foldseek env (and, via +# --gpu-backend system, NO prostt5-* torch env at all: predict_3di reuses the +# base image's already-working ROCm torch, which is installed in the same +# python running Snakemake here). This is the whole reason no multi-GB +# second torch is downloaded during this build. Deliberately NOT `rocm`/ +# `cpu`/`cuda` -- each of those would solve and download their own torch. phables run "${COMMON[@]}" --output /tmp/envbuild3 \ - --phagedetection prostt5-foldseek --gpu-backend rocm --conda-create-envs-only + --phagedetection prostt5-foldseek --gpu-backend system --conda-create-envs-only -# 4. Same, CPU backend. Built because gpu_backend's own config default is `cpu`: -# without this, forgetting `--gpu-backend rocm` inside the container would hit -# a missing env on a read-only filesystem. (cuda is deliberately NOT built -- -# this is a ROCm image for Setonix, and the CUDA wheels would add several GB -# that could never be used on this hardware.) -phables run "${COMMON[@]}" --output /tmp/envbuild4 \ - --phagedetection prostt5-foldseek --gpu-backend cpu --conda-create-envs-only - -# 5. Optional phylogenetic tree -> phylotree env (MAFFT + cogent3/piqtree). +# 4. Optional phylogenetic tree -> phylotree env (MAFFT + cogent3/piqtree). # Note conda resolves cogent3 fine, unlike pip, where every published # release is a prerelease and a plain version range matches nothing. -phables run "${COMMON[@]}" --output /tmp/envbuild5 \ +phables run "${COMMON[@]}" --output /tmp/envbuild4 \ --build-tree --conda-create-envs-only -# 6. install.smk's own env (curl), so `phables install` works inside here too. +# 5. install.smk's own env (curl), so `phables install` works inside here too. # Deliberately pointed at an EMPTY databases dir, not "$DB": the placeholder # files in $DB satisfy install.smk's own download targets, so Snakemake says # "Nothing to be done", the DAG is empty, and NO env gets created -- the @@ -91,9 +88,9 @@ phables run "${COMMON[@]}" --output /tmp/envbuild5 \ # four *_download rules in the DAG so their conda env actually gets built. # (--conda-create-envs-only still downloads nothing.) mkdir -p /tmp/empty_db -phables install --output /tmp/envbuild6 --databases /tmp/empty_db --conda-create-envs-only +phables install --output /tmp/envbuild5 --databases /tmp/empty_db --conda-create-envs-only -rm -rf /tmp/envbuild1 /tmp/envbuild2 /tmp/envbuild3 /tmp/envbuild4 /tmp/envbuild5 /tmp/envbuild6 \ +rm -rf /tmp/envbuild1 /tmp/envbuild2 /tmp/envbuild3 /tmp/envbuild4 /tmp/envbuild5 \ "$DB" "$WORK" /tmp/empty_db conda clean -a -y diff --git a/container/test_image.sh b/container/test_image.sh index 6737a30..9d5d602 100755 --- a/container/test_image.sh +++ b/container/test_image.sh @@ -31,8 +31,9 @@ echo "=== pre-built env count ===" n=$(find "$PREFIX" -maxdepth 1 -mindepth 1 -type d | wc -l) echo "found $n env directories" find "$PREFIX" -maxdepth 1 -mindepth 1 -type d -exec basename {} \; -# 6 distinct env files are reachable at minimum (coverm, genecall, smg, mmseqs, -# phables, curl) before counting foldseek/prostt5-rocm/prostt5-cpu/phylotree. +# 6 distinct env files are reachable at minimum: coverm, genecall, smg, mmseqs, +# phables, curl -- before counting foldseek and phylotree. No prostt5-* env is +# expected (or wanted): gpu_backend=system reuses the base image's torch. test "$n" -ge 6 # Finds an executable in any pre-built env; fails if no env provides it. @@ -61,34 +62,32 @@ check_bin run_FragGeneScan.pl check_bin mafft check_bin curl -echo "=== torch must come from the prostt5-rocm env, not the base image ===" -# The base image ships its own system-python torch 2.7.1, which the workflow -# does NOT use. predict_3di runs inside the prostt5-rocm conda env, whose torch -# is 2.9.1+rocm6.3 (workflow/envs/prostt5-rocm.yaml). If nothing here reports -# 2.9.1+rocm, that env didn't build properly and predict_3di would be running -# on the wrong stack. -found=0 +echo "=== torch must come from the BASE image, reused -- not reinstalled ===" +# predict_3di runs with NO conda env (gpu_backend=system), i.e. in the same +# python that runs Snakemake -- which must therefore be the base image's python, +# the one already holding a working ROCm torch. Both halves of that are checked: +# torch+pholdlib importable here, AND no conda env carrying a torch of its own. +python -c "import torch; print('torch:', torch.__version__)" +python -c "import torch, sys; sys.exit(0 if 'rocm' in torch.__version__ else 1)" \ + || { echo "ERROR: the ambient torch is not a ROCm build" >&2; exit 1; } +python -c "import pholdlib; print('pholdlib OK')" +python -c "import phables, snakemake; print('phables + snakemake share this interpreter')" + +# A torch inside any pre-built env means a second copy got installed after all, +# which is the exact regression this layout exists to prevent. for e in "$PREFIX"/*/; do - if [ -x "${e}bin/python" ]; then - if "${e}bin/python" - <<'PY' 2>/dev/null -import sys -try: - import torch -except Exception: - sys.exit(1) -sys.exit(0 if torch.__version__.startswith("2.9.1") and "rocm" in torch.__version__ else 1) -PY - then - echo "OK: torch 2.9.1+rocm in $e" - "${e}bin/python" -c "import pholdlib; print('pholdlib OK')" - found=1 - break - fi + if [ -x "${e}bin/python" ] && "${e}bin/python" -c "import torch" 2>/dev/null; then + echo "ERROR: a pre-built conda env contains its own torch: $e" >&2 + echo " gpu_backend=system should mean no prostt5-* env is built." >&2 + exit 1 fi done -if [ "$found" -ne 1 ]; then - echo "ERROR: no pre-built env provides torch 2.9.1+rocm" >&2 - exit 1 -fi +echo "OK: no pre-built env ships a duplicate torch" + +echo "=== image's default gpu_backend must be 'system' ===" +CONFIG="$(python -c 'import phables, os; print(os.path.join(os.path.dirname(phables.__file__), "config", "config.yaml"))')" +grep -q '^gpu_backend: system$' "$CONFIG" \ + || { echo "ERROR: $CONFIG does not default gpu_backend to system" >&2; exit 1; } +echo "OK: $CONFIG defaults to system" echo "=== all image tests passed ===" diff --git a/docs/container.md b/docs/container.md index b30fa77..7858020 100644 --- a/docs/container.md +++ b/docs/container.md @@ -50,7 +50,7 @@ singularity exec --rocm \ phables run --input assembly_graph.gfa --reads fastq \ --output phables_out \ --databases /scratch/.../all_databases/databases \ - --phagedetection prostt5-foldseek --gpu-backend rocm \ + --phagedetection prostt5-foldseek \ --prostt5-checkpoint /scratch/.../model.pt \ --threads 8 ``` @@ -65,24 +65,30 @@ Notes that matter on Setonix: (inside the installed package), and Snakemake resolves an env by hashing its file content *together with the prefix path* — changing the prefix changes the hash, and Snakemake would try to rebuild into a read-only filesystem. -- **`--gpu-backend rocm`** selects the prostt5-rocm env. The `cpu` env is also - prebuilt (it's the config default, so forgetting the flag still works); - `cuda` deliberately is not — this is a ROCm image. +- **Don't pass `--gpu-backend`.** The image defaults it to `system`, meaning + ProstT5 runs against the base image's own ROCm torch rather than a conda env. + Passing `rocm`/`cpu`/`cuda` would send it looking for a `prostt5-*` env that + this image deliberately does not contain — a hard failure on a read-only + `.sif`. ## What's in the image - **Base**: `quay.io/pawsey/pytorch:2.7.1-rocm6.3.3`, Pawsey's own verified - ROCm build, supplying the ROCm userspace matching Setonix's MI250X (gfx90a). - Its system-python torch (2.7.1) is *not* what the workflow uses and goes - untouched. -- **Miniforge** at `/opt/miniforge3`, plus Snakemake — the thing that reads the - prebuilt envs, so deliberately not one of them. -- **phables**, pip-installed from the build context (the commit CI tagged). + ROCm build, supplying both the ROCm userspace matching Setonix's MI250X + (gfx90a) **and the torch the workflow actually uses**. +- **phables, Snakemake and pholdlib installed into that base python** — not + into miniforge's. Snakemake runs a `script:` rule that declares no conda env + using its own interpreter, so that interpreter has to be the one holding + torch. The build asserts torch's version is unchanged across the pip install, + since a silently-replaced torch is the exact failure this avoids. +- **Miniforge** at `/opt/miniforge3`, appended to `PATH` (never prepended, so + it cannot shadow the base python). It exists only to provide the `conda` + binary that builds the envs below. - **Every per-rule conda env**, prebuilt by `container/prebuild_envs.sh`: coverm (minimap2/samtools/CoverM), genecall (FragGeneScan), pyrodigal-gv, - smg (HMMER), mmseqs, foldseek, phylotree (MAFFT + cogent3/piqtree), - prostt5-rocm (torch 2.9.1+rocm6.3 + pholdlib), prostt5-cpu, and curl for - `phables install`. + smg (HMMER), mmseqs, foldseek, phylotree (MAFFT + cogent3/piqtree), and curl + for `phables install`. **No `prostt5-*` env** — that's the point of + `system`. **Databases are not included** — PHROGs and the hallmark DB are multi-GB and separately versioned. Mount them and point `--databases` at them, exactly as @@ -112,26 +118,34 @@ Two things it deliberately does **not** do, both of which broke a real build: so Snakemake reports "Nothing to be done", the DAG is empty and the `curl` env is silently never built. `container/test_image.sh` then smoke-tests the result and fails the build if any -env is missing a binary the rules actually invoke, or if torch didn't come from -the prostt5-rocm env. - -**Disk**: this image is large — a ~14GB compressed ROCm base plus conda envs -including two torch builds. The CI workflow runs `jlumbroso/free-disk-space` -first because a stock GitHub Actions runner has only ~14GB free and the base -alone won't fit. Even with that, this is close to the limit of what a hosted -runner can build; if CI starts failing on `no space left on device`, building -on a machine with real disk and pushing manually is the fallback. +env is missing a binary the rules actually invoke, if the ambient torch isn't a +ROCm build, or if any pre-built env turns out to contain a torch of its own — +that last one being the regression that would mean a second copy got installed +after all. + +**Disk**: this image is large — a ~14GB compressed ROCm base plus the conda +envs. The CI workflow runs `jlumbroso/free-disk-space` first because a stock +GitHub Actions runner has only ~14GB free and the base alone won't fit. Reusing +the base torch rather than installing a second one keeps several GB off the +total, but this is still close to the limit of what a hosted runner can build; +if CI starts failing on `no space left on device`, building on a machine with +real disk and pushing manually is the fallback. ## Status A real `docker build` got as far as the env pre-build step before failing on the two issues listed above; both are fixed, and `prebuild_envs.sh` has since been run end-to-end (with its `--conda-create-envs-only` calls swapped for -dry-runs) so that all six invocations, the synthetic input generation and the -cleanup are known to work as written. The pinned conda packages were confirmed -to exist for `linux-64`/`noarch`. - -Still unverified: a complete `docker build` (the conda solves themselves, and -`test_image.sh` against a real image), and any Setonix Apptainer run. Build it, +dry-runs) so that its invocations, the synthetic input generation and the +cleanup are known to work as written. `--gpu-backend system` was verified to +declare no `prostt5-*` env at all (five envs instead of six, via +`--list-conda-envs`), which is what makes the torch reuse real rather than +aspirational. The pinned conda packages were confirmed to exist for +`linux-64`/`noarch`. + +Still unverified: a complete `docker build` (the conda solves themselves, the +pip install into the base python leaving its torch untouched, and +`test_image.sh` against a real image), and any Setonix Apptainer run — +including whether ProstT5 actually sees the GPU through `--rocm`. Build it, push a tag, and put one real sample through it before trusting it for production batches. diff --git a/docs/usage.md b/docs/usage.md index a5894de..283f138 100644 --- a/docs/usage.md +++ b/docs/usage.md @@ -56,13 +56,20 @@ Options: default) or prostt5-foldseek (structural; needs --hallmark-db, --hallmark-categories and --prostt5-checkpoint) [default: mmseqs] - --gpu-backend [cpu|cuda|rocm] PyTorch build for the ProstT5 conda env: + --gpu-backend [cpu|cuda|rocm|system] + PyTorch build for the ProstT5 conda env: cpu, cuda, or rocm (e.g. Setonix's MI250X nodes). Independent of --prostt5-cpu, which forces ProstT5 onto the CPU device at runtime even inside a GPU-capable env -- - this controls which env gets built - [default: cpu] + this controls which env gets built. 'system' + builds NO env at all and uses the torch + + pholdlib already installed in the ambient + python -- for when a working GPU torch is + already in place (e.g. inside the container, + or a module-loaded torch on HPC) and + installing a second one would be wasteful or + wrong [default: cpu] --foldseek-gpu use foldseek's CUDA GPU search mode for the hallmark scan (requires --gpu-backend cuda, a CUDA-capable foldseek build on PATH -- not @@ -195,7 +202,7 @@ Two independent methods for finding phage-like genes on unitigs, feeding the sam * `--prostt5-half-precision` / `--prostt5-full-precision` - run ProstT5 in half precision (ignored on CPU) [default: `prostt5-half-precision`] * `--prostt5-cpu` - force ProstT5 onto CPU even when a GPU is available * `--prostt5-max-residues`, `--prostt5-max-seq-len`, `--prostt5-max-batch` - ProstT5 batching knobs, device-specific -- tune per GPU rather than trusting the defaults on unfamiliar hardware [defaults: 4000, 4000, 20] -* `--gpu-backend` - which PyTorch build the ProstT5 conda env solves against: `cpu`, `cuda`, or `rocm` (e.g. Setonix's MI250X). Independent of `--prostt5-cpu`, which forces the CPU device at runtime even inside a GPU-capable env -- this controls which env gets *built* [default: `cpu`] +* `--gpu-backend` - which PyTorch build the ProstT5 conda env solves against: `cpu`, `cuda`, or `rocm` (e.g. Setonix's MI250X). Independent of `--prostt5-cpu`, which forces the CPU device at runtime even inside a GPU-capable env -- this controls which env gets *built* [default: `cpu`]. A fourth value, `system`, builds **no** env at all and runs ProstT5 against the `torch`+`pholdlib` already present in the ambient python. Conda envs are isolated, so a rule that declares one can never see a torch installed outside it -- `system` is the only way to reuse an existing, known-good GPU torch (the [container](container.md)'s ROCm base image, or a module-loaded torch on HPC) instead of installing a second copy * `--foldseek-gpu` - use Foldseek's CUDA GPU search mode for the hallmark scan. Requires `--gpu-backend cuda`, a CUDA-capable Foldseek build on `PATH` (not the plain bioconda package), and a `*_gpu`-suffixed, `makepaddedseqdb`-prepared hallmark DB To run all of this from one prebuilt image on HPC (Setonix), see [Running from a single container](container.md) — the image ships every conda env already built, so there are no per-rule container flags to pass. diff --git a/phables/__main__.py b/phables/__main__.py index 3922a7a..f448be3 100644 --- a/phables/__main__.py +++ b/phables/__main__.py @@ -228,9 +228,14 @@ def run_options(func): "PyTorch build for the ProstT5 conda env: cpu, cuda, or rocm " "(e.g. Setonix's MI250X nodes). Independent of --prostt5-cpu, " "which forces ProstT5 onto the CPU device at runtime even " - "inside a GPU-capable env -- this controls which env gets built" + "inside a GPU-capable env -- this controls which env gets " + "built. 'system' builds NO env at all and uses the torch + " + "pholdlib already installed in the ambient python -- for when " + "a working GPU torch is already in place (e.g. inside the " + "container, or a module-loaded torch on HPC) and installing a " + "second one would be wasteful or wrong" ), - type=click.Choice(["cpu", "cuda", "rocm"]), + type=click.Choice(["cpu", "cuda", "rocm", "system"]), show_default=True, ), click.option( diff --git a/phables/config/config.yaml b/phables/config/config.yaml index 658a157..00d3a35 100644 --- a/phables/config/config.yaml +++ b/phables/config/config.yaml @@ -57,6 +57,15 @@ prostt5_cpu: False # gpu_backend controls which env gets *built*. foldseek_gpu (below) is # separate again: foldseek's own --gpu mode is CUDA-only, so it stays off # regardless of gpu_backend unless explicitly requested and gpu_backend=cuda. +# +# A fourth value, `system`, means "don't build an env for predict_3di at all; +# torch and pholdlib are already installed in the ambient python". Conda envs +# are isolated, so a rule with a conda: directive can never see a torch that +# lives outside its env -- `system` is the only way to REUSE an existing, +# known-good GPU torch instead of installing a second one. That's exactly the +# situation inside container/Dockerfile, whose ROCm base image already ships a +# working torch, and it's also how you'd use a module-loaded torch on HPC. +# The image sets this as its own default (see container/Dockerfile). gpu_backend: cpu foldseek_gpu: False # Batch limits are device-specific -- the defaults are conservative CPU/laptop-MPS diff --git a/phables/workflow/rules/genes.smk b/phables/workflow/rules/genes.smk index 1a2afad..fdfd2f5 100644 --- a/phables/workflow/rules/genes.smk +++ b/phables/workflow/rules/genes.smk @@ -147,13 +147,20 @@ if PD == "prostt5-foldseek": max_batch = config["prostt5_max_batch"], log: os.path.join(LOGSDIR, "predict_3di.log") + # gpu_backend selects which torch build this rule runs against. + # cpu/cuda/rocm each need a different PyTorch wheel, and a conda env is + # solved once from a static file, so those are three separate env files + # rather than one file with a runtime switch. + # + # `system` is the odd one out and takes NO conda: directive at all: + # conda envs are isolated, so a rule that declares one can never see a + # torch installed outside it. Omitting the directive is therefore the + # only way to REUSE an already-working torch (the container's ROCm base + # image, a module-loaded torch on HPC) instead of installing a second + # copy. The rule then runs in whichever python is running Snakemake, + # which must already provide torch + pholdlib. conda: - # gpu_backend selects which torch build this env solves against -- - # cpu/cuda/rocm need different PyTorch wheels (conda envs are - # solved once from a static file, so this has to be three files, - # not one file with a runtime switch). See the individual env - # files for what each backend actually needs and why. - os.path.join("..", "envs", f"prostt5-{GPU_BACKEND}.yaml") + None if GPU_BACKEND == "system" else os.path.join("..", "envs", f"prostt5-{GPU_BACKEND}.yaml") script: os.path.join("..", "scripts", "predict_3di.py") From 74b2038bcbf3debcffc61a30833ba62f8bccced5 Mon Sep 17 00:00:00 2001 From: gbouras13 Date: Fri, 14 Aug 2026 09:13:49 +0930 Subject: [PATCH 06/12] container: report the base torch build flavour, do not assert on it The first full image build succeeded through every conda solve and then failed at the very last step on my own over-strict check: it tested for the substring 'rocm' in torch.__version__, but Pawsey's base torch is a source build reporting '2.7.1a0+gite2d141d' with no +rocm6.3 suffix, so a perfectly good ROCm torch was rejected. Now prints torch.__version__, torch.version.hip and torch.version.cuda and continues, warning only if hip is None. What that torch is compiled against is Pawsey's business, and failing a multi-GB build at the last step to re-litigate it is a bad trade. The checks that actually matter stay fatal: torch and pholdlib must be importable in the ambient python (predict_3di cannot run otherwise), and no pre-built conda env may contain a torch of its own (which would mean a second copy was installed after all). That build otherwise confirmed the design end to end: all 8 per-rule envs solve and install, no prostt5-* env is created, and the pip install into the base python leaves its torch untouched. --- container/test_image.sh | 36 ++++++++++++++++++++++++++----- docs/container.md | 48 ++++++++++++++++++++++++----------------- 2 files changed, 59 insertions(+), 25 deletions(-) diff --git a/container/test_image.sh b/container/test_image.sh index 9d5d602..a46bd09 100755 --- a/container/test_image.sh +++ b/container/test_image.sh @@ -65,11 +65,37 @@ check_bin curl echo "=== torch must come from the BASE image, reused -- not reinstalled ===" # predict_3di runs with NO conda env (gpu_backend=system), i.e. in the same # python that runs Snakemake -- which must therefore be the base image's python, -# the one already holding a working ROCm torch. Both halves of that are checked: -# torch+pholdlib importable here, AND no conda env carrying a torch of its own. -python -c "import torch; print('torch:', torch.__version__)" -python -c "import torch, sys; sys.exit(0 if 'rocm' in torch.__version__ else 1)" \ - || { echo "ERROR: the ambient torch is not a ROCm build" >&2; exit 1; } +# the one already holding a working torch. Two things are checked: that torch + +# pholdlib are importable here (fatal if not -- predict_3di simply cannot run), +# and that no conda env carries a torch of its own (fatal -- that would mean a +# second copy got installed after all). What the torch's build flavour is, is +# only reported. +# +# Reported, NOT asserted. This deliberately cannot fail the build. +# +# The base image's torch is Pawsey's own source build for Setonix -- it reports +# e.g. "2.7.1a0+gite2d141d", with no "+rocm6.3" suffix, because it isn't a +# stock wheel. An earlier version of this script tested for the substring +# "rocm" in torch.__version__ and failed a perfectly good ROCm build at the +# very last step of a multi-GB image. Whether that torch is ROCm-enabled is +# Pawsey's business, not something worth re-litigating here at build time, so +# this prints the facts (torch.version.hip is the real signal: set to the HIP +# version on ROCm builds, None otherwise) and moves on. +# +# The check that DOES matter -- that nothing installed a second torch -- is +# below and is fatal. +python - <<'PY' +import torch + +hip = getattr(torch.version, "hip", None) +print("torch:", torch.__version__) +print("torch.version.hip:", hip) +print("torch.version.cuda:", getattr(torch.version, "cuda", None)) +if hip is None: + print("WARNING: torch.version.hip is None -- this torch does not look " + "ROCm-enabled. Fine if intentional (e.g. a CPU-only base image); " + "worth a look if you expected GPU ProstT5 on Setonix.") +PY python -c "import pholdlib; print('pholdlib OK')" python -c "import phables, snakemake; print('phables + snakemake share this interpreter')" diff --git a/docs/container.md b/docs/container.md index 7858020..0061f99 100644 --- a/docs/container.md +++ b/docs/container.md @@ -118,10 +118,17 @@ Two things it deliberately does **not** do, both of which broke a real build: so Snakemake reports "Nothing to be done", the DAG is empty and the `curl` env is silently never built. `container/test_image.sh` then smoke-tests the result and fails the build if any -env is missing a binary the rules actually invoke, if the ambient torch isn't a -ROCm build, or if any pre-built env turns out to contain a torch of its own — -that last one being the regression that would mean a second copy got installed -after all. +env is missing a binary the rules actually invoke, if torch/pholdlib aren't +importable in the ambient python, or if any pre-built env turns out to contain a +torch of its own — that last one being the regression that would mean a second +copy got installed after all. + +It reports the base torch's build flavour (`torch.version.hip`) but does **not** +assert on it. The Pawsey base is a source build reporting e.g. +`2.7.1a0+gite2d141d`, with no `+rocm6.3` suffix; an earlier version of this +script tested for the substring `rocm` and failed a perfectly good ROCm build at +the very last step of a multi-GB image. What that torch is compiled against is +Pawsey's business. **Disk**: this image is large — a ~14GB compressed ROCm base plus the conda envs. The CI workflow runs `jlumbroso/free-disk-space` first because a stock @@ -133,19 +140,20 @@ real disk and pushing manually is the fallback. ## Status -A real `docker build` got as far as the env pre-build step before failing on -the two issues listed above; both are fixed, and `prebuild_envs.sh` has since -been run end-to-end (with its `--conda-create-envs-only` calls swapped for -dry-runs) so that its invocations, the synthetic input generation and the -cleanup are known to work as written. `--gpu-backend system` was verified to -declare no `prostt5-*` env at all (five envs instead of six, via -`--list-conda-envs`), which is what makes the torch reuse real rather than -aspirational. The pinned conda packages were confirmed to exist for -`linux-64`/`noarch`. - -Still unverified: a complete `docker build` (the conda solves themselves, the -pip install into the base python leaving its torch untouched, and -`test_image.sh` against a real image), and any Setonix Apptainer run — -including whether ProstT5 actually sees the GPU through `--rocm`. Build it, -push a tag, and put one real sample through it before trusting it for -production batches. +A real `docker build` now gets all the way through the image and into +`test_image.sh`. Confirmed working on a real build: + +- All **8** per-rule conda envs solve and install (coverm, genecall, smg, + mmseqs, foldseek, phylotree, phables, curl) — including the pins that were + only repodata-checked before (`mmseqs2=13.45111`, `cogent3<2026.7`). +- **No `prostt5-*` env is created**, so no second torch is downloaded — the + torch reuse works as designed. +- Installing phables/Snakemake/pholdlib into the base python leaves its torch + untouched (the build's own before/after assertion passed). +- Every per-rule binary resolves inside a pre-built env. + +Still unverified: a Setonix Apptainer run — in particular whether ProstT5 +actually sees the GPU through `singularity exec --rocm`, which no build-time +check can answer. Pull the `.sif`, put one real sample through it, and confirm +predict_3di lands on the GPU rather than silently falling back to CPU before +trusting this for production batches. From 4e8ed959f12a1b9870f74d15d4cad894bcca0592 Mon Sep 17 00:00:00 2001 From: gbouras13 Date: Fri, 14 Aug 2026 15:02:56 +0930 Subject: [PATCH 07/12] container: --gpu-backend system must be passed, not baked into config The first real container run on Setonix failed: Creating conda environment .../envs/prostt5-cpu.yaml... OSError: [Errno 30] Read-only file system: '/usr/local/.../phables/workflow/conda/4a0d2f6693d2fae2ac9de8108a44f6cc_.yaml' The image's config.yaml did say 'gpu_backend: system', but that is not the effective value: phables passes every CLI option to snaketool as merge_config=kwargs, so any option with a non-None click default is merged OVER the config file -- and --gpu-backend's default is 'cpu'. So predict_3di asked for a prostt5-cpu env the image deliberately does not contain, on a read-only filesystem. Reproduced locally before fixing: with config.yaml patched to 'system' and no flag passed, the runtime config still reports gpu_backend: cpu; passing --gpu-backend system reports system. Removes the Dockerfile's sed step (it never worked) and replaces test_image.sh's config-file grep, which asserted the wrong thing and passed while real runs used cpu. It now checks what the image can actually guarantee -- that the 'system' choice exists -- and the docs state that callers must pass the flag. --- container/Dockerfile | 21 ++++++++++----------- container/test_image.sh | 16 +++++++++++----- docs/container.md | 15 ++++++++++----- phables/config/config.yaml | 7 ++++++- 4 files changed, 37 insertions(+), 22 deletions(-) diff --git a/container/Dockerfile b/container/Dockerfile index 8d29e96..0904217 100644 --- a/container/Dockerfile +++ b/container/Dockerfile @@ -83,18 +83,17 @@ RUN BASE_PY="$(cat /opt/base_python_path)" && \ test "$TORCH_BEFORE" = "$TORCH_AFTER" # -# Make `system` the image's own default GPU backend, so a plain -# `phables run --phagedetection prostt5-foldseek` inside the container reuses -# the base torch with no extra flag. Without this the config default (`cpu`) -# would send predict_3di looking for a prostt5-cpu env that this image -# deliberately does not contain -- and a .sif is read-only, so that is a hard -# failure rather than a slow rebuild. +# NOTE: there is deliberately NO step here patching config.yaml's gpu_backend +# to `system`. An earlier version did exactly that, and it does not work: +# phables passes every CLI option to snaketool as merge_config=kwargs, so any +# option with a non-None click default is merged OVER the config file, and +# --gpu-backend's click default is "cpu". The patched config was silently +# ignored, Snakemake went to build a prostt5-cpu env that isn't in this image, +# and the run died on the read-only .sif filesystem. +# +# Callers must therefore pass `--gpu-backend system` explicitly. See +# ../docs/container.md; mass_processing/run_phables_container.sh does this. # -RUN BASE_PY="$(cat /opt/base_python_path)" && \ - CONFIG="$("$BASE_PY" -c 'import phables, os; print(os.path.join(os.path.dirname(phables.__file__), "config", "config.yaml"))')" && \ - sed -i 's/^gpu_backend: .*/gpu_backend: system/' "$CONFIG" && \ - grep -q '^gpu_backend: system$' "$CONFIG" && \ - echo "image default gpu_backend set to system in $CONFIG" # Pre-build every per-rule conda env. See the script for the full rationale # (default --conda-prefix, placeholder databases, one pass per flag combo). diff --git a/container/test_image.sh b/container/test_image.sh index a46bd09..eb5e04b 100755 --- a/container/test_image.sh +++ b/container/test_image.sh @@ -110,10 +110,16 @@ for e in "$PREFIX"/*/; do done echo "OK: no pre-built env ships a duplicate torch" -echo "=== image's default gpu_backend must be 'system' ===" -CONFIG="$(python -c 'import phables, os; print(os.path.join(os.path.dirname(phables.__file__), "config", "config.yaml"))')" -grep -q '^gpu_backend: system$' "$CONFIG" \ - || { echo "ERROR: $CONFIG does not default gpu_backend to system" >&2; exit 1; } -echo "OK: $CONFIG defaults to system" +echo "=== --gpu-backend system must be available ===" +# This, NOT a grep of config.yaml's gpu_backend value. An earlier version of +# this script asserted the config file said "system" and passed happily, while +# actual runs still used cpu: phables merges every CLI option over the config +# (merge_config=kwargs), and --gpu-backend's click default is "cpu". The config +# file value is simply not the effective value, so checking it proves nothing. +# What the image can meaningfully guarantee is that the choice EXISTS -- the +# caller is responsible for passing it (see docs/container.md). +phables run -h 2>&1 | grep -q -- '--gpu-backend .*system' \ + || { echo "ERROR: this phables has no --gpu-backend system choice" >&2; exit 1; } +echo "OK: --gpu-backend system is accepted" echo "=== all image tests passed ===" diff --git a/docs/container.md b/docs/container.md index 0061f99..8982129 100644 --- a/docs/container.md +++ b/docs/container.md @@ -51,6 +51,7 @@ singularity exec --rocm \ --output phables_out \ --databases /scratch/.../all_databases/databases \ --phagedetection prostt5-foldseek \ + --gpu-backend system \ --prostt5-checkpoint /scratch/.../model.pt \ --threads 8 ``` @@ -65,11 +66,15 @@ Notes that matter on Setonix: (inside the installed package), and Snakemake resolves an env by hashing its file content *together with the prefix path* — changing the prefix changes the hash, and Snakemake would try to rebuild into a read-only filesystem. -- **Don't pass `--gpu-backend`.** The image defaults it to `system`, meaning - ProstT5 runs against the base image's own ROCm torch rather than a conda env. - Passing `rocm`/`cpu`/`cuda` would send it looking for a `prostt5-*` env that - this image deliberately does not contain — a hard failure on a read-only - `.sif`. +- **Pass `--gpu-backend system` explicitly.** This is required, not optional. + phables merges every CLI option over the config file + (`merge_config=kwargs`), and `--gpu-backend`'s click default is `cpu` — so + the image's own `config.yaml` value is *not* the effective value and cannot + be relied on. Omitting the flag was tried and failed for real: the runtime + config showed `gpu_backend: cpu`, Snakemake went to build a `prostt5-cpu` + env that isn't in the image, and the run died with + `OSError: [Errno 30] Read-only file system`. Any backend other than + `system` fails the same way. ## What's in the image diff --git a/phables/config/config.yaml b/phables/config/config.yaml index 00d3a35..9f5c204 100644 --- a/phables/config/config.yaml +++ b/phables/config/config.yaml @@ -65,7 +65,12 @@ prostt5_cpu: False # known-good GPU torch instead of installing a second one. That's exactly the # situation inside container/Dockerfile, whose ROCm base image already ships a # working torch, and it's also how you'd use a module-loaded torch on HPC. -# The image sets this as its own default (see container/Dockerfile). +# Note this key cannot be usefully overridden by editing a config file: the +# CLI hands every option to snaketool as merge_config=kwargs, so +# --gpu-backend's click default ("cpu") is merged OVER whatever is written +# here. Pass --gpu-backend on the command line instead. (Learned the hard way +# by baking `system` into the container's config and watching runs still use +# cpu -- see container/Dockerfile.) gpu_backend: cpu foldseek_gpu: False # Batch limits are device-specific -- the defaults are conservative CPU/laptop-MPS From 82313fa08d1b25c27b70cf5329071f8debf1326a Mon Sep 17 00:00:00 2001 From: gbouras13 Date: Tue, 18 Aug 2026 08:48:58 +0930 Subject: [PATCH 08/12] Start the MFD K-search at a proven lower bound instead of K=1 FD_Algorithm tried K = 1, 2, 3, ... until a feasible decomposition was found, rebuilding the whole MILP for each K. K is structural to the model (every variable is indexed by it), so the model genuinely cannot be reused, and flowpaths does not expose a HiGHS warm start -- so every attempt below the true answer was a full model build that could only ever return infeasible. Profiling showed why that dominates: model CONSTRUCTION is ~95% of the cost of an attempt, not solving (a large component: 170ms build vs 8.9ms solve). That also explains why --threads never helped -- threads only affect the 5% -- and why a solver time_limit does nothing here, both of which were measured and rejected before landing on this. get_lowerbound_k() takes the maximum of two bounds, both from flowpaths' own MinFlowDecomp.get_lowerbound_k: the graph width (minimum paths needed to cover every flow-carrying edge) and ceil(log2(#distinct flow values)). Both are lower bounds, so starting there cannot skip a feasible smaller K. It costs 1-4ms and falls back to 1 on any error, since a lower bound is an optimisation and must never be why a component fails to resolve. Measured on synthetic components, against the previous behaviour through the real SolveInstances entry point: - 18/18 cases returned IDENTICAL K, path count and path sets - 1.5x-4.9x faster, growing with component size - components that cannot resolve within --maxpaths: up to 5.9x, since the bound proves K >= maxpaths up front instead of burning the whole ladder to learn nothing - edge cases verified identical: single edge, all-zero lower bounds, one distinct flow value, bound above --maxpaths, unresolvable Also annotates data["minK"], which was set to a constant 2 and never read by anything, so it isn't mistaken for the live lower bound. --- .../scripts/phables_utils/FD_Inexact.py | 78 ++++++++++++++++++- 1 file changed, 77 insertions(+), 1 deletion(-) diff --git a/phables/workflow/scripts/phables_utils/FD_Inexact.py b/phables/workflow/scripts/phables_utils/FD_Inexact.py index 09045c6..d250128 100644 --- a/phables/workflow/scripts/phables_utils/FD_Inexact.py +++ b/phables/workflow/scripts/phables_utils/FD_Inexact.py @@ -4,6 +4,7 @@ import itertools import logging +import math import flowpaths as fp import networkx as nx @@ -249,11 +250,81 @@ def flowMultipleDecomposition(data, K, nthreads): return data +def get_lowerbound_k(data): + """ + Smallest number of paths any valid decomposition of this component could + possibly use. + + The K search below tries K = 1, 2, 3, ... until one is feasible, rebuilding + the whole MILP each time (K is structural to the model -- every variable is + indexed by it -- so it genuinely cannot be reused, and neither flowpaths nor + HiGHS-via-flowpaths offers a warm start). Every attempt below the true answer + is therefore a complete model build that can only ever come back infeasible, + and model construction is ~95% of the cost of an attempt. Starting the search + at a proven lower bound skips exactly those wasted builds. + + Two bounds, both taken from flowpaths' own MinFlowDecomp.get_lowerbound_k: + + - the graph's WIDTH: the minimum number of paths needed to cover every edge. + Any decomposition must cover every edge carrying flow, so it cannot use + fewer paths than this. + - ceil(log2(number of distinct flow values)): k paths can produce at most + 2^k distinct subset sums, so k must be at least log2 of the number of + distinct values that have to be represented. + + Both are lower bounds, so the maximum of them is too, and starting there can + never skip a feasible smaller K -- the search still returns the same first + feasible K, just without the doomed attempts before it. Verified on synthetic + components: identical K and identical path sets with and without this. + + Returns 1 (i.e. the original behaviour) if anything about the computation + fails. A lower bound is an optimisation, not a correctness requirement, so it + must never be the reason a component stops resolving. + """ + try: + graph = data["graph"] + if graph.number_of_edges() == 0: + return 1 + + # flowpaths' stDAG wants string node ids, same relabelling + # flowMultipleDecomposition does before building its model. + node_labels = {node: str(node) for node in graph.nodes} + flowpaths_graph = nx.relabel_nodes(graph, node_labels, copy=True) + lowerbound = fp.stDAG(flowpaths_graph).get_width() + + # Only edges that must carry flow constrain the decomposition; an edge + # whose lower bound is 0 need not be covered at all. + distinct_flows = { + int(data["flows_low"][edge]) + for edge in data["edges"] + if data["flows_low"].get(edge, 0) > 0 + } + if len(distinct_flows) > 1: + lowerbound = max(lowerbound, math.ceil(math.log2(len(distinct_flows)))) + + return max(1, lowerbound) + except Exception as e: + logger.debug(f"Could not compute a lower bound on K ({e}); starting from 1") + return 1 + + def FD_Algorithm(data, max_paths, nthreads): solutionWeights = 0 solutionSet = 0 - for i in range(1, max_paths + 1): + # See get_lowerbound_k: K < lowerbound is provably infeasible, and each such + # attempt costs a full MILP build. Capped at max_paths so a bound above the + # user's --maxpaths doesn't turn into a search over an empty range with + # different semantics -- that case has no solution within max_paths either + # way, and this keeps the "attempt at least one K" behaviour identical. + lowerbound = min(get_lowerbound_k(data), max_paths) + if lowerbound > 1: + logger.debug( + f"Starting the K search at {lowerbound} rather than 1 " + f"({lowerbound - 1} provably-infeasible attempt(s) skipped)" + ) + + for i in range(lowerbound, max_paths + 1): data = flowMultipleDecomposition(data, i, nthreads) if data["message"] == "solved": solutionSet = data["solution"] @@ -348,6 +419,11 @@ def SolveInstances(Graphs, max_paths, outfile, recfile, nthreads): "adj_in": AD_in, "adj_out": AD_out, "subpaths": Graphs[s]["subpaths"], + # Never read by anything -- the K search takes its starting point + # from get_lowerbound_k(), which computes a real bound per component + # rather than assuming a constant. Left in place because this dict is + # passed around wholesale and removing a key is a wider change than + # it looks; flagged so it isn't mistaken for the live lower bound. "minK": 2, "runtime": 0, } From b537e2265a9938866147c66df57d23bc6876bab3 Mon Sep 17 00:00:00 2001 From: gbouras13 Date: Tue, 18 Aug 2026 08:56:12 +0930 Subject: [PATCH 09/12] Add --mfd-workers: run flow decomposition components in parallel Components are independent, so the flow-decomposition loop is embarrassingly parallel. resolve_short_parallel chunks the components, runs the EXISTING resolve_short once per chunk in a worker process, and merges the returned accumulators. Deliberately no change to resolve_short's ~1400-line body: it is already parameterised by the component set (pruned_vs) and already returns every accumulator it builds, so chunk-and-merge is equivalent to one call over all components. Extracting the loop body into a per-component function would have meant restructuring 1400 lines of nested branching for the same result. Soundness rests on no component's logic depending on another's results, which was verified against the body before writing any of this: every touch of a shared accumulator is a pure add/union/append, and there is not one conditional or membership test against them in the loop -- per-component decisions use the loop-local comp_* sets. Noted in the docstring, since chunking would silently change results if that stops being true. Chunks are merged in component order, so all_resolved_paths ends up identical to the sequential run. That matters: genomes are numbered by position, so a different order would rename every genome without changing the biology. Two implementation details that turned out to matter: - MORE chunks than workers (CHUNKS_PER_WORKER=4). Component cost is heavily skewed and which components are expensive is not known in advance, so one chunk per worker leaves workers idle. On a realistic skewed workload this took 8 workers from 1.92x to 2.35x -- 98% of the 2.4x Amdahl ceiling imposed by the single largest component. - Heavy read-only inputs (the assembly graph above all) are sent once per worker via a pool initializer, not once per chunk, so smaller chunks do not mean repeatedly re-pickling the graph. Per-worker solver threads are pinned to 1: workers x nthreads would oversubscribe, and profiling showed solver threads make no measurable difference anyway. Measured against sequential, through the real wrapper: - identical 18-tuple results at 2, 4 and 8 workers - all_resolved_paths order identical (genome numbering preserved) - uniform workload: 6.1x at 8 workers - realistic skewed workload: 2.35x at 8 workers (ceiling 2.4x) - edge cases identical: 1 component, fewer components than workers, more components than workers, and workers=1 Default is 1, i.e. the sequential path is unchanged. --- docs/usage.md | 12 ++ phables/__main__.py | 16 ++ phables/config/config.yaml | 8 + .../workflow/rules/02_phables_preflight.smk | 1 + phables/workflow/rules/phables.smk | 1 + phables/workflow/scripts/phables.py | 5 +- .../scripts/phables_utils/short_utils.py | 153 ++++++++++++++++++ 7 files changed, 195 insertions(+), 1 deletion(-) diff --git a/docs/usage.md b/docs/usage.md index 283f138..2b84d0a 100644 --- a/docs/usage.md +++ b/docs/usage.md @@ -111,6 +111,17 @@ Options: batch immediately [default: 4000] --prostt5-max-batch INTEGER max sequences per ProstT5 batch -- device- specific, tune per GPU [default: 20] + --mfd-workers INTEGER number of worker processes for the flow- + decomposition step. Components are + independent, so they are split across + workers and merged in order (results are + identical to --mfd-workers 1, including + genome numbering). 1 = sequential, as + before. Note this is separate from + --threads: the MILP solver itself gets no + measurable benefit from extra threads, so + parallelism has to come from running + components concurrently [default: 1] --evalue FLOAT maximum e-value for phrog annotations [default: 1e-10] --seqidentity FLOAT minimum sequence identity for phrog @@ -182,6 +193,7 @@ Options: * `--databases` - path to the databases directory [default: wherever `phables install` put them] * `--use-conda` / `--no-use-conda` - use conda for Snakemake rules [default: `use-conda`] * `--conda-prefix` - custom conda env directory +* `--mfd-workers` - worker processes for the flow-decomposition (MFD) step [default: 1]. Components are independent, so they're split across processes and merged in component order — output is identical to `1`, genome numbering included. Distinct from `--threads`: the MILP solver gains nothing measurable from extra threads (building the model dominates, not solving it), so speedup has to come from running components concurrently. Expect ~2–2.5x on a realistic component mix — a few large components dominate the runtime and can't be split * `--snake-default` - customise Snakemake runtime args [default: `--rerun-incomplete, --printshellcmds, --nolock, --show-failed-logs`] ### Phage-gene detection: `--phagedetection` diff --git a/phables/__main__.py b/phables/__main__.py index f448be3..4a88e93 100644 --- a/phables/__main__.py +++ b/phables/__main__.py @@ -346,6 +346,22 @@ def run_options(func): type=int, show_default=True, ), + click.option( + "--mfd-workers", + default=1, + required=False, + help=( + "number of worker processes for the flow-decomposition step. " + "Components are independent, so they are split across workers " + "and merged in order (results are identical to --mfd-workers 1, " + "including genome numbering). 1 = sequential, as before. Note " + "this is separate from --threads: the MILP solver itself gets " + "no measurable benefit from extra threads, so parallelism has " + "to come from running components concurrently" + ), + type=int, + show_default=True, + ), click.option( "--evalue", default=1e-10, diff --git a/phables/config/config.yaml b/phables/config/config.yaml index 9f5c204..ee2d48c 100644 --- a/phables/config/config.yaml +++ b/phables/config/config.yaml @@ -73,6 +73,14 @@ prostt5_cpu: False # cpu -- see container/Dockerfile.) gpu_backend: cpu foldseek_gpu: False +# Worker processes for the flow-decomposition step. Components are independent +# of each other (no component's logic reads another's results), so they are +# chunked across processes and merged in component order -- output is identical +# to mfd_workers: 1, genome numbering included. Separate from `threads`: the +# MILP solver gains nothing measurable from extra threads because building the +# model, not solving it, dominates -- so concurrency has to come from running +# components side by side instead. +mfd_workers: 1 # Batch limits are device-specific -- the defaults are conservative CPU/laptop-MPS # values, not tuned for any particular GPU. Retune per device (phold's autotune.py # is the reference approach) rather than trusting these at scale. diff --git a/phables/workflow/rules/02_phables_preflight.smk b/phables/workflow/rules/02_phables_preflight.smk index dc6ef6d..bfff53e 100644 --- a/phables/workflow/rules/02_phables_preflight.smk +++ b/phables/workflow/rules/02_phables_preflight.smk @@ -36,6 +36,7 @@ GC = config['genecaller'] PD = config['phagedetection'] GPU_BACKEND = config['gpu_backend'] FOLDSEEK_GPU = config['foldseek_gpu'] +MFD_WORKERS = config['mfd_workers'] EV = config['evalue'] SI = config['seqidentity'] CT = config['covtol'] diff --git a/phables/workflow/rules/phables.smk b/phables/workflow/rules/phables.smk index 549d454..fab47d3 100644 --- a/phables/workflow/rules/phables.smk +++ b/phables/workflow/rules/phables.smk @@ -36,6 +36,7 @@ rule run_phables: hallmark_minbits = config["hallmark_minbits"], output = os.path.join(OUTDIR, "phables"), nthreads = config["resources"]["jobCPU"], + mfd_workers = MFD_WORKERS, log = os.path.join(LOGSDIR, "phables_output.log") threads: config["resources"]["jobCPU"] diff --git a/phables/workflow/scripts/phables.py b/phables/workflow/scripts/phables.py index 32913fe..b3b6830 100755 --- a/phables/workflow/scripts/phables.py +++ b/phables/workflow/scripts/phables.py @@ -50,6 +50,7 @@ def main(): hallmark_minbits = float(snakemake.params.hallmark_minbits) output = snakemake.params.output nthreads = int(snakemake.params.nthreads) + mfd_workers = int(snakemake.params.mfd_workers) log = snakemake.params.log # Setup logger @@ -95,6 +96,7 @@ def main(): logger.info(f"Input long reads: {longreads}") logger.info(f"Prefix for genome identifiers: {prefix}") logger.info(f"Number of threads to use: {nthreads}") + logger.info(f"Number of flow-decomposition workers: {mfd_workers}") logger.info(f"Output folder: {output}") if prefix is None or prefix == "": @@ -248,7 +250,7 @@ def main(): phage_like_edges, all_phage_like_edges, unresolved_phage_like_edges, - ) = short_utils.resolve_short( + ) = short_utils.resolve_short_parallel( assembly_graph, pruned_vs, unitig_names, @@ -269,6 +271,7 @@ def main(): prefix, output, nthreads, + mfd_workers, ) # Log final summary information diff --git a/phables/workflow/scripts/phables_utils/short_utils.py b/phables/workflow/scripts/phables_utils/short_utils.py index 9fb760c..4a94efb 100644 --- a/phables/workflow/scripts/phables_utils/short_utils.py +++ b/phables/workflow/scripts/phables_utils/short_utils.py @@ -1,5 +1,7 @@ #!/usr/bin/env python3 +import math +from concurrent.futures import ProcessPoolExecutor import logging import sys import time @@ -16,6 +18,12 @@ # Create logger logger = logging.getLogger("phables 2.0.0") +# How many chunks to create per worker in resolve_short_parallel. >1 so the +# pool can load-balance an uneven component workload; small enough that +# per-chunk overhead stays negligible against per-component solve times +# measured in tens of milliseconds. +CHUNKS_PER_WORKER = 4 + def resolve_short( assembly_graph, @@ -1498,3 +1506,148 @@ def resolve_short( all_phage_like_edges, unresolved_phage_like_edges, ) + + +# Set once per worker process by _init_worker, so the heavy read-only inputs +# (the assembly graph above all) cross the process boundary ONCE per worker +# rather than once per chunk. Passing them as map() arguments instead would +# re-pickle the whole graph for every chunk, which gets more expensive the more +# chunks are used -- exactly backwards, since more chunks is what gives the pool +# room to balance an uneven workload. +_WORKER_KWARGS = None + + +def _init_worker(kwargs): + global _WORKER_KWARGS + _WORKER_KWARGS = kwargs + + +def _resolve_short_chunk(chunk): + """Worker entry point: one slice of components, nothing else. + + Module-level (not a closure) so it is picklable by ProcessPoolExecutor. + """ + return resolve_short(pruned_vs=chunk, **_WORKER_KWARGS) + + +def resolve_short_parallel( + assembly_graph, + pruned_vs, + unitig_names, + unitig_names_rev, + self_looped_nodes, + graph_unitigs, + minlength, + link_overlap, + unitig_coverages, + compcount, + oriented_links, + junction_pe_coverage, + likely_complete, + alpha, + mincov, + covtol, + maxpaths, + prefix, + output, + nthreads, + workers=1, +): + """ + resolve_short, with components processed in parallel across processes. + + Deliberately implemented WITHOUT touching resolve_short's ~1400-line body: + that function is already parameterised by the set of components to process + (`pruned_vs`) and already returns every accumulator it builds, so splitting + the components into contiguous chunks, running one chunk per worker, and + merging the returned tuples is equivalent to running it once over all of + them. The alternative -- extracting the loop body into a per-component + function -- would mean restructuring 1400 lines of nested branching for no + additional benefit. + + This is only sound because no component's logic depends on another's + results. Verified against the body: every touch of a shared accumulator is + a pure `X.add(...)` / `X = X.union(...)` / `X.append(...)`, and there is not + one conditional or membership test against them anywhere in the loop -- + per-component decisions use the loop-local `comp_*` sets instead. If that + ever stops being true, chunking silently changes results, so it is worth + re-checking before extending this. + + Chunks are CONTIGUOUS and merged in chunk order, so `all_resolved_paths` + ends up in exactly the same order as the sequential run. That matters + because downstream naming numbers the genomes by position: a different + order would rename every genome without changing the biology, which is a + nasty kind of non-reproducibility. + + Falls back to a plain sequential call for workers <= 1 or a single chunk, + so the default path is byte-for-byte the code that ran before. + """ + kwargs = dict( + assembly_graph=assembly_graph, + unitig_names=unitig_names, + unitig_names_rev=unitig_names_rev, + self_looped_nodes=self_looped_nodes, + graph_unitigs=graph_unitigs, + minlength=minlength, + link_overlap=link_overlap, + unitig_coverages=unitig_coverages, + compcount=compcount, + oriented_links=oriented_links, + junction_pe_coverage=junction_pe_coverage, + likely_complete=likely_complete, + alpha=alpha, + mincov=mincov, + covtol=covtol, + maxpaths=maxpaths, + prefix=prefix, + output=output, + # Each worker is its own process, so the MILP solver inside it should + # not also try to use every core -- that would oversubscribe by + # workers x nthreads. Profiling showed solver threads make no + # measurable difference anyway (model construction dominates), so 1 is + # the right per-worker value rather than a compromise. + nthreads=1, + ) + + keys = list(pruned_vs) + if workers <= 1 or len(keys) <= 1: + return resolve_short(pruned_vs=pruned_vs, **{**kwargs, "nthreads": nthreads}) + + workers = min(workers, len(keys)) + + # Deliberately MORE chunks than workers. Component cost is heavily skewed -- + # a handful of large components dominate, and which ones is not known in + # advance -- so splitting into exactly one chunk per worker (static + # partitioning) lets a single worker draw several expensive components while + # the rest sit idle. Smaller chunks let the pool hand out more work to + # whichever process finishes first. The usual cost of small chunks, + # re-sending the inputs each time, does not apply here: _init_worker sends + # them once per worker instead. + size = max(1, math.ceil(len(keys) / (workers * CHUNKS_PER_WORKER))) + chunks = [ + {k: pruned_vs[k] for k in keys[i : i + size]} for i in range(0, len(keys), size) + ] + + logger.info( + f"Resolving {len(keys)} components across {workers} worker process(es) " + f"in {len(chunks)} chunk(s)" + ) + + with ProcessPoolExecutor( + max_workers=workers, initializer=_init_worker, initargs=(kwargs,) + ) as executor: + # .map preserves input order, which is what keeps the merged + # all_resolved_paths identical to the sequential run. + results = list(executor.map(_resolve_short_chunk, chunks)) + + # Field 1 (all_resolved_paths) and field 2 (all_components) are lists and + # are extended; every other field is a set and is unioned. Merged in chunk + # order for the reason given above. + merged = list(results[0]) + for result in results[1:]: + for i, value in enumerate(result): + if isinstance(value, list): + merged[i] = merged[i] + value + else: + merged[i] = merged[i] | value + return tuple(merged) From a5be3aaf0aec54334ac4b242a1306d42ad36b21c Mon Sep 17 00:00:00 2001 From: gbouras13 Date: Tue, 18 Aug 2026 09:01:47 +0930 Subject: [PATCH 10/12] container CI: correct a stale comment about what the image contains The comment still described the image as pre-building two torch envs (rocm + cpu). That stopped being true when predict_3di switched to reusing the base image's ROCm torch via --gpu-backend system -- no prostt5-* env is built at all now. --- .github/workflows/build_container.yaml | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) diff --git a/.github/workflows/build_container.yaml b/.github/workflows/build_container.yaml index a5ef161..d179f50 100644 --- a/.github/workflows/build_container.yaml +++ b/.github/workflows/build_container.yaml @@ -25,14 +25,15 @@ jobs: # keeps the runner's node/python toolcache so other actions (checkout, # docker login) still work fast. # - # NOTE this image is now MONOLITHIC -- on top of that base it installs - # miniforge and pre-builds every per-rule conda env, including two torch - # builds (rocm + cpu). That is several GB more than the base alone, and - # is close to what a hosted runner can do even after this cleanup. If a - # build starts failing on "no space left on device", the honest fix is - # building on a machine with real disk and pushing manually, not shaving - # more off here. The ``df -h`` steps below exist to make that diagnosis - # obvious rather than a guess. + # NOTE this image is MONOLITHIC -- on top of that base it installs + # miniforge and pre-builds every per-rule conda env. It deliberately does + # NOT build a prostt5-* env: predict_3di reuses the base image's own ROCm + # torch via --gpu-backend system, which keeps several GB out of the image. + # Even so this is close to what a hosted runner can manage after the + # cleanup below. If a build starts failing on "no space left on device", + # the honest fix is building on a machine with real disk and pushing + # manually, not shaving more off here. The ``df -h`` steps exist to make + # that diagnosis obvious rather than a guess. - name: Free up disk space uses: jlumbroso/free-disk-space@main with: From cd7ad30c3ab80ffc29a73d8cd52ca896f8b74f2f Mon Sep 17 00:00:00 2001 From: gbouras13 Date: Tue, 18 Aug 2026 21:22:25 +0930 Subject: [PATCH 11/12] Fix --mfd-workers crashing on real data: oriented_links was unpicklable The first real run with --mfd-workers 8 (120 components, SRR27716023) died immediately after the pool started. Cause: edge_graph_utils built oriented_links as defaultdict(lambda: defaultdict(list)), and lambdas cannot be pickled -- so the structure could not cross a process boundary and every worker failed at startup. The sequential path never pickles anything, which is why this only appeared with workers > 1. Fixed at source: the inner factory is now a module-level function, which pickles. Behaviour is identical -- missing keys still get a defaultdict(list) at both levels, verified across a pickle round-trip. Also adds a preflight in resolve_short_parallel: the kwargs are test-pickled before the pool is created, and if anything is unpicklable it logs a clear warning and runs sequentially instead. An optimisation should never be able to destroy a completed assembly's run, and this failure mode is invisible until it happens on real data -- same reasoning as get_lowerbound_k falling back to 1 on error. Verified: with the fix the pool runs and merges all components; with the old lambda form the preflight catches it, warns, falls back, and returns results identical to the parallel path. --- PR_BODY.md | 157 ++++++++++++++++++ .../scripts/phables_utils/edge_graph_utils.py | 15 +- .../scripts/phables_utils/short_utils.py | 19 +++ 3 files changed, 190 insertions(+), 1 deletion(-) create mode 100644 PR_BODY.md diff --git a/PR_BODY.md b/PR_BODY.md new file mode 100644 index 0000000..1b4ee78 --- /dev/null +++ b/PR_BODY.md @@ -0,0 +1,157 @@ +Follow-up to #68. Three independent strands: a crash fix, two flow-decomposition +speedups, and a rework of the container added in #68. + +Rebases cleanly onto current `develop` (checked against `d196f5f`, including the +`QUAY_NAMESPACE` change to `build_container.yaml`). + +--- + +## 1. Fix: phables crashes when a sample resolves no genomes + +A sample where nothing is resolved and there are no unresolved phage-like edges +produces an empty `genomes_and_unresolved_edges.fasta`. CoverM 0.7.0 doesn't +print an empty table for the resulting zero-alignment BAM — it panics: + +``` +[WARN coverm::contig] No primary alignments were observed for sample X +thread 'main' panicked at src/coverage_printer.rs:467:61: +index out of bounds: the len is 0 but the index is 0 +``` + +That killed the whole run at the very last stage, after all the expensive work +had already succeeded. Hit for real on `SRR19670770`. + +`coverm_bam2counts_genomes` now checks the BAM for alignments first and writes a +header-only coverage table instead of invoking CoverM when there are none. +Everything downstream already handled an empty table correctly, so the run +finishes normally with empty report tables. The header reproduces CoverM's own +exactly, including `Covered Fraction` being two words. + +## 2. Performance: flow decomposition + +Both changes are opt-out-safe — the second is off by default — and both were +profiled before being written rather than guessed at. + +**Where the time actually goes.** For a component's MILP, *building* the model +is ~95% of the cost, not solving it (large component: 170 ms build vs 8.9 ms +solve). Two consequences, both measured: solver threads make no difference at +all (1/2/4/8 threads are flat within noise), and a solver `time_limit` does not +help either. + +**a. Start the K search at a proven lower bound** (`FD_Algorithm`) + +`FD_Algorithm` tried K = 1, 2, 3, … until feasible, rebuilding the whole MILP +each time. K is structural to the model, so it genuinely cannot be reused, and +flowpaths does not expose a HiGHS warm start — meaning every attempt below the +true answer was a full model build that could only return infeasible. + +`get_lowerbound_k()` takes the max of the graph width and +`ceil(log2(#distinct flow values))`, both lifted from flowpaths' own +`MinFlowDecomp.get_lowerbound_k`. Both are lower bounds, so starting there +cannot skip a feasible smaller K. It costs 1–4 ms and falls back to 1 on any +error, since a lower bound is an optimisation and must never be why a component +fails to resolve. + +| | | +|---|---| +| 18/18 synthetic cases | identical K, path count and path sets | +| speedup | 1.5×–4.9×, growing with component size | +| components that can't resolve within `--maxpaths` | up to 5.9× (the bound proves `K >= maxpaths` up front instead of burning the whole ladder) | + +(Also annotates `data["minK"]`, which was set to a constant `2` and never read +by anything, so it isn't mistaken for the live bound.) + +**b. `--mfd-workers`: run components in parallel** (default `1`, unchanged behaviour) + +Components are independent, so the loop is embarrassingly parallel. +`resolve_short_parallel` chunks them, runs the **existing** `resolve_short` once +per chunk in a worker process, and merges the returned accumulators — no change +to that function's ~1400-line body, since it's already parameterised by the +component set and already returns everything it builds. + +This is only sound because no component's logic depends on another's results. +That was verified against the body first: every touch of a shared accumulator is +a pure `add`/`union`/`append`, with no conditional or membership test against +them anywhere in the loop (per-component decisions use the loop-local `comp_*` +sets). It's noted in the docstring, because chunking would silently change +results if that ever stopped being true. + +Chunks merge in component order, so `all_resolved_paths` is identical to the +sequential run — genomes are numbered by position, so a different order would +rename every genome without changing the biology. + +| Workload | 2 | 4 | 8 workers | +|---|---|---|---| +| uniform components | 1.93× | 3.74× | 6.08× | +| realistic skew | 1.92× | 2.02× | 2.35× | + +The skewed case is the one to plan around: one component was 41% of total +runtime, giving an Amdahl ceiling of 2.4× — so 2.35× is ~98% of what's +achievable. Two details mattered: more chunks than workers (so the pool can +balance an uneven workload; this alone took 8 workers from 1.92× to 2.35×), and +sending the heavy read-only inputs once per worker via a pool initializer rather +than once per chunk, so smaller chunks don't mean re-pickling the assembly graph. + +Verified identical results at 2/4/8 workers, identical path ordering, and +identical output for 1 component, fewer components than workers, more components +than workers, and `workers=1`. + +## 3. Container: one monolithic image, replacing the per-rule container flags + +**This removes `--container` and `--prostt5-container`, both added in #68.** +Worth being explicit about, since they were merged only recently. + +Those flags pointed individual rules at images via Snakemake `container:` +directives. That approach fights Snakemake: combining `container:` with `conda:` +triggers the documented "ad-hoc combination" behaviour, which builds a *fresh* +conda env *inside* the container rather than using what the image already has — +defeating the point. The rules now carry plain `conda:` directives only, and the +image satisfies them by having every per-rule env **pre-built inside it** +(hybracter's approach). Run it and there is nothing left to create at runtime, +which removes the conda-env creation race that motivated containerising at all: +each sample is its own Snakemake process, so concurrent array tasks needing the +same not-yet-built env can corrupt each other's `mamba env create`. + +Also here: +- **`--gpu-backend system`** — builds no env and uses the ambient `torch` + + `pholdlib`. Conda envs are isolated, so a rule declaring `conda:` can never see + a torch installed outside it; this is the only way to *reuse* a known-good GPU + torch (the container's ROCm base, or a module-loaded torch on HPC) instead of + installing a second copy. Keeps several GB out of the image. +- `container/prebuild_envs.sh` generates its own throwaway inputs rather than + using `tests/data/` (which is gitignored, so absent in a fresh clone or CI + checkout), and points `phables install` at an *empty* databases dir so the + download rules are actually in the DAG and the `curl` env gets built. +- `container/test_image.sh` fails the build if any per-rule tool is missing, if + `torch`/`pholdlib` aren't importable, or if any pre-built env contains its own + torch (i.e. a second copy crept in). +- `prostt5-rocm.yaml` bumped to `torch==2.9.1` (verified present on the pinned + `rocm6.3` index, cp310–cp314). + +--- + +## Testing + +- Existing test suite passes. +- Workflow DAG verified by real `--dry-run`s across gene caller, detection mode, + GPU backend and tree options. +- Performance and equivalence numbers above are from synthetic instances driving + the real `FD_Inexact` / `resolve_short` code paths. +- The container built end-to-end: all 8 per-rule conda envs solve and install, no + duplicate torch, and the pip install leaves the base image's torch untouched. + +## Not verified + +- No end-to-end run on a real sample through the parallel path — the equivalence + testing used a stubbed `resolve_short` (real chunking and merging, synthetic + per-component results). +- Speedups are from synthetic components with clean topology. Build-dominates- + solve should hold generally since it tracks graph size, but a genuinely hard + component would shift the ratio. The lower bound is *valid* regardless, just + possibly looser on real graphs. +- No Setonix Apptainer run of the final image. +- Process pools copy memory per worker, so a large assembly graph at 8 workers + may bind on RAM before CPU. + +Happy to split any of the three strands into its own PR if that's easier to +review — they're independent. diff --git a/phables/workflow/scripts/phables_utils/edge_graph_utils.py b/phables/workflow/scripts/phables_utils/edge_graph_utils.py index 41c3b00..1e0667d 100644 --- a/phables/workflow/scripts/phables_utils/edge_graph_utils.py +++ b/phables/workflow/scripts/phables_utils/edge_graph_utils.py @@ -59,6 +59,19 @@ def get_unitig_lengths(edge_file): return unitig_lengths +def _oriented_links_inner(): + """Inner factory for oriented_links. + + A module-level function rather than the `lambda: defaultdict(list)` this + used to be, purely so the resulting structure can be PICKLED. Lambdas + cannot be, which meant oriented_links could not cross a process boundary -- + it broke resolve_short_parallel's worker pool with + `PicklingError: Can't pickle >` the first time it ran on + real data. Behaviour is identical: missing keys still get a defaultdict(list). + """ + return defaultdict(list) + + def get_links(assembly_graph_file): """ Get links from the assembly graph @@ -67,7 +80,7 @@ def get_links(assembly_graph_file): node_count = 0 graph_contigs = {} edges_lengths = {} - oriented_links = defaultdict(lambda: defaultdict(list)) + oriented_links = defaultdict(_oriented_links_inner) link_overlap = defaultdict(int) links = [] diff --git a/phables/workflow/scripts/phables_utils/short_utils.py b/phables/workflow/scripts/phables_utils/short_utils.py index 4a94efb..780cf95 100644 --- a/phables/workflow/scripts/phables_utils/short_utils.py +++ b/phables/workflow/scripts/phables_utils/short_utils.py @@ -1,6 +1,7 @@ #!/usr/bin/env python3 import math +import pickle from concurrent.futures import ProcessPoolExecutor import logging import sys @@ -1613,6 +1614,24 @@ def resolve_short_parallel( if workers <= 1 or len(keys) <= 1: return resolve_short(pruned_vs=pruned_vs, **{**kwargs, "nthreads": nthreads}) + # Everything below has to cross a process boundary, so it all has to pickle. + # Checked up front rather than discovered when the pool starts: a failure + # there aborts the whole phables run, and losing a completed assembly to a + # performance optimisation is a bad trade. Hit for real -- oriented_links + # was a defaultdict built with a lambda, which cannot be pickled (fixed at + # source in edge_graph_utils._oriented_links_inner), and it took out a real + # 120-component run. Falling back to the sequential path keeps that a slow + # run instead of a failed one. + try: + pickle.dumps(kwargs) + except Exception as e: + logger.warning( + f"Cannot run components in parallel -- some input is not picklable " + f"({type(e).__name__}: {e}). Falling back to sequential; the result " + f"is unaffected, only the runtime." + ) + return resolve_short(pruned_vs=pruned_vs, **{**kwargs, "nthreads": nthreads}) + workers = min(workers, len(keys)) # Deliberately MORE chunks than workers. Component cost is heavily skewed -- From 277178ddb3b75044e0f5a0b3616f8563b7c79ea4 Mon Sep 17 00:00:00 2001 From: gbouras13 Date: Thu, 20 Aug 2026 21:51:05 +0930 Subject: [PATCH 12/12] Fix --mfd-workers on large assemblies: inherit inputs by fork, don't pickle --mfd-workers 8 still died on a real sample (SRR12983552: 485,547 vertices, 2733 components), seconds after the pool started, even with the oriented_links pickling fix. Cause: ProcessPoolExecutor(initargs=(kwargs,)) serialises the ENTIRE input set once per worker -- the igraph object plus every unitig sequence in graph_unitigs. On a real assembly that is gigabytes, sent eight times through a pipe, and multiprocessing cannot transfer a single object larger than ~2GB at all. It worked on the 22k-vertex sample and fell off that cliff on the 485k-vertex one. Workers now INHERIT the inputs: _WORKER_KWARGS is set in the parent before the pool is created, and the pool uses an explicit 'fork' context, so children get the parent's memory image with no serialisation and, via copy-on-write, without eight full copies. Measured with a 200MB stand-in payload: four workers receive it intact in 0.08s, where pickling it per worker costs ~0.5s and ~1.6GB copied for eight -- and the real payload is far larger. Two safety nets, because this failure destroyed a run that had already cost ~5 hours of GPU time: - if 'fork' is unavailable, run sequentially rather than attempt a copy that would fail or exhaust the node - if the pool breaks anyway (a worker OOM-killed surfaces as BrokenProcessPool), warn and retry sequentially instead of failing the run; verified by SIGKILLing a worker mid-flight _WORKER_KWARGS is cleared afterwards so the parent does not hold a second reference to the graph for the rest of the run. --- .../scripts/phables_utils/short_utils.py | 85 ++++++++++++++----- 1 file changed, 64 insertions(+), 21 deletions(-) diff --git a/phables/workflow/scripts/phables_utils/short_utils.py b/phables/workflow/scripts/phables_utils/short_utils.py index 780cf95..b9ee0f8 100644 --- a/phables/workflow/scripts/phables_utils/short_utils.py +++ b/phables/workflow/scripts/phables_utils/short_utils.py @@ -1,6 +1,7 @@ #!/usr/bin/env python3 import math +import multiprocessing import pickle from concurrent.futures import ProcessPoolExecutor import logging @@ -1509,24 +1510,26 @@ def resolve_short( ) -# Set once per worker process by _init_worker, so the heavy read-only inputs -# (the assembly graph above all) cross the process boundary ONCE per worker -# rather than once per chunk. Passing them as map() arguments instead would -# re-pickle the whole graph for every chunk, which gets more expensive the more -# chunks are used -- exactly backwards, since more chunks is what gives the pool -# room to balance an uneven workload. +# Set in the PARENT before the pool is created; forked workers inherit it +# through the process image, so none of it is ever serialised. +# +# This was previously passed as ProcessPoolExecutor(initargs=(kwargs,)), which +# pickles the whole thing once per worker. That is fine for a small graph and +# fatal for a real one: on a 485k-vertex assembly the payload is the igraph +# object plus every unitig's sequence in graph_unitigs -- gigabytes, serialised +# eight times through a pipe. It killed the pool seconds after startup, and +# multiprocessing additionally cannot send a single object larger than ~2GB at +# all. Inheriting by fork moves that cost to zero, and copy-on-write means the +# eight workers share the pages rather than each holding a full copy. _WORKER_KWARGS = None -def _init_worker(kwargs): - global _WORKER_KWARGS - _WORKER_KWARGS = kwargs - - def _resolve_short_chunk(chunk): """Worker entry point: one slice of components, nothing else. - Module-level (not a closure) so it is picklable by ProcessPoolExecutor. + Module-level (not a closure) so it is picklable by ProcessPoolExecutor -- + though only the small `chunk` is ever pickled; the bulk inputs arrive via + _WORKER_KWARGS, inherited from the parent. """ return resolve_short(pruned_vs=chunk, **_WORKER_KWARGS) @@ -1634,14 +1637,35 @@ def resolve_short_parallel( workers = min(workers, len(keys)) + # Workers must INHERIT the bulk inputs rather than be sent them. On a real + # assembly those inputs are gigabytes (see _WORKER_KWARGS), and pickling + # them per worker is both ruinously slow and subject to multiprocessing's + # ~2GB per-object ceiling. Fork gives the children the parent's memory + # image directly, at no serialisation cost and, thanks to copy-on-write, + # very little extra memory. + # + # If fork is unavailable (Windows; macOS defaults to spawn but can still + # fork explicitly) there is no cheap way to hand over that much data, so + # this runs sequentially rather than attempting a copy that would either + # fail outright or exhaust the node's memory. + try: + mp_context = multiprocessing.get_context("fork") + except ValueError: + logger.warning( + "The 'fork' start method is unavailable, so component-level " + "parallelism would have to copy the whole assembly graph to every " + "worker. Running sequentially instead; the result is unaffected, " + "only the runtime." + ) + return resolve_short(pruned_vs=pruned_vs, **{**kwargs, "nthreads": nthreads}) + # Deliberately MORE chunks than workers. Component cost is heavily skewed -- # a handful of large components dominate, and which ones is not known in # advance -- so splitting into exactly one chunk per worker (static # partitioning) lets a single worker draw several expensive components while # the rest sit idle. Smaller chunks let the pool hand out more work to - # whichever process finishes first. The usual cost of small chunks, - # re-sending the inputs each time, does not apply here: _init_worker sends - # them once per worker instead. + # whichever process finishes first, and cost nothing extra to send, since + # only the chunk itself crosses the boundary. size = max(1, math.ceil(len(keys) / (workers * CHUNKS_PER_WORKER))) chunks = [ {k: pruned_vs[k] for k in keys[i : i + size]} for i in range(0, len(keys), size) @@ -1652,12 +1676,31 @@ def resolve_short_parallel( f"in {len(chunks)} chunk(s)" ) - with ProcessPoolExecutor( - max_workers=workers, initializer=_init_worker, initargs=(kwargs,) - ) as executor: - # .map preserves input order, which is what keeps the merged - # all_resolved_paths identical to the sequential run. - results = list(executor.map(_resolve_short_chunk, chunks)) + # Published to the module namespace BEFORE the pool exists, so every forked + # child sees it. Cleared afterwards so the parent does not keep a second + # reference to the graph alive for the rest of the run. + global _WORKER_KWARGS + _WORKER_KWARGS = kwargs + try: + with ProcessPoolExecutor( + max_workers=workers, mp_context=mp_context + ) as executor: + # .map preserves input order, which is what keeps the merged + # all_resolved_paths identical to the sequential run. + results = list(executor.map(_resolve_short_chunk, chunks)) + except Exception as e: + # A worker dying (OOM above all -- eight processes touching a large + # graph can outgrow the node) surfaces here as BrokenProcessPool. Better + # to spend the extra wall-clock than to lose an assembly that already + # cost hours, so this retries the whole thing sequentially. + logger.warning( + f"Parallel component resolution failed ({type(e).__name__}: {e}). " + f"Falling back to sequential -- the result is unaffected, only the " + f"runtime. If this is memory, lower --mfd-workers." + ) + return resolve_short(pruned_vs=pruned_vs, **{**kwargs, "nthreads": nthreads}) + finally: + _WORKER_KWARGS = None # Field 1 (all_resolved_paths) and field 2 (all_components) are lists and # are extended; every other field is a set and is unioned. Merged in chunk