Update phylo workflow

phylogenetic/config/hku1/genemap.gff

@@ -1,13 +1,13 @@
 ##sequence-region NC_006577.2 1 29926
 NC_006577.2	feature	source	1	29926	.	+	.	gene=nuc
-NC_006577.2	feature	gene	206	13600	.	+	.	gene=ORF1a
-NC_006577.2	feature	gene	13600	21753	.	+	.	gene=ORF1b
-NC_006577.2	feature	gene	21773	22933	.	+	.	gene=HE
-NC_006577.2	feature	gene	22942	27012	.	+	.	gene=Spike
-NC_006577.2	feature	gene	22978	25221	.	+	.	gene=S1
-NC_006577.2	feature	gene	27051	27380	.	+	.	gene=S2
-NC_006577.2	feature	gene	27051	27380	.	+	.	gene=ORF4
-NC_006577.2	feature	gene	27373	27621	.	+	.	gene=E
-NC_006577.2	feature	gene	27633	28304	.	+	.	gene=M
-NC_006577.2	feature	gene	28320	29645	.	+	.	gene=N
-NC_006577.2	feature	gene	28342	28959	.	+	.	gene=N2
+NC_006577.2	feature	gene	206	13600	.	+	.	gene_name=ORF1a
+NC_006577.2	feature	gene	13600	21753	.	+	.	gene_name=ORF1b
+NC_006577.2	feature	gene	21773	22933	.	+	.	gene_name=HE
+NC_006577.2	feature	gene	22942	27012	.	+	.	gene_name=Spike
+NC_006577.2	feature	gene	22978	25221	.	+	.	gene_name=S1
+NC_006577.2	feature	gene	27051	27380	.	+	.	gene_name=S2
+NC_006577.2	feature	gene	27051	27380	.	+	.	gene_name=ORF4
+NC_006577.2	feature	gene	27373	27621	.	+	.	gene_name=E
+NC_006577.2	feature	gene	27633	28304	.	+	.	gene_name=M
+NC_006577.2	feature	gene	28320	29645	.	+	.	gene_name=N
+NC_006577.2	feature	gene	28342	28959	.	+	.	gene_name=N2

phylogenetic/rules/prepare_sequences.smk

@@ -10,6 +10,7 @@ This part of the workflow expects a metadata TSV and FASTA file as inputs
 and will produce an aligned FASTA file of subsampled sequences as an output.
 """
 
+
 rule filter:
     message:
         """
@@ -18,15 +19,15 @@ rule filter:
           - minimum genome length of {params.min_length}
         """
     input:
-        sequences = lambda wildcards:config[wildcards.virus]["prepare_sequences"]["sequences"],
-        metadata = lambda wildcards:config[wildcards.virus]["metadata"],
-        exclude = "config/{virus}/dropped_strains.txt"
+        sequences=lambda wildcards: config[wildcards.virus]["prepare_sequences"]["sequences"],
+        metadata=lambda wildcards: config[wildcards.virus]["metadata"],
+        exclude="config/{virus}/dropped_strains.txt",
     output:
-        sequences = "results/{virus}/filtered.fasta"
+        sequences="results/{virus}/filtered.fasta",
     params:
-        group_by = lambda wildcards:config[wildcards.virus]["prepare_sequences"]["group_by"],
-        sequences_per_group = lambda wildcards:config[wildcards.virus]["prepare_sequences"]["sequences_per_group"],
-        min_length = lambda wildcards:config[wildcards.virus]["prepare_sequences"]["min_length"]
+        group_by=lambda wildcards: config[wildcards.virus]["prepare_sequences"]["group_by"],
+        sequences_per_group=lambda wildcards: config[wildcards.virus]["prepare_sequences"]["sequences_per_group"],
+        min_length=lambda wildcards: config[wildcards.virus]["prepare_sequences"]["min_length"],
     shell:
         """
         augur filter \
@@ -39,26 +40,27 @@ rule filter:
             --sequences-per-group {params.sequences_per_group} \
             --min-length {params.min_length}
         """
+
+
 rule align:
     message:
         """
         Aligning sequences to {input.reference}
           - filling gaps with N
         """
     input:
-        sequences = rules.filter.output.sequences,
-        reference= lambda wildcards:config[wildcards.virus]["reference"],
-        genemap= lambda wildcards:config[wildcards.virus]["genemap"],
+        sequences=rules.filter.output.sequences,
+        reference=lambda wildcards: config[wildcards.virus]["reference"],
+        genemap=lambda wildcards: config[wildcards.virus]["genemap"],
     output:
-        alignment= "results/{virus}/aligned.fasta",
-        insertions= "results/{virus}/insertions.fasta"
+        alignment="results/{virus}/aligned.fasta",
+        insertions="results/{virus}/insertions.tsv",
     shell:
         """
-        nextalign run \
-            --reference {input.reference} \
-            --genemap {input.genemap} \
-            --retry-reverse-complement \
+        nextclade run \
+            --input-ref {input.reference} \
+            --input-annotation {input.genemap} \
             --output-fasta - \
-            --output-insertions {output.insertions} \
+            --output-tsv {output.insertions} \
             {input.sequences} | seqkit seq -i > {output.alignment}
         """