From a4089e71ea0e704c1c20700e20132b3592817f34 Mon Sep 17 00:00:00 2001 From: Jover Lee Date: Mon, 16 Dec 2024 16:33:44 -0800 Subject: [PATCH 1/3] Switch accession links to UCSC file Adds a rule to fetch GISAID and GenBank accession links from UCSC and uses the output in the `transform_genbank` and `transform_gisaid` rules. The UCSC file is modified to keep the headers and format exactly the same as the current `source-data/accessions.tsv.gz` file so that it can be directly replaced. --- bin/fetch-accession-links | 9 +++++++++ workflow/snakemake_rules/curate.smk | 21 +++++++++++++++++++-- 2 files changed, 28 insertions(+), 2 deletions(-) create mode 100755 bin/fetch-accession-links diff --git a/bin/fetch-accession-links b/bin/fetch-accession-links new file mode 100755 index 00000000..7378ed4d --- /dev/null +++ b/bin/fetch-accession-links @@ -0,0 +1,9 @@ +#!/bin/bash +set -euo pipefail + +curl "https://hgwdev.gi.ucsc.edu/~angie/epiToPublicAndDate.latest" \ + --fail --silent --show-error \ + --header 'User-Agent: https://github.com/nextstrain/ncov-ingest (hello@nextstrain.org)' \ + | csvtk -t add-header --names gisaid_epi_isl,genbank_accession,strain,date \ + | csvtk -t cut --fields genbank_accession,gisaid_epi_isl \ + | gzip -c diff --git a/workflow/snakemake_rules/curate.smk b/workflow/snakemake_rules/curate.smk index 75b65f48..6fee03d2 100644 --- a/workflow/snakemake_rules/curate.smk +++ b/workflow/snakemake_rules/curate.smk @@ -23,6 +23,19 @@ Produces different output files for GISAID vs GenBank: """ +rule fetch_accession_links: + """ + Fetch the accession links between GISAID and GenBank + """ + output: + accessions="data/accessions.tsv.gz", + retries: 5 + shell: + """ + ./bin/fetch-accession-links > {output.accessions:q} + """ + + rule transform_rki_data: input: ndjson="data/rki.ndjson", @@ -60,7 +73,8 @@ rule transform_genbank_data: biosample = "data/genbank/biosample.tsv", ndjson = "data/genbank.ndjson", cog_uk_accessions = "data/cog_uk_accessions.tsv", - cog_uk_metadata = "data/cog_uk_metadata.csv.gz" + cog_uk_metadata = "data/cog_uk_metadata.csv.gz", + accessions = "data/accessions.tsv.gz", output: fasta = "data/genbank_sequences.fasta", metadata = "data/genbank_metadata_transformed.tsv", @@ -75,6 +89,7 @@ rule transform_genbank_data: --duplicate-biosample {output.duplicate_biosample} \ --cog-uk-accessions {input.cog_uk_accessions} \ --cog-uk-metadata {input.cog_uk_metadata} \ + --accessions {input.accessions} \ --output-metadata {output.metadata} \ --output-fasta {output.fasta} > {output.flagged_annotations} """ @@ -105,7 +120,8 @@ rule merge_open_data: rule transform_gisaid_data: input: - ndjson = "data/gisaid.ndjson" + ndjson = "data/gisaid.ndjson", + accessions = "data/accessions.tsv.gz", output: fasta = "data/gisaid/sequences.fasta", metadata = "data/gisaid/metadata_transformed.tsv", @@ -116,6 +132,7 @@ rule transform_gisaid_data: shell: """ ./bin/transform-gisaid {input.ndjson} \ + --accessions {input.accessions} \ --output-metadata {output.metadata} \ --output-fasta {output.fasta} \ --output-additional-info {output.additional_info} \ From 1e1f01489530ee205bd8326019dcb279835b6b15 Mon Sep 17 00:00:00 2001 From: Jover Lee Date: Tue, 17 Dec 2024 15:37:44 -0800 Subject: [PATCH 2/3] Concatenate accession link files MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Concatenates the downloaded accession links file to the existing `source-data/accessions.tsv.gz` since the downloaded file seems to be missing some of the existing links.¹ The transform scripts apply the accession links in order, so the last matching accession link is used in the final metadata. This allows us to default to the downloaded file which is the latest data. ¹ --- bin/fetch-accession-links | 3 +-- workflow/snakemake_rules/curate.smk | 21 ++++++++++++++++++--- 2 files changed, 19 insertions(+), 5 deletions(-) diff --git a/bin/fetch-accession-links b/bin/fetch-accession-links index 7378ed4d..42aa6e6b 100755 --- a/bin/fetch-accession-links +++ b/bin/fetch-accession-links @@ -5,5 +5,4 @@ curl "https://hgwdev.gi.ucsc.edu/~angie/epiToPublicAndDate.latest" \ --fail --silent --show-error \ --header 'User-Agent: https://github.com/nextstrain/ncov-ingest (hello@nextstrain.org)' \ | csvtk -t add-header --names gisaid_epi_isl,genbank_accession,strain,date \ - | csvtk -t cut --fields genbank_accession,gisaid_epi_isl \ - | gzip -c + | csvtk -t cut --fields genbank_accession,gisaid_epi_isl diff --git a/workflow/snakemake_rules/curate.smk b/workflow/snakemake_rules/curate.smk index 6fee03d2..f8c00765 100644 --- a/workflow/snakemake_rules/curate.smk +++ b/workflow/snakemake_rules/curate.smk @@ -28,7 +28,7 @@ rule fetch_accession_links: Fetch the accession links between GISAID and GenBank """ output: - accessions="data/accessions.tsv.gz", + accessions=temp("data/accessions.tsv"), retries: 5 shell: """ @@ -36,6 +36,21 @@ rule fetch_accession_links: """ +rule concat_accession_links: + input: + source_data="source-data/accessions.tsv.gz", + accessions="data/accessions.tsv", + output: + all_accessions="data/all_accessions.tsv.gz" + shell: + r""" + gunzip -kcfq {input.source_data:q} \ + | csvtk concat -t - {input.accessions:q} \ + | csvtk uniq -t -f genbank_accession,gisaid_epi_isl \ + | gzip -c > {output.all_accessions:q} + """ + + rule transform_rki_data: input: ndjson="data/rki.ndjson", @@ -74,7 +89,7 @@ rule transform_genbank_data: ndjson = "data/genbank.ndjson", cog_uk_accessions = "data/cog_uk_accessions.tsv", cog_uk_metadata = "data/cog_uk_metadata.csv.gz", - accessions = "data/accessions.tsv.gz", + accessions = "data/all_accessions.tsv.gz", output: fasta = "data/genbank_sequences.fasta", metadata = "data/genbank_metadata_transformed.tsv", @@ -121,7 +136,7 @@ rule merge_open_data: rule transform_gisaid_data: input: ndjson = "data/gisaid.ndjson", - accessions = "data/accessions.tsv.gz", + accessions = "data/all_accessions.tsv.gz", output: fasta = "data/gisaid/sequences.fasta", metadata = "data/gisaid/metadata_transformed.tsv", From fcb407c0072a3a58081e3098aab5d7f8d29b4ebf Mon Sep 17 00:00:00 2001 From: Jover Lee Date: Thu, 19 Dec 2024 11:54:30 -0800 Subject: [PATCH 3/3] transform: Move annotations to _after_ accessions Allows overrides of accessions if there is something that we want to manually correct in the accession links. --- bin/transform-genbank | 3 +-- bin/transform-gisaid | 2 +- 2 files changed, 2 insertions(+), 3 deletions(-) diff --git a/bin/transform-genbank b/bin/transform-genbank index 95c00f2d..5af323e4 100755 --- a/bin/transform-genbank +++ b/bin/transform-genbank @@ -207,8 +207,8 @@ if __name__ == '__main__': | ParseGeographicColumnsGenbank( base / 'source-data/us-state-codes.tsv' ) | AbbreviateAuthors() | ApplyUserGeoLocationSubstitutionRules(geoRules) - | MergeUserAnnotatedMetadata(annotations, idKey = 'genbank_accession' ) | MergeUserAnnotatedMetadata(accessions, idKey = 'genbank_accession_rev' ) + | MergeUserAnnotatedMetadata(annotations, idKey = 'genbank_accession' ) | FillDefaultLocationData() | patchUKData(args.cog_uk_accessions, args.cog_uk_metadata) | GenbankProblematicFilter( args.problem_data, @@ -301,4 +301,3 @@ if __name__ == '__main__': strain_name = updated_strain_names_by_line_no[entry[LINE_NUMBER_KEY]] print( '>' , strain_name , sep='' , file= fasta_OUT) print( entry['sequence'] , file= fasta_OUT) - diff --git a/bin/transform-gisaid b/bin/transform-gisaid index eb1bea0d..e8f95a9e 100755 --- a/bin/transform-gisaid +++ b/bin/transform-gisaid @@ -182,8 +182,8 @@ if __name__ == '__main__': pipeline = (pipeline | ApplyUserGeoLocationSubstitutionRules(geoRules) - | MergeUserAnnotatedMetadata(annotations) | MergeUserAnnotatedMetadata(accessions) + | MergeUserAnnotatedMetadata(annotations) | FillDefaultLocationData() )