Mercurial > repos > kkonganti > cfsan_cronology

annotate 0.2.0/bin/dl_pdg_metadata.py @ 17:b571995ddb51

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
planemo upload
author kkonganti
date Mon, 15 Jul 2024 19:01:29 -0400
parents a5f31c44f8c9
children
rev   line source
kkonganti@11 1 #!/usr/bin/env python3
kkonganti@11 2
kkonganti@11 3 # Kranti Konganti
kkonganti@11 4
kkonganti@11 5 import argparse
kkonganti@11 6 import inspect
kkonganti@11 7 import logging
kkonganti@11 8 import os
kkonganti@11 9 import re
kkonganti@11 10 import shutil
kkonganti@11 11 import tempfile
kkonganti@11 12 from html.parser import HTMLParser
kkonganti@11 13 from urllib.request import urlopen
kkonganti@11 14
kkonganti@11 15 # Set logging format.
kkonganti@11 16 logging.basicConfig(
kkonganti@11 17 format="\n" + "=" * 55 + "\n%(asctime)s - %(levelname)s\n" + "=" * 55 + "\n%(message)s\n",
kkonganti@11 18 level=logging.DEBUG,
kkonganti@11 19 )
kkonganti@11 20
kkonganti@11 21
kkonganti@11 22 # Multiple inheritence for pretty printing of help text.
kkonganti@11 23 class MultiArgFormatClasses(argparse.RawTextHelpFormatter, argparse.ArgumentDefaultsHelpFormatter):
kkonganti@11 24 pass
kkonganti@11 25
kkonganti@11 26
kkonganti@11 27 # HTMLParser override class to get PDG release and latest Cluster .tsv file
kkonganti@11 28 class NCBIPathogensHTMLParser(HTMLParser):
kkonganti@11 29 def __init__(self, *, convert_charrefs: bool = ...) -> None:
kkonganti@11 30 super().__init__(convert_charrefs=convert_charrefs)
kkonganti@11 31 self.reset()
kkonganti@11 32 self.href_data = list()
kkonganti@11 33
kkonganti@11 34 def handle_data(self, data):
kkonganti@11 35 self.href_data.append(data)
kkonganti@11 36
kkonganti@11 37
kkonganti@11 38 def dl_pdg(**kwargs) -> None:
kkonganti@11 39 """
kkonganti@11 40 Method to save the PDG metadata file and
kkonganti@11 41 return the latest PDG release.
kkonganti@11 42 """
kkonganti@11 43 db_path, url, regex, suffix, overwrite, release = [kwargs[k] for k in kwargs.keys()]
kkonganti@11 44
kkonganti@11 45 if (db_path or url) == None:
kkonganti@11 46 logging.error("Please provide absolute UNIX path\n" + "to store the result DB flat files.")
kkonganti@11 47 exit(1)
kkonganti@11 48
kkonganti@11 49 if re.match(r"^PDG\d+\.\d+$", release):
kkonganti@11 50 url = re.sub("latest_snps", release.strip(), url)
kkonganti@11 51
kkonganti@11 52 html_parser = NCBIPathogensHTMLParser()
kkonganti@11 53 logging.info(f"Finding latest NCBI PDG release at:\n{url}")
kkonganti@11 54
kkonganti@11 55 with urlopen(url) as response:
kkonganti@11 56 with tempfile.NamedTemporaryFile(delete=False) as tmp_html_file:
kkonganti@11 57 shutil.copyfileobj(response, tmp_html_file)
kkonganti@11 58
kkonganti@11 59 with open(tmp_html_file.name, "r") as html:
kkonganti@11 60 html_parser.feed("".join(html.readlines()))
kkonganti@11 61
kkonganti@11 62 pdg_filename = re.search(regex, "".join(html_parser.href_data)).group(0)
kkonganti@11 63 pdg_release = pdg_filename.rstrip(suffix)
kkonganti@11 64 pdg_metadata_url = "/".join([url, pdg_filename])
kkonganti@11 65 pdg_release = pdg_filename.rstrip(suffix)
kkonganti@11 66 dest_dir = os.path.join(db_path, pdg_release)
kkonganti@11 67
kkonganti@11 68 logging.info(f"Found NCBI PDG file:\n{pdg_metadata_url}")
kkonganti@11 69
kkonganti@11 70 if (
kkonganti@11 71 not overwrite
kkonganti@11 72 and re.match(r".+?\.metadata\.tsv$", pdg_filename)
kkonganti@11 73 and os.path.exists(dest_dir)
kkonganti@11 74 ):
kkonganti@11 75 logging.error(f"DB path\n{dest_dir}\nalready exists. Please use -f to overwrite.")
kkonganti@11 76 exit(1)
kkonganti@11 77 elif overwrite and not re.match(r".+?\.reference_target\.cluster_list\.tsv$", pdg_filename):
kkonganti@11 78 shutil.rmtree(dest_dir, ignore_errors=True) if os.path.exists(dest_dir) else None
kkonganti@11 79 os.makedirs(dest_dir)
kkonganti@11 80 elif (
kkonganti@11 81 not overwrite
kkonganti@11 82 and re.match(r".+?\.metadata\.tsv$", pdg_filename)
kkonganti@11 83 and not os.path.exists(dest_dir)
kkonganti@11 84 ):
kkonganti@11 85 os.makedirs(dest_dir)
kkonganti@11 86
kkonganti@11 87 tsv_at = os.path.join(dest_dir, pdg_filename)
kkonganti@11 88 logging.info(f"Saving to:\n{tsv_at}")
kkonganti@11 89
kkonganti@11 90 with urlopen(pdg_metadata_url) as response:
kkonganti@11 91 with open(tsv_at, "w") as tsv:
kkonganti@11 92 tsv.writelines(response.read().decode("utf-8"))
kkonganti@11 93
kkonganti@11 94 html.close()
kkonganti@11 95 tmp_html_file.close()
kkonganti@11 96 os.unlink(tmp_html_file.name)
kkonganti@11 97 tsv.close()
kkonganti@11 98 response.close()
kkonganti@11 99
kkonganti@11 100 return tsv_at, dest_dir
kkonganti@11 101
kkonganti@11 102
kkonganti@11 103 def main() -> None:
kkonganti@11 104 """
kkonganti@11 105 This script is part of the `cronology_db` Nextflow workflow and is only
kkonganti@11 106 tested on POSIX sytems.
kkonganti@11 107 It:
kkonganti@11 108 1. Downloads the latest NCBI Pathogens Release metadata file, which
kkonganti@11 109 looks like PDGXXXXXXXXXX.2504.metadata.csv and also the SNP cluster
kkonganti@11 110 information file which looks like PDGXXXXXXXXXX.2504.reference_target.cluster_list.tsv
kkonganti@11 111 2. Generates a new metadata file with only required information such as
kkonganti@11 112 computed_serotype, isolates GenBank or RefSeq downloadable genome FASTA
kkonganti@11 113 URL.
kkonganti@11 114 """
kkonganti@11 115
kkonganti@11 116 prog_name = os.path.basename(inspect.stack()[0].filename)
kkonganti@11 117
kkonganti@11 118 parser = argparse.ArgumentParser(
kkonganti@11 119 prog=prog_name, description=main.__doc__, formatter_class=MultiArgFormatClasses
kkonganti@11 120 )
kkonganti@11 121
kkonganti@11 122 required = parser.add_argument_group("required arguments")
kkonganti@11 123
kkonganti@11 124 parser.add_argument(
kkonganti@11 125 "-db",
kkonganti@11 126 dest="db_path",
kkonganti@11 127 default=os.getcwd(),
kkonganti@11 128 required=False,
kkonganti@11 129 help="Absolute UNIX path to a path where all results files are\nstored.",
kkonganti@11 130 )
kkonganti@11 131 parser.add_argument(
kkonganti@11 132 "-f",
kkonganti@11 133 dest="overwrite_db",
kkonganti@11 134 default=False,
kkonganti@11 135 required=False,
kkonganti@11 136 action="store_true",
kkonganti@11 137 help="Force overwrite a PDG release directory at DB path\nmentioned with -db.",
kkonganti@11 138 )
kkonganti@11 139 parser.add_argument(
kkonganti@11 140 "-org",
kkonganti@11 141 dest="organism",
kkonganti@11 142 default="Cronobacter",
kkonganti@11 143 required=False,
kkonganti@11 144 help="The organism to create the DB flat files\nfor.",
kkonganti@11 145 )
kkonganti@11 146 required.add_argument(
kkonganti@11 147 "-rel",
kkonganti@11 148 dest="release",
kkonganti@11 149 default=False,
kkonganti@11 150 required=False,
kkonganti@11 151 help="If you get a 404 error, try mentioning the actual release identifier.\n"
kkonganti@11 152 + "Ex: For Cronobacter, you can get the release identifier by going to:\n"
kkonganti@11 153 + " https://ftp.ncbi.nlm.nih.gov/pathogen/Results/Cronobacter\n"
kkonganti@11 154 + "Ex: If you want metadata beloginging to release PDG000000002.2507, then you\n"
kkonganti@11 155 + " would use this command-line option as:\n -rel PDG000000002.2507",
kkonganti@11 156 )
kkonganti@11 157
kkonganti@11 158 args = parser.parse_args()
kkonganti@11 159 db_path = args.db_path
kkonganti@11 160 org = args.organism
kkonganti@11 161 overwrite = args.overwrite_db
kkonganti@11 162 release = args.release
kkonganti@11 163 ncbi_pathogens_loc = "/".join(
kkonganti@11 164 ["https://ftp.ncbi.nlm.nih.gov/pathogen/Results", org, "latest_snps"]
kkonganti@11 165 )
kkonganti@11 166
kkonganti@11 167 if not db_path:
kkonganti@11 168 db_path = os.getcwd()
kkonganti@11 169
kkonganti@11 170 # Save metadata
kkonganti@11 171 file, dest_dir = dl_pdg(
kkonganti@11 172 db_path=db_path,
kkonganti@11 173 url="/".join([ncbi_pathogens_loc, "Metadata"]),
kkonganti@11 174 regex=re.compile(r"PDG\d+\.\d+\.metadata\.tsv"),
kkonganti@11 175 suffix=".metadata.tsv",
kkonganti@11 176 overwrite=overwrite,
kkonganti@11 177 release=release,
kkonganti@11 178 )
kkonganti@11 179
kkonganti@11 180 # Save cluster to target mapping
kkonganti@11 181 dl_pdg(
kkonganti@11 182 db_path=db_path,
kkonganti@11 183 url="/".join([ncbi_pathogens_loc, "Clusters"]),
kkonganti@11 184 regex=re.compile(r"PDG\d+\.\d+\.reference_target\.cluster_list\.tsv"),
kkonganti@11 185 suffix="reference_target\.cluster_list\.tsv",
kkonganti@11 186 overwrite=overwrite,
kkonganti@11 187 release=release,
kkonganti@11 188 )
kkonganti@11 189
kkonganti@11 190 # Create accs.txt for dataformat to fetch required ACC fields
kkonganti@11 191 accs_file = os.path.join(dest_dir, "accs_all.txt")
kkonganti@11 192 with open(file, "r") as pdg_metadata_fh:
kkonganti@11 193 with open(accs_file, "w") as accs_fh:
kkonganti@11 194 for line in pdg_metadata_fh.readlines():
kkonganti@11 195 if re.match(r"^#", line) or line in ["\n", "\n\r", "\r"]:
kkonganti@11 196 continue
kkonganti@11 197 cols = line.strip().split("\t")
kkonganti@11 198 asm_acc = cols[9]
kkonganti@11 199 accs_fh.write(f"{asm_acc}\n") if (asm_acc != "NULL") else None
kkonganti@11 200 accs_fh.close()
kkonganti@11 201 pdg_metadata_fh.close()
kkonganti@11 202
kkonganti@11 203 logging.info("Finished writing accessions for dataformat tool.")
kkonganti@11 204
kkonganti@11 205
kkonganti@11 206 if __name__ == "__main__":
kkonganti@11 207 main()