diff --git a/workflow/Snakefile b/workflow/Snakefile index 646e988..f037d9d 100644 --- a/workflow/Snakefile +++ b/workflow/Snakefile @@ -21,6 +21,7 @@ configfile: workflow.source_path("../config/container_config.yaml") gencode_or_ucsc = None gencode2ensembl_file = "resources/gencode2ensembl_human.tsv" ucsc_genome_build = "hg38" +ENSEMBL_ORGANISM = "homo_sapiens" validate(config, workflow.source_path("schemas/config.schema.yaml")) @@ -34,6 +35,8 @@ if config["organism"] == "mouse": gencode_or_ucsc = "gencode" if config["genome_build"] < "GRCm39" else "UCSC" ucsc_genome_build = "mm39" if config["genome_build"] == "GRCm39" else "mm10" gencode2ensembl_file = "resources/gencode2ensembl_mouse.tsv" + ENSEMBL_ORGANISM = "mus_musculus" + # translate the gencode version to ensembl version for ensembl specific resources gencode2ensembl = pd.read_csv( @@ -44,6 +47,13 @@ ENSEMBL_VERSION = gencode2ensembl[ gencode2ensembl.GENCODE_release == config["release"] ].Ensembl_release.values[0] +ENSEMBL_ASSEMBLY_BUILD = ( + gencode2ensembl[gencode2ensembl.GENCODE_release == config["release"]] + .Genome_assembly_version.values[0] + .removeprefix("GRCh") + .removeprefix("GRCm") +) + GENCODE2ENSEMBL_CHROM_MAPPING = workflow.source_path( "resources/GRCh38_ensembl2gencode.txt" ) diff --git a/workflow/rules/common.smk b/workflow/rules/common.smk index 3c3838b..e8cc80e 100644 --- a/workflow/rules/common.smk +++ b/workflow/rules/common.smk @@ -226,3 +226,21 @@ def get_bowtie2_prefix(index_files: list[str]): # Join the dirname back in and get everything up to the first '.'. return os.path.join(dirname, basename.split(".")[0]) + + +def uniprot_snapshot_url(wildcards): + """ + Function to build the URL for the uniprot snapshot tar.gz file based on checkpoint output. + """ + release_file = checkpoints.extract_uniprot_release_from_ensembl_external_data.get( + **wildcards + ).output.uniprot_release + with open(release_file, "r") as f: + release = f.read().strip().split("\t")[1] + return storage( + ( + config["uniprot_url"] + + f"/previous_releases/release-{release}/" + + f"knowledgebase/knowledgebase{release}.tar.gz" + ) + ) diff --git a/workflow/rules/pull_resources.smk b/workflow/rules/pull_resources.smk index fbbf69a..d2e2c87 100644 --- a/workflow/rules/pull_resources.smk +++ b/workflow/rules/pull_resources.smk @@ -934,3 +934,101 @@ Write the license information file for the pulled resources. exec &> "{log}" cp "{input.license_info}" "{output.license_info}" """ + + +rule download_ensembl_external_data: + """ +Download external_data table from ENSEMBL MySQL dump to extract uniprot release used in ENSEMBL pipeline. +""" + input: + ensembl_remote=storage( + f"https://ftp.ensembl.org/pub/release-{ENSEMBL_VERSION}/mysql/{ENSEMBL_ORGANISM}_core_{ENSEMBL_VERSION}_{ENSEMBL_ASSEMBLY_BUILD}/external_db.txt.gz" + ), + output: + external_data=temp("resources/uniprot/external_data.txt.gz"), + log: + "logs/pull_resources/ensembl_external_data.log", + benchmark: + "benchmarks/pull_resources/ensembl_external_data.txt" + conda: + "../envs/shellutils.yaml" + container: + config["container"].get("shell_utils") + shell: + """ + exec &> "{log}" + cp "{input.ensembl_remote}" "{output.external_data}" + """ + + +checkpoint extract_uniprot_release_from_ensembl_external_data: + """ +Extract the UniProt release version from the ENSEMBL external_data MySQL table. +""" + input: + external_data="resources/uniprot/external_data.txt.gz", + script=workflow.source_path("../scripts/extract_uniprot_release.py"), + output: + uniprot_release="resources/uniprot/uniprot_release.txt", + log: + "logs/pull_resources/extract_uniprot_release_from_ensembl_external_data.log", + benchmark: + "benchmarks/pull_resources/extract_uniprot_release_from_ensembl_external_data.txt" + conda: + "../envs/pull_uniprot.yaml" + container: + config["container"].get("notebook") + shell: + """ + exec &> "{log}" + python {input.script} \ + --external_data "{input.external_data}" \ + --outfile "{output.uniprot_release}" + """ + + +rule download_uniprot_snapshot: + """ +Download the UniProt snapshot release from the UniProt FTP server. +""" + input: + uniprot_release=uniprot_snapshot_url, + output: + uniprot_snapshot=temp("resources/uniprot/uniprot_snapshot.tar.gz"), + log: + "logs/pull_resources/download_uniprot_snapshot.log", + benchmark: + "benchmarks/pull_resources/download_uniprot_snapshot.txt" + conda: + "../envs/shellutils.yaml" + container: + config["container"].get("shell_utils") + shell: + """ + exec &> "{log}" + cp "{input.uniprot_release}" "{output.uniprot_snapshot}" + """ + + +rule extract_uniprot_snapshot: + """ +Extract the UniProt snapshot release to obtain SwissProt and TrEMBL data. +""" + input: + uniprot_snapshot="resources/uniprot/uniprot_snapshot.tar.gz", + output: + swissprot="resources/uniprot/swisprot.dat.gz", + trembl="resources/uniprot/trembl.dat.gz", + log: + "logs/pull_resources/extract_uniprot_snapshot.log", + benchmark: + "benchmarks/pull_resources/extract_uniprot_snapshot.txt" + conda: + "../envs/shellutils.yaml" + container: + config["container"].get("shell_utils") + shell: + """ + exec &> "{log}" + tar xzf "{input.uniprot_snapshot}" -C resources/uniprot uniprot_trembl.dat.gz uniprot_sprot.dat.gz + """ diff --git a/workflow/schemas/config.schema.yaml b/workflow/schemas/config.schema.yaml index 14379e7..a243ca9 100644 --- a/workflow/schemas/config.schema.yaml +++ b/workflow/schemas/config.schema.yaml @@ -96,6 +96,11 @@ properties: needed when organism is "human". default: https://storage.googleapis.com/gcp-public-data--gnomad/release type: string + uniprot_url: + description: >- + URL from which the UniProt resources may be downloaded. + default: https://ftp.uniprot.org/pub/databases/uniprot + type: string chrom_filter: description: >- List of chromosome names to be included in the genomelib. Only needed when diff --git a/workflow/scripts/extract_uniprot_release.py b/workflow/scripts/extract_uniprot_release.py new file mode 100644 index 0000000..fc53453 --- /dev/null +++ b/workflow/scripts/extract_uniprot_release.py @@ -0,0 +1,55 @@ +import pathlib +import re +import argparse +import gzip +import sys + +RELEASE_PATTERN = re.compile(r"UniProtKB/Swiss-Prot\s+Release\s+(\d{4}_\d{2})") + + +def parse_uniprot_release(external_data_file: pathlib.Path) -> str: + """ + Parse the uniprot release from the ENSEMBL external_data.txt.gz file. + + Args: + external_data_file: Path to the external_data.txt.gz file. + """ + with gzip.open(external_data_file, "rt", encoding="utf-8") as file: + for line in file: + match = RELEASE_PATTERN.search(line) + if match: + release = match.group(1) + return release + + print("Uniprot release not found in external_data.txt.gz file.") + return + + +# Main + +parser = argparse.ArgumentParser( + description="Extract uniprot release from ENSEMBL external_data.txt.gz file" +) +parser.add_argument( + "--outfile", + type=str, + help="Output TSV file", + required=True, +) +parser.add_argument( + "--external_data", + type=str, + help="Path external_data.txt.gz file", + required=True, +) + + +args = parser.parse_args() + +release = parse_uniprot_release(args.external_data) + +if release is None: + sys.exit(1) + +with open(args.outfile, "w") as file_hande: + file_hande.write(f"Release\t{release}\n")