From 41991ff7fa03cf92e2a6c657687f74f65a6cb425 Mon Sep 17 00:00:00 2001 From: Victor Lin <13424970+victorlin@users.noreply.github.com> Date: Fri, 6 Dec 2024 16:15:11 -0800 Subject: [PATCH 1/3] Add columns in a separate rule Narrows the scope of the append_usvi rule to just merging data files. --- phylogenetic/rules/merge_sequences_usvi.smk | 39 +++++++++++++-------- 1 file changed, 25 insertions(+), 14 deletions(-) diff --git a/phylogenetic/rules/merge_sequences_usvi.smk b/phylogenetic/rules/merge_sequences_usvi.smk index 2dd36dd..2067d3b 100644 --- a/phylogenetic/rules/merge_sequences_usvi.smk +++ b/phylogenetic/rules/merge_sequences_usvi.smk @@ -21,26 +21,20 @@ This part of the workflow usually includes the following steps: """ -rule append_usvi: - """Appending USVI sequences +rule add_metadata_columns: + """Add columns to metadata Notable columns: - - accession: Either the GenBank accession or USVI accession. - - genbank_accession: GenBank accession for Auspice to generate a URL to the NCBI GenBank record. Empty for USVI sequences. - - url: URL used in Auspice, to either link to the USVI github repo (https://github.com/blab/zika-usvi/) or link to the NCBI GenBank record ('https://www.ncbi.nlm.nih.gov/nuccore/*') + - genbank_accession: GenBank accession for Auspice to generate a URL to the NCBI GenBank record. + - [NEW] accession: The GenBank accession. Added to go alongside USVI accession. + - [NEW] url: URL linking to the NCBI GenBank record ('https://www.ncbi.nlm.nih.gov/nuccore/*'). Added to go alongside USVI url. """ input: - sequences = "data/sequences.fasta", - metadata = "data/metadata.tsv", - usvi_sequences = "data/sequences_usvi.fasta", - usvi_metadata = "data/metadata_usvi.tsv" + metadata = "data/metadata.tsv" output: - sequences = "data/sequences_all.fasta", - metadata = "data/metadata_all.tsv" + metadata = "data/metadata_modified.tsv" shell: """ - cat {input.sequences} {input.usvi_sequences} > {output.sequences} - csvtk mutate2 -tl \ -n url \ -e '"https://www.ncbi.nlm.nih.gov/nuccore/" + $genbank_accession' \ @@ -48,7 +42,24 @@ rule append_usvi: | csvtk mutate2 -tl \ -n accession \ -e '$genbank_accession' \ - | csvtk concat -tl - {input.usvi_metadata} \ + > {output.metadata} + """ + +rule append_usvi: + """Appending USVI sequences""" + input: + sequences = "data/sequences.fasta", + metadata = "data/metadata_modified.tsv", + usvi_sequences = "data/sequences_usvi.fasta", + usvi_metadata = "data/metadata_usvi.tsv" + output: + sequences = "data/sequences_all.fasta", + metadata = "data/metadata_all.tsv" + shell: + """ + cat {input.sequences} {input.usvi_sequences} > {output.sequences} + + csvtk concat -tl {input.metadata} {input.usvi_metadata} \ | tsv-select -H -f accession --rest last \ > {output.metadata} """ \ No newline at end of file From 14b89631ac357c2353707c32112f790fbd37aec5 Mon Sep 17 00:00:00 2001 From: Victor Lin <13424970+victorlin@users.noreply.github.com> Date: Fri, 6 Dec 2024 16:39:04 -0800 Subject: [PATCH 2/3] Replace csvtk concat with augur merge Use our own tooling designed for this type of operation. There is a slight change in behavior: augur merge adds quotes around values with spaces whereas csvtk doesn't, but I think that's ok. --- phylogenetic/rules/merge_sequences_usvi.smk | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/phylogenetic/rules/merge_sequences_usvi.smk b/phylogenetic/rules/merge_sequences_usvi.smk index 2067d3b..fe8296e 100644 --- a/phylogenetic/rules/merge_sequences_usvi.smk +++ b/phylogenetic/rules/merge_sequences_usvi.smk @@ -59,7 +59,8 @@ rule append_usvi: """ cat {input.sequences} {input.usvi_sequences} > {output.sequences} - csvtk concat -tl {input.metadata} {input.usvi_metadata} \ - | tsv-select -H -f accession --rest last \ - > {output.metadata} + augur merge \ + --metadata ingest={input.metadata} usvi={input.usvi_metadata} \ + --metadata-id-columns accession \ + --output-metadata {output.metadata} """ \ No newline at end of file From b92c80ce34c4ba2c0c27dbe9ee2df1d44642d0ef Mon Sep 17 00:00:00 2001 From: Victor Lin <13424970+victorlin@users.noreply.github.com> Date: Fri, 6 Dec 2024 17:16:43 -0800 Subject: [PATCH 3/3] Replace cat with seqkit rmdup This better reflects the duplicate handling used by augur merge. In case of duplicates, augur merge will keep the last one from the input list. seqkit rmdup keeps the first one from the input list so the order is reversed. There are some behavior changes, none of which should have any impact on downstream usage: 1. The order of sequences in the file is reversed on a file level. 2. seqkit's default behavior wraps lines at 60 characters. --- phylogenetic/rules/merge_sequences_usvi.smk | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/phylogenetic/rules/merge_sequences_usvi.smk b/phylogenetic/rules/merge_sequences_usvi.smk index fe8296e..ffc7a50 100644 --- a/phylogenetic/rules/merge_sequences_usvi.smk +++ b/phylogenetic/rules/merge_sequences_usvi.smk @@ -57,7 +57,7 @@ rule append_usvi: metadata = "data/metadata_all.tsv" shell: """ - cat {input.sequences} {input.usvi_sequences} > {output.sequences} + seqkit rmdup {input.usvi_sequences} {input.sequences} > {output.sequences} augur merge \ --metadata ingest={input.metadata} usvi={input.usvi_metadata} \