Skip to content

Lessons Learned

Every failure encountered during pipeline development (Mar 2026), documented so they don't happen again.

Docker Image Issues

AnnotSV: Official image doesn't exist

  • Failed: bioinfochrustrasbourg/annotsv:3.4.4 — no such image on Docker Hub
  • Fix (then): Used getwilds/annotsv:3.4.4 (Fred Hutch maintained). That image is retired here; the pipeline now pins the Bioconda build quay.io/biocontainers/annotsv (see ANNOTSV_IMAGE in versions.env).

SnpEff/SnpSift: Combined package

  • Failed: quay.io/biocontainers/snpsift:5.2--hdfd78af_1 — no such manifest
  • Fix: Use quay.io/biocontainers/snpeff:5.2--hdfd78af_1 — SnpSift is bundled inside the SnpEff package

SnpEff Database: Azure blob storage outage

  • Failed: snpEff download GRCh38.105 and all database names — Azure blob storage returned 0-byte files for ALL URLs including GRCh38.105, GRCh38.mane.1.2.ensembl
  • Root cause: SnpEff databases are hosted on Azure blob storage which was down
  • Fix: Pivot to Ensembl VEP as alternative functional annotation tool

ExpansionHunter: Different image formats

  • Failed: quay.io/biocontainers/expansionhunter:5.0.0--hd03093a_1 — manifest not found. Also mgibio/expansionhunter:latest — not found
  • Fix: Use weisburd/expansionhunter:latest — binary at /ExpansionHunter/bin/ExpansionHunter, variant catalogs at /pathogenic_repeats/GRCh38/
  • Note: This is ExpansionHunter v2.5.5 which uses --repeat-specs (directory), not --variant-catalog (single JSON)

ExpansionHunter: Missing required --log parameter

  • Failed: Container exits immediately with "the option '--log' is required but missing"
  • Fix: Always include --log /output/sample_eh.log in the command

StellarPGx: Empty Docker image

  • Failed: twesigomwedavid/stellarpgx:latest — image exists but contains no StellarPGx binaries
  • Status: UNSOLVED. No working Docker image found for StellarPGx as of Mar 2026.

Tool-Specific Issues

TelomereHunter: Permission denied

  • Failed: OSError: [Errno 13] Permission denied: '/output/<sample>' when writing output (the sample name in the error path will vary)
  • Fix: Add --user root flag to docker run

TelomereHunter: pip install on host fails

  • Failed: pip install telomerehunter on the host gives UnicodeDecodeError — Python environment issues
  • Fix: Use the Docker image (lgalarno/telomerehunter, digest-pinned in versions.env) instead of native install

HLA-LA: Graph not serialized (3 failures)

  • Failed attempt 1: zlskidmore/hla-la:latest — graph at src/additionalReferences/PRG_MHC_GRCh38_withIMGT/ exists but is NOT serialized. HLA-LA exits silently.
  • Failed attempt 2: Copied graph files to graphs/PRG_MHC_GRCh38_withIMGT/ — "graph not complete"
  • Failed attempt 3: Ran --action prepareGraph in detached container — container exited during prep, graph still not serialized
  • Root cause: HLA-LA requires a pre-serialized graph (~40GB RAM to prepare, takes hours). No standard Docker image includes it.
  • Fix: Use jiachenzdocker/hla-la:latest (27.5GB image with pre-built graph) OR switch to T1K

HLA-LA: Binary crash with pre-built graph image

  • Failed: jiachenzdocker/hla-la:latest — read extraction succeeds but HLA-LA C++ binary crashes during graph alignment even with 32GB RAM and 8 threads. Error: "HLA-LA execution not successful."
  • Root cause: Likely an incompatibility between the pre-built binary and the BAM data format, or an unmet memory requirement (the graph deserialization may need >32GB)
  • Status: UNSOLVED. HLA-LA from WGS BAMs is unreliable in Docker. Alternative: a clinical lab result, if you have one, or arcas-hla or T1K with partial coordinates

T1K: Coordinate file with wrong values

  • Failed: t1k-build.pl -d hla.dat -g reference.fasta.fai produced coordinate file with chr19 -1 -1 + for all HLA genes
  • Root cause: The -g parameter expects the actual reference FASTA (3.1GB), not the FAI index (158KB). The AddGeneCoord.pl script needs to align HLA sequences against the genome to find coordinates.
  • Fix: Use -g Homo_sapiens_assembly38.fasta (the full FASTA, not the .fai)

T1K: BAM extraction produces 0-byte FASTQ

  • Failed: bam-extractor runs but produces empty _candidate_1.fq and _candidate_2.fq
  • Root cause: Coordinate file had -1 -1 coordinates (see above), so no genomic region was extracted
  • Fix: Regenerate coordinate file with proper reference FASTA

Cyrius CYP2D6: Inconclusive on short-read WGS

  • Result: Cyrius returned None for CYP2D6 star alleles
  • Root cause: CYP2D6 has extensive pseudogene homology (CYP2D7, CYP2D8) making short-read WGS unreliable
  • Mitigation: Use a clinical lab result, if you have one, for CYP2D6; consider long-read sequencing in future

bcftools/htslib Issues

"The index file is older than the data file" on published BAMs and VCFs

  • Observed: samtools warned about the published <sample>_sorted.bam.bai, and htslib about the VEP VCF's .tbi, on every region query after a Nextflow run, although each task wrote its index after its data file.
  • Root cause: publishDir (mode copy) copies a task's output files in parallel. The small index finishes first and the large BAM or VCF gets the later time.
  • Why it matters: the warning is noise there, and it hides the real case: an old index left beside a BAM or VCF that was rewritten, which serves the wrong regions without an error.
  • Fix: main.nf's completion handler runs after the last copy. It sets each published index that this run wrote, and that is older than its data file, to the data file's time. An older index from before the run, or a hard link into work/, is left as it is and named in a warning. validate-setup.sh warns about any index under the sample folder older than its data file: index it again when the data changed, touch it when the file was only copied.

bgzip/tabix not in bcftools image PATH

  • Failed: staphb/bcftools:1.21 does not include bgzip or tabix in $PATH
  • Fix: Use bcftools view -Oz -o output.vcf.gz (native bgzip output) and bcftools index -t (native tabix index) instead of piping to bgzip/tabix
  • Alternative: Use quay.io/biocontainers/samtools:1.21 which includes all htslib tools

MIS VCF conversion: 0-byte output files

  • Failed: All 22 chr files in mis_ready/ were 0 bytes
  • Root cause: Script used bgzip -c > output.vcf.gz which failed silently because bgzip wasn't available
  • Fix: See above — use bcftools view -Oz -o

VEP Cache Issues

VEP INSTALL.pl: Permission denied on temp directory

  • Failed: Cannot open Local file /opt/vep/.vep/tmp/homo_sapiens_vep_112_GRCh38.tar.gz
  • Fix: Run with --user root and pre-create the temp dir: mkdir -p /opt/vep/.vep/tmp && chmod 777 /opt/vep/.vep/tmp

VEP INSTALL.pl: Silent download failure

  • Failed: Container exits after "downloading..." with no cache files extracted
  • Root cause: The 17GB download timed out or was interrupted; INSTALL.pl doesn't retry
  • Fix: Download manually with wget -c (supports resume), then extract with tar xzf, then run VEP with --cache --dir_cache /path/to/cache

General Docker Tips

Always use resource limits

docker run --cpus 4 --memory 8g ...
Without limits, tools like DeepVariant or minimap2 will consume ALL available RAM and crash the host.

Build for amd64 from macOS

docker build --platform linux/amd64 ...
macOS is arm64; most Linux servers are amd64. Images built on Mac without --platform won't run on amd64 servers.

Use --rm for one-shot containers

Always use --rm for analysis containers to avoid accumulating stopped containers. Use -d (detached) for long-running jobs.

Always use --user root for write access

Most bioinformatics containers run as non-root users. If writing to bind-mounted volumes, add --user root to avoid permission issues.

CI / Workflow Issues

ShellCheck warnings still fail the GitHub Action

  • Observed: ludeeus/action-shellcheck@master exits non-zero even when configured with severity: warning
  • Impact: Once the check is required on main, "warning-only" findings still block merges
  • Fix: Clear ShellCheck warnings before enabling required status checks, or explicitly relax the workflow instead of assuming warnings are advisory only

Protect main only after CI is green

  • Observed: Required status checks become a trap if you enable branch protection while the default branch or active PR branch is still red
  • Fix: Get Lint and Smoke Tests green first, then enable required checks, block force-pushes, and block deletion

MToolBox Issues

MToolBox: No working Docker image exists

  • Failed: robertopreste/mtoolbox:latest — "repository does not exist or may require docker login"
  • Also checked: No image on quay.io/biocontainers, ghcr.io, or Docker Hub
  • Root cause: MToolBox was never officially containerized. GitHub issue #107 (Mar 2022) confirms: "Not at the moment."
  • Fix: Use GATK Mutect2 in mitochondrial mode instead (broadinstitute/gatk:4.6.2.0, the version the scripts pin). Mutect2 handles mitochondrial heteroplasmy detection natively and is well-maintained.

CNVnator Issues

CNVnator: Biocontainer tag with wrong build hash

  • Failed: quay.io/biocontainers/cnvnator:0.4.1--py312hc02a2a2_7 — manifest not found
  • Fix: Use quay.io/biocontainers/cnvnator:0.4.1--py312h99c8fb2_11. Biocontainer hashes encode the conda build hash and change between builds. Always verify at quay.io/repository/biocontainers/cnvnator.

Delly Issues

Delly: Biocontainer tag doesn't exist

  • Failed: quay.io/biocontainers/delly:1.2.9--ha41ced6_0 — manifest not found
  • Fix: Use quay.io/biocontainers/delly:1.7.3--hd6466ae_0 (latest as of Mar 2026). Biocontainer tags are version-specific and change frequently — always verify at quay.io/repository/biocontainers/delly.

Delly: SV annotation phase takes 2-3 hours

  • Observed: Delly's "SV annotation" step runs for 2-3 hours at 100% CPU on a 30X WGS genome. No new log output appears during this time, which can look like the process is stuck.
  • This is normal. Delly genotypes every candidate SV site against the reference, which is CPU-intensive. Total runtime for 30X WGS: ~3-4 hours.
  • Tip: Use docker stats to confirm the container is still using CPU. If CPU is at 0%, the process may actually be stuck.

Delly: Output is BCF format, not VCF

  • Gotcha: Delly writes BCF (binary VCF), not VCF. The output file has a .bcf extension.
  • Fix: Convert with bcftools view input.bcf -Oz -o output.vcf.gz and index with bcftools index -t output.vcf.gz. The pipeline script handles this automatically.

CNVnator: ROOT file appears empty (266 bytes) during tree extraction

  • Observed: During the -tree step, the .root file stays at 266 bytes (just the ROOT header) until the entire BAM is parsed.
  • This is normal. For a 30X WGS (~80GB BAM), the tree step takes ~5-10 minutes. The ROOT file grows to ~900MB-1.2GB only at the very end when the tree is flushed to disk.
  • If the container exits and the file is still 266 bytes: Check if a corrupt .root file from a previous failed run is blocking it. Delete and retry.

CPSR/PCGR Issues

CPSR: --pcgr_dir path confusion (PCGR 1.x, historical)

  • Failed: cpsr --pcgr_dir /genome/pcgr_data/data → "Data directory (/genome/pcgr_data/data/data) does not exist"
  • Root cause: CPSR 1.x internally appends /data to whatever --pcgr_dir you pass. If you point to the data/ directory inside the extracted bundle, it looks for data/data/.
  • Fix (1.x): Point --pcgr_dir to the parent of the data/ directory: --pcgr_dir /genome/pcgr_data (not /genome/pcgr_data/data)
  • Superseded by PCGR 2.x — the --pcgr_dir flag no longer exists. See migration notes below.

CPSR: Docker image is inside PCGR

  • Failed: sigven/cpsr:2.0.0 does not exist on Docker Hub
  • Fix: Use sigven/pcgr:2.2.5 which bundles both pcgr and cpsr binaries at /usr/local/bin/

PCGR 2.x Migration (1.4.1 to 2.2.5)

  • CLI completely changed: The --pcgr_dir flag is gone. Replaced by --refdata_dir (for the ref data bundle) and --vep_dir (for the VEP cache). These are separate mount points inside the container.
  • Docker volumes changed: PCGR 1.x used a single -v ${GENOME_DIR}:/genome mount. PCGR 2.x requires four separate mounts: VEP cache (:/mnt/.vep), ref data bundle (:/mnt/bundle), input VCFs (:/mnt/inputs), and outputs (:/mnt/outputs).
  • Data bundle is smaller and different: The old monolithic ~21 GB bundle (pcgr.databundle.grch38.20220203.tgz) that included VEP cache is replaced by a smaller ~5 GB ref data bundle (pcgr_ref_data.20250314.grch38.tgz). VEP cache is now provided separately (reuse the same cache from step 13).
  • Bundle extraction requires extra step: After tar xzf, the extracted data/ directory must be moved into a version-stamped directory: mkdir -p 20250314 && mv data/ 20250314/. The --refdata_dir mount points to this version directory.
  • Docker tag: sigven/pcgr:1.4.1 → sigven/pcgr:2.2.5. The image still bundles both pcgr and cpsr binaries.
  • Old data bundle is incompatible: If you have the 1.x bundle, you must download the 2.x bundle fresh. The directory structure and expected paths are completely different.

Michigan Imputation Server Notes

Minimum 20 samples per job

  • MIS is designed for genotyping array data from cohort studies
  • Single-sample WGS submissions may be rejected (20-sample minimum)
  • For individual WGS: imputation adds minimal value (you already have 90%+ variant coverage)
  • Main benefit for WGS would be phasing, not imputation

Registration required

  • Must create account at imputationserver.sph.umich.edu
  • API tokens expire after 30 days
  • Results auto-deleted after 7 days

TOPMed panel is best for Europeans

  • TOPMed Freeze 8 (r2): 132K samples, 705M variants, GRCh38 native
  • Uses chr prefix (which GRCh38 BAMs already have)
  • HRC r1.1 (32K samples) is European-centric but hg19 only

Alternative Variant Caller Issues

FreeBayes 1.3.7: SIGILL crash (exit code 132)

  • Failed: quay.io/biocontainers/freebayes:1.3.7--h1870644_0 — freebayes --version works, but actual variant calling triggers SIGILL (illegal instruction, exit code 132)
  • Tested on: Intel i5-14500
  • Root cause: Likely a build-time CPU optimization mismatch in the 1.3.7 biocontainer binary
  • Fix: Use quay.io/biocontainers/freebayes:1.3.6--hbfe0e7f_2 which works correctly

FreeBayes: Memory grows to ~13 GB on full genome

  • Observed: FreeBayes memory usage grows unpredictably during full-genome runs: 463MB at 30 min, 6.4GB at 60 min, 12.7GB at 90 min, then stabilizes ~12GB
  • Original limit: --memory 16g was too tight — would have OOM-killed at 80% usage
  • Fix: Use --memory 32g for full-genome runs. Peak observed was 12.8GB but growth is non-linear and region-dependent.

FreeBayes: Single-threaded, no parallelism

  • Observed: FreeBayes has no -t or --threads flag. Full 30X WGS takes 8-12 hours.
  • Workaround: Use --region chr22 (or INTERVALS=chr22) for quick testing (~20-40 min)
  • For production: Consider GNU parallel with per-chromosome regions, then merge VCFs

GATK HaplotypeCaller: bcftools index fails on existing .tbi

  • Failed: bcftools index -t fails with "index file exists" after GATK already creates its own .tbi
  • Fix: Use bcftools index -ft (with -f force flag) to overwrite the GATK-generated index

GATK HaplotypeCaller: Requires .dict file

  • Failed: GATK HaplotypeCaller fails if Homo_sapiens_assembly38.dict is missing
  • Fix: Generate once with gatk CreateSequenceDictionary -R /genome/reference/Homo_sapiens_assembly38.fasta

bcftools isec: -R vs -r for region strings

  • Failed: bcftools isec -R chr22 treats -R (uppercase) as a BED file path, fails with "file not found"
  • Fix: Use -r chr22 (lowercase) for region strings. -R expects a file.

TIDDIT >=3.9: Requires BWA index for local assembly

  • Failed: tiddit --sv exits with "The reference must be indexed using bwa index; run bwa index, or skip local assembly (--skip_assembly)"
  • Root cause: TIDDIT 3.9+ uses local assembly for breakpoint refinement, which requires BWA index files alongside the reference
  • Fix: Use --skip_assembly when no classic BWA index is available. TIDDIT's assembly step calls classic bwa, which needs the classic index files (.amb .ann .bwt .pac .sa). BWA-MEM2's index (.bwt.2bit.64, .0123) is a different format and does not work, even after aligning with BWA-MEM2.

TIDDIT: Image tag 3.7.0 doesn't exist on quay.io

  • Failed: quay.io/biocontainers/tiddit:3.7.0--py312h24f4cff_1 — manifest unknown
  • Fix: Use quay.io/biocontainers/tiddit:3.9.5--py312h6e8b409_0. Always verify tags at quay.io/repository/biocontainers/tiddit.

Strelka2: --callRegions needs bgzipped + tabixed BED

  • Failed: --callRegions reference.fasta.fai → "Can't find expected call-regions bed index file"
  • Fix: Create a proper bgzipped BED file with tabix index. Use GATK container for bgzip/tabix (not in bcftools or samtools staphb images).

bgzip/tabix not in staphb/samtools or staphb/bcftools images

  • Observed: Neither staphb/samtools:1.20 nor staphb/bcftools:1.21 include bgzip or tabix in PATH
  • Fix: Use broadinstitute/gatk:4.6.1.0 which has both at /usr/bin/bgzip and /usr/bin/tabix. Or use bcftools view -Oz as a bgzip alternative.

FreeBayes chr22 variant count (3x more than DeepVariant)

  • Observed: FreeBayes calls ~247K variants on chr22 vs DeepVariant ~91K and GATK ~69K
  • Interpretation: The ~200K FreeBayes-unique variants are mostly false positives. FreeBayes maximizes sensitivity at the cost of precision.
  • Recommendation: Always quality-filter FreeBayes output with bcftools filter or vcffilter before use.

Chip Data Conversion (Genotyping Arrays → VCF)

  • What failed: plink --23file (1.9) to import + plink2 --ref-from-fa force to fix REF/ALT
  • Why: For single-sample data, ALL homozygous positions are monomorphic. plink's .bim stores only one allele for these. --ref-from-fa cannot create a proper ALT because there's no second allele slot. Homozygous ALT genotypes silently become homozygous REF.
  • Verified: on one homozygous-alternate chip position (chip genotype AA where the reference base is T): plink wrote REF=A, ALT=., GT=0/0 (WRONG); bcftools wrote REF=T, ALT=A, GT=1/1 (CORRECT). ~66K positions (11%) corrupted.
  • Fix: Use bcftools convert --tsv2vcf -f <reference.fa>. Single command, no intermediate binary format.
  • --allow-extra-chr cannot be used with --23file
  • Female samples with Y calls (MyHeritage GSA PAR region) error with sex=F
  • Sex inference defaults to male unless explicitly set

MyHeritage CSV must be converted to TSV

  • Quoted CSV with "RSID","CHROMOSOME","POSITION","RESULT" columns
  • Strip ## comments, header, quotes; convert commas to tabs

bcftools hg19 VCF needs chr prefix before liftover

  • hg19 reference uses numeric chromosomes; chain file expects chr prefix
  • bcftools annotate --rename-chrs between conversion and liftover

PharmCAT chip vs WGS results (MyHeritage GSA, verified 2026-03-31)

  • Correct: CYP2B6, CYP4F2, DPYD and NUDT15 matched the WGS diplotypes
  • Missed: CYP2C19 (25 missing), VKORC1 (1 missing)
  • Wrong: CYP3A5 came out as the reference diplotype because 4 of its defining positions are not on the chip; WGS called a different diplotype
  • Fix: treat a chip-only call as unresolved when any of its defining positions are missing from the chip, and confirm it from WGS before relying on it
  • Total: 888 missing PGx positions from the GSA chip

ROH and PRS need special flags for chip data

  • ROH: -G30 required (no FORMAT/PL in chip VCF)
  • PRS: no-mean-imputation required (single sample lacks allele frequencies)
  • PRS matching: chip ~12% of large scoring files vs WGS ~28%

v0.3.0 Tool Additions (Apr 2026)

ExpansionHunter v5.0.0: Completely different CLI from v2.5.5

  • Old (v2.5.5): --bam, --ref-fasta, --repeat-specs (directory), --vcf, --json, --log (all required)
  • New (v5.0.0): --reads, --reference, --variant-catalog (single JSON), --output-prefix, --threads
  • The --log flag is gone in v5.0.0. Do NOT pass it or the command will fail.
  • Biocontainer image: quay.io/biocontainers/expansionhunter:5.0.0--hc26b3af_5. Binary is ExpansionHunter (on PATH, not at /ExpansionHunter/bin/).
  • Bundled catalog: /usr/local/share/ExpansionHunter/variant_catalog/grch38/variant_catalog.json (31 loci). No need to download separately.

GRIDSS: Requires BWA index (not minimap2)

  • Failed: GRIDSS exits with "BWA index not found" when using default minimap2 alignment
  • Fix: Build the classic BWA index (bwa index on the reference, about an hour). GRIDSS runs classic bwa internally, so BWA-MEM2's index files (.bwt.2bit.64) do not help, and aligning with BWA-MEM2 first does not provide the index it needs. The 04b-gridss.sh script checks for the classic files and prints the command.
  • Note: GRIDSS outputs ALL SVs as BND (breakend) notation. Standard DEL/DUP/INV types require post-processing conversion for SURVIVOR merge compatibility.

GRIDSS: 32 GB memory requirement

  • Observed: GRIDSS assembly-based SV calling needs ~28 GB JVM heap for 30X WGS
  • Fix: Container runs with --memory 32g and -Xmx28g JVM argument. Will fail silently on machines with < 32 GB RAM.

GRIDSS: ENCODE blacklist download

  • Observed: GRIDSS benefits from an ENCODE blacklist to suppress known artifact regions
  • Fix: Script auto-downloads ENCFF356LFX.bed.gz (hg38 blacklist) from ENCODE on first run. If download fails (offline), GRIDSS runs without it (lower precision but still functional).

fastp: Maximum 16 threads despite --workers flag

  • Observed: fastp accepts -w (workers) up to 16. Values above 16 are clamped to 16. The -w flag controls I/O worker threads; actual adapter detection is single-threaded.
  • Recommendation: Use -w matching THREADS up to 16. For most WGS runs, -w 4 is sufficient.

fastp: BGI/MGI adapter auto-detection

  • Observed: fastp's --detect_adapter_for_pe works for Illumina, BGI, and MGI adapters without specifying adapter sequences. BGI/MGI adapters are compiled into fastp's knownadapters.h.
  • No action needed: The --detect_adapter_for_pe flag handles all common sequencing platforms.

mosdepth: --fast-mode skips per-base output

  • Observed: --fast-mode uses a simpler, faster counting method and does NOT write the per-base .per-base.bed.gz file. This saves ~2 GB of output and cuts runtime by ~40%.
  • Fix: Always use --fast-mode unless per-base resolution is specifically needed.

MultiQC: Auto-discovers fastp JSON by content, not filename

  • Observed: MultiQC identifies fastp output by looking for "before_filtering": { in JSON files, not by filename pattern. Files must end in .json.
  • Tip: Ensure fastp's -j output uses .json extension.

Octopus: No issues observed

  • Docker image dancooke/octopus:0.7.4 works out of the box for germline calling
  • Supports --threads for parallelism (unlike FreeBayes)
  • Typical memory: 8-12 GB for 30X WGS (much less than FreeBayes peak of 13 GB)

PharmCAT 3.x Migration (2.15.5 to 3.2.0)

Preprocessor script renamed (no .py extension)

  • Old (2.15.5): python3 /pharmcat/pharmcat_vcf_preprocessor.py
  • New (3.2.0): python3 /pharmcat/pharmcat_vcf_preprocessor
  • Impact: Step 7 preprocessor command must drop the .py suffix or the container exits with "No such file"

Reporter flags: must be explicit for both formats

  • Old (2.15.5): -reporterJson produced JSON; HTML was always generated by default
  • New (3.2.0): If ANY format flag is specified, ONLY those formats are produced. To get both HTML and JSON, you must pass -reporterJson -reporterHtml
  • Impact: Step 7 now passes both flags explicitly. Without -reporterHtml, the HTML report (used for manual review) would silently stop being generated.

JSON property rename: wildtypeAllele to referenceAllele

  • Old (2.15.5): wildtypeAllele property in gene result JSON objects
  • New (3.2.0): Renamed to referenceAllele
  • Impact: Step 27's JSON parser does not use this property directly, so no code change was needed. Any downstream scripts or notebooks that parse wildtypeAllele must be updated.

New features in 3.2.0

  • NAT2 calling: PharmCAT 3.x includes improved NAT2 acetylator status calling
  • BCF support: Preprocessor now accepts BCF input files directly (no conversion needed)
  • Single-gene calling: New -g flag allows running PharmCAT on a single gene (useful for targeted re-analysis)

Nextflow version compatibility (2026-06)

The pipeline runs on Nextflow 25.10.x and later; 24.x fails at parse time

  • Observed: A full run requires the validated Nextflow release (NEXTFLOW_VERSION in versions.env; 25.10.4 when this was found, 25.10.8 until October 2026, 26.04.7 until 2026-10-14, 26.10.x since). Other versions fail before any process executes:
  • 26.04.4 — the strict config parser rejects top-level def/variable declarations in nextflow.config ("Variable declarations cannot be mixed with config statements"), and then the def check_max(...) function in conf/base.config ("Unexpected input: '('").
  • 24.04.4 — the DSL2 module parser flags the optional annotation inputs in modules/local/vcfanno/main.nf as "Variable already defined in the process scope" (cadd_snv/cadd_indel/spliceai_*/revel/alphamissense, referenced inside the def has_nochr/def has_chr expressions). 25.10.4 tolerates this; 24.04.4 does not.
  • Fix status: nextflow.config is strict-parser-clean — the execution-report timestamp is inlined into each report path (no top-level def; see #30/#31), which also preserves per-run report history. conf/base.config has since replaced check_max() with process.resourceLimits, which clears the second 26.x parse error above. The last 26.x error was VCF_PRECHECK setting FILTER_STATUS and FILTER_COUNTS without declaring them as env() outputs; it declares them now. CI (nextflow.yml) runs nextflow lint and the all-tools stub run on both the pinned release and the newest 26.10.x release; the real-data runs (e2e.yml) use the pinned release. 25.10.x still passes the manifest floor, but CI no longer runs it. The vcfanno input scope is not refactored, so 24.x still fails.
  • Tip: NXF_VER=26.10.5 nextflow run main.nf ... runs the release the real-data tests use; run-all.sh sets NXF_VER from versions.env unless you set it yourself. The manifest.nextflowVersion floor is !>=25.10.4, so 24.x is rejected up front; there is no cap, and 26.10.x passes the lint and stub checks.

CYP2D6 copy number: no single caller settles it (withdrawn lesson)

  • Withdrawn: an earlier version of this entry told readers to trust pypgx on the BAM over Cyrius and PharmCAT for CYP2D6 deletion and duplication alleles, because pypgx made a deletion call where Cyrius returned None/None and PharmCAT No Result. That conclusion does not hold.
  • Why: pypgx calls copy number from read depth. On a reference with ALT contigs and an aligner that is not run ALT-aware, reads from the CYP2D locus can split between the primary copy and an ALT copy, depth on the primary drops, and a depth-based caller can report a whole-gene deletion that is not there.
  • Rule: on a reference with ALT contigs, compare CYP2D6 depth with its flanks before trusting any copy-number call from pypgx or Cyrius, and report CYP2D6 only when two callers agree. A Cyrius None/None means "no call", not "no deletion" and not "deletion". If you have a clinical lab result for CYP2D6, it outranks all of these.

A reference with ALT contigs and an aligner that is not run ALT-aware thin the depth at paralogous loci

  • Failed: setup.sh installed the Broad Homo_sapiens_assembly38.fasta (3,366 sequences, with ALT, HLA and decoy contigs), and step 02 aligned with plain minimap2 -x sr. A read that matches a primary locus and its ALT copy equally well got MAPQ 0, and callers skip MAPQ 0 reads, so depth thinned at CYP2D6, the MHC and KIR. Measured on the HG002 fixture (ALT depth A/B workflow): depth at MAPQ >= 1 was 0.0x at CYP2D6 and HLA-A where the no-ALT analysis set gives 16.8x and 26.9x, and 85% of the MHC reads were aligned to ALT or HLA contigs. See realignment.
  • Root cause: ALT contigs only help an aligner that is run ALT-aware (BWA-MEM with the .alt file and its post-processing, or DRAGEN's graph reference); the pipeline runs none. The extra contigs also forced per-contig workarounds in CNVpytor and Delly.
  • Fix: the default reference is NCBI's GRCh38 no-ALT analysis set, reference/GRCh38_no_alt_analysis_set.fasta (195 sequences), under a new name so an old file is never read by accident. validate-setup.sh fails on a reference with _alt or HLA- contigs (unless ALLOW_ALT_REFERENCE=true) and on a BAM whose @SQ names and lengths differ from the .fai, naming the first difference.
  • Rule: a change of reference means realigning every sample, from FASTQ or from the old BAM (realignment). Never run a step on a BAM whose header does not match the reference; validate-setup.sh checks it.

CNVpytor migration (2026-07)

CNVpytor 1.3.2 biocontainer ships without GC/mask data and its downloader is broken

  • Failed: cnvpytor -his aborts with Some reference genome resource files are missing. Run 'cnvpytor -download' — the quay.io/biocontainers/cnvpytor:1.3.2--pyhdfd78af_0 image's cnvpytor/data/ dir contains only an empty readme.txt.
  • Failed: cnvpytor -download itself crashes in 1.3.2 (AttributeError: 'PosixPath' object has no attribute 'split' in genome.py), so it cannot self-heal.
  • Fix: Pre-download the pinned v1.3.2 GC/mask files and bind-mount them onto the container's package data dir (/usr/local/lib/python3.12/site-packages/cnvpytor/data). Runs then work fully offline (verified with --network none). See docs/00-reference-setup.md.

CNVpytor's resource check requires every genome's files to exist, not just hg38

  • Observed: genome.py:check_resources() iterates every bundled reference genome (hg19, hg38, chm13v2.0, chm13v1.1, kn99) and os.path.exists()-checks each gc_file/mask_file. Mounting only gc_hg38.pytor+mask_hg38.pytor still fails the check.
  • Fix: Provide all seven files. Only gc_hg38.pytor/mask_hg38.pytor are actually read for an hg38 BAM; the others just need to exist (the check is existence-only). All seven total ~90 MB.

CNVpytor container has no bcftools/bgzip

  • Observed: Unlike the old CNVnator biocontainer, cnvpytor:1.3.2 bundles no bcftools/bgzip/tabix/samtools.
  • Fix: The bash step runs a separate staphb/bcftools container for VCF normalization; the Nextflow module splits into two processes (CNVPYTOR calls → CNVPYTOR_VCF reheader/sort/index). cnvpytor -view emits a proper VCFv4.2 (SVTYPE/END/SVLEN, ALT DEL/DUP/LOH, GT/CN) but only carries ##contig lines for processed chromosomes — reheader from the reference .fai before merging with Manta/Delly in step 22.

Real-data lessons from full-WGS runs (2026-07)

CNVpytor chokes on GRCh38 ALT/HLA/decoy contigs — restrict -rd to canonical chromosomes

  • Failed: On a real full-reference GRCh38 BAM (hundreds of ALT/HLA/decoy contigs), cnvpytor -rd (no -chrom) crashes/stalls during read-depth import and produces no calls (only a tiny stub .pytor). The bundled GC-correction data covers only the main chromosomes. A chr20-only validation did not surface this.
  • Fix: Pass -chrom chr1 … chr22 chrX chrY to the -rd step (scripts/18-cnvpytor.sh + modules/local/cnvpytor). Validated on two real 30× genomes (≈2200 and ≈2400 canonical CNV calls, ~99% concordant with the prior CNVnator counts). Lesson: validate depth-based callers on a full BAM, not a single chromosome.

pypgx image must match the bundle version — 0.27.0 mismatched 0.26.0 and shipped pandas 3.0.3

  • Failed: pypgx:0.27.0--pyh106432d_0 failed all genes at runtime with pandas.errors.LossySetitemError / TypeError: Invalid value — pypgx assigns floats into int columns, which pandas ≥2.1 (the 0.27.0 image bundles 3.0.3) rejects. It was also mismatched against the pipeline's 0.26.0 pypgx-bundle. The smoke test only runs pypgx --version, so CI could not catch it.
  • Fix (v0.8.2): pin the image to pypgx:0.26.0--pyh7e72e81_0 to match the bundle (pandas 2.3.3, no error). Verified end-to-end on a real 30× genome: 23/23 genes succeed. Lesson: the pypgx image version MUST track the pypgx-bundle branch, and bumping one without the other (and without a real gene-call test) silently breaks step 32.

Stranger over-flags RFC1 (CANVAS) from short reads — do not read it as a diagnosis

  • Observed: Stranger can report RFC1 STR_STATUS=full_mutation for a modest expansion (e.g. 51/73 of the degenerate AARRG motif). CANVAS requires the AAGGG motif specifically, biallelic, at ~400–2000+ repeats — short-read ExpansionHunter cannot resolve AAGGG vs the benign AAAAG, and the catalog's STR_PATHOLOGIC_MIN=12 is not the clinical threshold.
  • Interpretation: Treat an RFC1 flag as uninterpretable from short-read WGS — confirm with motif-aware/flanking-PCR testing only if clinically indicated. (Documented in docs/09b-stranger.md.)
  • Since then: ExpansionHunter's count still cannot tell the motifs apart, but when Stranger flags RFC1, step 9b counts the motif in the reads over the repeat (bin/rfc1_motif.py) and gives a verdict; without the sample's BAM and its index the verdict is not_checked (STR Annotation, Notes).

Vendor VCF intake (2026-10)

PharmCAT's Java step stops on a backslash in a valid header line

  • Failed: pharmcat.jar exits with Error parsing metadata: character to be escaped is missing when a ##FILTER or ##INFO line holds a backslash. bcftools writes such a line itself for a soft filter with a quoted string (bcftools filter -s LowDP -e 'FORMAT/DP<10 && GT!="0/0"'), and the VCF 4.3 spec asks for the escape. The preprocessor copies the line through; the bundled vcf-parser 0.3.1 refuses it. PharmCAT 3.4.0 bundles the same parser.
  • Fix: step 07 and the PHARMCAT module rewrite the header of PharmCAT's own copy (\" to ', any other \ to /, on ## lines) before the Java step. The calls are the same as without the line. Remove the rewrite once a PharmCAT release bundles a newer vcf-parser; bumping PharmCAT is not the fix.

Ensembl contig names fail silently

  • Failed: with contigs named 1, MT, bcftools view -r chrM returns no record with exit 0, so haplogrep3 writes a header-only file, and the ROH summary's chrX|chrY filter lets X through. Only the ClinVar screen failed loudly, and only when selected.
  • Fix: VCF_PRECHECK stops the run before any analysis when no contig holding records is chr-named, and prints the bcftools annotate --rename-chrs command. Lesson: a step that selects a region by name must not read zero records as a clean result.

bcftools norm -f stops at the first contig the reference lacks

  • Failed: exit 255 (The sequence "NT_113889.1" was not found) on the full chr-renamed ClinVar file, and on a sample record on a scaffold the reference lacks (a half-done rename, or the no-alt reference).
  • Fix: the ClinVar screen reads both files only on the contigs they share with each other and the reference, with a targets file whose end is a constant (chrom 1 2147483647): a header without contig lengths would otherwise give an empty region and zero records.

workflow.onComplete saw workflow as null

  • Failed: Cannot get property 'success' on null object on every run under Nextflow 25.10: the handler runs with the script binding's variable map as its delegate, and a map answers null for a name it lacks.
  • Fix: the handler reads local variables (run_info, run_log, outdir) set in the workflow body; closures resolve those where they are written.

PGx consensus and the step backlog (2026-10)

The planned CYP2D6 depth rule passed the broken BAM and flagged the clean one

  • Failed: the plan for step 32's CYP2D6 depth check was a ratio rule: unreliable when CYP2D6 has less than 0.6 of the depth of its flanks. On the e2e fixture's CYP2D slice, mapped to the no-ALT reference, CYP2D6 has 0.62 of the flank depth with all reads and 0.57 at MAPQ >= 1, its normal value; mapped to the Broad hg38 FASTA with ALT contigs, CYP2D6 keeps 0.28 of the flank depth with all reads. On those measured depths the ratio cut at 0.6 flagged the clean BAM and passed the broken one (CYP2D6 depth check).
  • Why: CYP2D6 sits beside CYP2D7 and CYP2D8, so its depth is below the flanks on a correct alignment too. What an ALT contig changes is where the reads go: a third of them move to chr22_KI270879v1_alt, and the ones left get MAPQ 0.
  • Fix: the shipped rule (bin/cyp2d6_depth_check.py) is the share of reads with MAPQ >= 1: unreliable when more than 40% of the reads at CYP2D6 (when it has at least 15% of the flank depth) or in the flanks have MAPQ 0. The share is 0.91 without ALT contigs and 0.005 with them. The ALT depth A/B workflow runs the check on both mappings and fails unless it passes the first and flags the second.
  • Rule: measure a threshold on both the good and the broken input before writing it down; a ratio that sounds right can sit on the wrong side of the normal value.

SURVIVOR names its columns after the input samples and exits 0 on a missing file

  • Failed (found while replacing step 22's position binning): SURVIVOR merge takes each output column's name from the input VCF's sample. Three callers of one sample give three columns with one name, which bcftools then refuses. And when it cannot open an input it prints a message and exits 0.
  • Fix: step 22 and SURVIVOR_PREP rename each input's sample to the caller's name, and the step checks that SURVIVOR wrote a VCF instead of trusting the exit code.

Yleaf downloads the whole hg38 FASTA unless its config names one

  • Failed: Yleaf 3.2.1 checks for a full reference before anything else and downloads hg38 when its config file (inside the image, read-only) names none, also for a BAM, whose pileup never needs the sequence. Under --network none that fails.
  • Fix: step 37 and Y_HAPLOGROUP set Yleaf's reference constant to the pipeline's FASTA in a one-line Python launcher before calling Yleaf's main().

The Yleaf image has no samtools, and Yleaf then hangs instead of failing

  • Failed: the image tests of Yleaf ran into their one-hour limit. The Bioconda recipe of Yleaf 3.2.1 lists no samtools, which Yleaf calls for a BAM (samtools idxstats, samtools mpileup). Yleaf runs each BAM in a multiprocessing.Pool worker and stops on the failed call with SystemExit; that kills the worker, and Pool.map waits for a result that never comes. With a serial pool the same run stops at once with samtools: command not found.
  • Also missing: the package installs Yleaf's Python code only, without its data folder (marker positions, haplogroup tree): setup.py lists no package data. setup.sh --yleaf-data installs that folder from the release archive, checked by sha256, and the launcher points Yleaf at it.
  • Fix: bin/yleaf_run.py: step 37 and the pipeline make the idxstats and the pileup at Yleaf's marker positions in SAMTOOLS_IMAGE (with Yleaf's own flags), and run Yleaf on them with its samtools calls served from those files, its reference constant set and its pools replaced by a serial map. The image test checks that the functions and calls the launcher replaces are still in the image.

Dependency bumps (2026-10)

The first -resume on Nextflow 26.04 caches nothing from a 25.10 work directory

  • Failed: after the Nextflow 26.04.7 bump, a -resume over a work directory written by 25.10 reran every process. 26.04 hashes a map input by its keys and values where 25.10 hashed it differently, and every process here takes the meta map, so every task hash changed.
  • Fix: none needed in the code; the first run after the upgrade starts a fresh work directory, and later resumes behave as before. A BAM that an earlier run already produced goes in through the BAM entry of the samplesheet instead of a resume.

Renovate's commitMessagePrefix is used as written

  • Failed: "commitMessagePrefix": "deps" gave PR titles like deps Update python Docker tag to v3.14.8: no colon, so not a conventional commit, while this page and updating.md said deps:.
  • Fix: semanticCommits enabled with semanticCommitType deps and an empty scope (ci for the lint tools), which writes deps: update ... itself. The five PRs of 2026-10-09 were retitled by hand before the squash.

config:recommended turned the deps: titles back into chore:

  • Failed: Renovate's PR titles started with chore: although semanticCommitType said deps. config:recommended includes the :semanticPrefixFixDepsChoreOthers preset, whose package rules set chore for every package and fix for runtime dependencies. Package rules win over the top-level setting, and a dependency found by a custom regex manager has no dependency type, so every PR got chore:.
  • Fix: :semanticCommitTypeAll(deps) after config:recommended in extends: its rule sets deps for every package and comes later, so it wins. The ci rule for the lint tools and the Renovate CLI comes later still. The dry run lists the new pins; the first Renovate PR after the change shows the title.

Nextflow's local executor (2026-10)

The local executor can lose track of free memory, and a task that needs all of it then never starts

  • Failed: the Nextflow CI step that fails six kinds of stub task once on purpose (exit 137, as a kill for memory would) and reads the doubled memory of each retry went silent until the 45-minute job timeout, twice on real commits. Its config told the local executor it had 16 CPUs and 64 GB, and three of the retries ask for 64 GB each. A loop of that one step on Nextflow 26.04.7, on GitHub's 4-CPU runners, stopped the same way in 22 of 600 runs, and in 34 of 600 more runs with trace logging on.
  • Why: Nextflow 26.04.7 keeps the local executor's free CPUs and memory in two plain fields (availCpus and availMemory in LocalPollingMonitor). One thread subtracts a task's request when it starts the task, another adds it back when the task ends, and nothing locks the two. When both happen at the same moment one of the changes is lost. In a hung run captured with trace logging, all 30 submitted tasks had completed, .nextflow.log said No more task to compute, and the executor still counted availCpus: 16 and availMemory: 62 GB of 64. Two GB had gone missing from the count. The two retries waiting for 64 GB could not start, and no running task was left to free anything. In 33 of the 34 hangs captured that way the count of free memory was short. In the other one memory was right and the count of free CPUs was 4 of 16. In four of the 34 the executor counted more free CPUs than the 16 it was given, which is the same race the other way round: a lost subtraction lets it start more than it was told it has.
  • Fix: nothing in this repo can repair the count, so the CI step no longer depends on it. It tells the executor it has 512 CPUs and 2048 GB, more than the 248 CPUs and 856 GB that all 46 task attempts of the run ask for added together, and fails when that stops being true. With that config the loop ran 600 times with no hang. The step also has a 10-minute limit of its own, and a failed job prints the tasks that were running or waiting and uploads .nextflow.log.
  • In a real run: the same code counts the machine's CPUs and RAM, so a real run can stop the same way; we have not seen one do it. The loss needs one task to end in the instant another starts. Stub tasks that fail within milliseconds do that often; how often real tools do is not measured. What to look for and what to do is in Troubleshooting. We reported the race upstream as nextflow-io/nextflow#7816.
  • Fix status: Nextflow 26.10.0 keeps both counts in atomic counters (ResourceTracker, from nextflow-io/nextflow#7783, made for another purpose), so a start and an end at the same moment no longer lose a change. The releases before it have the race (we read the source of 25.10.4, 26.04.0 and 26.04.7). The pipeline moved to 26.10.0 on 2026-10-14. The CI step keeps its large executor, which also covers the newest-patch leg and anyone who runs it on an older release.

Polygenic scores on the default list (2026-10)

pgsc_calc's own retry gives its first task more threads under the same memory cap

  • Failed: in one run on a ~30x sample, the Nextflow PRS path stopped at its last tool. pgsc_calc v2.3.0's FORMAT_SCOREFILES ran on 2 CPUs in an 8 GB container and was killed (exit 137); pgsc_calc retried it on 4 CPUs in the same 8 GB, where it died with BrokenProcessPool (exit 1). pgsc_calc then exited 1, which our retry list does not cover, so the run ended. In another run PRS_PREPARE was killed at 2 GB and again at 4 GB.
  • Why: pgsc_calc's process_medium label asks for 2 CPUs and 16 GB times the attempt, capped at --max_memory, and pgscatalog-format holds one scoring file per thread in memory. Under a fixed cap the retry raises the need and keeps the limit. Measured on a GitHub runner on the default list: about 14.2 GB on one thread, 15.7 GB on two, 23.8 GB on four.
  • Fix: both entry points give FORMAT_SCOREFILES one thread, the whole of pgsc_calc's memory and no retry inside pgsc_calc. PRS asks for 24 GB and turns a memory kill inside pgsc_calc into exit 137, so Nextflow's retry doubles the cap; PRS_PREPARE asks for 16 GB. tests/test_prs_memory.sh checks the requests against the recorded figures on every pull request, and measures them again on dispatch.

maxRetries = 0 still allows one retry

  • Failed: to stop pgsc_calc retrying FORMAT_SCOREFILES under the same memory cap, the PRS steps first gave it maxRetries = 0. The Nextflow CI stand-in for pgsc_calc, which exits 2 on any attempt after the first, showed pgsc_calc retried it anyway.
  • Why: Nextflow 26.04.7, and still 26.10.0, reads maxRetries as result ? result as int : 1 (TaskConfig.getMaxRetries). Zero is false in Groovy, so maxRetries = 0 means one retry.
  • Fix: the PRS steps set errorStrategy = 'finish' on FORMAT_SCOREFILES, which never retries. To turn retries off for a process, set its errorStrategy, not maxRetries.

A PRS summary with a Variants_Used column comes from the old scoring path

  • Failed: before 2026-10-01 step 25 scored with plink2 directly and wrote <sample>_prs_summary.tsv with the columns Condition, PGS_ID, Score, Variants_Used, Variants_Total. Score held plink2's SCORE1_AVG (an average, not the sum), Variants_Used held NAMED_ALLELE_DOSAGE_SUM (a dosage sum, not a count of matched variants), and the condition labels came from a hard-coded list that named PGS000738, a vitiligo score, as schizophrenia.
  • Why: the script read the last two columns of plink2's .sscore by position and took the labels from its own list instead of the PGS Catalog.
  • Fix: the summary now has Score_SUM and Variants_Matched and labels each score with the catalog's trait. Such old files can still sit in a sample folder. Do not reuse them or compare them with new ones: a header with Variants_Used marks one. Rerun step 25 instead (Keeping versions up to date).

Fine-grained ancestry (2026-10)

ADMIXTURE's bootstrap crashes when it is combined with projection

  • Failed: admixture -P -B100 sample.bed 6 printed Bootstrap replicate #0 and died with a segmentation fault, on the 1.3.0 BioContainer, on every input we tried: a planted person, a real one, 235 SNVs and 1,003 SNVs.
  • Why: not investigated beyond the crash; 1.3.0 is the only build BioContainers publishes, so a newer ADMIXTURE would need an image built here.
  • Fix: step 26b runs its own block bootstrap. bin/ancestry_fine_report.py sample-bed writes run_admixture.sh, which projects the sample once per level, then for each replicate draws 100 blocks of neighbouring SNVs with replacement (awk, seeded by the replicate number), writes them as a PLINK .ped with the matching rows of the frequency file, and projects that. 100 replicates of 1,003 SNVs take about a second in total.

ADMIXTURE reads a .ped the other way round from a .bed

  • Failed: the same genotypes projected as a .bed and as a .ped gave different fractions (0.96 against 0.46 for the largest group).
  • Why: in PLINK's "12" coding, allele 2 is the .bim's A1 (what plink2 wrote as the ALT allele here), not allele 1. Writing homozygous-ALT as 1 1 mirrors every genotype.
  • Fix: PED_ALLELES in bin/ancestry_fine_report.py maps the .bed codes to 2 2, 1 2 and 1 1, and tests/test_ancestry_fine_report.py checks both encodings of the same three genotypes.

ADMIXTURE stops on a SNV with no genotype, which for one person is every missing site

  • Failed: admixture -P exited 255 with detected that all genotypes are missing for a SNP locus.
  • Why: ADMIXTURE refuses a SNV whose genotypes are all missing. With one person in the input, every site where that person has no genotype is such a SNV.
  • Fix: the sample's input holds only the SNVs it has a genotype for, and each level's .P.in is cut to the same rows, so the frequencies still line up.

RFMix dies with a segmentation fault on a short chromosome

  • Failed: RFMix 2.03 printed mm::allocate() ... requesting zero or negative size (0) at setup_ref_haplotypes() and segfaulted on the fixture's slices (0.1 to 2.7 cM each), while the same command on chr6:29.9-33.1 Mb and on chr22:20-35 Mb ran to the end.
  • Why: not investigated beyond its message: on those slices setup_ref_haplotypes() asks for a buffer of size zero, and the longer ranges do not hit it.
  • Fix: the step paints a chromosome only with 100 genotyped SNVs or more spanning 5 cM or more (bin/ancestry_fine_report.py chroms, from the panel's genetic map). Every autosome of a whole genome spans more than 40 cM; on a slice panel the step says it painted nothing instead of crashing.

htslib writes a remote index into the working directory

  • Failed: bcftools view https://.../hgdp1kgp_chr1...bcf inside the container could not create hgdp1kgp_chr1....bcf.csi in its working directory, so it could not read the index and failed.
  • Why: htslib saves the index of a remote file in the current directory, and the container's working directory is not writable.
  • Fix: the install runs that command with the writable panel folder as its working directory and deletes the downloaded .csi files when it ends.

E2E parity checks (2026-10)

A tabix index differs between two runs with the same records

  • Failed: family-joint-1 checks that a single-sample run publishes the same files as main. It failed on three runs in a row (PRs 149 and 151, unrelated changes), each time on HG002/clinvar/HG002_pass.vcf.gz.tbi alone, with a different md5 every time.
  • Why: the check hashes each file's text without its date and command lines, but an index is binary. It holds byte offsets into the bgzip file, and the VCF header carries bcftools' dated command line, so the first block compresses to a different size from run to run and every offset moves.
  • Fix: the check lists .tbi and .csi files by name only, so a missing index still fails it. The records behind them are compared through the VCF itself.

The raw-data export's site table (2026-10)

grep -Fwf with 620,000 rsIDs never finished reading dbSNP

  • Failed: the first full build of step 40's site table (setup.sh --consumer-sites) read 0.1% of dbSNP's 30 GB GRCh38 VCF in two hours, with grep at one full core. On a 200,000-line slice (30 MB) the same command printed nothing in 20 minutes. CI never builds the full table, so nothing had run it at full size.
  • Why: we did not profile grep. What we measured: with the same 620,000 patterns, GNU grep 3.11 is fast without -w and stalls with it. The likely cause is the word check: many rsIDs are prefixes of others (rs1 of rs10, rs100, ...), so nearly every line has hits that fail the check and send grep back to retry the other patterns there.
  • Fix: each rsID is wrapped in the tabs that surround dbSNP's ID column and matched with grep -F alone. On a 5-million-line slice it took 4.6 s, kept the same lines as an exact lookup in awk, and gave a table identical to Python reading the slice unfiltered. The full table builds in about 25 minutes.

A rerun after a failed build downloaded dbSNP again

  • Failed: the build's failure message said the dbSNP download was kept and to run the command again, but the rerun started the 30 GB download from zero.
  • Why: fetch always downloads into <name>.part and renames it when the checksum matches. It never looks at a finished file under the final name, so the kept file was not used.
  • Fix: install_consumer_sites skips the download when the finished file exists. That file is whole: only fetch writes it, after the md5 matched.

Pharmacogenomics from the gVCF (2026-10)

PharmCAT sets aside a one-base record at a position longer than one base

  • Failed: after step 7 expanded the gVCF's reference blocks into one record per base, the alleles PharmCAT defines by an indel (DPYD, UGT1A1 *28, TPMT, CYP3A4 and others) stayed uncallable, even where the gVCF showed the reference over the whole span.
  • Why: PharmCAT 3.4.0's preprocessor does not read a one-base 0/0 record as a call at a position whose REF or ALT is longer than one base. It sets the record aside with FILTER PCATxINDEL (or PCATxREF for a two-base substitution) and lists the position as missing; where PharmCAT's REF is longer than one base, the matcher then discards it.
  • Fix: bin/pgx_refblock_fill.py writes PharmCAT's own record (its ID, REF and ALT) as 0/0 at each such position when every base of the span lies in reference blocks with GQ 20 or more and nothing else overlaps it. The details are in Pharmacogenomics, Input the preprocessor or PharmCAT refuses.

PharmCAT reads MT-RNR1 only as an outside call

  • Failed: MT-RNR1 (aminoglycoside-induced hearing loss) was No Result on every sample, although step 20 calls chrM with Mutect2 at high depth.
  • Why: PharmCAT 3.4.0 has no MT-RNR1 allele definition for its matcher, so chrM records in its VCF input are never matched for the gene. It lists MT-RNR1 among the genes it prefers as an outside call, and its own MT-RNR1 test passes the gene that way.
  • Fix: bin/pgx_mito_merge.py turns step 20's chrM calls into one MT-RNR1 outside call (a homoplasmic named allele, or Reference only when the BAM covers every named position), added to step 36's outside calls. Step 27 checks only CYP2D6, HLA-A and HLA-B against step 36's table: MT-RNR1 never appears there, and checking it marked the new call stale. The details are in Pharmacogenomics, Input the preprocessor or PharmCAT refuses.

PharmCAT's log file needs a writable working directory

  • Failed: every step 07 run printed a logback status dump with a pharmcat.log (Permission denied) stack trace. The calls were not affected.
  • Why: the jar's built-in logback.xml writes pharmcat.log by a relative path, into the working directory. In the image that is /pharmcat, which the calling user cannot write. No -D property moves the file.
  • Fix: step 07 runs the jar with -w /tmp. A Nextflow task already runs in its own writable work directory, so the module needed no change.

Nextflow's task wrapper header (2026-10)

Nextflow 26.10.3 quotes task names with double quotes

  • Failed: the nextflow (26.10) job, which runs the newest 26.10.x release, counted 0 T1K_BUILD, HLA_TYPING, INDEXCOV and DEEPVARIANT tasks on every branch from 2026-10-20, although the run's log shows them completed.
  • Why: each task's .command.run starts with a YAML header. Up to 26.10.2 Nextflow wrote ### name: 'X'; 26.10.3 writes ### name: "X" (BashWrapperBuilder now renders the values as YAML double-quoted strings). The counts matched single quotes only.
  • Fix: every reader of that header (four steps in .github/workflows/nextflow.yml and tests/e2e/couple-carrier-1.sh) accepts either quote, so the pinned release and the newest one both count.

MUC1 VNTR with VNtyper (2026-10)

VNtyper's image runs only its own command words

  • Failed: docker run VNTYPER_IMAGE bash -c ... and Nextflow's task wrapper started VNtyper's entrypoint, which accepts vntyper, celery, beat, health and bash, and ignores the arguments of bash.
  • Why: the image's /usr/local/bin/entrypoint.sh dispatches on its first word. VNtyper also reads its references through paths relative to /opt/vntyper.
  • Fix: step 52 passes --entrypoint "" and runs cd /opt/vntyper && vntyper pipeline ...; nextflow.config's docker profile gives the VNTYPER process the same --entrypoint "".

VNtyper refuses an output folder inside the BAM's folder

  • Failed: the VNTYPER process read "not assessed" on a BAM the bash step called positive. Its log said Pipeline log file is inside an operator-owned input tree.
  • Why: VNtyper 2.0.47 will not write its log under the folder of --bam, the link's folder or the file's. Nextflow stages the BAM as a link in the task folder, which is where run/ goes. The bash step writes to vntyper/, beside aligned/, so it never hit this.
  • Fix: the process links the BAM and its index into in/ and passes in/<bam>; run/ is then outside both folders. The e2e case vntyper-1 compares the process's table with the bash step's byte for byte.

Nextflow's trace file (2026-10)

The trace file can list a finished task twice

  • Failed: the nextflow (pinned) job's retry-memory step failed at random, once with FAIL MANTA: 5 first and 4 second attempts in the trace and once with FAIL DEEPVARIANT: 3 first and 2 second attempts in the trace. A rerun of the job passed both times. No retry was missing: the run's log shows every task submitted once and retried once.
  • Why: the extra row was the same failed first attempt, written a second time as the last line of the trace. Nextflow's trace observer keeps a map of running tasks: a submit or start event puts a task in it, the complete event removes it and writes its row, and at the end of the run every task still in the map is written again. A stub task that exits 137 at once can have its submit event handled after its complete event (the failing run's log prints the retry note before the task's own "Submitted process" line), so the task goes back into the map and gets a second row. The code is the same on 26.10.2 and on Nextflow's main branch.
  • Fix: the step adds task_id to the trace fields and counts each task once. A repeated row must match the first one (only its exit status may read -), and each task name must have exactly one first attempt and one retry. The two steps that read back the memory retries of ALIGN_MINIMAP2 and PRS key their rows on task_id the same way. Any other trace reader that counts rows should too.

Mobile element insertions with MEGAnE (2026-10)

MEGAnE stops on a BAM that lacks one of its GRCh38 contigs

  • Failed: call_genotype_38 on a BAM of chr20 alone stopped with "Chromosomes ... were NOT found" before reading a read.
  • Why: the wrapper passes -mainchr with its own list of every GRCh38 contig, ALT contigs included, and MEGAnE refuses a BAM that misses any of them. A BAM aligned to the no-ALT analysis set has none of the ALT contigs.
  • Fix: step 50 writes the list from the BAM's header (chr1 to chr22, chrX, chrY as the BAM names them) and passes it as -mainchr after the wrapper's own, so it wins.

MEGAnE stops when chrX is missing, unless told the sample has no sex chromosomes

  • Failed: with that list, a chr20-only BAM stopped with "Sex chromosome {'chrX', 'X'} was NOT found".
  • Why: MEGAnE measures the depth of chrX to genotype it against the autosomes.
  • Fix: when the BAM has no chrX, step 50 passes -no_sex_chr, which reads every chromosome as an autosome.

MEGAnE refuses a .fai that is not newer than its FASTA

  • Failed: MEGAnE stops with "Fasta index ... was older than fasta" when the .fai is not strictly newer than the FASTA, which a copied or unpacked reference can be (MEGAnE's scripts/setup.py).
  • Why: MEGAnE checks that the index was written after the FASTA. A Nextflow task stages both as links, which carry the times of the files they point to.
  • Fix: step 50 and the MEI process give MEGAnE a link to the FASTA with a fresh copy of the .fai beside it.

On one chromosome MEGAnE writes MEI_final_failed_genotyped.vcf, not the percentile VCF

  • Failed: on HG002's chr20 the step found no MEI_final_percentile_genotyped.vcf, the file MEGAnE's log recommends for one sample, although the run ended with exit 0.
  • Why: MEGAnE fits the read support of its candidates to set its filter thresholds and needs at least three clean candidates for the fit. One chromosome gave fewer, so it logged "Cannot automatically determine thresholds" and wrote its calls under fixed thresholds to MEI_final_failed_genotyped.vcf.
  • Fix: step 50 reads the percentile VCF and falls back to that file with a warning. A whole genome gives the percentile VCF.

MEGAnE crashes on a few hundred kilobases of reads

  • Failed: on the image test's 500 kb slice, call_genotype_38 stopped with an IndexError in process_mapped_seq.py (check_true_pos).
  • Why: with so few reads one of its read lists is empty, and the code indexes it anyway.
  • Fix: the image test builds the k-mer set and checks the options; the call is tested on all of chr20 in the e2e case mei-1.

MEGAnE's image has no ps, and every Nextflow task in it stops

  • Failed: the MEI process exited 1 within seconds with "Command 'ps' required by nextflow to collect task metrics cannot be found"; the bash step on the same BAM ran.
  • Why: with a trace, report or timeline on, Nextflow's task wrapper lists the task's processes with ps -e -o pid= -o ppid= to measure its memory, and checks for ps first. MEGAnE's image does not have it, and it is the only image there is.
  • Fix: modules/local/mei/resources/usr/bin/ps prints those two columns from /proc. It is a module binary (nextflow.enable.moduleBinaries = true), so it is on the PATH of the three tasks of modules/local/mei only.

The staphb samtools image cannot read a BAM over https

  • Failed: streaming the GIAB BAM's chr20 with SAMTOOLS_IMAGE failed with libcurl error 60.
  • Why: the image has no CA certificate bundle.
  • Fix: the e2e case mounts the runner's /etc/ssl/certs and sets CURL_CA_BUNDLE, as the LPA case does.

MEGAnE's image leaves the e2e runner short of disk

  • Failed: with steps 50, 51 and 52 all on main, the e2e job stopped with "No space left on device" in the cases after mei-1 (mocha-1, later the nextflow-* cases). Removing only the case's files was not enough: the free disk went from 7.6 GB to 8.9 GB, and the job still ran out.
  • Why: the e2e job pulls every image before the first case, and MEGAnE's image is about 6.4 GB to download and much more unpacked. It stayed on the disk after mei-1, the only case that uses it.
  • Fix: tests/e2e/mei-1.sh removes its whole sample directory and the image by its id when it ends, and logs the free disk before and after. With the image gone every case passed. The image is pinned by digest, so the e2e image cache never saves it anyway.
  • Image sizes step: the job's "Image sizes" step then stopped with numfmt: invalid number: ''. For an image that is gone, docker image inspect -f prints an empty line before the || echo 0 fallback, so the step's own "not here" branch never ran. The fix in .github/workflows/e2e.yml reads the size into a variable first (size=$(docker image inspect ...) || size=0) and prints ${size:-0}.

Downloads (2026-11)

One bad transfer of a large file failed the whole setup

  • Failed: setup.sh --archaic stopped on a wrong md5 for the 1000 Genomes phase 3 sites VCF (about 1.5 GB) in the e2e case archaic-ancestry-1. Other runs on the same day downloaded the same file with the right md5.
  • Why: fetch tried a download again only when the transfer itself failed. A transfer that finished with the wrong bytes failed its check, and fetch deleted it and returned 1 at once, so one bad transfer ended the run.
  • Fix: a download that fails its checksum, gzip -t or non-empty check is deleted and downloaded again from the start, up to FETCH_TRIES tries FETCH_WAIT apart, the same budget as a broken transfer. A broken transfer still keeps its .part and the next try resumes it. tests/test_fetch_retry.sh serves a bad file once and then the right one, and a bad file every time, and checks both against a mutant that returns at the first failed check.

CI flakes (2026-11)

AnnotSV's merge run in setup.sh failed in two seconds

  • Failed: tests/test_annotsv_benign_mane.sh failed now and then with setup.sh --annotsv-data exiting 1, and the rerun passed. Its log showed AnnotSV's merge run ending with exit 125 two seconds after it started, and the gnomAD-SV file left unmerged.
  • Why: the merge run was the first use of the AnnotSV image on the runner, so docker run pulled it. One pull broke off (short read: expected 603618355 bytes but got 8388608: unexpected EOF), and docker run does not try a pull again: it exits 125 before AnnotSV starts.
  • Fix: annotsv_merge_benign pulls the image first with pull_image, which tries FETCH_TRIES times. The Repeat a test workflow (.github/workflows/repeat-test.yml) runs the test 30 times on fresh runners to show it no longer flakes.

A slow download used up an e2e case's hour

  • Failed: archaic-ancestry-1 ran into the 3600-second case timeout on one e2e run (run 40933547961, attempt 1) and pushed the job past its own limit; the rerun passed.
  • Why: setup.sh --archaic downloads the 1000 Genomes phase 3 sites VCF (1.36 GB) from EBI. On that run EBI sent it at about 440 KB/s, and the download alone took 53 minutes, so the case's Nextflow half started with minutes left. The download never stopped: it was slow, not stalled. MoChA's chr20 panel in mocha-1 comes from the same server.
  • Fix: the e2e job caches both files (keyed on their pins) and scripts/ci/e2e-run.sh hands them to the two cases, so a run with a cache hit downloads neither. Separately, _get now bounds a real stall, which nothing did before: curl gets --connect-timeout 30 --speed-limit 1024 --speed-time ${FETCH_STALL:-120}, so a transfer under 1 KB/s for two minutes fails and fetch tries again, resuming the .part. wget gets --connect-timeout=30 --read-timeout=${FETCH_STALL:-120}, which fails only when no data arrives for that long. A slow but moving transfer like the one above is left alone.