Skip to content

feat: add germline variant calling job bundle using bioconda - #275

Merged
mwiebe merged 1 commit into
aws-deadline:mainlinefrom
mwiebe:feat/variant-calling-bioconda-sample
Aug 6, 2026
Merged

feat: add germline variant calling job bundle using bioconda#275
mwiebe merged 1 commit into
aws-deadline:mainlinefrom
mwiebe:feat/variant-calling-bioconda-sample

Conversation

@mwiebe

@mwiebe mwiebe commented Aug 6, 2026

Copy link
Copy Markdown
Contributor

What was the problem/requirement? (What/Why)

The repository had no sequence-bioinformatics sample. The existing life-sciences
samples cover protein structure prediction (esmfold_predict),
molecular dynamics (gromacs_md), and molecular docking
(virtual_screening_vina/), but nothing that
starts from DNA sequencing reads. Variant calling is the workhorse analysis in genomics,
and it suits Deadline Cloud well because both expensive stages divide cleanly: samples
align independently, and separate regions of the genome are called independently.

No sample sourced software from the bioconda channel
either. Nearly all open source bioinformatics tooling is published there, so a customer
evaluating Deadline Cloud for genomics had no worked example of the channel setup that
this domain requires.

What was the solution? (How)

Adds job_bundles/variant_calling_bwa/, a germline variant calling pipeline that goes
from paired-end FASTQ reads to a merged VCF using bwa, samtools, and bcftools.

Five steps across two independent axes of parallelism:

QualityControl  (1 task/sample)     BuildIndex  (1 task)
  fastqc on each read pair            samtools faidx + bwa index
           |                                  |
           |                        AlignReads  (1 task/sample)
           |                          bwa mem | samtools sort
           |                                  |
           |                        CallVariants  (1 task/region)  <- the scatter
           |                          bcftools mpileup | call | filter
           +----------------+-----------------+
                            |
                  MergeVariants  (1 task)                          <- the gather
                    bcftools concat | norm, stats, MultiQC

AlignReads fans out over samples and CallVariants over genome regions, using
different task parameters in different steps. With the shipped defaults that comes to
10 tasks.

The tool sequence follows the
AWS HealthOmics WDL variant-calling tutorial pipeline,
reimplemented for Open Job Description. Having the same pipeline expressed in both
specifications makes their designs comparable. The WDL version runs each tool in its own
BioContainers image; because those images are themselves built from bioconda packages,
the tool builds are equivalent while the delivery mechanism collapses to a single
CondaPackages line. This bundle also adds the second scatter over regions, which the
original lacks, and takes paired-end reads where the original takes one FASTQ per sample.

Design decisions worth calling out, each recorded in the commit message and the sample
README:

  • CondaChannels defaults to conda-forge bioconda, in that order. bioconda depends
    on conda-forge for its runtime libraries and documents that priority. Reversing the
    order causes resolution failures that are hard to read.
  • The reference is a ReferenceDir directory plus a ReferenceFastaName filename
    rather than one FILE-typed parameter. Job attachments stages exactly the path a
    FILE parameter names and does not sweep in siblings, so a .fai or bwa index built
    beside the FASTA reaches the worker only if the whole directory is staged.
  • The shell lives in scripts/ rather than inline in the template, following stage 2
    of the job development progression. At roughly 400
    lines the bundle is past the point where a self-contained template stays readable.
    scripts/common.sh holds the per-region VCF filename derivation that CallVariants
    and MergeVariants must agree on.
  • Every path is derived from a job parameter rather than passed between steps, because
    each step runs in its own session directory, and each step declares every step whose
    output it reads
    , because dependency entitlements do not chain.

Test data is downloaded by sample_inputs/fetch_test_data.py (about 950 KB from the
public aws-genomics-static-us-east-1 bucket, originating from nf-core/test-datasets,
MIT licensed) rather than committed.

The change also fixes an unrelated repository validation failure: check_external_links.py
now reports link locations with as_posix(), because str() on a relative path yields
docs\guide.md on a Windows checkout and the test expected forward slashes.

What is the impact of this change?

Additive. One new job bundle, one row in the job_bundles/README.md category table, and
a one-line fix in scripts/check_external_links.py. No existing sample changes behavior.

The bundle requires a Linux fleet and a conda queue environment whose channels include
conda-forge and bioconda. When a tool is missing, the BioToolchain job environment
fails on entry and names the exact packages and channels to configure, rather than
failing partway into a step.

Running the job incurs Deadline Cloud worker, storage, and data transfer charges. The
README documents that, and notes that genomic sequence data is often subject to consent,
privacy, and jurisdictional restrictions, so access controls deserve attention before the
bundle is pointed at real data.

How was this change tested?

Submitted to a Deadline Cloud farm with the shipped defaults. All 10 tasks succeeded
on a Linux service-managed fleet in about 100 seconds of wall clock. The merged VCF held
variants for both samples:

=== Variant Calling Summary ===
Regions called: 4
Total variants: 17
Samples:        tiny_n tiny_t

--- variants per contig ---
  1            17

--- bcftools stats ---
number of samples:	2
number of records:	17
number of SNPs:	17
number of indels:	0

The alignment is genuinely good rather than vacuous. samtools flagstat for tiny_n:

1174 + 0 in total (QC-passed reads + QC-failed reads)
1164 + 0 mapped (99.15% : N/A)
1144 + 0 properly paired (97.44% : N/A)

Every expected artifact landed: variants.vcf.gz and its .tbi, four per-region VCFs
each with an index, two sorted BAMs with .bai and flagstat summaries, four FastQC
reports, and a MultiQC report. MultiQC's source manifest confirms it parsed the real
FastQC archives and flagstat files.

The bioconda path is confirmed in the worker log, which shows the job's channel override
taking effect over the queue environment's default:

CondaChannels: conda-forge bioconda
CondaPackages: bwa samtools bcftools fastqc multiqc

Repository validation. python3 scripts/validate_repository.py passes (35 unit tests,
static local-link checks). python3 scripts/check_external_links.py passes, which the
as_posix() fix was needed for on a Windows checkout.

Template and script checks. The template decodes and builds a job through the
openjd model library, producing the documented task counts (BuildIndex 1,
QualityControl 2, AlignReads 2, CallVariants 4, MergeVariants 1). bash -n is
clean on all seven scripts.

Shared helper behavior, exercised directly under set -euo pipefail, the mode the
scripts run in:

Input to parse_list Result
"tiny_n, tiny_t" [tiny_n] [tiny_t]
" lead , both , trail " [lead] [both] [trail]
"1:131000-133999, 1:134000-136999" region strings preserved intact

Trimming only the surrounding whitespace matters: deleting interior spaces too would
rewrite a name like tiny n to tinyn before validate_sample_names saw it, so the job
would fail on a filename the user never typed. An empty list parameter is rejected up
front, because zero entries would otherwise be reported by every per-index task as past
the end of the list while the gather found nothing missing, and the failure would surface
inside bcftools handed no BAM arguments.

Not covered by this testing: the run reused a cached named conda environment, so it
did not exercise a cold bioconda solve. An earlier end-to-end run did cover that path.

Was this change documented?

Yes.

  • job_bundles/variant_calling_bwa/README.md
    covers what the sample demonstrates, prerequisites, the conda channel order and why it
    matters, setup, local iteration with the Open Job Description CLI, a full parameter and
    output table, how steps pass files to each other, why the ...Range parameters exist,
    security and cost, and troubleshooting. It also states what the sample leaves out that
    a production analysis would add (duplicate marking, base quality score recalibration,
    adapter trimming, per-contig ploidy), so nobody mistakes it for a clinical-grade
    pipeline.
  • job_bundles/variant_calling_bwa/sample_inputs/README.md
    documents the test data provenance and licensing.
  • The new row in job_bundles/README.md follows the existing
    category table format.
  • Root navigation is unchanged, because no recommended starting point changed.

By submitting this pull request, I confirm that you can use, modify, copy, and redistribute this contribution, under the terms of your choice.

@mwiebe
mwiebe requested a review from a team as a code owner August 6, 2026 05:50
@github-actions github-actions Bot added the waiting-on-maintainers Waiting on the maintainers to review. label Aug 6, 2026
@mwiebe
mwiebe force-pushed the feat/variant-calling-bioconda-sample branch from 834337c to 56f4740 Compare August 6, 2026 05:52
Adds a sequence-bioinformatics sample that goes from paired-end FASTQ reads to a
merged VCF using bwa, samtools, and bcftools installed from the bioconda channel.
This is the first sample in the repository to source software from bioconda rather
than conda-forge alone.

The pipeline aligns each sample in parallel, calls variants in parallel across
reference regions (the scatter), and concatenates the per-region calls into a
single VCF (the gather). It reimplements the AWS HealthOmics WDL variant-calling
tutorial pipeline, whose per-tool BioContainers images are themselves built from
bioconda, so the tool builds are equivalent while the delivery mechanism collapses
to a single CondaPackages line.

Notes on the design:

- CondaChannels defaults to "conda-forge bioconda" in that order. bioconda depends
  on conda-forge for its runtime libraries and documents this priority; reversing
  it causes resolution failures.
- The reference is a ReferenceDir directory plus a ReferenceFastaName filename
  within it, rather than one FILE-typed path. Job attachments stages exactly the
  path a FILE parameter names and does not sweep in siblings, so a .fai or bwa
  index built beside the FASTA only reaches the worker if the whole directory is
  staged. BuildIndex reuses a staged .fai and bwa index and builds whichever are
  missing, requiring all five bwa index files before reusing that set so a
  partial one is rebuilt rather than copied for bwa mem to fail on later.
- The shell for each step lives in scripts/ rather than embedded in the template,
  following stage 2 of the job development progression. Values reach the scripts
  as named command line arguments, since a staged script cannot use {{Param}}
  substitution. scripts/common.sh holds the per-region VCF filename derivation
  that CallVariants and MergeVariants both have to agree on, plus the shared
  comma-separated list parsing.
- parse_list trims only the whitespace surrounding each entry. Deleting interior
  spaces as well would rewrite a name like 'tiny n' to 'tinyn' before
  validate_sample_names ever saw it, failing on a filename the user never typed
  rather than on the name they did.
- Every list parameter is checked for at least one entry. An empty list parses to
  zero elements, which each per-index task reports as past the end of the list
  while the gather's completeness checks find nothing missing, so without the
  check the job would fail inside bcftools handed no BAM arguments instead of
  naming the parameter that was left empty.
- Samples and Regions are each paired with a SampleRange/RegionRange task
  parameter range rather than a maximum index, matching how the other samples
  parameterize a range. A range also expresses a sparse selection (0,2) or a
  single index, which is what re-running one failed region needs.
- Every path is derived from a job parameter rather than passed between steps,
  because each step runs in its own session directory.
- Steps declare every step whose output they read, since dependency entitlements
  do not chain.
- Regions defaults to four windows of contig 1 covering where the sample reads
  align, so each scatter task calls real variants. The list stays a parameter so
  the whole scatter is visible in the submission command, which caps it at around
  100 entries against the 1024-character parameter limit; the README sketches
  moving the list to a staged data file when a real whole-genome scatter needs
  more.

Test data is downloaded by sample_inputs/fetch_test_data.py from the public
nf-core/test-datasets mirror (MIT licensed) rather than committed.

Also reports external link locations with as_posix(), fixing a repository
validation test that failed on Windows checkouts because str() on a relative
path yields the native separator.

Verified end to end on a Linux service-managed fleet: all tasks succeeded and the
merged VCF held variants for both samples.

Signed-off-by: Mark <399551+mwiebe@users.noreply.github.com>
@mwiebe
mwiebe force-pushed the feat/variant-calling-bioconda-sample branch from 56f4740 to 3146d0e Compare August 6, 2026 06:03

@crowecawcaw crowecawcaw left a comment

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Cool! Does this workflow generate any interesting visuals that we could use as a hero image for the sample README?

@mwiebe

mwiebe commented Aug 6, 2026

Copy link
Copy Markdown
Contributor Author

Cool! Does this workflow generate any interesting visuals that we could use as a hero image for the sample README?

Not that I'm aware of, but I didn't dig into what options the tools it uses have for visualizing results. Would be a nice follow-up!

@mwiebe
mwiebe merged commit f846f1a into aws-deadline:mainline Aug 6, 2026
11 checks passed
@mwiebe
mwiebe deleted the feat/variant-calling-bioconda-sample branch August 6, 2026 20:22
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

waiting-on-maintainers Waiting on the maintainers to review.

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants