From f1540198876bb7ceeb391f201cafb06ec691bdc1 Mon Sep 17 00:00:00 2001 From: Victor Lin Date: Tue, 11 Aug 2026 11:21:53 -0700 Subject: [PATCH 1/5] Configuration schema and validation Adds a YAML schema file for the structure of the current 'config' variable and uses it for validation on the variable at the start of the workflow. Date values in the default config file are explicitly quoted to avoid being parsed into Python as date objects by Snakemake's YAML parser. --- CHANGELOG.md | 4 + config.schema.yaml | 226 ++++++++++++++++++++++++++++ config/configfile.yaml | 4 +- workflow/snakemake_rules/config.smk | 40 ++++- 4 files changed, 271 insertions(+), 3 deletions(-) create mode 100644 config.schema.yaml diff --git a/CHANGELOG.md b/CHANGELOG.md index d4d7144..10b32e9 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,3 +1,7 @@ # CHANGELOG We use this CHANGELOG to document breaking changes, new features, bug fixes, and config value changes that may affect both the usage of the workflows and the outputs of the workflows. + +## 2026 + +* TBD: Phylogenetic workflow configuration is now validated against a strict schema. The workflow will error if your configuration has extraneous entries that were previously ignored. diff --git a/config.schema.yaml b/config.schema.yaml new file mode 100644 index 0000000..774060b --- /dev/null +++ b/config.schema.yaml @@ -0,0 +1,226 @@ +$schema: "https://json-schema.org/draft/2020-12/schema" +title: RSV Phylogenetic Workflow Configuration +description: >- + This is the schema for the Nextstrain rsv phylogenetic workflow's + configuration file. + +$defs: + per_subtype_map: &per_subtype_map + type: object + additionalProperties: false + propertyNames: + title: Subtype + description: Subtype ('a' or 'b') + + per_build_map: &per_build_map + type: object + additionalProperties: false + propertyNames: + title: Build name + description: Build name ('genome', 'G', 'F', or 'F-antibody-escape') + + per_resolution_map: &per_resolution_map + type: object + additionalProperties: false + propertyNames: + title: Resolution name + description: Resolution name ('all-time', '6y', or '3y') + + input_item: + type: object + additionalProperties: false + required: + - name + anyOf: + - required: [metadata] + - required: [sequences] + properties: + name: + type: string + metadata: + type: string + sequences: + type: string + +type: object +additionalProperties: false +properties: + conda_environment: + type: string + genesforglycosylation: + type: array + items: + type: string + builds_to_run: + type: array + items: + type: string + resolutions_to_run: + type: array + items: + type: string + subtypes: + type: array + items: + type: string + inputs: + type: array + items: + $ref: "#/$defs/input_item" + additional_inputs: + type: array + items: + $ref: "#/$defs/input_item" + exclude: + type: string + description: + type: string + strain_id_field: + type: string + display_strain_field: + type: string + filter: + type: object + additionalProperties: false + properties: + group_by: + type: string + min_coverage: + <<: *per_build_map + patternProperties: + "^.*$": + type: number + min_length: + <<: *per_build_map + patternProperties: + "^.*$": + type: integer + resolutions: + <<: *per_resolution_map + patternProperties: + "^.*$": + type: object + additionalProperties: false + properties: + min_date: + type: string + background_min_date: + type: string + subsample_max_sequences: + <<: *per_build_map + patternProperties: + "^.*$": + type: integer + exclude_where: + type: object + additionalProperties: false + properties: + recent: + type: array + items: + type: string + background: + type: array + items: + type: string + missing_data_threshold: + type: integer + files: + type: object + additionalProperties: false + properties: + auspice_config: + type: string + auspice_config_additional_colorings: + type: string + auspice_config_f_antibody_escape: + type: string + auspice_config_non-genome_builds: + type: string + refine: + type: object + additionalProperties: false + properties: + coalescent: + type: string + date_inference: + type: string + clock_filter_iqd: + type: number + divergence_units: + type: string + ancestral: + type: object + additionalProperties: false + properties: + inference: + type: string + cds: + <<: *per_build_map + patternProperties: + "^.*$": + type: string + traits: + type: object + additionalProperties: false + properties: + columns: + type: [string, array] + items: + type: string + frequencies: + type: object + additionalProperties: false + properties: + resolutions: + <<: *per_resolution_map + patternProperties: + "^.*$": + type: object + additionalProperties: false + properties: + min_date: + type: string + nextclade_attributes: + <<: *per_subtype_map + patternProperties: + "^.*$": + type: object + additionalProperties: false + properties: + name: + type: string + reference_name: + type: string + accession: + type: string + f_dms_data: + type: string + f_dms_antibodies: + type: array + items: + type: string + dms_only_positive_escape: + type: boolean + enrich_antibody_escape: + <<: *per_build_map + patternProperties: + "^.*$": + type: object + additionalProperties: false + properties: + nseqs_per_antibody_scoretype: + type: integer + group_by: + type: array + items: + type: string + max_identical_f_prot_muts: + type: integer + max_identical_max_escape_mut: + type: integer + custom_rules: + type: array + description: Custom Snakemake rule files to include. If used, this will disable config schema validation. + items: + type: string diff --git a/config/configfile.yaml b/config/configfile.yaml index 3aa617b..02f67c3 100644 --- a/config/configfile.yaml +++ b/config/configfile.yaml @@ -44,7 +44,7 @@ filter: F-antibody-escape: 1200 resolutions: all-time: - min_date: 1975-01-01 + min_date: "1975-01-01" 6y: min_date: 6Y background_min_date: 12Y @@ -91,7 +91,7 @@ traits: frequencies: resolutions: all-time: - min_date: 1975-01-01 + min_date: "1975-01-01" 6y: min_date: 6Y 3y: diff --git a/workflow/snakemake_rules/config.smk b/workflow/snakemake_rules/config.smk index 199b4a9..ab284ff 100644 --- a/workflow/snakemake_rules/config.smk +++ b/workflow/snakemake_rules/config.smk @@ -5,4 +5,42 @@ OUTPUTS: results/run_config.yaml """ -write_config("results/run_config.yaml") +import sys +from augur.validate import load_json_schema_locally, validate_json, ValidateError +from pathlib import Path + +def main(): + dump_and_validate( + "results/run_config.yaml", + Path(workflow.basedir) / "config.schema.yaml" + ) + + +# TODO: move this to nextstrain/shared +# Copied from measles with some cleanups +# +def dump_and_validate(dump_path, schema_path): + """ + Write Snakemake's 'config' variable to a file, then validate it against the + schema. Do both in the same function so that the validation output can + easily reference the path of the dumped config for inspection. + """ + global config + + write_config(dump_path) + + if "custom_rules" in config: + print("WARNING: Skipping config schema validation because custom rules are defined.", file=sys.stderr) + return + + try: + validator = load_json_schema_locally(schema_path) + validate_json(config, validator, dump_path) + except ValidateError as e: + raise InvalidConfigError(str(e)) from e + +try: + main() +except InvalidConfigError as e: + print(f"ERROR: {e}", file=sys.stderr) + exit(1) From ac456a6e1ac04a7ada579dd7ee57c3f90fedf71a Mon Sep 17 00:00:00 2001 From: Victor Lin Date: Wed, 5 Aug 2026 17:24:47 -0700 Subject: [PATCH 2/5] Add separate config section for filter_for_pre_subsample_alignment Similar to "Add separate frequencies config" (0b221851), the filter_for_pre_subsample_alignment rule shouldn't rely on config from another rule. --- CHANGELOG.md | 1 + config.schema.yaml | 27 +++++++++++++++++++++++++++ config/configfile.yaml | 19 +++++++++++++++++++ workflow/snakemake_rules/core.smk | 6 +++--- 4 files changed, 50 insertions(+), 3 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 10b32e9..754e587 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -4,4 +4,5 @@ We use this CHANGELOG to document breaking changes, new features, bug fixes, and ## 2026 +* TBD: Phylogenetic workflow configuration for initial quality filtering has moved to its own section, `filter_for_pre_subsample_alignment`. * TBD: Phylogenetic workflow configuration is now validated against a strict schema. The workflow will error if your configuration has extraneous entries that were previously ignored. diff --git a/config.schema.yaml b/config.schema.yaml index 774060b..0b1e5b2 100644 --- a/config.schema.yaml +++ b/config.schema.yaml @@ -194,6 +194,33 @@ properties: type: string accession: type: string + filter_for_pre_subsample_alignment: + type: object + additionalProperties: false + properties: + group_by: + type: string + min_length: + <<: *per_build_map + patternProperties: + "^.*$": + type: integer + min_coverage: + <<: *per_build_map + patternProperties: + "^.*$": + type: number + resolutions: + <<: *per_resolution_map + patternProperties: + "^.*$": + type: object + additionalProperties: false + properties: + min_date: + type: string + background_min_date: + type: string f_dms_data: type: string f_dms_antibodies: diff --git a/config/configfile.yaml b/config/configfile.yaml index 02f67c3..48467a2 100644 --- a/config/configfile.yaml +++ b/config/configfile.yaml @@ -107,6 +107,25 @@ nextclade_attributes: reference_name: "hRSV/B/Australia/VIC-RCH056/2019" accession: "EPI_ISL_1653999" +filter_for_pre_subsample_alignment: + min_length: + genome: 10000 + G: 600 + F: 1200 + F-antibody-escape: 1200 + min_coverage: + genome: 0.3 + G: 0.3 + F: 0.3 + F-antibody-escape: 0.75 + resolutions: + all-time: + min_date: "1975-01-01" + 6y: + min_date: 6Y + 3y: + min_date: 3Y + # configuration specific to the F deep mutational scanning antibody escape data f_dms_data: dms-data/all_antibodies.csv f_dms_antibodies: # columns in `f_dms_data` with per-mutation escape diff --git a/workflow/snakemake_rules/core.smk b/workflow/snakemake_rules/core.smk index 96f2fe9..750f117 100644 --- a/workflow/snakemake_rules/core.smk +++ b/workflow/snakemake_rules/core.smk @@ -235,10 +235,10 @@ rule filter_for_pre_subsample_alignment: benchmark: "benchmarks/filter_for_pre_subsample_alignment_{a_or_b}_{build_name}_{resolution}.txt" params: - min_coverage=lambda w: f'{w.build_name.split("-")[0]}_coverage>{config["filter"]["min_coverage"][w.build_name]}', - min_length=lambda w: config["filter"]["min_length"][w.build_name], + min_coverage=lambda w: f'{w.build_name.split("-")[0]}_coverage>{config["filter_for_pre_subsample_alignment"]["min_coverage"][w.build_name]}', + min_length=lambda w: config["filter_for_pre_subsample_alignment"]["min_length"][w.build_name], strain_id=config["strain_id_field"], - min_date=lambda w: config["filter"]["resolutions"][w.resolution]["min_date"], + min_date=lambda w: config["filter_for_pre_subsample_alignment"]["resolutions"][w.resolution]["min_date"], shell: r""" exec &> >(tee {log:q}) From 494f760c50e77bb4bc8a0ddb209d03df715e2e57 Mon Sep 17 00:00:00 2001 From: Victor Lin Date: Wed, 5 Aug 2026 17:08:53 -0700 Subject: [PATCH 3/5] Use augur subsample The previous subsampling implementation was fixed to a two-sample recent+background split with some hardcoded parameters. Replacing it with augur subsample allows for more flexible configuration. In Snakemake, implementation is mostly copied from pathogen repos that have switched over to augur subsample. One notable difference is that the combine_samples rule must stay to handle output from the enrich_antibody_escape rule. In the config YAML, the subsampling configuration is much more verbose as a byproduct of increased flexibility. It was generated using a script, which I'll add in another commit since it makes additional changes. --- CHANGELOG.md | 3 +- config.schema.yaml | 66 +--- config/configfile.yaml | 572 ++++++++++++++++++++++++++-- workflow/snakemake_rules/config.smk | 16 + workflow/snakemake_rules/core.smk | 105 +---- 5 files changed, 586 insertions(+), 176 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 754e587..9d7bc6d 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -4,5 +4,6 @@ We use this CHANGELOG to document breaking changes, new features, bug fixes, and ## 2026 -* TBD: Phylogenetic workflow configuration for initial quality filtering has moved to its own section, `filter_for_pre_subsample_alignment`. +* TBD: The `filter` section in phylogenetic workflow configuration has been replaced by `subsample`/`custom_subsample` for subsampling, and `filter_for_pre_subsample_alignment` for initial quality filtering. **This is a breaking change**. + * NOTE: The workflow does not yet support proximal samples. * TBD: Phylogenetic workflow configuration is now validated against a strict schema. The workflow will error if your configuration has extraneous entries that were previously ignored. diff --git a/config.schema.yaml b/config.schema.yaml index 0b1e5b2..84f0d63 100644 --- a/config.schema.yaml +++ b/config.schema.yaml @@ -5,6 +5,13 @@ description: >- configuration file. $defs: + per_full_build_map: &per_full_build_map + type: object + additionalProperties: false + propertyNames: + title: Full build name + description: Full build name (e.g. 'a/genome/all-time') + per_subtype_map: &per_subtype_map type: object additionalProperties: false @@ -79,52 +86,19 @@ properties: type: string display_strain_field: type: string - filter: - type: object - additionalProperties: false - properties: - group_by: - type: string - min_coverage: - <<: *per_build_map - patternProperties: - "^.*$": - type: number - min_length: - <<: *per_build_map - patternProperties: - "^.*$": - type: integer - resolutions: - <<: *per_resolution_map - patternProperties: - "^.*$": - type: object - additionalProperties: false - properties: - min_date: - type: string - background_min_date: - type: string - subsample_max_sequences: - <<: *per_build_map - patternProperties: - "^.*$": - type: integer - exclude_where: - type: object - additionalProperties: false - properties: - recent: - type: array - items: - type: string - background: - type: array - items: - type: string - missing_data_threshold: - type: integer + subsample: &subsample_config + <<: *per_full_build_map + description: >- + Subsampling configuration. When using --configfile, it is recommended to + use 'custom_subsample' instead to ignore default subsampling configuration. + patternProperties: + "^.*$": + $ref: "https://nextstrain.org/schemas/augur/subsample-config/v1#/$defs/schemaForUnalignedSequences" + custom_subsample: + <<: *subsample_config + description: >- + Custom subsampling configuration. When using --configfile, this is + recommended over 'subsample' to ignore default subsampling configuration. files: type: object additionalProperties: false diff --git a/config/configfile.yaml b/config/configfile.yaml index 48467a2..5939779 100644 --- a/config/configfile.yaml +++ b/config/configfile.yaml @@ -29,41 +29,543 @@ strain_id_field: "accession" display_strain_field: "strain" -filter: - group_by: "year country" - min_coverage: - genome: 0.3 - G: 0.3 - F: 0.3 - F-antibody-escape: 0.75 - - min_length: - genome: 10000 - G: 600 - F: 1200 - F-antibody-escape: 1200 - resolutions: - all-time: - min_date: "1975-01-01" - 6y: - min_date: 6Y - background_min_date: 12Y - 3y: - min_date: 3Y - background_min_date: 12Y - - subsample_max_sequences: - genome: 3000 - G: 3000 - F: 3000 - F-antibody-escape: 2000 - - exclude_where: - recent: ["qc.overallStatus=bad"] - background: ["qc.overallStatus=bad", "qc.overallStatus=mediocre"] - - missing_data_threshold: 1000 - +subsample: + a/genome/all-time: + samples: + sample: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: '1975-01-01' + min_length: 10000 + query: genome_coverage>0.3 & missing_data<1000 + a/genome/6y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: 6Y + min_length: 10000 + query: genome_coverage>0.3 & missing_data<1000 + background: + include: config/include_a.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 300 + min_date: 12Y + max_date: 6Y + min_length: 10000 + query: genome_coverage>0.3 & missing_data<1000 & clade.str.startswith("A.D", na=False) + a/genome/3y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: 3Y + min_length: 10000 + query: genome_coverage>0.3 & missing_data<1000 + background: + include: config/include_a.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 300 + min_date: 12Y + max_date: 3Y + min_length: 10000 + query: genome_coverage>0.3 & missing_data<1000 & clade.str.startswith("A.D", na=False) + a/G/all-time: + samples: + sample: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: '1975-01-01' + min_length: 600 + query: G_coverage>0.3 & missing_data<1000 + a/G/6y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: 6Y + min_length: 600 + query: G_coverage>0.3 & missing_data<1000 + background: + include: config/include_a.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 300 + min_date: 12Y + max_date: 6Y + min_length: 600 + query: G_coverage>0.3 & missing_data<1000 & clade.str.startswith("A.D", na=False) + a/G/3y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: 3Y + min_length: 600 + query: G_coverage>0.3 & missing_data<1000 + background: + include: config/include_a.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 300 + min_date: 12Y + max_date: 3Y + min_length: 600 + query: G_coverage>0.3 & missing_data<1000 & clade.str.startswith("A.D", na=False) + a/F/all-time: + samples: + sample: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: '1975-01-01' + min_length: 1200 + query: F_coverage>0.3 & missing_data<1000 + a/F/6y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: 6Y + min_length: 1200 + query: F_coverage>0.3 & missing_data<1000 + background: + include: config/include_a.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 300 + min_date: 12Y + max_date: 6Y + min_length: 1200 + query: F_coverage>0.3 & missing_data<1000 & clade.str.startswith("A.D", na=False) + a/F/3y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: 3Y + min_length: 1200 + query: F_coverage>0.3 & missing_data<1000 + background: + include: config/include_a.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 300 + min_date: 12Y + max_date: 3Y + min_length: 1200 + query: F_coverage>0.3 & missing_data<1000 & clade.str.startswith("A.D", na=False) + a/F-antibody-escape/all-time: + samples: + sample: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 2000 + min_date: '1975-01-01' + min_length: 1200 + query: F_coverage>0.75 & missing_data<1000 + a/F-antibody-escape/6y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 2000 + min_date: 6Y + min_length: 1200 + query: F_coverage>0.75 & missing_data<1000 + background: + include: config/include_a.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 200 + min_date: 12Y + max_date: 6Y + min_length: 1200 + query: F_coverage>0.75 & missing_data<1000 & clade.str.startswith("A.D", na=False) + a/F-antibody-escape/3y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 2000 + min_date: 3Y + min_length: 1200 + query: F_coverage>0.75 & missing_data<1000 + background: + include: config/include_a.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 200 + min_date: 12Y + max_date: 3Y + min_length: 1200 + query: F_coverage>0.75 & missing_data<1000 & clade.str.startswith("A.D", na=False) + b/genome/all-time: + samples: + sample: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: '1975-01-01' + min_length: 10000 + query: genome_coverage>0.3 & missing_data<1000 + b/genome/6y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: 6Y + min_length: 10000 + query: genome_coverage>0.3 & missing_data<1000 + background: + include: config/include_b.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 300 + min_date: 12Y + max_date: 6Y + min_length: 10000 + query: genome_coverage>0.3 & missing_data<1000 & clade.str.startswith("B.D", na=False) + b/genome/3y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: 3Y + min_length: 10000 + query: genome_coverage>0.3 & missing_data<1000 + background: + include: config/include_b.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 300 + min_date: 12Y + max_date: 3Y + min_length: 10000 + query: genome_coverage>0.3 & missing_data<1000 & clade.str.startswith("B.D", na=False) + b/G/all-time: + samples: + sample: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: '1975-01-01' + min_length: 600 + query: G_coverage>0.3 & missing_data<1000 + b/G/6y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: 6Y + min_length: 600 + query: G_coverage>0.3 & missing_data<1000 + background: + include: config/include_b.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 300 + min_date: 12Y + max_date: 6Y + min_length: 600 + query: G_coverage>0.3 & missing_data<1000 & clade.str.startswith("B.D", na=False) + b/G/3y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: 3Y + min_length: 600 + query: G_coverage>0.3 & missing_data<1000 + background: + include: config/include_b.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 300 + min_date: 12Y + max_date: 3Y + min_length: 600 + query: G_coverage>0.3 & missing_data<1000 & clade.str.startswith("B.D", na=False) + b/F/all-time: + samples: + sample: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: '1975-01-01' + min_length: 1200 + query: F_coverage>0.3 & missing_data<1000 + b/F/6y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: 6Y + min_length: 1200 + query: F_coverage>0.3 & missing_data<1000 + background: + include: config/include_b.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 300 + min_date: 12Y + max_date: 6Y + min_length: 1200 + query: F_coverage>0.3 & missing_data<1000 & clade.str.startswith("B.D", na=False) + b/F/3y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 3000 + min_date: 3Y + min_length: 1200 + query: F_coverage>0.3 & missing_data<1000 + background: + include: config/include_b.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 300 + min_date: 12Y + max_date: 3Y + min_length: 1200 + query: F_coverage>0.3 & missing_data<1000 & clade.str.startswith("B.D", na=False) + b/F-antibody-escape/all-time: + samples: + sample: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 2000 + min_date: '1975-01-01' + min_length: 1200 + query: F_coverage>0.75 & missing_data<1000 + b/F-antibody-escape/6y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 2000 + min_date: 6Y + min_length: 1200 + query: F_coverage>0.75 & missing_data<1000 + background: + include: config/include_b.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 200 + min_date: 12Y + max_date: 6Y + min_length: 1200 + query: F_coverage>0.75 & missing_data<1000 & clade.str.startswith("B.D", na=False) + b/F-antibody-escape/3y: + samples: + recent: + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + group_by: + - year + - country + max_sequences: 2000 + min_date: 3Y + min_length: 1200 + query: F_coverage>0.75 & missing_data<1000 + background: + include: config/include_b.txt + exclude: config/outliers_ppx.txt + exclude_where: + - qc.overallStatus=bad + - qc.overallStatus=mediocre + group_by: + - year + - country + max_sequences: 200 + min_date: 12Y + max_date: 3Y + min_length: 1200 + query: F_coverage>0.75 & missing_data<1000 & clade.str.startswith("B.D", na=False) files: auspice_config: "config/auspice_config.json" auspice_config_additional_colorings: "config/auspice_config_additional_colorings.json" diff --git a/workflow/snakemake_rules/config.smk b/workflow/snakemake_rules/config.smk index ab284ff..1b537b6 100644 --- a/workflow/snakemake_rules/config.smk +++ b/workflow/snakemake_rules/config.smk @@ -14,6 +14,7 @@ def main(): "results/run_config.yaml", Path(workflow.basedir) / "config.schema.yaml" ) + write_subsample_config() # TODO: move this to nextstrain/shared @@ -39,6 +40,21 @@ def dump_and_validate(dump_path, schema_path): except ValidateError as e: raise InvalidConfigError(str(e)) from e + +def write_subsample_config(): + # TODO: Support custom build names in the workflow and infer from + # config["builds"]. + for a_or_b in ["a", "b"]: + for build_name in ["genome", "G", "F", "F-antibody-escape"]: + for resolution in ["all-time", "6y", "3y"]: + build = f"{a_or_b}/{build_name}/{resolution}" + if "custom_subsample" in config: + section = ["custom_subsample", build] + else: + section = ["subsample", build] + write_config(f"results/{build}/subsample_config.yaml", section=section) + + try: main() except InvalidConfigError as e: diff --git a/workflow/snakemake_rules/core.smk b/workflow/snakemake_rules/core.smk index 750f117..b528426 100644 --- a/workflow/snakemake_rules/core.smk +++ b/workflow/snakemake_rules/core.smk @@ -3,7 +3,7 @@ This part of the workflow expects input files sequences = "data/sequences.fasta" metadata = "data/metadata.tsv" """ - +from augur.subsample import get_referenced_files rule index_sequences: @@ -56,121 +56,38 @@ rule newreference: """ -rule filter_recent: - """ - filtering sequences - """ +rule subsample: input: sequences="results/{a_or_b}/sequences.fasta", metadata="results/{a_or_b}/metadata.tsv", sequence_index=rules.index_sequences.output, - exclude=config["exclude"], + config="results/{a_or_b}/{build_name}/{resolution}/subsample_config.yaml", + referenced_files=lambda w: get_referenced_files(f"results/{w.a_or_b}/{w.build_name}/{w.resolution}/subsample_config.yaml"), output: - sequences=build_dir - + "/{a_or_b}/{build_name}/{resolution}/filtered_recent.fasta", + sequences=build_dir + "/{a_or_b}/{build_name}/{resolution}/subsampled.fasta", log: - "logs/filter_recent_{a_or_b}_{build_name}_{resolution}.txt" + "logs/subsample_{a_or_b}_{build_name}_{resolution}.txt", benchmark: - "benchmarks/filter_recent_{a_or_b}_{build_name}_{resolution}.txt" + "benchmarks/subsample_{a_or_b}_{build_name}_{resolution}.txt", params: - group_by=config["filter"]["group_by"], - min_coverage=lambda w: f'{w.build_name.split("-")[0]}_coverage>{config["filter"]["min_coverage"][w.build_name]}', - min_length=lambda w: config["filter"]["min_length"][w.build_name], - subsample_max_sequences=lambda w: config["filter"][ - "subsample_max_sequences" - ][w.build_name], strain_id=config["strain_id_field"], - min_date=lambda w: config["filter"]["resolutions"][w.resolution]["min_date"], - exclude_where=config["filter"]["exclude_where"]["recent"], - missing_data_threshold=config["filter"]["missing_data_threshold"], shell: r""" exec &> >(tee {log:q}) - augur filter \ + augur subsample \ --sequences {input.sequences} \ --sequence-index {input.sequence_index} \ --metadata {input.metadata} \ --metadata-id-columns {params.strain_id} \ - --exclude {input.exclude} \ - --exclude-where {params.exclude_where:q} \ - --min-date {params.min_date} \ - --min-length {params.min_length} \ - --output {output.sequences} \ - --group-by {params.group_by} \ - --subsample-max-sequences {params.subsample_max_sequences} \ - --query '({params.min_coverage}) & missing_data<{params.missing_data_threshold}' - """ - - -rule filter_background: - """ - filtering sequences - """ - input: - sequences="results/{a_or_b}/sequences.fasta", - metadata="results/{a_or_b}/metadata.tsv", - sequence_index=rules.index_sequences.output, - include="config/include_{a_or_b}.txt", - exclude=config["exclude"], - output: - sequences=build_dir - + "/{a_or_b}/{build_name}/{resolution}/filtered_background.fasta", - metadata=build_dir - + "/{a_or_b}/{build_name}/{resolution}/filtered_background_metadata.tsv", - log: - "logs/filter_background_{a_or_b}_{build_name}_{resolution}.txt" - benchmark: - "benchmarks/filter_background_{a_or_b}_{build_name}_{resolution}.txt" - params: - group_by=config["filter"]["group_by"], - min_coverage=lambda w: f'{w.build_name.split("-")[0]}_coverage>{config["filter"]["min_coverage"][w.build_name]}', - min_length=lambda w: config["filter"]["min_length"][w.build_name], - subsample_max_sequences=lambda w: int( - config["filter"]["subsample_max_sequences"][w.build_name], - ) - // 10, - strain_id=config["strain_id_field"], - max_date=lambda w: config["filter"]["resolutions"][w.resolution]["min_date"], - min_date=lambda w: config["filter"]["resolutions"][w.resolution][ - "background_min_date" - ], - exclude_where=config["filter"]["exclude_where"]["background"], - missing_data_threshold=config["filter"]["missing_data_threshold"], - clade_prefix=lambda w: f"{w.a_or_b.upper()}.D", - shell: - r""" - exec &> >(tee {log:q}) - - augur filter \ - --sequences {input.sequences} \ - --sequence-index {input.sequence_index} \ - --metadata {input.metadata} \ - --metadata-id-columns {params.strain_id} \ - --include {input.include} \ - --exclude {input.exclude} \ - --exclude-where {params.exclude_where:q} \ - --min-date {params.min_date} \ - --max-date {params.max_date} \ - --min-length {params.min_length} \ - --output-sequences {output.sequences} \ - --output-metadata {output.metadata} \ - --group-by {params.group_by} \ - --subsample-max-sequences {params.subsample_max_sequences} \ - --query '({params.min_coverage}) & missing_data<{params.missing_data_threshold} & clade.str.startswith("{params.clade_prefix}", na=False)' + --config {input.config} \ + --output-sequences {output.sequences} """ rule combine_samples: input: subsamples=lambda w: ( - ( - [ - rules.filter_recent.output.sequences, - rules.filter_background.output.sequences, - ] - if "background_min_date" in config["filter"]["resolutions"][w.resolution] - else [rules.filter_recent.output.sequences] - ) + [rules.subsample.output.sequences] # potentially add sequences sampled to include maximum escape sequences + ( [ From 14ed05dd0fd8593069296cf86f9bf266e31bcb73 Mon Sep 17 00:00:00 2001 From: Victor Lin Date: Wed, 5 Aug 2026 17:22:55 -0700 Subject: [PATCH 4/5] Generate default config YAML This makes it easier to make changes to the subsample config. One downside is that the generated file is less readable with a strict YAML style and no comments. Comments have been moved to the script, but ideally they'd live in a schema which is used to generate user-facing docs. --- .gitattributes | 5 + README.md | 9 +- config/configfile.yaml | 135 +++++++------- scripts/generate_default_config.py | 272 +++++++++++++++++++++++++++++ 4 files changed, 347 insertions(+), 74 deletions(-) create mode 100644 scripts/generate_default_config.py diff --git a/.gitattributes b/.gitattributes index 46610df..7c5d1ec 100644 --- a/.gitattributes +++ b/.gitattributes @@ -1,3 +1,8 @@ # Allow Git to decide if file is text or binary # Always use LF line endings even on Windows. * text=auto eol=lf + +# This is a large generated file that, while text, it is not useful to +# routinely show the diff of. A diff can be forced as needed, e.g. with `git +# diff --text`. +/config/configfile.yaml -diff diff --git a/README.md b/README.md index fe171cc..13bc190 100644 --- a/README.md +++ b/README.md @@ -46,7 +46,14 @@ Once you've run the build, you can view the results with: ## Configuration -The default configuration is in [`config/configfile.yaml`](./config/configfile.yaml). +The default configuration is generated by +[scripts/generate_default_config.py](./scripts/generate_default_config.py). It +can be run with Nextstrain CLI: + +```sh +nextstrain shell . -c 'python scripts/generate_default_config.py' +``` + The workflow is contained in the [Snakefile](Snakefile) with included [rules](workflow/snakemake_rules/). Each rule specifies its file inputs and outputs and pulls its parameters from the config. There is little redirection and each diff --git a/config/configfile.yaml b/config/configfile.yaml index 5939779..38d15c0 100644 --- a/config/configfile.yaml +++ b/config/configfile.yaml @@ -1,34 +1,31 @@ -conda_environment: "workflow/envs/nextstrain.yaml" - -genesforglycosylation: ["G", "F"] - -# if the build name has a "-" in it, then the pipeline assumes -# the gene is the part before the first "-", so "F-antibody-escape" -# is parsed as the gene "F". -builds_to_run: ["genome", "G", "F", "F-antibody-escape"] - -resolutions_to_run: ["all-time", "6y", "3y"] - -subtypes: ['a', 'b'] - -# Both files must have a {a_or_b} expandable field to be replaced by "a" or "b" -# depending on if they are specified in the `subtypes` param above +# [DO NOT EDIT] This file was generated by scripts/generate_default_config.py. +conda_environment: workflow/envs/nextstrain.yaml +genesforglycosylation: +- G +- F +builds_to_run: +- genome +- G +- F +- F-antibody-escape +resolutions_to_run: +- all-time +- 6y +- 3y +subtypes: +- a +- b inputs: - - name: ppx_open - metadata: "https://data.nextstrain.org/files/workflows/rsv/{a_or_b}/metadata.tsv.gz" - sequences: "https://data.nextstrain.org/files/workflows/rsv/{a_or_b}/sequences.fasta.xz" - - name: ppx_restricted - metadata: "https://data.nextstrain.org/files/workflows/rsv/{a_or_b}/metadata_restricted.tsv.gz" - sequences: "https://data.nextstrain.org/files/workflows/rsv/{a_or_b}/sequences_restricted.fasta.xz" - -exclude: "config/outliers_ppx.txt" - -description: "config/description.md" - -strain_id_field: "accession" -display_strain_field: "strain" - - +- name: ppx_open + metadata: https://data.nextstrain.org/files/workflows/rsv/{a_or_b}/metadata.tsv.gz + sequences: https://data.nextstrain.org/files/workflows/rsv/{a_or_b}/sequences.fasta.xz +- name: ppx_restricted + metadata: https://data.nextstrain.org/files/workflows/rsv/{a_or_b}/metadata_restricted.tsv.gz + sequences: https://data.nextstrain.org/files/workflows/rsv/{a_or_b}/sequences_restricted.fasta.xz +exclude: config/outliers_ppx.txt +description: config/description.md +strain_id_field: accession +display_strain_field: strain subsample: a/genome/all-time: samples: @@ -567,48 +564,41 @@ subsample: min_length: 1200 query: F_coverage>0.75 & missing_data<1000 & clade.str.startswith("B.D", na=False) files: - auspice_config: "config/auspice_config.json" - auspice_config_additional_colorings: "config/auspice_config_additional_colorings.json" - auspice_config_f_antibody_escape: "config/auspice_config_dms-defaults.json" - auspice_config_non-genome_builds: "config/auspice_config_non-genome.json" - + auspice_config: config/auspice_config.json + auspice_config_additional_colorings: config/auspice_config_additional_colorings.json + auspice_config_f_antibody_escape: config/auspice_config_dms-defaults.json + auspice_config_non-genome_builds: config/auspice_config_non-genome.json refine: - coalescent: "opt" - date_inference: "marginal" + coalescent: opt + date_inference: marginal clock_filter_iqd: 4 - divergence_units: "mutations-per-site" - + divergence_units: mutations-per-site ancestral: - inference: "joint" - + inference: joint cds: - F: "F" - G: "G" - genome: "F" - F-antibody-escape: "F" - + F: F + G: G + genome: F + F-antibody-escape: F traits: - columns: "country region" - + columns: country region frequencies: resolutions: all-time: - min_date: "1975-01-01" + min_date: '1975-01-01' 6y: min_date: 6Y 3y: min_date: 3Y - nextclade_attributes: a: - name: "RSV-A NextClade using real-time tree" - reference_name: "hRSV/A/England/397/2017" - accession: "EPI_ISL_412866" + name: RSV-A NextClade using real-time tree + reference_name: hRSV/A/England/397/2017 + accession: EPI_ISL_412866 b: - name: "RSV-B NextClade using real-time tree" - reference_name: "hRSV/B/Australia/VIC-RCH056/2019" - accession: "EPI_ISL_1653999" - + name: RSV-B NextClade using real-time tree + reference_name: hRSV/B/Australia/VIC-RCH056/2019 + accession: EPI_ISL_1653999 filter_for_pre_subsample_alignment: min_length: genome: 10000 @@ -622,24 +612,23 @@ filter_for_pre_subsample_alignment: F-antibody-escape: 0.75 resolutions: all-time: - min_date: "1975-01-01" + min_date: '1975-01-01' 6y: min_date: 6Y 3y: min_date: 3Y - -# configuration specific to the F deep mutational scanning antibody escape data f_dms_data: dms-data/all_antibodies.csv -f_dms_antibodies: # columns in `f_dms_data` with per-mutation escape - - Clesrovimab-Fab - - Clesrovimab-IgG - - Nirsevimab-Fab - - Nirsevimab-IgG -dms_only_positive_escape: true # for DMS escape values, set any values < 0 to 0 -enrich_antibody_escape: # additional filtering for antibody escape sequences added to tree - F-antibody-escape: # for this build, enrich by these criteria - nseqs_per_antibody_scoretype: 500 # add this many sequences for each antibody for total and max escape - group_by: [country, year] # group by these variables - max_identical_f_prot_muts: 2 # for each group, no more than this many w identical F protein mutations - max_identical_max_escape_mut: 6 # for each group, no more than this many w same top F escape mutation - +f_dms_antibodies: +- Clesrovimab-Fab +- Clesrovimab-IgG +- Nirsevimab-Fab +- Nirsevimab-IgG +dms_only_positive_escape: true +enrich_antibody_escape: + F-antibody-escape: + nseqs_per_antibody_scoretype: 500 + group_by: + - country + - year + max_identical_f_prot_muts: 2 + max_identical_max_escape_mut: 6 diff --git a/scripts/generate_default_config.py b/scripts/generate_default_config.py new file mode 100644 index 0000000..0b92e63 --- /dev/null +++ b/scripts/generate_default_config.py @@ -0,0 +1,272 @@ +"""Generate the default configfile for the RSV workflow.""" + +from pathlib import Path +import yaml + + +class NoAliasDumper(yaml.SafeDumper): + def ignore_aliases(self, data): + return True + + +def main(): + config = generate_config() + + path = Path(__file__).resolve().parent.parent / "config" / "configfile.yaml" + + with open(path, "w") as f: + print("# [DO NOT EDIT] This file was generated by scripts/generate_default_config.py.", file=f) + yaml.dump(config, f, Dumper=NoAliasDumper, sort_keys=False, width=1000) + + +def generate_config(): + return { + "conda_environment": "workflow/envs/nextstrain.yaml", + + "genesforglycosylation": ["G", "F"], + # if the build name has a "-" in it, then the pipeline assumes + # the gene is the part before the first "-", so "F-antibody-escape" + # is parsed as the gene "F". + "builds_to_run": ["genome", "G", "F", "F-antibody-escape"], + + "resolutions_to_run": ["all-time", "6y", "3y"], + + "subtypes": ["a", "b"], + + # Both files must have a {a_or_b} expandable field to be replaced by "a" or "b" + # depending on if they are specified in the `subtypes` param above + "inputs": [ + { + "name": "ppx_open", + "metadata": "https://data.nextstrain.org/files/workflows/rsv/{a_or_b}/metadata.tsv.gz", + "sequences": "https://data.nextstrain.org/files/workflows/rsv/{a_or_b}/sequences.fasta.xz", + }, + { + "name": "ppx_restricted", + "metadata": "https://data.nextstrain.org/files/workflows/rsv/{a_or_b}/metadata_restricted.tsv.gz", + "sequences": "https://data.nextstrain.org/files/workflows/rsv/{a_or_b}/sequences_restricted.fasta.xz", + }, + ], + + "exclude": "config/outliers_ppx.txt", + + "description": "config/description.md", + + "strain_id_field": "accession", + "display_strain_field": "strain", + + "subsample": generate_subsample_config(), + + "files": { + "auspice_config": "config/auspice_config.json", + "auspice_config_additional_colorings": "config/auspice_config_additional_colorings.json", + "auspice_config_f_antibody_escape": "config/auspice_config_dms-defaults.json", + "auspice_config_non-genome_builds": "config/auspice_config_non-genome.json", + }, + + "refine": { + "coalescent": "opt", + "date_inference": "marginal", + "clock_filter_iqd": 4, + "divergence_units": "mutations-per-site", + }, + + "ancestral": { + "inference": "joint", + }, + + "cds": { + "F": "F", + "G": "G", + "genome": "F", + "F-antibody-escape": "F", + }, + + "traits": { + "columns": "country region", + }, + + "frequencies": { + "resolutions": { + "all-time": {"min_date": "1975-01-01"}, + "6y": {"min_date": "6Y"}, + "3y": {"min_date": "3Y"}, + }, + }, + + "nextclade_attributes": { + "a": { + "name": "RSV-A NextClade using real-time tree", + "reference_name": "hRSV/A/England/397/2017", + "accession": "EPI_ISL_412866", + }, + "b": { + "name": "RSV-B NextClade using real-time tree", + "reference_name": "hRSV/B/Australia/VIC-RCH056/2019", + "accession": "EPI_ISL_1653999", + }, + }, + + "filter_for_pre_subsample_alignment": { + "min_length": { + "genome": 10000, + "G": 600, + "F": 1200, + "F-antibody-escape": 1200, + }, + "min_coverage": { + "genome": 0.3, + "G": 0.3, + "F": 0.3, + "F-antibody-escape": 0.75, + }, + "resolutions": { + "all-time": {"min_date": "1975-01-01"}, + "6y": {"min_date": "6Y"}, + "3y": {"min_date": "3Y"}, + }, + }, + + # configuration specific to the F deep mutational scanning antibody escape data + "f_dms_data": "dms-data/all_antibodies.csv", + "f_dms_antibodies": [ # columns in `f_dms_data` with per-mutation escape + "Clesrovimab-Fab", + "Clesrovimab-IgG", + "Nirsevimab-Fab", + "Nirsevimab-IgG", + ], + "dms_only_positive_escape": True, # for DMS escape values, set any values < 0 to 0 + "enrich_antibody_escape": { # additional filtering for antibody escape sequences added to tree + "F-antibody-escape": { # for this build, enrich by these criteria + "nseqs_per_antibody_scoretype": 500, # add this many sequences for each antibody for total and max escape + "group_by": ["country", "year"], # group by these variables + "max_identical_f_prot_muts": 2, # for each group, no more than this many w identical F protein mutations + "max_identical_max_escape_mut": 6, # for each group, no more than this many w same top F escape mutation + }, + }, + } + + +def generate_subsample_config(): + """Generate the expanded subsample configuration for all builds.""" + subtypes = ["a", "b"] + + builds = { + "genome": { + "gene_coverage": "genome_coverage", + "min_coverage": 0.3, + "min_length": 10000, + "recent_max_seqs": 3000, + "bg_max_seqs": 300, + }, + "G": { + "gene_coverage": "G_coverage", + "min_coverage": 0.3, + "min_length": 600, + "recent_max_seqs": 3000, + "bg_max_seqs": 300, + }, + "F": { + "gene_coverage": "F_coverage", + "min_coverage": 0.3, + "min_length": 1200, + "recent_max_seqs": 3000, + "bg_max_seqs": 300, + }, + "F-antibody-escape": { + "gene_coverage": "F_coverage", + "min_coverage": 0.75, + "min_length": 1200, + "recent_max_seqs": 2000, + "bg_max_seqs": 200, + }, + } + + subsample = {} + for a_or_b in subtypes: + clade_prefix = f"{a_or_b.upper()}.D" + include_file = f"config/include_{a_or_b}.txt" + + for build_name, build_info in builds.items(): + coverage_var = build_info["gene_coverage"] + min_cov = build_info["min_coverage"] + min_len = build_info["min_length"] + recent_seqs = build_info["recent_max_seqs"] + bg_seqs = build_info["bg_max_seqs"] + + recent_query = f"{coverage_var}>{min_cov} & missing_data<1000" + bg_query = f"{coverage_var}>{min_cov} & missing_data<1000 & clade.str.startswith(\"{clade_prefix}\", na=False)" + + resolutions = { + "all-time": { + "sample": { + "exclude": "config/outliers_ppx.txt", + "exclude_where": ["qc.overallStatus=bad"], + "group_by": ["year", "country"], + "max_sequences": recent_seqs, + "min_date": "1975-01-01", + "min_length": min_len, + "query": recent_query, + }, + }, + "6y": { + "recent": { + "exclude": "config/outliers_ppx.txt", + "exclude_where": ["qc.overallStatus=bad"], + "group_by": ["year", "country"], + "max_sequences": recent_seqs, + "min_date": "6Y", + "min_length": min_len, + "query": recent_query, + }, + "background": { + "include": include_file, + "exclude": "config/outliers_ppx.txt", + "exclude_where": [ + "qc.overallStatus=bad", + "qc.overallStatus=mediocre", + ], + "group_by": ["year", "country"], + "max_sequences": bg_seqs, + "min_date": "12Y", + "max_date": "6Y", + "min_length": min_len, + "query": bg_query, + }, + }, + "3y": { + "recent": { + "exclude": "config/outliers_ppx.txt", + "exclude_where": ["qc.overallStatus=bad"], + "group_by": ["year", "country"], + "max_sequences": recent_seqs, + "min_date": "3Y", + "min_length": min_len, + "query": recent_query, + }, + "background": { + "include": include_file, + "exclude": "config/outliers_ppx.txt", + "exclude_where": [ + "qc.overallStatus=bad", + "qc.overallStatus=mediocre", + ], + "group_by": ["year", "country"], + "max_sequences": bg_seqs, + "min_date": "12Y", + "max_date": "3Y", + "min_length": min_len, + "query": bg_query, + }, + }, + } + + for resolution, samples in resolutions.items(): + key = f"{a_or_b}/{build_name}/{resolution}" + subsample[key] = {"samples": samples} + + return subsample + + +if __name__ == "__main__": + main() From 12efdc88f85581d1d51670042e74a9cbe280631c Mon Sep 17 00:00:00 2001 From: Victor Lin Date: Mon, 10 Aug 2026 15:18:13 -0700 Subject: [PATCH 5/5] Rename rule to filter_for_f_antibody_escape The new name makes it more obvious that this rule is only used when build_name=F-antibody-escape. --- CHANGELOG.md | 2 +- config.schema.yaml | 2 +- config/configfile.yaml | 2 +- scripts/generate_default_config.py | 2 +- workflow/snakemake_rules/core.smk | 14 +++++++------- 5 files changed, 11 insertions(+), 11 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 9d7bc6d..f72125a 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -4,6 +4,6 @@ We use this CHANGELOG to document breaking changes, new features, bug fixes, and ## 2026 -* TBD: The `filter` section in phylogenetic workflow configuration has been replaced by `subsample`/`custom_subsample` for subsampling, and `filter_for_pre_subsample_alignment` for initial quality filtering. **This is a breaking change**. +* TBD: The `filter` section in phylogenetic workflow configuration has been replaced by `subsample`/`custom_subsample` for subsampling, and `filter_for_f_antibody_escape` for initial quality filtering. **This is a breaking change**. * NOTE: The workflow does not yet support proximal samples. * TBD: Phylogenetic workflow configuration is now validated against a strict schema. The workflow will error if your configuration has extraneous entries that were previously ignored. diff --git a/config.schema.yaml b/config.schema.yaml index 84f0d63..80ea2a4 100644 --- a/config.schema.yaml +++ b/config.schema.yaml @@ -168,7 +168,7 @@ properties: type: string accession: type: string - filter_for_pre_subsample_alignment: + filter_for_f_antibody_escape: type: object additionalProperties: false properties: diff --git a/config/configfile.yaml b/config/configfile.yaml index 38d15c0..f9430e4 100644 --- a/config/configfile.yaml +++ b/config/configfile.yaml @@ -599,7 +599,7 @@ nextclade_attributes: name: RSV-B NextClade using real-time tree reference_name: hRSV/B/Australia/VIC-RCH056/2019 accession: EPI_ISL_1653999 -filter_for_pre_subsample_alignment: +filter_for_f_antibody_escape: min_length: genome: 10000 G: 600 diff --git a/scripts/generate_default_config.py b/scripts/generate_default_config.py index 0b92e63..dabdf24 100644 --- a/scripts/generate_default_config.py +++ b/scripts/generate_default_config.py @@ -107,7 +107,7 @@ def generate_config(): }, }, - "filter_for_pre_subsample_alignment": { + "filter_for_f_antibody_escape": { "min_length": { "genome": 10000, "G": 600, diff --git a/workflow/snakemake_rules/core.smk b/workflow/snakemake_rules/core.smk index b528426..e4a9af4 100644 --- a/workflow/snakemake_rules/core.smk +++ b/workflow/snakemake_rules/core.smk @@ -137,7 +137,7 @@ rule get_nextclade_dataset: """ -rule filter_for_pre_subsample_alignment: +rule filter_for_f_antibody_escape: """ Do the quality filtering applied to each sequence set before subsampling """ @@ -148,14 +148,14 @@ rule filter_for_pre_subsample_alignment: output: sequences=build_dir + "/{a_or_b}/{build_name}/{resolution}/pre_subsample/filtered_for_alignment.fasta", log: - "logs/filter_for_pre_subsample_alignment_{a_or_b}_{build_name}_{resolution}.txt" + "logs/filter_for_f_antibody_escape_{a_or_b}_{build_name}_{resolution}.txt" benchmark: - "benchmarks/filter_for_pre_subsample_alignment_{a_or_b}_{build_name}_{resolution}.txt" + "benchmarks/filter_for_f_antibody_escape_{a_or_b}_{build_name}_{resolution}.txt" params: - min_coverage=lambda w: f'{w.build_name.split("-")[0]}_coverage>{config["filter_for_pre_subsample_alignment"]["min_coverage"][w.build_name]}', - min_length=lambda w: config["filter_for_pre_subsample_alignment"]["min_length"][w.build_name], + min_coverage=lambda w: f'{w.build_name.split("-")[0]}_coverage>{config["filter_for_f_antibody_escape"]["min_coverage"][w.build_name]}', + min_length=lambda w: config["filter_for_f_antibody_escape"]["min_length"][w.build_name], strain_id=config["strain_id_field"], - min_date=lambda w: config["filter_for_pre_subsample_alignment"]["resolutions"][w.resolution]["min_date"], + min_date=lambda w: config["filter_for_f_antibody_escape"]["resolutions"][w.resolution]["min_date"], shell: r""" exec &> >(tee {log:q}) @@ -178,7 +178,7 @@ rule align_pre_subsample_sequences: Aligning all pre-subsampled quality-filtered sequences """ input: - sequences=rules.filter_for_pre_subsample_alignment.output.sequences, + sequences=rules.filter_for_f_antibody_escape.output.sequences, dataset=rules.get_nextclade_dataset.output.dataset, output: alignment=build_dir + "/{a_or_b}/{build_name}/{resolution}/pre_subsample/sequences.aligned.fasta",