From 3a13913a4654b4dc7564bafe55399374ab5507dc Mon Sep 17 00:00:00 2001 From: GeiserX <9169332+GeiserX@users.noreply.github.com> Date: Sat, 7 Nov 2026 17:14:35 +0100 Subject: [PATCH 01/10] fix(report): read a header-only Y tree path table as an empty path y_tree_place.py writes the header alone when no branch is derived; the report marked the whole Y tree section unreadable. Require the Branch column in the header, not a row. Bead: pgp-5bja.25 --- bin/collect_summary.py | 11 ++++++++--- tests/test_collect_summary.py | 12 ++++++++++-- 2 files changed, 18 insertions(+), 5 deletions(-) diff --git a/bin/collect_summary.py b/bin/collect_summary.py index 0caf4f4b..6dba8cea 100755 --- a/bin/collect_summary.py +++ b/bin/collect_summary.py @@ -2638,9 +2638,14 @@ def sec_y_tree(d, s): "explore": haplogroup_explore("y", hg)} lvl = first_existing(d, [f"y_haplogroup/{s}_y_tree_path.tsv"]) if lvl: - prows = read_tsv(lvl) - if not prows or "Branch" not in prows[0]: - raise Unreadable(f"{os.path.basename(lvl)} has no Branch column or no rows") + # bin/y_tree_place.py writes the header alone when no branch is derived: + # an empty path, not an unreadable file. + with open_text(lvl) as f: + rd = csv.DictReader(f, delimiter="\t") + prows = list(rd) + header = rd.fieldnames or [] + if "Branch" not in header: + raise Unreadable(f"{os.path.basename(lvl)} has no Branch column") out["path"] = [{k: (row.get(k) or "") for k in ("Level", "Branch", "SNPs", "On_tables", "Derived", "Ancestral", "No_call", "Ancestral_SNPs")} for row in prows] diff --git a/tests/test_collect_summary.py b/tests/test_collect_summary.py index 501aa928..fb800a62 100644 --- a/tests/test_collect_summary.py +++ b/tests/test_collect_summary.py @@ -2416,11 +2416,19 @@ def main(): "Path_no_call\tAncestral_in_terminal\tSNPs_called\tMin_depth\tTree\tPositions\t" "Outvoted_ancestral\nS\tNA\t0\t0\t0\t0\t.\t12\t10\tYFull YTree 14.06.0\t" "Yleaf 4.1.6 tables\t.\n") + # bin/y_tree_place.py writes the path table's header alone when no branch is derived. + put(f"{h}/y_haplogroup/S_y_tree_path.tsv", "Level\tBranch\tSNPs\tOn_tables\tDerived\tAncestral\t" + "No_call\tAncestral_SNPs\n") summ_na, txt_na, _ = render(h) check("Y tree: Hg NA reads as no derived call, and nothing ancestral is printed", - summ_na["sections"]["y_tree"]["data"].get("haplogroup") == "no branch has a derived call" + summ_na["sections"]["y_tree"]["state"] == "ok" + and summ_na["sections"]["y_tree"]["data"].get("haplogroup") == "no branch has a derived call" + and summ_na["sections"]["y_tree"]["data"].get("path") == [] and "Ancestral inside the terminal branch" not in txt_na, - summ_na["sections"]["y_tree"]["data"]) + summ_na["sections"]["y_tree"]) + put(f"{h}/y_haplogroup/S_y_tree_path.tsv", "Level\tSNPs\n") + check("Y tree: a path table with no Branch column is unreadable", + render(h)[0]["sections"]["y_tree"]["state"] == "unreadable") os.remove(f"{h}/y_haplogroup/S_y_tree.txt") os.remove(f"{h}/y_haplogroup/S_y_tree_path.tsv") From 5c2776d13b36d9b2f2b44938dbc5e4b0c4b85173 Mon Sep 17 00:00:00 2001 From: GeiserX <9169332+GeiserX@users.noreply.github.com> Date: Sat, 7 Nov 2026 17:16:24 +0100 Subject: [PATCH 02/10] fix(report): read no copy from a phased GT with an uncalled allele phase_haps gave {1} for '.|1', so a deletion '.|1' and a hit '0|1' in one phase set were called in cis and dropped as no second hit. A phased GT with a missing allele now gives no phase. Every e2e case that runs bin/cnvpytor_to_vcf.py now lists bin/family_report.py, which it imports. Bead: pgp-5bja.24 --- bin/collect_summary.py | 9 +++++++-- scripts/ci/e2e-case-paths.tsv | 18 +++++++++--------- tests/test_collect_summary.py | 10 ++++++++++ 3 files changed, 26 insertions(+), 11 deletions(-) diff --git a/bin/collect_summary.py b/bin/collect_summary.py index 6dba8cea..6df0cd61 100755 --- a/bin/collect_summary.py +++ b/bin/collect_summary.py @@ -3901,10 +3901,15 @@ def over(chrom, a, b): def phase_haps(gt): """The copies (haplotype indexes) a phased GT puts the alt allele on, or - None when the GT is not phased: '0|1' -> {1}, '1|0' -> {0}, '0/1' -> None.""" + None when the GT is not phased or a copy is not called: '0|1' -> {1}, + '1|0' -> {0}, '0/1' -> None, '.|1' -> None (the missing copy could carry + the other allele, so no cis or trans can be read from it).""" if not gt or "|" not in gt: return None - return frozenset(i for i, a in enumerate(gt.split("|")) if a not in ("0", ".")) + alleles = gt.split("|") + if "." in alleles: + return None + return frozenset(i for i, a in enumerate(alleles) if a != "0") def pair_phase(dl, v): diff --git a/scripts/ci/e2e-case-paths.tsv b/scripts/ci/e2e-case-paths.tsv index 5dc19dc9..7021aaec 100644 --- a/scripts/ci/e2e-case-paths.tsv +++ b/scripts/ci/e2e-case-paths.tsv @@ -48,9 +48,9 @@ alignment-markdup-gvcf-sex-3-homref assets/pgs_scores.tsv bin/collect_summary.py alignment-markdup-gvcf-sex-4-bwamem2 scripts/02a-alignment-bwamem2.sh scripts/ci/fake-docker/lib.sh ancestry-fine-grained-1 assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/mito_confirmed.tsv assets/pgs_scores.tsv assets/stub/sample.vcf.gz bin/ancestry_fine_report.py bin/collect_summary.py bin/render_report.py bin/results_index.py modules/local/ancestry_fine/* modules/local/html_report/* modules/local/prs/* modules/local/vcf_precheck/* scripts/24-html-report.sh scripts/25-prs.sh scripts/26-ancestry.sh scripts/26b-ancestry-fine.sh scripts/chip-to-vcf.sh scripts/ci/fake-docker/lib.sh scripts/setup.sh scripts/validate-setup.sh ancient-similarity-1 assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/mito_confirmed.tsv assets/stub/sample.vcf.gz bin/ancient_report.py bin/collect_summary.py bin/render_report.py bin/results_index.py modules/local/ancient/* modules/local/html_report/* modules/local/vcf_precheck/* scripts/24-html-report.sh scripts/44-ancient.sh scripts/ci/fake-docker/lib.sh scripts/setup.sh scripts/validate-setup.sh -annotsv-delly-depth-1-chain assets/par_grch38.bed assets/sv_rule_regions.tsv bin/annotsv_postfilter.py bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/sv_svlen_fix.py scripts/05-annotsv.sh scripts/15-duphold.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/ci/fake-docker/lib.sh tests/fixtures/annotsv/class45.sh -annotsv-report-section-1-chain assets/par_grch38.bed assets/sv_rule_regions.tsv bin/annotsv_postfilter.py bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/sv_svlen_fix.py scripts/05-annotsv.sh scripts/15-duphold.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/22-survivor-merge.sh scripts/24-html-report.sh scripts/ci/fake-docker/lib.sh scripts/generate-report.sh tests/fixtures/annotsv/class45.sh -annotsv-sv-inputs-1-manta assets/par_grch38.bed assets/sv_rule_regions.tsv bin/annotsv_postfilter.py bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/manta_inv_filter.py bin/sv_svlen_fix.py scripts/04-manta.sh scripts/05-annotsv.sh scripts/15-duphold.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/ci/fake-docker/lib.sh +annotsv-delly-depth-1-chain assets/par_grch38.bed assets/sv_rule_regions.tsv bin/annotsv_postfilter.py bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/family_report.py bin/sv_svlen_fix.py scripts/05-annotsv.sh scripts/15-duphold.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/ci/fake-docker/lib.sh tests/fixtures/annotsv/class45.sh +annotsv-report-section-1-chain assets/par_grch38.bed assets/sv_rule_regions.tsv bin/annotsv_postfilter.py bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/family_report.py bin/sv_svlen_fix.py scripts/05-annotsv.sh scripts/15-duphold.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/22-survivor-merge.sh scripts/24-html-report.sh scripts/ci/fake-docker/lib.sh scripts/generate-report.sh tests/fixtures/annotsv/class45.sh +annotsv-sv-inputs-1-manta assets/par_grch38.bed assets/sv_rule_regions.tsv bin/annotsv_postfilter.py bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/family_report.py bin/manta_inv_filter.py bin/sv_svlen_fix.py scripts/04-manta.sh scripts/05-annotsv.sh scripts/15-duphold.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/ci/fake-docker/lib.sh archaic-ancestry-1 assets/pgs_scores.tsv bin/archaic_report.py bin/collect_summary.py modules/local/archaic/* modules/local/prs/* modules/local/vcf_precheck/* scripts/25-prs.sh scripts/26-ancestry.sh scripts/43-archaic.sh scripts/ci/fake-docker/lib.sh scripts/generate-report.sh scripts/setup.sh scripts/validate-setup.sh bash-step-08-hla - bash-step-09-str assets/str_extra_loci.json assets/stub/db/catalog.json bin/rfc1_motif.py bin/strchive_convert.py scripts/09-expansion-hunter.sh scripts/09b-stranger.sh scripts/ci/fake-docker/lib.sh @@ -100,18 +100,18 @@ hla-second-caller-2-nextflow assets/acmg_sf.tsv assets/callability_genes.tsv ass ibd-segments-1 assets/families.example.tsv assets/par_grch38.bed assets/pedigree.example.tsv assets/stub/pedigree.tsv bin/family_report.py bin/family_segments.py bin/ibd_report.py scripts/38-family.sh scripts/ci/fake-docker/lib.sh scripts/setup.sh live-server-1-weblog assets/live/index.html assets/live/live.css assets/live/login.html assets/viewer/index.html bin/owngenome_live.py bin/render_report.py bin/run_state.py modules/local/roh/* modules/local/vcf_precheck/* scripts/11-roh-analysis.sh scripts/ci/fake-docker/lib.sh scripts/validate-setup.sh lpa-1 assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/mito_confirmed.tsv assets/par_grch38.bed assets/stub/sample.vcf.gz bin/cn_loci.tsv bin/collect_summary.py bin/copy_number.py bin/hba_copy_number.py bin/lpa_kiv2.py bin/render_report.py bin/results_index.py modules/local/deepvariant/* modules/local/html_report/* modules/local/lpa/* modules/local/vcf_precheck/* scripts/24-html-report.sh scripts/46-lpa.sh scripts/ci/fake-docker/lib.sh scripts/validate-setup.sh -mei-1 assets/mei_sites.tsv assets/par_grch38.bed assets/sv_rule_regions.tsv bin/annotsv_postfilter.py bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/mei_report.py bin/sv_svlen_fix.py modules/local/align_minimap2/* modules/local/annotsv/* modules/local/deepvariant/* modules/local/fastp/* modules/local/mei/* modules/local/vcf_precheck/* scripts/01b-fastp-qc.sh scripts/05-annotsv.sh scripts/15-duphold.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/50-mei.sh scripts/ci/check-images.sh scripts/ci/fake-docker/lib.sh scripts/setup.sh scripts/validate-setup.sh tests/fixtures/annotsv/class45.sh +mei-1 assets/mei_sites.tsv assets/par_grch38.bed assets/sv_rule_regions.tsv bin/annotsv_postfilter.py bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/family_report.py bin/mei_report.py bin/sv_svlen_fix.py modules/local/align_minimap2/* modules/local/annotsv/* modules/local/deepvariant/* modules/local/fastp/* modules/local/mei/* modules/local/vcf_precheck/* scripts/01b-fastp-qc.sh scripts/05-annotsv.sh scripts/15-duphold.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/50-mei.sh scripts/ci/check-images.sh scripts/ci/fake-docker/lib.sh scripts/setup.sh scripts/validate-setup.sh tests/fixtures/annotsv/class45.sh mito-extended-1-deletions assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/construct_variants.tsv assets/deep_intronic_sites.tsv assets/mito_confirmed.tsv assets/stub/db/clinvar.vcf.gz assets/stub/db/clinvar.vcf.gz.tbi assets/stub/sample.bam assets/stub/sample.bam.bai assets/stub/sample.vcf.gz bin/cdna_scan.py bin/clinvar_hits.awk bin/clinvar_join_adjacent.awk bin/clip_check.py bin/collect_summary.py bin/render_report.py bin/results_index.py modules/local/clinvar_screen/* modules/local/html_report/* modules/local/mito_variants/* modules/local/vcf_precheck/* scripts/06-clinvar-screen.sh scripts/16b-mosdepth.sh scripts/20-mtoolbox.sh scripts/24-html-report.sh scripts/ci/build-mito-confirmed.sh scripts/ci/fake-docker/lib.sh scripts/validate-setup.sh mocha-1 assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/mito_confirmed.tsv assets/stub/sample.vcf.gz bin/collect_summary.py bin/mocha_report.py bin/render_report.py bin/results_index.py modules/local/html_report/* modules/local/mocha/* modules/local/vcf_precheck/* scripts/24-html-report.sh scripts/51-mocha.sh scripts/ci/fake-docker/lib.sh scripts/setup.sh scripts/validate-setup.sh mosaic-loy-lox-1 assets/callability_genes.tsv bin/collect_summary.py bin/mosaic_sex_loss.py scripts/16b-mosdepth.sh scripts/33-sample-qc.sh scripts/ci/fake-docker/lib.sh scripts/setup.sh nextflow-from-fastq-1-bash assets/callability_genes.tsv assets/construct_variants.tsv assets/deep_intronic_sites.tsv assets/pgs_scores.tsv assets/pgx_genes.tsv assets/stub/db/hla.dat assets/stub/db/kir.dat bin/cdna_scan.py bin/clinvar_hits.awk bin/clinvar_join_adjacent.awk bin/clip_check.py bin/cn_loci.tsv bin/collect_summary.py bin/copy_number.py bin/cyp2d6_depth_check.py bin/hba_copy_number.py bin/lpa_kiv2.py bin/pgx_mito_merge.py bin/pgx_outside_calls.py bin/pgx_refblock_fill.py scripts/01b-fastp-qc.sh scripts/06-clinvar-screen.sh scripts/07-pharmacogenomics.sh scripts/08-hla-typing.sh scripts/11-roh-analysis.sh scripts/16b-mosdepth.sh scripts/20-mtoolbox.sh scripts/25-prs.sh scripts/32-pypgx.sh scripts/36-pgx-consensus.sh scripts/ci/fake-docker/lib.sh scripts/ci/parity-diff.sh scripts/validate-setup.sh nextflow-from-fastq-2-nextflow assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/construct_variants.tsv assets/deep_intronic_sites.tsv assets/live/index.html assets/mito_confirmed.tsv assets/par_grch38.bed assets/pgs_scores.tsv assets/pgx_contraindications.tsv assets/pgx_genes.tsv assets/pgx_inhibitors.tsv assets/stub/db/hla.dat assets/stub/db/kir.dat assets/stub/sample.bam assets/stub/sample.bam.bai assets/stub/sample.vcf.gz assets/viewer/index.html bin/cdna_scan.py bin/clinvar_hits.awk bin/clinvar_join_adjacent.awk bin/clip_check.py bin/cn_loci.tsv bin/collect_summary.py bin/copy_number.py bin/cyp2d6_depth_check.py bin/hba_copy_number.py bin/lpa_kiv2.py bin/pgx_mito_merge.py bin/pgx_outside_calls.py bin/pgx_parse.py bin/pgx_refblock_fill.py bin/render_report.py bin/results_index.py modules/local/align_minimap2/* modules/local/clinvar_screen/* modules/local/cpic_lookup/* modules/local/cyrius/* modules/local/deepvariant/* modules/local/fastp/* modules/local/hla_typing/* modules/local/html_report/* modules/local/indexcov/* modules/local/mito_haplogroup/* modules/local/mito_variants/* modules/local/mosdepth/* modules/local/multiqc/* modules/local/pharmcat/* modules/local/prs/* modules/local/roh/* modules/local/t1k_build/* modules/local/telomere_hunter/* modules/local/vcf_precheck/* modules/local/vcfanno/* scripts/01b-fastp-qc.sh scripts/06-clinvar-screen.sh scripts/07-pharmacogenomics.sh scripts/08-hla-typing.sh scripts/10-telomere-hunter.sh scripts/11-roh-analysis.sh scripts/12-mito-haplogroup.sh scripts/16-indexcov.sh scripts/16b-mosdepth.sh scripts/20-mtoolbox.sh scripts/21-cyrius.sh scripts/24-html-report.sh scripts/25-prs.sh scripts/26-ancestry.sh scripts/27-cpic-lookup.sh scripts/28-multiqc.sh scripts/30-vcfanno.sh scripts/32-pypgx.sh scripts/36-pgx-consensus.sh scripts/ci/check-images.sh scripts/ci/fake-docker/lib.sh scripts/cyrius-constraints.txt scripts/validate-setup.sh nextflow-from-fastq-3-sexcheck assets/callability_genes.tsv assets/live/index.html assets/viewer/index.html bin/collect_summary.py modules/local/indexcov/* modules/local/mosdepth/* modules/local/vcf_precheck/* scripts/16-indexcov.sh scripts/16b-mosdepth.sh scripts/ci/fake-docker/lib.sh scripts/validate-setup.sh -nextflow-hardening-1-manta assets/par_grch38.bed assets/sv_rule_regions.tsv bin/annotsv_postfilter.py bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/manta_inv_filter.py bin/sv_svlen_fix.py scripts/04-manta.sh scripts/05-annotsv.sh scripts/15-duphold.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/ci/fake-docker/lib.sh -nextflow-hardening-2-nf assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/construct_variants.tsv assets/deep_intronic_sites.tsv assets/live/index.html assets/mito_confirmed.tsv assets/par_grch38.bed assets/pgx_genes.tsv assets/stub/sample.bam assets/stub/sample.bam.bai assets/stub/sample.vcf.gz assets/sv_rule_regions.tsv assets/viewer/index.html bin/annotsv_postfilter.py bin/cdna_scan.py bin/clinvar_hits.awk bin/clinvar_join_adjacent.awk bin/clip_check.py bin/cn_loci.tsv bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/copy_number.py bin/cyp2d6_depth_check.py bin/hba_copy_number.py bin/lpa_kiv2.py bin/manta_inv_filter.py bin/pgx_mito_merge.py bin/pgx_outside_calls.py bin/pgx_refblock_fill.py bin/render_report.py bin/results_index.py bin/sv_svlen_fix.py bin/y_tree_place.py bin/yleaf_run.py modules/local/clinvar_screen/* modules/local/delly/* modules/local/html_report/* modules/local/indexcov/* modules/local/manta/* modules/local/mito_haplogroup/* modules/local/mosdepth/* modules/local/pgx_consensus/* modules/local/pharmcat/* modules/local/pypgx/* modules/local/survivor_merge/* modules/local/telomere_hunter/* modules/local/vcf_precheck/* modules/local/vcfanno/* modules/local/y_haplogroup/* scripts/04-manta.sh scripts/05-annotsv.sh scripts/06-clinvar-screen.sh scripts/07-pharmacogenomics.sh scripts/10-telomere-hunter.sh scripts/12-mito-haplogroup.sh scripts/15-duphold.sh scripts/16-indexcov.sh scripts/16b-mosdepth.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/20-mtoolbox.sh scripts/21-cyrius.sh scripts/22-survivor-merge.sh scripts/24-html-report.sh scripts/30-vcfanno.sh scripts/32-pypgx.sh scripts/36-pgx-consensus.sh scripts/37-y-haplogroup.sh scripts/ci/fake-docker/lib.sh scripts/ci/parity-diff.sh scripts/cyrius-constraints.txt scripts/setup.sh scripts/validate-setup.sh -nextflow-report-parity-1-sections assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/construct_variants.tsv assets/deep_intronic_sites.tsv assets/mito_confirmed.tsv assets/par_grch38.bed assets/pgx_genes.tsv assets/str_extra_loci.json assets/stub/db/catalog.json assets/stub/sample.bam assets/stub/sample.bam.bai assets/stub/sample.vcf.gz assets/sv_rule_regions.tsv bin/annotsv_postfilter.py bin/cdna_scan.py bin/clinvar_hits.awk bin/clinvar_join_adjacent.awk bin/clip_check.py bin/cn_loci.tsv bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/comphet_phase.py bin/constraint_join.awk bin/copy_number.py bin/csq_rank.awk bin/cyp2d6_depth_check.py bin/export_consumer.py bin/hba_copy_number.py bin/lpa_kiv2.py bin/pgx_outside_calls.py bin/render_report.py bin/results_index.py bin/rfc1_motif.py bin/strchive_convert.py bin/sv_svlen_fix.py modules/local/annotsv/* modules/local/clinvar_screen/* modules/local/expansion_hunter/* modules/local/export_consumer/* modules/local/html_report/* modules/local/mito_variants/* modules/local/pypgx/* modules/local/slivar/* modules/local/stranger/* modules/local/vcf_precheck/* modules/local/vep/* scripts/05-annotsv.sh scripts/06-clinvar-screen.sh scripts/09-expansion-hunter.sh scripts/09b-stranger.sh scripts/13-vep-annotation.sh scripts/15-duphold.sh scripts/16b-mosdepth.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/20-mtoolbox.sh scripts/24-html-report.sh scripts/31-slivar.sh scripts/32-pypgx.sh scripts/36-pgx-consensus.sh scripts/40-export-consumer.sh scripts/ci/fake-docker/lib.sh scripts/validate-setup.sh tests/fixtures/consumer_export/sites.tsv.gz +nextflow-hardening-1-manta assets/par_grch38.bed assets/sv_rule_regions.tsv bin/annotsv_postfilter.py bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/family_report.py bin/manta_inv_filter.py bin/sv_svlen_fix.py scripts/04-manta.sh scripts/05-annotsv.sh scripts/15-duphold.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/ci/fake-docker/lib.sh +nextflow-hardening-2-nf assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/construct_variants.tsv assets/deep_intronic_sites.tsv assets/live/index.html assets/mito_confirmed.tsv assets/par_grch38.bed assets/pgx_genes.tsv assets/stub/sample.bam assets/stub/sample.bam.bai assets/stub/sample.vcf.gz assets/sv_rule_regions.tsv assets/viewer/index.html bin/annotsv_postfilter.py bin/cdna_scan.py bin/clinvar_hits.awk bin/clinvar_join_adjacent.awk bin/clip_check.py bin/cn_loci.tsv bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/copy_number.py bin/cyp2d6_depth_check.py bin/family_report.py bin/hba_copy_number.py bin/lpa_kiv2.py bin/manta_inv_filter.py bin/pgx_mito_merge.py bin/pgx_outside_calls.py bin/pgx_refblock_fill.py bin/render_report.py bin/results_index.py bin/sv_svlen_fix.py bin/y_tree_place.py bin/yleaf_run.py modules/local/clinvar_screen/* modules/local/delly/* modules/local/html_report/* modules/local/indexcov/* modules/local/manta/* modules/local/mito_haplogroup/* modules/local/mosdepth/* modules/local/pgx_consensus/* modules/local/pharmcat/* modules/local/pypgx/* modules/local/survivor_merge/* modules/local/telomere_hunter/* modules/local/vcf_precheck/* modules/local/vcfanno/* modules/local/y_haplogroup/* scripts/04-manta.sh scripts/05-annotsv.sh scripts/06-clinvar-screen.sh scripts/07-pharmacogenomics.sh scripts/10-telomere-hunter.sh scripts/12-mito-haplogroup.sh scripts/15-duphold.sh scripts/16-indexcov.sh scripts/16b-mosdepth.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/20-mtoolbox.sh scripts/21-cyrius.sh scripts/22-survivor-merge.sh scripts/24-html-report.sh scripts/30-vcfanno.sh scripts/32-pypgx.sh scripts/36-pgx-consensus.sh scripts/37-y-haplogroup.sh scripts/ci/fake-docker/lib.sh scripts/ci/parity-diff.sh scripts/cyrius-constraints.txt scripts/setup.sh scripts/validate-setup.sh +nextflow-report-parity-1-sections assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/construct_variants.tsv assets/deep_intronic_sites.tsv assets/mito_confirmed.tsv assets/par_grch38.bed assets/pgx_genes.tsv assets/str_extra_loci.json assets/stub/db/catalog.json assets/stub/sample.bam assets/stub/sample.bam.bai assets/stub/sample.vcf.gz assets/sv_rule_regions.tsv bin/annotsv_postfilter.py bin/cdna_scan.py bin/clinvar_hits.awk bin/clinvar_join_adjacent.awk bin/clip_check.py bin/cn_loci.tsv bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/comphet_phase.py bin/constraint_join.awk bin/copy_number.py bin/csq_rank.awk bin/cyp2d6_depth_check.py bin/export_consumer.py bin/family_report.py bin/hba_copy_number.py bin/lpa_kiv2.py bin/pgx_outside_calls.py bin/render_report.py bin/results_index.py bin/rfc1_motif.py bin/strchive_convert.py bin/sv_svlen_fix.py modules/local/annotsv/* modules/local/clinvar_screen/* modules/local/expansion_hunter/* modules/local/export_consumer/* modules/local/html_report/* modules/local/mito_variants/* modules/local/pypgx/* modules/local/slivar/* modules/local/stranger/* modules/local/vcf_precheck/* modules/local/vep/* scripts/05-annotsv.sh scripts/06-clinvar-screen.sh scripts/09-expansion-hunter.sh scripts/09b-stranger.sh scripts/13-vep-annotation.sh scripts/15-duphold.sh scripts/16b-mosdepth.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/20-mtoolbox.sh scripts/24-html-report.sh scripts/31-slivar.sh scripts/32-pypgx.sh scripts/36-pgx-consensus.sh scripts/40-export-consumer.sh scripts/ci/fake-docker/lib.sh scripts/validate-setup.sh tests/fixtures/consumer_export/sites.tsv.gz nf-wiring-1-report-tables assets/acmg_sf.tsv assets/test_samplesheet.csv assets/stub/* bin/clip_check.py bin/collect_summary.py bin/cpsr_sample_id bin/render_report.py bin/results_index.py conf/test_all.config conf/test.config modules/local/align_minimap2/* modules/local/clinical_filter/* modules/local/clinvar_screen/* modules/local/cpsr/* modules/local/html_report/* modules/local/slivar/* modules/local/vep/* scripts/24-html-report.sh tests/fixtures/consumer_export/sites.tsv.gz -nf-wiring-2-sv-str-tools bin/cnvpytor_to_vcf.py bin/manta_inv_filter.py bin/strchive_convert.py bin/sv_svlen_fix.py modules/local/cnvpytor/* modules/local/duphold/* modules/local/expansion_hunter/* modules/local/indexcov/* modules/local/manta/* modules/local/vcf_precheck/* +nf-wiring-2-sv-str-tools bin/cnvpytor_to_vcf.py bin/family_report.py bin/manta_inv_filter.py bin/strchive_convert.py bin/sv_svlen_fix.py modules/local/cnvpytor/* modules/local/duphold/* modules/local/expansion_hunter/* modules/local/indexcov/* modules/local/manta/* modules/local/vcf_precheck/* output-contract-1-hg002 assets/* bin/* modules/* scripts/* paralog-catalog-1-gba assets/paralog_loci.tsv assets/parascopy_models_v1.7.txt bin/cn_loci.tsv bin/copy_number.py bin/hba_copy_number.py bin/lpa_kiv2.py modules/local/parascopy/* modules/local/vcf_precheck/* scripts/35-paralogs.sh scripts/ci/fake-docker/lib.sh scripts/setup.sh scripts/validate-setup.sh paralog-checks-and-doc-fixes-1-hba assets/par_grch38.bed assets/paralog_loci.tsv bin/cn_loci.tsv bin/collect_summary.py bin/copy_number.py bin/hba_copy_number.py bin/lpa_kiv2.py bin/paraphase_report.py modules/local/align_minimap2/* modules/local/deepvariant/* modules/local/fastp/* modules/local/paraphase/* modules/local/parascopy/* modules/local/vcf_precheck/* scripts/01b-fastp-qc.sh scripts/35-paralogs.sh scripts/35b-paraphase.sh scripts/ci/check-images.sh scripts/ci/fake-docker/lib.sh scripts/setup.sh scripts/validate-setup.sh @@ -139,7 +139,7 @@ qc-1-sample-qc bin/collect_summary.py bin/mosaic_sex_loss.py scripts/28-multiqc. qc-2-nextflow-sample-qc assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/families.example.tsv assets/live/index.html assets/mito_confirmed.tsv assets/par_grch38.bed assets/pedigree.example.tsv assets/pgx_contraindications.tsv assets/pgx_genes.tsv assets/pgx_inhibitors.tsv assets/stub/db/hla.dat assets/stub/db/kir.dat assets/stub/pedigree.tsv assets/stub/sample.vcf.gz assets/viewer/index.html bin/cn_loci.tsv bin/collect_summary.py bin/copy_number.py bin/cyp2d6_depth_check.py bin/family_report.py bin/family_segments.py bin/hba_copy_number.py bin/ibd_report.py bin/lpa_kiv2.py bin/mosaic_sex_loss.py bin/pgx_mito_merge.py bin/pgx_outside_calls.py bin/pgx_parse.py bin/pgx_refblock_fill.py bin/render_report.py bin/results_index.py modules/local/align_minimap2/* modules/local/cpic_lookup/* modules/local/deepvariant/* modules/local/fastp/* modules/local/hla_typing/* modules/local/html_report/* modules/local/indexcov/* modules/local/mito_haplogroup/* modules/local/mito_variants/* modules/local/mosdepth/* modules/local/pharmcat/* modules/local/pypgx/* modules/local/roh/* modules/local/somalier/* modules/local/t1k_build/* modules/local/telomere_hunter/* modules/local/vcf_precheck/* modules/local/verifybamid2/* scripts/01b-fastp-qc.sh scripts/07-pharmacogenomics.sh scripts/08-hla-typing.sh scripts/10-telomere-hunter.sh scripts/11-roh-analysis.sh scripts/12-mito-haplogroup.sh scripts/16-indexcov.sh scripts/16b-mosdepth.sh scripts/20-mtoolbox.sh scripts/24-html-report.sh scripts/27-cpic-lookup.sh scripts/32-pypgx.sh scripts/33-sample-qc.sh scripts/36-pgx-consensus.sh scripts/38-family.sh scripts/ci/check-images.sh scripts/ci/fake-docker/lib.sh scripts/validate-setup.sh qc-3-cram assets/callability_genes.tsv bin/collect_summary.py modules/local/cram_archive/* modules/local/mosdepth/* modules/local/vcf_precheck/* scripts/16b-mosdepth.sh scripts/34-cram-archive.sh scripts/ci/fake-docker/lib.sh scripts/validate-setup.sh reference-noalt scripts/ci/fake-docker/lib.sh scripts/validate-setup.sh -report-followups-1-sections assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/genotype_rules.tsv assets/mito_confirmed.tsv assets/par_grch38.bed assets/str_extra_loci.json assets/stub/db/catalog.json assets/stub/sample.vcf.gz assets/sv_rule_regions.tsv assets/traits.tsv bin/annotsv_postfilter.py bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/export_consumer.py bin/manta_inv_filter.py bin/render_report.py bin/results_index.py bin/rfc1_motif.py bin/strchive_convert.py bin/sv_svlen_fix.py bin/traits_report.py modules/local/annotsv/* modules/local/expansion_hunter/* modules/local/export_consumer/* modules/local/html_report/* modules/local/manta/* modules/local/traits/* modules/local/vcf_precheck/* scripts/04-manta.sh scripts/05-annotsv.sh scripts/09-expansion-hunter.sh scripts/15-duphold.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/24-html-report.sh scripts/40-export-consumer.sh scripts/42-traits.sh scripts/ci/fake-docker/lib.sh scripts/generate-report.sh scripts/validate-setup.sh tests/fixtures/consumer_export/sites.tsv.gz +report-followups-1-sections assets/acmg_sf.tsv assets/callability_genes.tsv assets/carrier_genes.tsv assets/genotype_rules.tsv assets/mito_confirmed.tsv assets/par_grch38.bed assets/str_extra_loci.json assets/stub/db/catalog.json assets/stub/sample.vcf.gz assets/sv_rule_regions.tsv assets/traits.tsv bin/annotsv_postfilter.py bin/cnvpytor_to_vcf.py bin/collect_summary.py bin/export_consumer.py bin/family_report.py bin/manta_inv_filter.py bin/render_report.py bin/results_index.py bin/rfc1_motif.py bin/strchive_convert.py bin/sv_svlen_fix.py bin/traits_report.py modules/local/annotsv/* modules/local/expansion_hunter/* modules/local/export_consumer/* modules/local/html_report/* modules/local/manta/* modules/local/traits/* modules/local/vcf_precheck/* scripts/04-manta.sh scripts/05-annotsv.sh scripts/09-expansion-hunter.sh scripts/15-duphold.sh scripts/18-cnvpytor.sh scripts/19-delly.sh scripts/24-html-report.sh scripts/40-export-consumer.sh scripts/42-traits.sh scripts/ci/fake-docker/lib.sh scripts/generate-report.sh scripts/validate-setup.sh tests/fixtures/consumer_export/sites.tsv.gz report-sections-mito-roh-telomere-cnv-1-mito assets/callability_genes.tsv bin/collect_summary.py scripts/10-telomere-hunter.sh scripts/16b-mosdepth.sh scripts/20-mtoolbox.sh scripts/ci/fake-docker/lib.sh scripts/generate-report.sh scripts/validate-setup.sh reports-and-pgx-correctness-clinical assets/construct_variants.tsv bin/clip_check.py bin/comphet_phase.py bin/constraint_join.awk bin/csq_rank.awk scripts/23-clinical-filter.sh scripts/31-slivar.sh scripts/ci/fake-docker/lib.sh reports-and-pgx-correctness-pgx assets/pgx_contraindications.tsv assets/pgx_genes.tsv assets/pgx_inhibitors.tsv bin/cn_loci.tsv bin/collect_summary.py bin/copy_number.py bin/cyp2d6_depth_check.py bin/hba_copy_number.py bin/lpa_kiv2.py bin/pgx_outside_calls.py bin/pgx_parse.py scripts/27-cpic-lookup.sh scripts/32-pypgx.sh scripts/36-pgx-consensus.sh scripts/ci/fake-docker/lib.sh scripts/validate-setup.sh diff --git a/tests/test_collect_summary.py b/tests/test_collect_summary.py index fb800a62..4079a649 100644 --- a/tests/test_collect_summary.py +++ b/tests/test_collect_summary.py @@ -1215,6 +1215,16 @@ def run(name, gt): _, t = run("homref", "0/0") check("sv phase: a deletion genotyped 0/0 pairs with nothing", t == [], t) + # A phased GT with one copy not called says nothing about cis or trans. + check("sv phase: '.|1' gives no copy, '0|1' gives copy 1", + (collect_summary.phase_haps(".|1"), collect_summary.phase_haps("0|1")) == (None, frozenset({1})), + (collect_summary.phase_haps(".|1"), collect_summary.phase_haps("0|1"))) + dl = {"inside": False, "chrom": "chr1", "start": 100, "end": 200, "gt": ".|1", "ps": "90"} + v = {"chrom": "chr1", "pos": 300, "gt": "0|1", "ps": "90", "copies": 1, "ploidy": 2} + verdict, why = collect_summary.pair_phase(dl, v) + check("sv phase: a deletion '.|1' and a hit '0|1' in one phase set are not called in cis", + verdict == "possibly biallelic" and "in cis" not in why, (verdict, why)) + def report_followups_sections(work, schema): """26. step 42's genotype rules, step 9b's motif table of the other From 9f19f150e8c333f523466ae0648c46bf73aaceb9 Mon Sep 17 00:00:00 2001 From: GeiserX <9169332+GeiserX@users.noreply.github.com> Date: Sat, 7 Nov 2026 17:24:26 +0100 Subject: [PATCH 03/10] fix(live): pin the nf-weblog plugin in versions.env with a Renovate rule The plugin pin sat in run-all.sh, the one image or plugin version kept outside versions.env, and no Renovate rule saw it. It is now NF_WEBLOG_PLUGIN in versions.env; a custom manager bumps it there and in the example config of docs/live-site.md in one PR; the e2e case reads it from versions.env. Bead: pgp-5bja.23 --- docs/live-site.md | 2 +- renovate.json | 14 ++++++++++++++ scripts/run-all.sh | 3 +-- tests/e2e/live-server-1-weblog.sh | 6 +++--- versions.env | 4 ++++ 5 files changed, 23 insertions(+), 6 deletions(-) diff --git a/docs/live-site.md b/docs/live-site.md index 26018b78..fe03b28e 100644 --- a/docs/live-site.md +++ b/docs/live-site.md @@ -29,7 +29,7 @@ OWNGENOME_LIVE_URL=http://127.0.0.1:8788 ./scripts/run-all.sh your_name male # `./scripts/owngenome-live.sh your_name` does the same and prints the address of that run's page. `hash-password` writes `GENOME_DIR/.owngenome-live/password` with mode 600. The file holds a PBKDF2-SHA256 hash, never the password. Open http://127.0.0.1:8787, log in, and pick the run. -With `OWNGENOME_LIVE_URL` set, `run-all.sh` writes `/nextflow/live.config` and passes it to Nextflow with `-c`. That config loads the nf-weblog plugin (pinned in `run-all.sh`; Nextflow downloads it on the first live run) and points it at `OWNGENOME_LIVE_URL/ingest/`, so every task's progress is posted. The numbered script steps that run after the pipeline (GRIDSS, imputation prep, the somatic caller, the extra callers) post when they start and how they end (the benchmark, `benchmark-variants.sh`, has no step number and posts nothing), and `run-all.sh` posts once more as it ends, with its exit code. A site that is down or slow costs one warning; the run never waits for it and never fails for it. `OWNGENOME_LIVE_URL` must be a plain `http://host:port` URL. +With `OWNGENOME_LIVE_URL` set, `run-all.sh` writes `/nextflow/live.config` and passes it to Nextflow with `-c`. That config loads the nf-weblog plugin (pinned as `NF_WEBLOG_PLUGIN` in `versions.env`; Nextflow downloads it on the first live run) and points it at `OWNGENOME_LIVE_URL/ingest/`, so every task's progress is posted. The numbered script steps that run after the pipeline (GRIDSS, imputation prep, the somatic caller, the extra callers) post when they start and how they end (the benchmark, `benchmark-variants.sh`, has no step number and posts nothing), and `run-all.sh` posts once more as it ends, with its exit code. A site that is down or slow costs one warning; the run never waits for it and never fails for it. `OWNGENOME_LIVE_URL` must be a plain `http://host:port` URL. A bare `nextflow run main.nf` can report too. Write a config like this one, with the sample id at the end of the URL (the page counts only the tasks tagged with it), and pass it with `-c`: diff --git a/renovate.json b/renovate.json index 4cf6b687..a91850c8 100644 --- a/renovate.json +++ b/renovate.json @@ -63,6 +63,20 @@ "extractVersionTemplate": "^v(?.+)$", "versioningTemplate": "regex:^(?\\d+)\\.(?\\d+)\\.(?\\d+)$" }, + { + "customType": "regex", + "description": "nf-weblog, the Nextflow plugin the live site's config loads: NF_WEBLOG_PLUGIN in versions.env and the example config in docs/live-site.md, in one PR. The e2e case live-server-1-weblog reads the pin from versions.env, so it holds no copy to match.", + "managerFilePatterns": [ + "/^versions\\.env$/", + "/^docs/live-site\\.md$/" + ], + "matchStrings": [ + "NF_WEBLOG_PLUGIN=\"nf-weblog@(?\\d+\\.\\d+\\.\\d+)\"", + "id 'nf-weblog@(?\\d+\\.\\d+\\.\\d+)'" + ], + "depNameTemplate": "nextflow-io/nf-weblog", + "datasourceTemplate": "github-releases" + }, { "customType": "regex", "description": "actionlint, downloaded by version in the lint workflow.", diff --git a/scripts/run-all.sh b/scripts/run-all.sh index ec353a40..5ec9d955 100755 --- a/scripts/run-all.sh +++ b/scripts/run-all.sh @@ -296,8 +296,7 @@ else echo " benchmark-variants skipped (only one caller VCF: set EXTRA_CALLERS, export NXF_VER=${NXF_VER:-$NEXTFLOW_VERSION} # The live site (docs/live-site.md). Nextflow's weblog is the nf-weblog plugin, which Nextflow 26.10 loads only when # a config names it, so the block lives in a config of this run's own, passed only when OWNGENOME_LIVE_URL is set: a run -# without it neither downloads the plugin nor warns about an unknown weblog scope. Bump the pin with its release. -NF_WEBLOG_PLUGIN=nf-weblog@1.2.0 +# without it neither downloads the plugin nor warns about an unknown weblog scope. The pin is NF_WEBLOG_PLUGIN in versions.env. if [ -n "${OWNGENOME_LIVE_URL:-}" ]; then # The ingest listener speaks plain HTTP and serves /ingest/ at its root, so https or a path would lose every post. [[ "$OWNGENOME_LIVE_URL" =~ ^http://[A-Za-z0-9.:_-]+/?$ ]] \ diff --git a/tests/e2e/live-server-1-weblog.sh b/tests/e2e/live-server-1-weblog.sh index c98b4e40..75035997 100644 --- a/tests/e2e/live-server-1-weblog.sh +++ b/tests/e2e/live-server-1-weblog.sh @@ -1,15 +1,15 @@ #!/usr/bin/env bash # The live site (bin/owngenome_live.py) on a real Nextflow run: the pinned Nextflow # loads the nf-weblog plugin from the config run-all.sh writes when -# OWNGENOME_LIVE_URL is set (the same text, the pin read from run-all.sh), posts each +# OWNGENOME_LIVE_URL is set (the same text, the pin read from versions.env), posts each # task to the site's ingest port, and the site stores the messages and shows # step 11 (ROH) done in state.json. A control reads state.json before the run: # step 11 is 'not planned' then, so the state check can fail. . "$(dirname "$0")/lib.sh" command -v nextflow >/dev/null || { fail "nextflow is not on PATH"; finish; } -PLUGIN=$(sed -n 's/^NF_WEBLOG_PLUGIN=\([^ ]*\)$/\1/p' "${REPO}/scripts/run-all.sh") -check "run-all.sh pins nf-weblog with a version" has '^nf-weblog@[0-9]+\.[0-9]+\.[0-9]+$' "$PLUGIN" +PLUGIN=$(sed -n 's/^NF_WEBLOG_PLUGIN="\([^"]*\)"$/\1/p' "${REPO}/versions.env") +check "versions.env pins nf-weblog with a version" has '^nf-weblog@[0-9]+\.[0-9]+\.[0-9]+$' "$PLUGIN" G="$GENOME_DIR" ROOT="${CASE_TMP}/root" STATE="${CASE_TMP}/state" diff --git a/versions.env b/versions.env index fb7d3876..b05b4bff 100644 --- a/versions.env +++ b/versions.env @@ -224,6 +224,10 @@ GRIDSS_BLACKLIST_COMMIT="bd7daff8fc32f349752ff713e1fd776361a57d7e" # times larger, has been read and the two steps adapted. # Nextflow release the pipeline is validated with. NEXTFLOW_VERSION="26.10.5" +# The nf-weblog plugin (github.com/nextflow-io/nf-weblog) the live site's +# config loads: run-all.sh writes it into /nextflow/live.config when +# OWNGENOME_LIVE_URL is set, and Nextflow downloads it on the first live run. +NF_WEBLOG_PLUGIN="nf-weblog@1.2.0" # Release of pgsc_calc (github.com/PGScatalog/pgsc_calc) steps 25 and 26 run; # the PGSC_*_IMAGE lines above belong to it. setup.sh pulls it with Nextflow. PGSC_CALC_VERSION="v2.3.0" From c7a136a435a9d03abe181f075089cbdc851728ae Mon Sep 17 00:00:00 2001 From: GeiserX <9169332+GeiserX@users.noreply.github.com> Date: Sat, 7 Nov 2026 17:24:54 +0100 Subject: [PATCH 04/10] fix(report): read the CFTR and UGT1A1 repeats without overclaiming - The setup hint is now the command that runs: STR_EXTRA_LOCI=yes ./scripts/setup.sh --strchive-data . - No 5T copy no longer says the tract changes nothing: the poly-T still sets how R117H reads, so the line reads R117H and names its copy. - A TA length with no named allele gets no function claim and no PharmCAT comparison. - The readings use the traits and CPIC section states: a stale table is read and said to be so, a failed or unreadable one is not read. Bead: pgp-5bja.28 --- bin/collect_summary.py | 103 ++++++++++++++++++++++++---------- docs/09-str-expansions.md | 6 +- tests/test_collect_summary.py | 47 ++++++++++++++-- 3 files changed, 118 insertions(+), 38 deletions(-) diff --git a/bin/collect_summary.py b/bin/collect_summary.py index 6df0cd61..953f7e0d 100755 --- a/bin/collect_summary.py +++ b/bin/collect_summary.py @@ -2099,6 +2099,10 @@ def str_not_graded(d, s): # The repeats STRchive leaves out (assets/str_extra_loci.json), read from # their counts when step 09 sized them (STR_EXTRA_LOCI=yes; docs/09-str-expansions.md). +STR_EXTRA_SETUP = "STR_EXTRA_LOCI=yes ./scripts/setup.sh --strchive-data " +# What a section state other than ok says about a table another step wrote. +OTHER_STEP_STATE = {"stale": "from an earlier run", "failed": "its step failed in the last run", + "unreadable": "unreadable"} # CFTR intron 9 (legacy IVS8): ExpansionHunter writes the TG and the poly-T # tract as two records, so which TG count sits on the copy of which T count # is known only when one of the two is the same on both copies. R117H @@ -2139,13 +2143,25 @@ def trait_genotype(d, s, rsid): return (row.get("status") or "", row.get("genotype") or "") if row else None -def cftr_reading(by, d, s): +def r117h_genotype(d, s, state): + """(genotype or None, why not read, note) of R117H in step 42's traits + table, given the traits section's state: a table of a failed or + unreadable step is not read, one from an earlier run is read and said so.""" + if state in ("failed", "unreadable"): + return None, f"step 42's traits table is {OTHER_STEP_STATE[state]}", "" + r = trait_genotype(d, s, CFTR_R117H) + if not r or r[0] != "called": + return None, "no step 42 traits row for rs78655421 was called", "" + return r[1], "", (" (R117H from step 42's traits table of an earlier run)" if state == "stale" else "") + + +def cftr_reading(by, d, s, traits_state="ok"): tg, t = by.get("CFTR_TG"), by.get("CFTR_T") out = {"locus": "CFTR", "repeat": "intron 9 (legacy IVS8) TG and poly-T tract", "source": CFTR_SOURCE} if not tg or not t: return dict(out, status="not_sized", genotype=".", phase=".", - reading="not in step 09's output: STR_EXTRA_LOCI=yes setup.sh --strchive-data adds the tract " - "to the catalog (docs/09-str-expansions.md).") + reading=f"not in step 09's output: {STR_EXTRA_SETUP} adds the tract to the catalog " + "(docs/09-str-expansions.md).") a, b = str_counts(tg), str_counts(t) if not a or not b: return dict(out, status="not_sized", genotype=".", phase=".", @@ -2161,10 +2177,21 @@ def cftr_reading(by, d, s): other = sorted({n for n in b if n not in (5, 7, 9)}) odd = (f" A poly-T length of {', '.join(f'{n}T' for n in other)} is none of the usual 5T, 7T and 9T." if other else "") + r117, why, note = r117h_genotype(d, s, traits_state) if not five: - return dict(out, five_t=False, phase="not_needed", - reading=f"No 5T copy ({b[0]}T/{b[1]}T), so the tract does not change how a CFTR finding is " - f"read.{odd}") + # No 5T risk of its own, but the poly-T still sets how R117H reads: + # with no 5T copy, an R117H copy sits on a 7T or 9T (or other) copy. + line = (f"No 5T copy ({b[0]}T/{b[1]}T): the tract adds no 5T allele of its own. The poly-T length still " + "matters beside R117H, which is milder and of variable penetrance on a 7T or 9T copy and CF-causing " + "with 5T in cis.") + if r117 is None: + line += f" R117H was not read ({why})." + elif r117.split("/").count("A"): + copies = f"{b[0]}T" if b[0] == b[1] else f"{b[0]}T or {b[1]}T" + line += f" R117H ({r117}) is read: with no 5T copy it sits on a {copies} copy.{note}" + else: + line += f" No R117H copy.{note}" + return dict(out, five_t=False, phase="not_needed", r117h=r117 or ".", reading=line + odd) if len(pairs) == 1: # One TG count on both copies, or 5T on both copies with their two TG counts. tgs = [a[0]] if a[0] == a[1] else a @@ -2174,36 +2201,34 @@ def cftr_reading(by, d, s): line = (f"5T: {five} copy{'s' if five > 1 else ''} with the 5T poly-T tract, {tg_text}. 5T lowers the splicing " "of exon 10 and is a variable-penetrance allele: with TG12 or TG13 on the same copy it more often " "causes CBAVD or another CFTR-related disorder when the other copy carries a CF-causing variant.") - r = trait_genotype(d, s, CFTR_R117H) - if not r or r[0] != "called": + if r117 is None: phase = "not_read" - r117 = (" R117H was not read (no step 42 traits row for rs78655421 was called), so whether it sits on the " - "5T copy is not known.") + text = f" R117H was not read ({why}), so whether it sits on the 5T copy is not known." else: - n = r[1].split("/").count("A") + n = r117.split("/").count("A") if n == 0: phase = "not_needed" - r117 = " No R117H copy: the 5T is read on its own." + text = " No R117H copy: the 5T is read on its own." elif n == 2 or five == 2: phase = "in_cis" - r117 = (f" R117H ({r[1]}) and 5T share a copy: R117H with 5T in cis is a CF-causing allele of " + text = (f" R117H ({r117}) and 5T share a copy: R117H with 5T in cis is a CF-causing allele of " "variable penetrance.") else: phase = "unknown" - r117 = (" One R117H copy and one 5T copy: phase unknown. R117H lies 17.6 kb from the tract, farther " + text = (" One R117H copy and one 5T copy: phase unknown. R117H lies 17.6 kb from the tract, farther " "than a read pair spans, so whether they share a copy (in cis, a CF-causing allele of variable " "penetrance) or not (in trans) is not read; testing a parent tells them apart.") - return dict(out, five_t=True, phase=phase, r117h=r[1] if r and r[0] == "called" else ".", - reading=line + r117 + odd) + text += note + return dict(out, five_t=True, phase=phase, r117h=r117 or ".", reading=line + text + odd) -def ugt1a1_reading(by, d, s): +def ugt1a1_reading(by, d, s, cpic_state="ok"): item = by.get("UGT1A1") out = {"locus": "UGT1A1", "repeat": "promoter TA repeat", "source": UGT1A1_SOURCE} if not item: return dict(out, status="not_sized", genotype=".", pharmcat=".", agreement="not_compared", - reading="not in step 09's output: STR_EXTRA_LOCI=yes setup.sh --strchive-data adds the repeat " - "to the catalog (docs/09-str-expansions.md). PharmCAT (step 7) still reads UGT1A1 *28.") + reading=f"not in step 09's output: {STR_EXTRA_SETUP} adds the repeat to the catalog " + "(docs/09-str-expansions.md). PharmCAT (step 7) still reads UGT1A1 *28.") c = str_counts(item) if not c: return dict(out, status="not_sized", genotype=".", pharmcat=".", agreement="not_compared", @@ -2213,15 +2238,20 @@ def ugt1a1_reading(by, d, s): genotype = "/".join(names) dec = sum(n in UGT1A1_DECREASED for n in names) ta = "/".join(f"TA{n - 1}" for n in c) - reading = (f"{genotype} ({ta}, {c[0]}/{c[1]} TA units as ExpansionHunter counts them): " - + ("no decreased-function TA allele." if not dec else - f"{dec} decreased-function allele{'s' if dec > 1 else ''} (*28 and *37 are decreased function).") - + ("" if all(n in UGT1A1_ALLELES for n in c) else " A TA length outside TA5 to TA8 is not a named allele here.")) out.update(status="read", genotype=genotype, ta_units=c) + reading = f"{genotype} ({ta}, {c[0]}/{c[1]} TA units as ExpansionHunter counts them): " + if not all(n in UGT1A1_ALLELES for n in c): + # A length with no named allele has no function to read, and so no + # count of decreased-function copies to set beside PharmCAT's. + return dict(out, pharmcat=".", agreement="not_compared", + reading=reading + "a TA length outside TA5 to TA8 is not a named allele here, so its function " + "is not read and it is not compared with PharmCAT (step 7).") + reading += ("no decreased-function TA allele." if not dec else + f"{dec} decreased-function allele{'s' if dec > 1 else ''} (*28 and *37 are decreased function).") # PharmCAT's UGT1A1 call (step 27's table, with the reference-block fill # of bin/pgx_refblock_fill.py): its count of *28 and *37 copies. pc, agreement = ".", "not_compared" - p = first_existing(d, [f"cpic/{s}_phenotypes.tsv"]) + p = first_existing(d, [f"cpic/{s}_phenotypes.tsv"]) if cpic_state not in ("failed", "unreadable") else None row = None if p: try: @@ -2235,20 +2265,26 @@ def ugt1a1_reading(by, d, s): # A combination allele reads "*80+*28" or "[*80 + *28]". n = sum(any(x in re.split(r"[+ \[\]]", a) for x in UGT1A1_DECREASED) for a in alleles) agreement = "agrees" if n == dec else "differs" + old = " (step 27's table is from an earlier run)" if cpic_state == "stale" and agreement != "not_compared" else "" if agreement == "agrees": - reading += f" PharmCAT (step 7) called {pc}: the same number of *28 and *37 copies." + reading += f" PharmCAT (step 7) called {pc}: the same number of *28 and *37 copies{old}." elif agreement == "differs": - reading += (f" PharmCAT (step 7) called {pc}, a different number of *28 and *37 copies: read the reads at " - "chr2:233760234 before using either.") + reading += (f" PharmCAT (step 7) called {pc}, a different number of *28 and *37 copies{old}: read the reads " + "at chr2:233760234 before using either.") + elif cpic_state in ("failed", "unreadable"): + reading += f" Not compared with PharmCAT (step 27's table is {OTHER_STEP_STATE[cpic_state]})." else: reading += " Not compared with PharmCAT (no unambiguous UGT1A1 call from step 27)." return dict(out, pharmcat=pc, agreement=agreement, reading=reading) -def str_readings(by, d, s): +def str_readings(by, d, s, states=None): """What the counts of the extra repeats mean, one entry per locus, and the - 5-HTTLPR line that no step can test.""" - return [cftr_reading(by, d, s), ugt1a1_reading(by, d, s), + 5-HTTLPR line that no step can test. states: the traits and cpic section + states, once collect() knows them (a stale, failed or unreadable table is + said so, not read as current).""" + states = states or {} + return [cftr_reading(by, d, s, states.get("traits", "ok")), ugt1a1_reading(by, d, s, states.get("cpic", "ok")), {"locus": "SLC6A4", "repeat": "5-HTTLPR promoter length polymorphism", "status": "not_tested", "genotype": ".", "reading": HTTLPR_READING, "source": "Bousman et al. 2023, PMID 37032427"}] @@ -4795,6 +4831,13 @@ def collect(sample, sample_dir, declared_sex=None): sections[key] = sec cyp2d6 = cyp2d6_block(sample, sample_dir, sections) + # The repeat readings lean on the traits (R117H) and CPIC (UGT1A1) + # tables, whose sections come after the expansions one: read them again + # with those sections' states. + ex = sections["expansions"] + if ex["state"] in ("ok", "stale"): + ex["data"]["readings"] = str_readings({l["locus"]: l for l in ex["data"]["loci"]}, sample_dir, sample, + {k: sections[k]["state"] for k in ("traits", "cpic")}) if sections["locityper"]["state"] in ("ok", "stale"): primary = locityper_primary(sections) for x in sections["locityper"]["data"]["loci"]: diff --git a/docs/09-str-expansions.md b/docs/09-str-expansions.md index 4e2cce94..2261c0a1 100644 --- a/docs/09-str-expansions.md +++ b/docs/09-str-expansions.md @@ -49,9 +49,9 @@ Each entry has a `NoThresholds` reason and no Stranger thresholds: step 9b gives The reports read the counts (`bin/collect_summary.py`, under the repeat expansions section, `readings` in `summary.json`): -- **CFTR**: the TG and poly-T counts of both copies, such as `TG11-7T/TG11-7T`. ExpansionHunter writes the two tracts as two records, so when both the TG and the poly-T counts differ between the copies, which TG sits with which poly-T is not known: the line says `TG11/TG12 5T/7T` and names both pairings. A 5T copy gets the 5T line: 5T is a variable-penetrance allele, more often disease-causing with TG12 or TG13 on the same copy, and it matters when the other copy carries a CF-causing variant. The line then reads R117H (rs78655421), a row of step 42's traits table: no R117H copy; R117H and 5T in cis when one of them is on both copies; otherwise phase unknown, because R117H lies 17.6 kb from the tract, farther than a read pair spans. No 5T copy gives no 5T line. -- **UGT1A1**: the TA count of each copy as its allele (`*36`, `*1`, `*28`, `*37`), and the number of decreased-function copies (`*28`, `*37`) compared with PharmCAT's UGT1A1 call of step 27, which reads `*28` from the small-variant calls with the reference-block fill of `bin/pgx_refblock_fill.py`. The line says whether the two agree; when they differ, read the reads before using either. PharmCAT's call is still the one step 27's drug guidance uses. -- Without the opt-in, both lines say the repeat was not sized and how to add it. +- **CFTR**: the TG and poly-T counts of both copies, such as `TG11-7T/TG11-7T`. ExpansionHunter writes the two tracts as two records, so when both the TG and the poly-T counts differ between the copies, which TG sits with which poly-T is not known: the line says `TG11/TG12 5T/7T` and names both pairings. A 5T copy gets the 5T line: 5T is a variable-penetrance allele, more often disease-causing with TG12 or TG13 on the same copy, and it matters when the other copy carries a CF-causing variant. The line then reads R117H (rs78655421), a row of step 42's traits table: no R117H copy; R117H and 5T in cis when one of them is on both copies; otherwise phase unknown, because R117H lies 17.6 kb from the tract, farther than a read pair spans. No 5T copy gives no 5T line, but the poly-T still matters beside R117H (milder on a 7T or 9T copy, CF-causing with 5T in cis), so the line still reads R117H and says which poly-T copy it sits on. A traits table from an earlier run is read and said to be so; one of a failed or unreadable step 42 is not read. +- **UGT1A1**: the TA count of each copy as its allele (`*36`, `*1`, `*28`, `*37`), and the number of decreased-function copies (`*28`, `*37`) compared with PharmCAT's UGT1A1 call of step 27, which reads `*28` from the small-variant calls with the reference-block fill of `bin/pgx_refblock_fill.py`. The line says whether the two agree; when they differ, read the reads before using either. A TA length outside TA5 to TA8 has no named allele, so the line makes no function claim and does not compare it. A step 27 table from an earlier run is compared and said to be so; one of a failed or unreadable step is not compared. PharmCAT's call is still the one step 27's drug guidance uses. +- Without the opt-in, both lines say the repeat was not sized and give the command that adds it: `STR_EXTRA_LOCI=yes ./scripts/setup.sh --strchive-data `. `tests/test_str_catalog.py` checks each region against the GRCh38 sequence, `tests/test_collect_summary.py` checks the readings on made-up counts (a TG11-7T/TG11-7T sample is the red control of the 5T line), the e2e case `genotype-rule-rows-1-str-extra` sizes both tracts on made-up reads with ExpansionHunter 5.0.0, and the e2e case `str-extra-repeats-1-readings` runs step 09 and both reports on made-up reads and checks the lines. diff --git a/tests/test_collect_summary.py b/tests/test_collect_summary.py index 4079a649..5f0f6e4a 100644 --- a/tests/test_collect_summary.py +++ b/tests/test_collect_summary.py @@ -1328,7 +1328,8 @@ def str_extra_readings(work, schema): import html as htmlmod fmt = "GT:REPCN" - def made_up(name, tg, t, ugt, r117h=None, ugt_pc=None, stranger=False): + def made_up(name, tg, t, ugt, r117h=None, ugt_pc=None, stranger=False, traits_text=None, side_age=0, + run_status=None): x = os.path.join(work, "str-extra", name) eh = (f"chr7\t117548606\t.\tT\t\t.\tPASS\tEND=117548628;REPID=CFTR;VARID=CFTR_TG\t{fmt}\t0/1:{tg}\n" f"chr7\t117548628\t.\tT\t\t.\tPASS\tEND=117548635;REPID=CFTR;VARID=CFTR_T\t{fmt}\t0/1:{t}\n" @@ -1342,12 +1343,15 @@ def made_up(name, tg, t, ugt, r117h=None, ugt_pc=None, stranger=False): put(f"{x}/expansion_hunter/{name}_eh_catalog.json", json.dumps([ {"LocusId": "CFTR", "NoThresholds": "the CFTR tract"}, {"LocusId": "UGT1A1", "NoThresholds": "the UGT1A1 repeat"}])) - if r117h: - put(f"{x}/traits/{name}_traits.tsv", "\t".join(collect_summary.traits_report.OUT_COLUMNS) + "\n" - + f"rs78655421\tCFTR\tchr7\t117530975\tCFTR R117H\t{r117h}\tcalled\tr\tno\tn\t32404922\n") + if r117h or traits_text: + put(f"{x}/traits/{name}_traits.tsv", traits_text or "\t".join(collect_summary.traits_report.OUT_COLUMNS) + + "\n" + f"rs78655421\tCFTR\tchr7\t117530975\tCFTR R117H\t{r117h}\tcalled\tr\tno\tn\t32404922\n", + age_days=side_age) if ugt_pc: put(f"{x}/cpic/{name}_phenotypes.tsv", "Gene\tDiplotype\tPhenotype\tStatus\n" - f"UGT1A1\t{ugt_pc}\tsome function\tnon-normal\n") + f"UGT1A1\t{ugt_pc}\tsome function\tnon-normal\n", age_days=side_age) + if run_status: + put(f"{x}/logs/run_status.tsv", run_status) summ, txt, html_ = render(x, name) ed = summ["sections"]["expansions"]["data"] return summ, ed, {r["locus"]: r for r in ed.get("readings") or []}, txt, htmlmod.unescape(html_) @@ -1413,6 +1417,39 @@ def made_up(name, tg, t, ugt, r117h=None, ugt_pc=None, stranger=False): [rd[k]["status"] for k in ("CFTR", "UGT1A1", "SLC6A4")] == ["not_sized", "not_sized", "not_tested"] and "STR_EXTRA_LOCI=yes" in rd["CFTR"]["reading"] and "CFTR intron 9 (legacy IVS8) TG and poly-T tract: " "not sized" in txt, rd) + check("extra repeats: the setup hint is a command that runs, with the script path and the genome folder", + all("STR_EXTRA_LOCI=yes ./scripts/setup.sh --strchive-data " in rd[k]["reading"] + for k in ("CFTR", "UGT1A1")), rd) + + # No 5T copy: the poly-T still matters beside R117H, so the line must not say it changes nothing. + _, _, rd, _, _ = made_up("EXH", "11/11", "7/9", "7/7", r117h="G/A") + c = rd["CFTR"] + check("extra repeats: no 5T copy with one R117H copy says R117H sits on the 7T or 9T copy", + c.get("five_t") is False and "R117H (G/A) is read" in c["reading"] and "7T or 9T copy" in c["reading"] + and "does not change how a CFTR finding is read" not in c["reading"] and c.get("r117h") == "G/A", c) + _, _, rd, _, _ = made_up("EXI", "11/11", "7/7", "7/7") + check("extra repeats: no 5T copy and no traits table says R117H was not read", + "R117H was not read" in rd["CFTR"]["reading"], rd["CFTR"]) + + # A TA length with no named allele: no function claim and no PharmCAT comparison. + _, _, rd, _, _ = made_up("EXJ", None, None, "7/11", ugt_pc="*1/*1") + u = rd["UGT1A1"] + check("extra repeats: UGT1A1 TA6/TA10 makes no function claim and is not compared with PharmCAT", + u.get("genotype") == "*1/TA10" and u.get("agreement") == "not_compared" + and "decreased-function" not in u["reading"] and "not compared with PharmCAT" in u["reading"], u) + + # The traits and CPIC tables are read with their sections' states. + bad_traits = "rsid\tgenotype\tstatus\nrs78655421\tG/A\tcalled\n" + _, _, rd, _, _ = made_up("EXK", "12/12", "5/7", "7/8", traits_text=bad_traits) + check("extra repeats: an unreadable traits table gives no R117H, and says why", + rd["CFTR"].get("phase") == "not_read" and "traits table is unreadable" in rd["CFTR"]["reading"], rd["CFTR"]) + status = f"meta\tstarted_epoch\t{time.time() - DAY}\nstep\t09\tok\nstep\t27\tskipped\nstep\t42\tskipped\n" + _, _, rd, _, _ = made_up("EXL", "12/12", "5/7", "7/8", r117h="G/A", ugt_pc="*1/*28", side_age=3, + run_status=status) + check("extra repeats: R117H and PharmCAT from tables of an earlier run are said to be so", + "R117H from step 42's traits table of an earlier run" in rd["CFTR"]["reading"] + and rd["UGT1A1"]["agreement"] == "agrees" and "step 27's table is from an earlier run" in rd["UGT1A1"]["reading"], + rd) SMN_COLS = ("chrom\tstart\tend\tlocus\tagCN_filter\tagCN\tagCN_qual\tpsCN_filter\tpsCN\tpsCN_qual\t" From 70c074227012ff5ca03828ff4c566beabf870179 Mon Sep 17 00:00:00 2001 From: GeiserX <9169332+GeiserX@users.noreply.github.com> Date: Sat, 7 Nov 2026 17:25:16 +0100 Subject: [PATCH 05/10] docs: say what the report search, the Nextflow report and the tests really cover - 24-html-report: a section's locus name (repeat, Locityper, copy-number locus) is indexed too, not only gene symbols; the Nextflow report lacks only slivar's compound-het counts. - testing: nf-wiring-2-sv-str-tools runs the Nextflow duphold, CNVpytor and ExpansionHunter modules with their real images. - 53-locityper: the image test repeats the Locityper commands of the database build and the step; it does not run setup.sh or step 53. Bead: pgp-5bja.27 Bead: pgp-5bja.28 --- docs/24-html-report.md | 4 ++-- docs/53-locityper.md | 2 +- docs/testing.md | 6 +++--- 3 files changed, 6 insertions(+), 6 deletions(-) diff --git a/docs/24-html-report.md b/docs/24-html-report.md index 5efddee2..31ab1503 100644 --- a/docs/24-html-report.md +++ b/docs/24-html-report.md @@ -68,7 +68,7 @@ After them, "What this report leaves out" lists the steps not run and what no st - **Links to the step's own files**, relative to the report, when the summary says they exist: the PharmCAT HTML report, the CPIC recommendations per gene, the CPSR HTML report (or the `cpsr/` folder that holds it) and the Exomiser ranking table. The report sits in the sample folder beside each step's folder, in the bash layout and in a Nextflow `--outdir`. MultiQC and the couple report are not linked: the summary does not record them, and the couple report sits in the family folder. - **Nothing is fetched.** A Content-Security-Policy (`default-src 'none'; img-src data:; style-src 'unsafe-inline'`) lets the browser load nothing but the file itself and its inline pictures, and there is no JavaScript. The only addresses off this computer are the haplogroup links (FamilyTreeDNA, YFull, the docs) and the project link in the footer, each marked "opens a website"; the FamilyTreeDNA and YFull ones add that the address carries the haplogroup. The genome viewer links go to `127.0.0.1`, this computer. - **A Viewer column on every findings table.** The ClinVar hits, the secondary-findings genes, the mitochondrial variants, the known deep intronic sites, the cDNA scan, the AnnotSV class 4 and 5 SVs and the calls moved to class 3, slivar's pairs in trans, the mosaic alterations and Exomiser's ranked variants each end with a Viewer column: one link per row, to that row's position (a range for an SV, a region or a pair). The link opens the reads there in the [genome viewer](viewer.md) once `scripts/owngenome-view.sh ` runs; a note under the table says so. The column is new, so every other cell keeps its text; it replaces the single line of links that sat under three of these tables. A row with no position has an empty Viewer cell. On the [live site](live-site.md) the served report points the same links at the site's own viewer of the run. -- **Genes and rsIDs in these results.** The closing list is built by `bin/results_index.py`. It holds every gene symbol and rsID that a section in state complete or stale names in `summary.json`, by the same rule `scripts/owngenome-lookup.sh` uses for "in the reports" (`names_of` in `bin/genome_lookup.py`). It adds the genes of the clinical filter (step 23) and variant prioritization (step 31), whose cards only count them: from their summary tables on the bash path, from the VCFs on the Nextflow path. The two give the same genes when a variant's most severe consequence is on the transcript the report ranks first (MANE, else canonical); the step 31 table keeps the gene of the most severe consequence, so where an overlapping gene's transcript carries the worse one, the two forms name different genes. A ClinVar hit's rsIDs come from step 13's VEP output (`Existing_variation`), because the hits file's ID column holds ClinVar's own id; without step 13 a hit has none. The Nextflow `HTML_REPORT` process does not get step 13's VEP output, so its index has no ClinVar rsIDs; step 24, which `run-all.sh` runs over the output folder afterwards, adds them. A failed, missing or unreadable section adds nothing. The browser's find (Ctrl+F) searches the list; there is still no script. A position is not in it: ask `scripts/owngenome-lookup.sh` for a genotype ([genotype lookup](lookup.md)). +- **Genes and rsIDs in these results.** The closing list is built by `bin/results_index.py`. It holds every gene symbol and rsID that a section in state complete or stale names in `summary.json`, by the same rule `scripts/owngenome-lookup.sh` uses for "in the reports" (`names_of` in `bin/genome_lookup.py`). A section's locus name counts as a term too, when it holds an upper-case letter: a repeat expansion (HTT, FMR1, CFTR_TG), a Locityper locus (HLA-DRB1, FCGR) and a copy-number locus are indexed by their `locus` field, so a term is not always a gene symbol. It adds the genes of the clinical filter (step 23) and variant prioritization (step 31), whose cards only count them: from their summary tables on the bash path, from the VCFs on the Nextflow path. The two give the same genes when a variant's most severe consequence is on the transcript the report ranks first (MANE, else canonical); the step 31 table keeps the gene of the most severe consequence, so where an overlapping gene's transcript carries the worse one, the two forms name different genes. A ClinVar hit's rsIDs come from step 13's VEP output (`Existing_variation`), because the hits file's ID column holds ClinVar's own id; without step 13 a hit has none. The Nextflow `HTML_REPORT` process does not get step 13's VEP output, so its index has no ClinVar rsIDs; step 24, which `run-all.sh` runs over the output folder afterwards, adds them. A failed, missing or unreadable section adds nothing. The browser's find (Ctrl+F) searches the list; there is still no script. A position is not in it: ask `scripts/owngenome-lookup.sh` for a genotype ([genotype lookup](lookup.md)). - **Dark on screen, light on paper.** - **Print one section.** Open a section from the contents (or its "link to this section"), then print: only that section prints, with the header (sample, generation date, pipeline commit, ClinVar file date, ClinVar in CPSR's bundle, HLA database), the list of what is not assessed and the disclaimer. A chapter prints the same way from its title in the contents. Without an anchor in the address the whole report prints, as it does in a browser without CSS `:has()`. @@ -177,4 +177,4 @@ start ${GENOME_DIR}/${SAMPLE}/${SAMPLE}_report.html - Steps that were not run show "Not run" and are listed under Steps Not Run — this is expected. A step recorded as failed shows "Failed" instead - The report contains health findings: the ClinVar table lists your pathogenic and likely pathogenic variants, and other sections summarise pharmacogenomics, CPSR and the other steps. It holds no raw reads or full VCF, but share it only as you would share a medical record - Re-run this script anytime to update the report after running additional steps -- The Nextflow `HTML_REPORT` module runs the same `bin/render_report.py` on the outputs of its run's QC, ClinVar, PharmCAT, CPIC, CPSR, clinical filter, slivar, ROH and mitochondrial haplogroup steps, so both reports print the same numbers for the same file. Its clinical filter and slivar cards show counts only. For every section on a Nextflow run, run this script with `GENOME_DIR` set to the `--outdir` +- The Nextflow `HTML_REPORT` module runs the same `bin/render_report.py` on the outputs of its run's QC, ClinVar, PharmCAT, CPIC, CPSR, clinical filter, slivar, ROH and mitochondrial haplogroup steps, so both reports print the same numbers for the same file. Its slivar card has no compound-het variant and gene counts: the Nextflow slivar module writes no `_compound_hets.tsv`. For every section on a Nextflow run, run this script with `GENOME_DIR` set to the `--outdir` diff --git a/docs/53-locityper.md b/docs/53-locityper.md index 5b4f415d..9c89183c 100644 --- a/docs/53-locityper.md +++ b/docs/53-locityper.md @@ -88,7 +88,7 @@ Both reports carry a section "Locityper Second Opinion at Hard Loci" in the medi ## How it was checked -The image test (`tests/smoke/locityper.sh`, the `LOCITYPER_IMAGE` row of `tests/smoke/commands.tsv`) runs the database build and the step on the fixture's GIAB HG002 reads at HLA-DRB1, HLA-DQB1 and CYP2D6. The panel holds HG002's own two assembled haplotypes, which are the truth. Both called haplotypes at each locus must match them at QV 20 or more (one difference per 100 bases or fewer), aligned with `locityper align`. That QV measures how close each called haplotype is to HG002's true one; it is not the `quality` column, Locityper's own genotype quality. +The image test (`tests/smoke/locityper.sh`, the `LOCITYPER_IMAGE` row of `tests/smoke/commands.tsv`) repeats, inside the image, the Locityper commands that `setup.sh --locityper-data` and step 53 run (the database build, `preproc` and `genotype`) on the fixture's GIAB HG002 reads at HLA-DRB1, HLA-DQB1 and CYP2D6. It does not run `setup.sh` or `scripts/53-locityper.sh` themselves. The panel holds HG002's own two assembled haplotypes, which are the truth. Both called haplotypes at each locus must match them at QV 20 or more (one difference per 100 bases or fewer), aligned with `locityper align`. That QV measures how close each called haplotype is to HG002's true one; it is not the `quality` column, Locityper's own genotype quality. The red control takes out of the panel HG002's haplotypes and every other haplotype within QV 20 of them, rebuilds the database and genotypes again; the same check must then fail at every locus. Leaving HG002 out with `--leave-out` is not enough: other people in the panel carry near copies of HG002's haplotypes, and the call stays close. diff --git a/docs/testing.md b/docs/testing.md index fdabfd79..43b14519 100644 --- a/docs/testing.md +++ b/docs/testing.md @@ -168,15 +168,15 @@ This is every step script that no e2e case runs, plus the full-size index build | 02b (long-read alignment) | long reads; the e2e fixture has none | the minimap2 image's smoke test, which maps synthetic HiFi reads with `-x map-hifi`; nothing runs the script | | 03e (Clair3) | long reads; the e2e fixture has none | the image's smoke test, with a recall check on synthetic HiFi reads and a check that the two model directories the script names exist; the fake-docker suite runs the script | | 04c (Sniffles2) | long reads; the e2e fixture has none | the image's smoke test on the synthetic HiFi alignment; `tests/fake-docker/bash-step-partial-outputs.sh` runs the script | -| 15 (duphold) | no case yet | the image's smoke test on a planted deletion; `tests/fake-docker/bash-step-sv.sh` runs the script; the Nextflow module is stub-run | -| 18 (CNVpytor) | no case yet; it needs its own resource files | the image's smoke test; the SV merge case reads a planted CNVpytor VCF, not one CNVpytor wrote | +| 15 (duphold) | the bash step: no case yet | case `nf-wiring-2-sv-str-tools` runs the Nextflow `DUPHOLD` and `DUPHOLD_FILTER` with the real image on the fixture's Manta calls; the image's smoke test on a planted deletion; `tests/fake-docker/bash-step-sv.sh` runs the script | +| 18 (CNVpytor) | the bash step: no case yet | case `nf-wiring-2-sv-str-tools` runs the Nextflow `CNVPYTOR` and `CNVPYTOR_VCF` with the real image and the pinned GC and mask files; the image's smoke test; the SV merge case reads a planted CNVpytor VCF, not one CNVpytor wrote | | 29 (somatic Mutect2) | no case yet | `tests/fake-docker/bash-step-somatic-mito.sh` checks the Mutect2 command line | | 03c (Strelka2) | no case yet | the image's smoke test, with a recall check against the truth slice | | 03d (Octopus) | no case yet | the image's smoke test, with a recall check against the truth slice; case `bash-step-bash42` runs the script under bash 4.2 against the fake docker | | 04a (TIDDIT) | no case yet | the image's smoke test; settle-doubts question 4 runs the script with only a BWA-MEM2 index; case `bash-step-bash42` and `tests/fake-docker/bash-step-sv.sh` run it against the fake docker | | chip-to-vcf (array data to VCF) on a vendor download | case `consumer-export-1` runs it only on the 23andMe and AncestryDNA files step 40 wrote | `tests/test_chip_to_vcf.sh` in the unit-tests job runs the script with the real bcftools and Picard images on two ten-row vendor files | -The Nextflow ExpansionHunter module is only stub-run (`-profile test_all -stub` in `nextflow.yml`); the bash step 09 runs in case `bash-step-09-str`. +Case `nf-wiring-2-sv-str-tools` runs the Nextflow `EXPANSION_HUNTER` module with the real image on the fixture's BAM (the bundled catalog cut to the fixture's contigs); the bash step 09 runs in case `bash-step-09-str`. ## The settle-doubts job From 7a4c99e1fe449835e619c89a83001c10232e2997 Mon Sep 17 00:00:00 2001 From: GeiserX <9169332+GeiserX@users.noreply.github.com> Date: Sat, 7 Nov 2026 17:27:45 +0100 Subject: [PATCH 06/10] fix(schema): name the cram,crai and gvcf,gvcf_index samplesheet columns in --input main.nf reads both pairs, but the schema's --input description named only the FASTQ, BAM and VCF columns. Bead: pgp-5bja.29 --- nextflow_schema.json | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/nextflow_schema.json b/nextflow_schema.json index 6795fd58..d22fecf3 100644 --- a/nextflow_schema.json +++ b/nextflow_schema.json @@ -15,7 +15,7 @@ "type": "string", "format": "file-path", "mimetype": "text/csv", - "description": "Path to comma-separated samplesheet: sample, then fastq_1,fastq_2 (aligned and called), or bam,bam_index (called), or vcf,vcf_index (a VCF with one sample) with an optional bam,bam_index; sex (male or female) is required on every row that is called and when expansion_hunter is enabled.", + "description": "Path to comma-separated samplesheet: sample, then fastq_1,fastq_2 (aligned and called), or bam,bam_index (called), or vcf,vcf_index (a VCF with one sample) with an optional bam,bam_index; sex (male or female) is required on every row that is called and when expansion_hunter is enabled. A row can also give cram,crai instead of the BAM, and gvcf,gvcf_index beside the VCF.", "fa_icon": "fas fa-file-csv" }, "outdir": { From 7d310a41b5c6cd73d8c9504e14bc2c36b413c673 Mon Sep 17 00:00:00 2001 From: GeiserX <9169332+GeiserX@users.noreply.github.com> Date: Sat, 7 Nov 2026 17:28:00 +0100 Subject: [PATCH 07/10] test(docs): fail when a schema parameter has no row in docs/parameters.md The page is written by hand from nextflow_schema.json, so a new parameter could ship without a row. The test also fails on a row for a parameter the schema lacks (main.nf's --pedigree is allowed by name), and runs a dropped row and a planted row as controls. Bead: pgp-5bja.30 --- tests/test_parameters_doc.py | 76 ++++++++++++++++++++++++++++++++++++ 1 file changed, 76 insertions(+) create mode 100644 tests/test_parameters_doc.py diff --git a/tests/test_parameters_doc.py b/tests/test_parameters_doc.py new file mode 100644 index 00000000..a04274e7 --- /dev/null +++ b/tests/test_parameters_doc.py @@ -0,0 +1,76 @@ +#!/usr/bin/env python3 +"""docs/parameters.md lists every parameter of nextflow_schema.json. + +Checks: + 1. every parameter of the schema (each group's properties, and any + top-level ones) has a row "| `--` |" in docs/parameters.md; + 2. every row of the page names a schema parameter, or a parameter main.nf + reads that the schema does not hold yet (DOC_ONLY, each checked against + main.nf), so a removed or renamed parameter cannot keep its row; + 3. controls: the page with one row taken out fails check 1 and names that + parameter; the page with a planted row fails check 2. +""" +import json +import os +import re +import sys + +ROOT = os.path.dirname(os.path.dirname(os.path.abspath(__file__))) +SCHEMA = os.path.join(ROOT, "nextflow_schema.json") +PAGE = os.path.join(ROOT, "docs", "parameters.md") +MAIN = os.path.join(ROOT, "main.nf") +# Read by main.nf, not (yet) in the schema; the page writes their rows by hand. +DOC_ONLY = {"pedigree"} +ROW = re.compile(r"^\| `--([A-Za-z0-9_]+)` \|", re.M) + +failures = 0 + + +def check(desc, ok, detail=""): + global failures + print(f"[{'PASS' if ok else 'FAIL'}] {desc}" + ("" if ok else f" -- {detail}")) + if not ok: + failures += 1 + + +def schema_params(schema): + names = set(schema.get("properties", {})) + for group in schema.get("definitions", schema.get("$defs", {})).values(): + names |= set(group.get("properties", {})) + return names + + +def problems(params, page_text, main_text): + rows = ROW.findall(page_text) + missing = sorted(params - set(rows)) + extra = sorted(r for r in set(rows) - params + if r not in DOC_ONLY or not re.search(rf"\bparams\.{r}\b|'{r}'", main_text)) + return missing, extra + + +with open(SCHEMA, encoding="utf-8") as f: + params = schema_params(json.load(f)) +with open(PAGE, encoding="utf-8") as f: + page = f.read() +with open(MAIN, encoding="utf-8") as f: + main_nf = f.read() + +check("the schema has parameters to compare", len(params) > 20, len(params)) +missing, extra = problems(params, page, main_nf) +check("every schema parameter has a row in docs/parameters.md", not missing, f"no row for: {', '.join(missing)}") +check("every row of docs/parameters.md is a schema parameter or one main.nf reads", not extra, + f"not in the schema: {', '.join(extra)}") +check("each parameter has one row", len(ROW.findall(page)) == len(set(ROW.findall(page))), + sorted(r for r in set(ROW.findall(page)) if ROW.findall(page).count(r) > 1)) + +# Controls: the same comparison must catch a dropped row and a planted one. +victim = sorted(params)[0] +dropped = re.sub(rf"^\| `--{victim}` \|.*\n", "", page, flags=re.M) +check(f"control: the page without the --{victim} row fails, naming it", + problems(params, dropped, main_nf)[0] == [victim], problems(params, dropped, main_nf)) +planted = page + "| `--no_such_parameter` | none | planted |\n" +check("control: a planted row for a parameter the schema lacks fails", + problems(params, planted, main_nf)[1] == ["no_such_parameter"], problems(params, planted, main_nf)) + +print(f"\nRESULT: {'ALL PASS' if not failures else f'{failures} FAILED'}") +sys.exit(1 if failures else 0) From d85e9d5759bf05d3d6aeec7dfb54cd3a3b8d2067 Mon Sep 17 00:00:00 2001 From: GeiserX <9169332+GeiserX@users.noreply.github.com> Date: Sat, 7 Nov 2026 18:07:56 +0100 Subject: [PATCH 08/10] fix(report): say a homozygous R117H with no 5T copy sits on both copies --- bin/collect_summary.py | 3 +++ tests/test_collect_summary.py | 9 +++++++++ 2 files changed, 12 insertions(+) diff --git a/bin/collect_summary.py b/bin/collect_summary.py index 953f7e0d..25f2e285 100755 --- a/bin/collect_summary.py +++ b/bin/collect_summary.py @@ -2186,6 +2186,9 @@ def cftr_reading(by, d, s, traits_state="ok"): "with 5T in cis.") if r117 is None: line += f" R117H was not read ({why})." + elif r117.split("/").count("A") == 2: + copies = f"both {b[0]}T copies" if b[0] == b[1] else f"both copies ({b[0]}T and {b[1]}T)" + line += f" R117H ({r117}) is read: it sits on {copies}.{note}" elif r117.split("/").count("A"): copies = f"{b[0]}T" if b[0] == b[1] else f"{b[0]}T or {b[1]}T" line += f" R117H ({r117}) is read: with no 5T copy it sits on a {copies} copy.{note}" diff --git a/tests/test_collect_summary.py b/tests/test_collect_summary.py index 5f0f6e4a..5ef81ccc 100644 --- a/tests/test_collect_summary.py +++ b/tests/test_collect_summary.py @@ -1427,6 +1427,15 @@ def made_up(name, tg, t, ugt, r117h=None, ugt_pc=None, stranger=False, traits_te check("extra repeats: no 5T copy with one R117H copy says R117H sits on the 7T or 9T copy", c.get("five_t") is False and "R117H (G/A) is read" in c["reading"] and "7T or 9T copy" in c["reading"] and "does not change how a CFTR finding is read" not in c["reading"] and c.get("r117h") == "G/A", c) + # Homozygous R117H with no 5T copy sits on both copies, not on one or the other. + _, _, rd, _, _ = made_up("EXK", "11/11", "7/9", "7/7", r117h="A/A") + c = rd["CFTR"] + check("extra repeats: no 5T copy with two R117H copies says R117H sits on both copies (7T and 9T)", + c.get("five_t") is False and "R117H (A/A) is read: it sits on both copies (7T and 9T)" in c["reading"] + and "it sits on a 7T or 9T copy" not in c["reading"], c) + _, _, rd, _, _ = made_up("EXL", "11/11", "7/7", "7/7", r117h="A/A") + check("extra repeats: two R117H copies on two 7T copies says both 7T copies", + "R117H (A/A) is read: it sits on both 7T copies" in rd["CFTR"]["reading"], rd["CFTR"]) _, _, rd, _, _ = made_up("EXI", "11/11", "7/7", "7/7") check("extra repeats: no 5T copy and no traits table says R117H was not read", "R117H was not read" in rd["CFTR"]["reading"], rd["CFTR"]) From 6d8c717987198d936debef029da04083a4d199b3 Mon Sep 17 00:00:00 2001 From: GeiserX <9169332+GeiserX@users.noreply.github.com> Date: Sat, 7 Nov 2026 18:33:17 +0100 Subject: [PATCH 09/10] docs(report): say a chrX, chrY or chrM expansion position can be a search term --- docs/24-html-report.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/24-html-report.md b/docs/24-html-report.md index 31ab1503..ee465b71 100644 --- a/docs/24-html-report.md +++ b/docs/24-html-report.md @@ -68,7 +68,7 @@ After them, "What this report leaves out" lists the steps not run and what no st - **Links to the step's own files**, relative to the report, when the summary says they exist: the PharmCAT HTML report, the CPIC recommendations per gene, the CPSR HTML report (or the `cpsr/` folder that holds it) and the Exomiser ranking table. The report sits in the sample folder beside each step's folder, in the bash layout and in a Nextflow `--outdir`. MultiQC and the couple report are not linked: the summary does not record them, and the couple report sits in the family folder. - **Nothing is fetched.** A Content-Security-Policy (`default-src 'none'; img-src data:; style-src 'unsafe-inline'`) lets the browser load nothing but the file itself and its inline pictures, and there is no JavaScript. The only addresses off this computer are the haplogroup links (FamilyTreeDNA, YFull, the docs) and the project link in the footer, each marked "opens a website"; the FamilyTreeDNA and YFull ones add that the address carries the haplogroup. The genome viewer links go to `127.0.0.1`, this computer. - **A Viewer column on every findings table.** The ClinVar hits, the secondary-findings genes, the mitochondrial variants, the known deep intronic sites, the cDNA scan, the AnnotSV class 4 and 5 SVs and the calls moved to class 3, slivar's pairs in trans, the mosaic alterations and Exomiser's ranked variants each end with a Viewer column: one link per row, to that row's position (a range for an SV, a region or a pair). The link opens the reads there in the [genome viewer](viewer.md) once `scripts/owngenome-view.sh ` runs; a note under the table says so. The column is new, so every other cell keeps its text; it replaces the single line of links that sat under three of these tables. A row with no position has an empty Viewer cell. On the [live site](live-site.md) the served report points the same links at the site's own viewer of the run. -- **Genes and rsIDs in these results.** The closing list is built by `bin/results_index.py`. It holds every gene symbol and rsID that a section in state complete or stale names in `summary.json`, by the same rule `scripts/owngenome-lookup.sh` uses for "in the reports" (`names_of` in `bin/genome_lookup.py`). A section's locus name counts as a term too, when it holds an upper-case letter: a repeat expansion (HTT, FMR1, CFTR_TG), a Locityper locus (HLA-DRB1, FCGR) and a copy-number locus are indexed by their `locus` field, so a term is not always a gene symbol. It adds the genes of the clinical filter (step 23) and variant prioritization (step 31), whose cards only count them: from their summary tables on the bash path, from the VCFs on the Nextflow path. The two give the same genes when a variant's most severe consequence is on the transcript the report ranks first (MANE, else canonical); the step 31 table keeps the gene of the most severe consequence, so where an overlapping gene's transcript carries the worse one, the two forms name different genes. A ClinVar hit's rsIDs come from step 13's VEP output (`Existing_variation`), because the hits file's ID column holds ClinVar's own id; without step 13 a hit has none. The Nextflow `HTML_REPORT` process does not get step 13's VEP output, so its index has no ClinVar rsIDs; step 24, which `run-all.sh` runs over the output folder afterwards, adds them. A failed, missing or unreadable section adds nothing. The browser's find (Ctrl+F) searches the list; there is still no script. A position is not in it: ask `scripts/owngenome-lookup.sh` for a genotype ([genotype lookup](lookup.md)). +- **Genes and rsIDs in these results.** The closing list is built by `bin/results_index.py`. It holds every gene symbol and rsID that a section in state complete or stale names in `summary.json`, by the same rule `scripts/owngenome-lookup.sh` uses for "in the reports" (`names_of` in `bin/genome_lookup.py`). A section's locus name counts as a term too, when it holds an upper-case letter: a repeat expansion (HTT, FMR1, CFTR_TG), a Locityper locus (HLA-DRB1, FCGR) and a copy-number locus are indexed by their `locus` field, so a term is not always a gene symbol. An expansion record with no REPID or VARID takes its position as its locus (`chrX:147912050`), and that counts as a term only when it holds an upper-case letter (a chrX, chrY or chrM position). It adds the genes of the clinical filter (step 23) and variant prioritization (step 31), whose cards only count them: from their summary tables on the bash path, from the VCFs on the Nextflow path. The two give the same genes when a variant's most severe consequence is on the transcript the report ranks first (MANE, else canonical); the step 31 table keeps the gene of the most severe consequence, so where an overlapping gene's transcript carries the worse one, the two forms name different genes. A ClinVar hit's rsIDs come from step 13's VEP output (`Existing_variation`), because the hits file's ID column holds ClinVar's own id; without step 13 a hit has none. The Nextflow `HTML_REPORT` process does not get step 13's VEP output, so its index has no ClinVar rsIDs; step 24, which `run-all.sh` runs over the output folder afterwards, adds them. A failed, missing or unreadable section adds nothing. The browser's find (Ctrl+F) searches the list; there is still no script. Apart from those locus names a position is not in it, and the search does not answer a position: ask `scripts/owngenome-lookup.sh` for a genotype ([genotype lookup](lookup.md)). - **Dark on screen, light on paper.** - **Print one section.** Open a section from the contents (or its "link to this section"), then print: only that section prints, with the header (sample, generation date, pipeline commit, ClinVar file date, ClinVar in CPSR's bundle, HLA database), the list of what is not assessed and the disclaimer. A chapter prints the same way from its title in the contents. Without an anchor in the address the whole report prints, as it does in a browser without CSS `:has()`. From c5c465c396f710786c604595b7108d5bfc66a4a8 Mon Sep 17 00:00:00 2001 From: GeiserX <9169332+GeiserX@users.noreply.github.com> Date: Sun, 8 Nov 2026 02:46:52 +0100 Subject: [PATCH 10/10] docs(export): say why the round trip leaves chrM out now that chip-to-vcf moves 23andMe MT sites --- docs/40-export-consumer.md | 2 +- tests/e2e/consumer-export-1.sh | 6 ++++-- 2 files changed, 5 insertions(+), 3 deletions(-) diff --git a/docs/40-export-consumer.md b/docs/40-export-consumer.md index c2d1bbc5..c281cf47 100644 --- a/docs/40-export-consumer.md +++ b/docs/40-export-consumer.md @@ -81,4 +81,4 @@ On the CI fixture (the gVCF of the fixture's 15 regions, 7,623 sites of the CI t - The AncestryDNA file has the 23andMe v5 sites, not the sites of an AncestryDNA chip; the two chips share part of their sites. A service that matches on the sites it expects uses the ones the files have in common. - The 23andMe file has no internal `i` sites and none of the dropped sites, so it has fewer rows than a file from 23andMe. - In a duplicated region the list's GRCh37 position and dbSNP's GRCh38 position can sit on two different copies. The file then carries the list's position with the genotype read at dbSNP's. On the CI table this is one site of 7,338 off chrM (in SMN1/SMN2), found by the round trip below. -- `scripts/chip-to-vcf.sh` reads both files back. The e2e case `consumer-export-1` does that on the CI fixture, lifts the result to GRCh38 and checks it against the fixture's calls; chrM is left out there, because the hg19 chain maps hg19's chrM, not the rCRS positions the file uses. +- `scripts/chip-to-vcf.sh` reads both files back. The e2e case `consumer-export-1` does that on the CI fixture, lifts the result to GRCh38 and checks it against the fixture's calls; chrM is left out of both files there. `chip-to-vcf.sh` treats the MT positions of a 23andMe-layout file as 23andMe numbers them, on hg19's chrM, and moves each one to the rCRS through the chain's chrM block. This file already writes rCRS positions, so its MT sites would land one or two bases off. The AncestryDNA file's MT positions are read as rCRS positions and keep their place, but the case checks neither. diff --git a/tests/e2e/consumer-export-1.sh b/tests/e2e/consumer-export-1.sh index b858f3ea..238d3aa9 100755 --- a/tests/e2e/consumer-export-1.sh +++ b/tests/e2e/consumer-export-1.sh @@ -15,8 +15,10 @@ # table's GRCh38 position and has the alleles of case 21's calls there # (hom-ref where case 21 called nothing), but for at most 0.1% of the # sites, where the list's GRCh37 position and dbSNP's GRCh38 position -# fall on two copies of a duplicated region. chrM is left out: the hg19 -# chain maps hg19's chrM, not the rCRS positions the file uses; +# fall on two copies of a duplicated region. chrM is left out of both +# files: chip-to-vcf.sh moves a 23andMe-layout file's MT positions from +# hg19's chrM numbering to the rCRS, and step 40 already writes rCRS +# positions, so those sites would land one or two bases off; # 5. red controls: the table shifted by one GRCh38 position stops the # export (reference bases differ); shifted by one GRCh37 position, the # export runs and the round trip fails;