Phylogenetic: Add frequencies panel #89

Author: j23414

CHANGELOG.md

@@ -4,6 +4,7 @@ We use this CHANGELOG to document breaking changes, new features, bug fixes, and
 
 ## 2025
 
+* 24 March 2025: Add frequencies panel @j23414
 * 21 March 2025: Migrated WA specific config to NW-PaGe [PR#90] @j23414, @DOH-LMT2303, @DOH-PNT0303
 * 18 March 2025: Rename "global" to "all-lineages, update various auspice_config settings @j23414
 * 14 March 2025: Pushed from staging to live @j23414

phylogenetic/Snakefile

@@ -21,7 +21,8 @@ wildcard_constraints:
 
 rule all:
     input:
-        auspice = expand("auspice/WNV_{build}.json", build=builds)
+        auspice = expand("auspice/WNV_{build}.json", build=builds),
+        tip_frequencies_json= expand("auspice/WNV_{build}_tip-frequencies.json", build=builds)
 
 # These rules are imported in the order that they are expected to run.
 # Each Snakefile will have documented inputs and outputs that should be kept as

phylogenetic/defaults/all-lineages/auspice_config.json

@@ -40,6 +40,12 @@
     "host_genus",
     "host_type"
   ],
+  "panels": [
+    "tree",
+    "map",
+    "entropy",
+    "frequencies"
+  ],
   "display_defaults": {
     "color_by": "lineage",
     "map_triplicate": true,

phylogenetic/defaults/all-lineages/reference.gb

@@ -1,6 +1,6 @@
-LOCUS       AF260968               11029 bp    RNA     linear   VRL 27-AUG-2000
+LOCUS       AF260968_REF               11029 bp    RNA     linear   VRL 27-AUG-2000
 DEFINITION  West Nile virus strain Eg101, complete genome.
-ACCESSION   AF260968
+ACCESSION   AF260968_REF
 VERSION     AF260968.1
 KEYWORDS    .
 SOURCE      West Nile virus (WNV)

phylogenetic/defaults/config.yaml

@@ -78,6 +78,14 @@ traits:
     'lineage',
   ]
 
+tip_frequencies:
+  # 2000 since there is an increase in WNV at that time
+  min_date: "2000-01-01"
+  max_date: "6M"
+  # Quarterly narrow_bandwidth or every 3 months (3 /12.0 = 0.25)
+  narrow_bandwidth: 0.25
+  proportion_wide: 0.0
+
 export:
   description: "defaults/description.md"
   auspice_config: "defaults/all-lineages/auspice_config.json"

phylogenetic/defaults/include.txt

@@ -1,4 +1,4 @@
-# AF260968 # Egypt 1951 will be used as augur align reference
+AF260968 # Egypt 1951 all-lineages reference
 NC_001563 # Lineage 2 reference
 NC_009942 # Lineage 1 reference
 HM051416 # Isreal 1953

phylogenetic/defaults/lineage-1A/auspice_config.json

@@ -40,6 +40,12 @@
     "host_genus",
     "host_type"
   ],
+  "panels": [
+    "tree",
+    "map",
+    "entropy",
+    "frequencies"
+  ],
   "display_defaults": {
     "color_by": "lineage",
     "map_triplicate": true,

phylogenetic/defaults/lineage-2/auspice_config.json

@@ -40,6 +40,12 @@
     "host_genus",
     "host_type"
   ],
+  "panels": [
+    "tree",
+    "map",
+    "entropy",
+    "frequencies"
+  ],
   "display_defaults": {
     "color_by": "lineage",
     "map_triplicate": true,

phylogenetic/rules/export.smk

@@ -54,3 +54,32 @@ rule export:
             --include-root-sequence-inline \
             --output {output.auspice} 2>&1 | tee {log}
         """
+
+rule tip_frequencies:
+    """
+    Estimating KDE frequencies for tips
+    """
+    input:
+        tree = "results/{build}/tree.nwk",
+        metadata = "results/{build}/metadata_filtered.tsv",
+    output:
+        tip_freq = "auspice/WNV_{build}_tip-frequencies.json"
+    params:
+        strain_id = config["strain_id_field"],
+        min_date = config["tip_frequencies"]["min_date"],
+        max_date = config["tip_frequencies"]["max_date"],
+        narrow_bandwidth = config["tip_frequencies"]["narrow_bandwidth"],
+        proportion_wide = config["tip_frequencies"]["proportion_wide"]
+    shell:
+        r"""
+        augur frequencies \
+            --method kde \
+            --tree {input.tree} \
+            --metadata {input.metadata} \
+            --metadata-id-columns {params.strain_id} \
+            --min-date {params.min_date} \
+            --max-date {params.max_date} \
+            --narrow-bandwidth {params.narrow_bandwidth} \
+            --proportion-wide {params.proportion_wide} \
+            --output {output.tip_freq}
+        """

phylogenetic/rules/prepare_sequences.smk

@@ -39,5 +39,6 @@ rule align:
             --output {output.alignment} \
             --fill-gaps \
             --reference-sequence {input.reference} \
+            --remove-reference \
             --nthreads {threads:q} 2>&1 | tee {log}
         """