aboutsummaryrefslogtreecommitdiff
path: root/scripts
diff options
context:
space:
mode:
Diffstat (limited to 'scripts')
-rw-r--r--scripts/cli/__init__.py3
-rw-r--r--scripts/cli/logging.py18
-rw-r--r--scripts/cli/options.py56
-rw-r--r--scripts/cli/validators.py10
-rw-r--r--scripts/cli_parser.py23
-rw-r--r--scripts/compute_phenotype_means.py101
-rw-r--r--scripts/insert_data.py6
-rw-r--r--scripts/insert_samples.py40
-rw-r--r--scripts/load_phenotypes_to_db.py551
-rw-r--r--scripts/phenotypes/__init__.py1
-rw-r--r--scripts/phenotypes/delete_phenotypes.py173
-rw-r--r--scripts/process_rqtl2_bundle.py4
-rw-r--r--scripts/qc_on_rqtl2_bundle.py13
-rw-r--r--scripts/redis_logger.py2
-rw-r--r--scripts/rqtl2/entry.py30
-rw-r--r--scripts/rqtl2/install_genotypes.py6
-rw-r--r--scripts/rqtl2/install_phenos.py7
-rw-r--r--scripts/rqtl2/phenotypes_qc.py64
-rw-r--r--scripts/run_qtlreaper.py238
-rw-r--r--scripts/worker.py2
20 files changed, 1256 insertions, 92 deletions
diff --git a/scripts/cli/__init__.py b/scripts/cli/__init__.py
new file mode 100644
index 0000000..45bbda9
--- /dev/null
+++ b/scripts/cli/__init__.py
@@ -0,0 +1,3 @@
+"""Package to hold CLI-specific utilities."""
+
+from . import options
diff --git a/scripts/cli/logging.py b/scripts/cli/logging.py
new file mode 100644
index 0000000..30ecf17
--- /dev/null
+++ b/scripts/cli/logging.py
@@ -0,0 +1,18 @@
+"""Logging for scripts."""
+import logging
+
+def setup_logging(
+ script_logger: logging.Logger,
+ loglevel: str,
+ modules: tuple[str, ...] = tuple()
+):
+ """Setup module-level loggers to the same log-level as the application."""
+ logging.basicConfig(
+ encoding="utf-8",
+ format=("%(asctime)s — %(filename)s:%(lineno)s — %(levelname)s: "
+ "%(message)s"),
+ level=logging.INFO)
+ script_logger.setLevel(getattr(logging, loglevel.upper()))
+ effective_loglevel = logging.getLevelName(script_logger.getEffectiveLevel())
+ for module in modules:
+ logging.getLogger(module).setLevel(effective_loglevel)
diff --git a/scripts/cli/options.py b/scripts/cli/options.py
new file mode 100644
index 0000000..58d3df4
--- /dev/null
+++ b/scripts/cli/options.py
@@ -0,0 +1,56 @@
+"""General options to be added to CLI scripts."""
+from argparse import ArgumentParser
+
+
+def add_logging(parser: ArgumentParser) -> ArgumentParser:
+ """Add optional log-level option"""
+ loglevels = ("debug", "info", "warning", "error", "critical")
+ parser.add_argument(
+ "--log_level",
+ "--log-level",
+ "--loglevel",
+ metavar="LOG-LEVEL",
+ type=str,
+ default="INFO",
+ choices=loglevels,
+ help=("Controls the severity of events to log. Valid values are: " +
+ ", ".join(f"'{level}'" for level in loglevels)))
+ return parser
+
+
+def add_mariadb_uri(parser: ArgumentParser) -> ArgumentParser:
+ """Add the MySQL/MariaDB URI argument."""
+ parser.add_argument("db_uri",
+ metavar="DB-URI",
+ type=str,
+ help="MariaDB/MySQL connection URL")
+ return parser
+
+
+def add_species_id(parser: ArgumentParser) -> ArgumentParser:
+ """Add species-id as a mandatory argument."""
+ parser.add_argument("species_id",
+ metavar="SPECIES-ID",
+ type=int,
+ help="The species to operate on.")
+ return parser
+
+
+def add_population_id(parser: ArgumentParser) -> ArgumentParser:
+ """Add population-id as a mandatory argument."""
+ parser = add_species_id(parser)
+ parser.add_argument("population_id",
+ metavar="POPULATION-ID",
+ type=int,
+ help="The ID for the population to operate on.")
+ return parser
+
+
+def add_dataset_id(parser: ArgumentParser) -> ArgumentParser:
+ """Add dataset-id as a mandatory argument."""
+ parser = add_population_id(parser)
+ parser.add_argument("dataset_id",
+ metavar="DATASET-ID",
+ type=int,
+ help="The ID for the dataset to operate on.")
+ return parser
diff --git a/scripts/cli/validators.py b/scripts/cli/validators.py
new file mode 100644
index 0000000..6d16e4c
--- /dev/null
+++ b/scripts/cli/validators.py
@@ -0,0 +1,10 @@
+"""CLI options validators."""
+from pathlib import Path
+
+
+def directory_exists(val: str) -> Path:
+ """Check that directory path specified actually exists."""
+ _dir = Path(val).absolute()
+ if _dir.is_dir() and _dir.exists():
+ return _dir
+ raise FileNotFoundError(f"The path '{_dir}' MUST exist and be a directory.")
diff --git a/scripts/cli_parser.py b/scripts/cli_parser.py
index d42ae66..bf39731 100644
--- a/scripts/cli_parser.py
+++ b/scripts/cli_parser.py
@@ -3,6 +3,20 @@ from uuid import UUID
from typing import Optional
from argparse import ArgumentParser
+
+def add_logging_option(parser: ArgumentParser) -> ArgumentParser:
+ """Add optional log-level option"""
+ parser.add_argument(
+ "--log-level",
+ "--loglevel",
+ type=str,
+ default="INFO",
+ choices=["DEBUG", "INFO", "WARNING", "ERROR", "CRITICAL",
+ "debug", "info", "warning", "error", "critical"],
+ help="The severity of events to track with the logger.")
+ return parser
+
+
def init_cli_parser(program: str, description: Optional[str] = None) -> ArgumentParser:
"""Initialise the CLI arguments parser."""
parser = ArgumentParser(prog=program, description=description)
@@ -19,13 +33,8 @@ def init_cli_parser(program: str, description: Optional[str] = None) -> Argument
type=int,
default=86400,
help="How long to keep any redis keys around.")
- parser.add_argument(
- "--loglevel",
- type=str,
- default="INFO",
- choices=["DEBUG", "INFO", "WARNING", "ERROR", "CRITICAL"],
- help="The severity of events to track with the logger.")
- return parser
+ return add_logging_option(parser)
+
def add_global_data_arguments(parser: ArgumentParser) -> ArgumentParser:
"""Add the global (present in nearly ALL scripts) CLI arguments."""
diff --git a/scripts/compute_phenotype_means.py b/scripts/compute_phenotype_means.py
new file mode 100644
index 0000000..6d39ace
--- /dev/null
+++ b/scripts/compute_phenotype_means.py
@@ -0,0 +1,101 @@
+"""Compute phenotype means."""
+import sys
+import logging
+from pathlib import Path
+from typing import TypeVar
+from argparse import Namespace, ArgumentParser
+
+import MySQLdb
+
+from gn_libs import mysqldb
+from uploader import setup_modules_logging
+
+from .cli_parser import add_logging_option
+from .load_phenotypes_to_db import update_means
+
+logger = logging.getLogger(__name__)
+logging.basicConfig(
+ encoding="utf-8",
+ format="%(asctime)s - %(name)s - %(levelname)s — %(message)s",
+ level=logging.INFO)
+
+
+def fetch_xref_id(conn: mysqldb.Connection, population_id: int) -> tuple[int, ...]:
+ """Fetch a population's cross-reference IDs."""
+ logger.debug("Fetching the xref IDs.")
+ with conn.cursor(cursorclass=MySQLdb.cursors.DictCursor) as cursor:
+ query = "SELECT Id FROM PublishXRef WHERE InbredSetId=%(population_id)s"
+ cursor.execute(query, {"population_id": population_id})
+ return tuple(int(row["Id"]) for row in cursor.fetchall())
+
+
+def run(args) -> int:
+ """Run the script."""
+ logger.debug("Running the script!")
+ with mysqldb.database_connection(args.db_uri) as mariadb_conn:
+ xref_ids = args.cross_ref_ids or fetch_xref_id(mariadb_conn, args.population_id)
+ if len(xref_ids):
+ update_means(mariadb_conn,
+ args.population_id,
+ xref_ids)
+ logger.debug("Successfully computed means for %02d phenotypes.",
+ len(xref_ids))
+ return 0
+ _reasons = (
+ f"no population exists with the ID {args.population_id}",
+ "the population exists but it has no phenotypes linked to it yet")
+ logger.error(
+ "No cross-reference IDs to run against. Likely causes are: %s",
+ " OR ".join(_reasons) + ".")
+ return 1
+
+
+T = TypeVar("T")
+def comma_separated_list(val: str, itemstype: type = str) -> tuple[T, ...]:
+ """Convert val into a list of items of type 'itemstype'."""
+ return tuple(itemstype(item.strip()) for item in val.split(","))
+
+
+def comma_separated_list_of_integers(val: str) -> tuple[int, ...]:
+ """Convert 'val' into list of items of type 'int'."""
+ return comma_separated_list(val, int)
+
+
+if __name__ == "__main__":
+ def parse_args() -> Namespace:
+ """Define and parse the CLI parsers accepted by this script."""
+ parser = ArgumentParser(
+ "compute-phenotype-means",
+ description="Compute/Recompute the phenotype means.")
+ parser.add_argument("db_uri",
+ metavar="db-uri",
+ type=str,
+ help="MariaDB/MySQL connection URL")
+ parser.add_argument("jobs_db_path",
+ metavar="jobs-db-path",
+ type=Path,
+ help="Path to jobs' SQLite database.")
+ parser.add_argument("population_id",
+ metavar="population-id",
+ type=int,
+ help=("Identifier for the InbredSet group/"
+ "population to run means against."))
+ ## Optional arguments
+ parser = add_logging_option(parser)
+ parser.add_argument(
+ "--cross-ref-ids",
+ type=comma_separated_list_of_integers,
+ help=("Provide cross-reference IDs to narrow the number of "
+ "phenotypes that the means are computed against."),
+ default=[])
+
+ return parser.parse_args()
+
+ def main() -> int:
+ """compute-phenotype-means: Entry-point function."""
+ args = parse_args()
+ logger.setLevel(getattr(logging, args.log_level.upper()))
+ setup_modules_logging(logger, ("scripts.load_phenotypes_to_db",))
+ return run(args)
+
+ sys.exit(main())
diff --git a/scripts/insert_data.py b/scripts/insert_data.py
index 67038f8..aec0251 100644
--- a/scripts/insert_data.py
+++ b/scripts/insert_data.py
@@ -197,7 +197,7 @@ def probeset_ids(dbconn: mdb.Connection,
break
yield row
-def insert_means(# pylint: disable=[too-many-locals, too-many-arguments]
+def insert_means(# pylint: disable=[too-many-locals, too-many-arguments, too-many-positional-arguments]
filepath: str, speciesid: int, platform_id: int, datasetid: int,
dbconn: mdb.Connection, rconn: Redis) -> int: # pylint: disable=[unused-argument]
"Insert the means/averages data into the database"
@@ -232,7 +232,7 @@ def insert_means(# pylint: disable=[too-many-locals, too-many-arguments]
item for sublist in
read_datavalues(filepath, headings, strains).values()
for item in sublist),
- start=(last_data_id(dbconn)+1)))
+ start=last_data_id(dbconn)+1))
with dbconn.cursor(cursorclass=DictCursor) as cursor:
while True:
means = tuple(take(the_means, 10000))
@@ -245,7 +245,7 @@ def insert_means(# pylint: disable=[too-many-locals, too-many-arguments]
cursor.executemany(xref_query, means)
return 0
-def insert_se(# pylint: disable = [too-many-arguments,too-many-locals]
+def insert_se(# pylint: disable = [too-many-arguments,too-many-locals, too-many-positional-arguments]
filepath: str, speciesid: int, platformid: int, datasetid: int,
dbconn: mdb.Connection, rconn: Redis) -> int: # pylint: disable=[unused-argument]
"Insert the standard-error data into the database"
diff --git a/scripts/insert_samples.py b/scripts/insert_samples.py
index 1b0a052..96ae8e2 100644
--- a/scripts/insert_samples.py
+++ b/scripts/insert_samples.py
@@ -3,12 +3,13 @@ import sys
import logging
import pathlib
import argparse
+import traceback
import MySQLdb as mdb
-from redis import Redis
+
from gn_libs.mysqldb import database_connection
-from uploader.check_connections import check_db, check_redis
+from uploader.check_connections import check_db
from uploader.species.models import species_by_id
from uploader.population.models import population_by_id
from uploader.samples.models import (
@@ -33,8 +34,7 @@ class SeparatorAction(argparse.Action):
"""Process the value passed in."""
setattr(namespace, self.dest, (chr(9) if values == "\\t" else values))
-def insert_samples(conn: mdb.Connection,# pylint: disable=[too-many-arguments]
- rconn: Redis,# pylint: disable=[unused-argument]
+def insert_samples(conn: mdb.Connection,# pylint: disable=[too-many-arguments, too-many-positional-arguments]
speciesid: int,
populationid: int,
samplesfile: pathlib.Path,
@@ -73,6 +73,7 @@ def insert_samples(conn: mdb.Connection,# pylint: disable=[too-many-arguments]
print("Samples upload successfully completed.")
return 0
+
if __name__ == "__main__":
def cli_args():
@@ -117,33 +118,30 @@ if __name__ == "__main__":
help=("The character used to delimit (surround?) the value in "
"each column."))
- # == Script-specific extras ==
- parser.add_argument("--redisuri",
- help="URL to initialise connection to redis",
- default="redis:///")
-
args = parser.parse_args()
return args
def main():
"""Run script to insert samples into the database."""
-
+ status_code = 1 # Exit with an Exception
args = cli_args()
check_db(args.databaseuri)
- check_redis(args.redisuri)
if not args.samplesfile.exists():
logging.error("File not found: '%s'.", args.samplesfile)
return 2
- with (Redis.from_url(args.redisuri, decode_responses=True) as rconn,
- database_connection(args.databaseuri) as dbconn):
- return insert_samples(dbconn,
- rconn,
- args.speciesid,
- args.populationid,
- args.samplesfile,
- args.separator,
- args.firstlineheading,
- args.quotechar)
+ with database_connection(args.databaseuri) as dbconn:
+ try:
+ status_code = insert_samples(dbconn,
+ args.speciesid,
+ args.populationid,
+ args.samplesfile,
+ args.separator,
+ args.firstlineheading,
+ args.quotechar)
+ except Exception as _exc:# pylint: disable=[broad-exception-caught]
+ print(traceback.format_exc(), file=sys.stderr)
+
+ return status_code
sys.exit(main())
diff --git a/scripts/load_phenotypes_to_db.py b/scripts/load_phenotypes_to_db.py
new file mode 100644
index 0000000..31eb715
--- /dev/null
+++ b/scripts/load_phenotypes_to_db.py
@@ -0,0 +1,551 @@
+"""Load phenotypes and their data provided in files into the database."""
+import sys
+import uuid
+import json
+import time
+import logging
+import argparse
+from pathlib import Path
+from zipfile import ZipFile
+from datetime import datetime
+from typing import Any, Iterable
+from urllib.parse import urljoin
+from functools import reduce, partial
+
+from MySQLdb.cursors import DictCursor
+
+from gn_libs import jobs, mysqldb, sqlite3, monadic_requests as mrequests
+
+from r_qtl import r_qtl2 as rqtl2
+from uploader.species.models import species_by_id
+from uploader.population.models import population_by_species_and_id
+from uploader.samples.models import samples_by_species_and_population
+from uploader.phenotypes.models import (
+ dataset_by_id,
+ save_phenotypes_data,
+ create_new_phenotypes,
+ quick_save_phenotypes_data)
+from uploader.publications.models import fetch_publication_by_id
+
+from scripts.rqtl2.bundleutils import build_line_joiner, build_line_splitter
+
+from functional_tools import take
+
+logging.basicConfig(
+ format="%(asctime)s — %(filename)s:%(lineno)s — %(levelname)s: %(message)s")
+logger = logging.getLogger(__name__)
+
+
+
+def __replace_na_strings__(line, na_strings):
+ return ((None if value in na_strings else value) for value in line)
+
+
+def save_phenotypes(
+ conn: mysqldb.Connection,
+ control_data: dict[str, Any],
+ population_id,
+ publication_id,
+ filesdir: Path
+) -> tuple[dict, ...]:
+ """Read `phenofiles` and save the phenotypes therein."""
+ phenofiles = tuple(filesdir.joinpath(_file) for _file in control_data["phenocovar"])
+ if len(phenofiles) <= 0:
+ return tuple()
+
+ if control_data["phenocovar_transposed"]:
+ logger.info("Undoing transposition of the files rows and columns.")
+ phenofiles = tuple(
+ rqtl2.transpose_csv_with_rename(
+ _file,
+ build_line_splitter(control_data),
+ build_line_joiner(control_data))
+ for _file in phenofiles)
+
+ _headers = rqtl2.read_csv_file_headers(phenofiles[0],
+ control_data["phenocovar_transposed"],
+ control_data["sep"],
+ control_data["comment.char"])
+ return create_new_phenotypes(
+ conn,
+ population_id,
+ publication_id,
+ (dict(zip(_headers,
+ __replace_na_strings__(line, control_data["na.strings"])))
+ for filecontent
+ in (rqtl2.read_csv_file(path,
+ separator=control_data["sep"],
+ comment_char=control_data["comment.char"])
+ for path in phenofiles)
+ for idx, line in enumerate(filecontent)
+ if idx != 0))
+
+
+def __row_to_dataitems__(
+ sample_row: dict,
+ dataidmap: dict,
+ pheno_name2id: dict[str, int],
+ samples: dict
+) -> Iterable[dict]:
+ samplename = sample_row["id"]
+
+ return ({
+ "phenotype_id": dataidmap[pheno_name2id[phenoname]]["phenotype_id"],
+ "data_id": dataidmap[pheno_name2id[phenoname]]["data_id"],
+ "sample_name": samplename,
+ "sample_id": samples[samplename]["Id"],
+ "value": phenovalue
+ } for phenoname, phenovalue in sample_row.items() if phenoname != "id")
+
+
+def __build_dataitems__(
+ phenofiles,
+ control_data,
+ samples,
+ dataidmap,
+ pheno_name2id
+):
+ _headers = rqtl2.read_csv_file_headers(
+ phenofiles[0],
+ False, # Any transposed files have been un-transposed by this point
+ control_data["sep"],
+ control_data["comment.char"])
+ _filescontents = (
+ rqtl2.read_csv_file(path,
+ separator=control_data["sep"],
+ comment_char=control_data["comment.char"])
+ for path in phenofiles)
+ _linescontents = (
+ __row_to_dataitems__(
+ dict(zip(("id",) + _headers[1:],
+ __replace_na_strings__(line, control_data["na.strings"]))),
+ dataidmap,
+ pheno_name2id,
+ samples)
+ for linenum, line in (enumline for filecontent in _filescontents
+ for enumline in enumerate(filecontent))
+ if linenum > 0)
+ return (item for items in _linescontents
+ for item in items
+ if item["value"] is not None)
+
+
+def save_numeric_data(# pylint: disable=[too-many-positional-arguments,too-many-arguments]
+ conn: mysqldb.Connection,
+ dataidmap: dict,
+ pheno_name2id: dict[str, int],
+ samples: dict,
+ control_data: dict,
+ filesdir: Path,
+ filetype: str,
+ table: str
+):
+ """Read data from files and save to the database."""
+ phenofiles = tuple(
+ filesdir.joinpath(_file) for _file in control_data[filetype])
+ if len(phenofiles) <= 0:
+ return tuple()
+
+ if control_data[f"{filetype}_transposed"]:
+ logger.info("Undoing transposition of the files rows and columns.")
+ phenofiles = tuple(
+ rqtl2.transpose_csv_with_rename(
+ _file,
+ build_line_splitter(control_data),
+ build_line_joiner(control_data))
+ for _file in phenofiles)
+
+ try:
+ logger.debug("Attempt quick save with `LOAD … INFILE`.")
+ return quick_save_phenotypes_data(
+ conn,
+ table,
+ __build_dataitems__(
+ phenofiles,
+ control_data,
+ samples,
+ dataidmap,
+ pheno_name2id),
+ filesdir)
+ except Exception as _exc:# pylint: disable=[broad-exception-caught]
+ logger.debug("Could not use `LOAD … INFILE`, using raw query",
+ exc_info=True)
+ time.sleep(60)
+ return save_phenotypes_data(
+ conn,
+ table,
+ __build_dataitems__(
+ phenofiles,
+ control_data,
+ samples,
+ dataidmap,
+ pheno_name2id))
+
+
+save_pheno_data = partial(save_numeric_data,
+ filetype="pheno",
+ table="PublishData")
+
+
+save_phenotypes_se = partial(save_numeric_data,
+ filetype="phenose",
+ table="PublishSE")
+
+
+save_phenotypes_n = partial(save_numeric_data,
+ filetype="phenonum",
+ table="NStrain")
+
+
+def update_auth(# pylint: disable=[too-many-locals,too-many-positional-arguments,too-many-arguments]
+ auth_details,
+ resource_details,
+ species,
+ population,
+ dataset,
+ xrefdata):
+ """Grant the user access to their data."""
+ logger.info("Updating authorisation for the data.")
+ logger.debug("Resource details for the authorisation: %s", resource_details)
+ authserver, token = auth_details
+ _tries = 0
+ _delay = 1
+ headers = {
+ "Authorization": f"Bearer {token}",
+ "Content-Type": "application/json"
+ }
+ def authserveruri(endpoint):
+ return urljoin(authserver, endpoint)
+
+ def __fetch_user_details__():
+ logger.info("… Fetching user details")
+ return mrequests.get(
+ authserveruri("/auth/user/"),
+ headers=headers
+ )
+
+ def __link_data__(user):
+ logger.info("… linking uploaded data to user's group")
+ return mrequests.post(
+ authserveruri("/auth/data/link/phenotype"),
+ headers=headers,
+ json={
+ "species_name": species["Name"],
+ "group_id": user["group"]["group_id"],
+ "selected": [
+ {
+ "SpeciesId": species["SpeciesId"],
+ "InbredSetId": population["Id"],
+ "PublishFreezeId": dataset["Id"],
+ "dataset_name": dataset["Name"],
+ "dataset_fullname": dataset["FullName"],
+ "dataset_shortname": dataset["ShortName"],
+ "PublishXRefId": item["xref_id"]
+ }
+ for item in xrefdata
+ ],
+ "using-raw-ids": "on"
+ }).then(lambda ld_results: (user, ld_results))
+
+ def __fetch_phenotype_category_details__(user, linkeddata):
+ logger.info("… fetching phenotype category details")
+ return mrequests.get(
+ authserveruri("/auth/resource/categories"),
+ headers=headers
+ ).then(
+ lambda categories: (
+ user,
+ linkeddata,
+ next(category for category in categories
+ if category["resource_category_key"] == "phenotype"))
+ )
+
+ def __create_resource__(user, linkeddata, category):
+ logger.info("… creating authorisation resource object")
+ return mrequests.post(
+ authserveruri("/auth/resource/create"),
+ headers=headers,
+ json={
+ **resource_details,
+ "resource_category": category["resource_category_id"],
+ "public": "off"
+ }).then(lambda cr_results: (user, linkeddata, cr_results))
+
+ def __attach_data_to_resource__(user, linkeddata, resource):
+ logger.info("… attaching data to authorisation resource object")
+ return mrequests.post(
+ authserveruri("/auth/resource/data/link"),
+ headers=headers,
+ json={
+ "dataset_type": "phenotype",
+ "resource_id": resource["resource_id"],
+ "data_link_ids": [
+ item["data_link_id"] for item in linkeddata["traits"]]
+ }).then(lambda attc: (user, linkeddata, resource, attc))
+
+ def __handle_error__(resp):
+ error = resp.json()
+ if error.get("error") == "IntegrityError":
+ # This is hacky. If the auth already exists, something went wrong
+ # somewhere.
+ # This needs investigation to recover correctly.
+ logger.error(
+ "Error: The authorisation for the data was already set up.")
+ return 0
+ logger.error("ERROR: Updating the authorisation for the data failed.")
+ logger.debug(
+ "ERROR: The response from the authorisation server was:\n\t%s",
+ error)
+ return 1
+
+ def __handle_success__(_val):
+ logger.info(
+ "The authorisation for the data has been updated successfully.")
+ return 0
+
+ return __fetch_user_details__().then(__link_data__).then(
+ lambda result: __fetch_phenotype_category_details__(*result)
+ ).then(
+ lambda result: __create_resource__(*result)
+ ).then(
+ lambda result: __attach_data_to_resource__(*result)
+ ).either(__handle_error__, __handle_success__)
+
+
+def load_data(# pylint: disable=[too-many-locals]
+ conn: mysqldb.Connection, job: dict
+) -> tuple[dict, dict, dict, tuple[int, ...]]:
+ """Load the data attached in the given job."""
+ _job_metadata = job["metadata"]
+ # Steps
+ # 0. Read data from the files: can be multiple files per type
+ #
+ _species = species_by_id(conn, int(_job_metadata["species_id"]))
+ _population = population_by_species_and_id(
+ conn,
+ _species["SpeciesId"],
+ int(_job_metadata["population_id"]))
+ _dataset = dataset_by_id(
+ conn,
+ _species["SpeciesId"],
+ _population["Id"],
+ int(_job_metadata["dataset_id"]))
+ # 1. Just retrive the publication: Don't create publications for now.
+ _publication = fetch_publication_by_id(
+ conn, int(_job_metadata.get("publication_id", "0"))) or {"Id": 0}
+ # 2. Save all new phenotypes:
+ # -> return phenotype IDs
+ bundle = Path(_job_metadata["bundle_file"])
+ _control_data = rqtl2.control_data(bundle)
+ logger.info("Extracting the zipped bundle of files.")
+ _outdir = Path(bundle.parent, f"bundle_{bundle.stem}")
+ with ZipFile(str(bundle), "r") as zfile:
+ _files = rqtl2.extract(zfile, _outdir)
+ logger.info("Saving new phenotypes.")
+ _phenos = save_phenotypes(conn,
+ _control_data,
+ _population["Id"],
+ _publication["Id"],
+ _outdir)
+
+ def __build_phenos_maps__(accumulator, row):
+ return ({
+ **accumulator[0],
+ row["phenotype_id"]: {
+ "population_id": _population["Id"],
+ "phenotype_id": row["phenotype_id"],
+ "data_id": row["data_id"],
+ "publication_id": row["publication_id"],
+ }
+ }, {
+ **accumulator[1],
+ row["pre_publication_abbreviation"]: row["phenotype_id"]
+ }, (
+ accumulator[2] + ({
+ "xref_id": row["xref_id"],
+ "population_id": row["population_id"],
+ "phenotype_id": row["phenotype_id"],
+ "publication_id": row["publication_id"],
+ "data_id": row["data_id"]
+ },)))
+ dataidmap, pheno_name2id, _xrefs = reduce(# type: ignore[var-annotated]
+ __build_phenos_maps__, _phenos, ({},{}, tuple()))
+ # 3. a. Fetch the strain names and IDS: create name->ID map
+ samples = {
+ row["Name"]: row
+ for row in samples_by_species_and_population(
+ conn, _species["SpeciesId"], _population["Id"])}
+ # b. Save all the data items (DataIds are vibes), return new IDs
+ logger.info("Saving new phenotypes data.")
+ _num_data_rows = save_pheno_data(conn=conn,
+ dataidmap=dataidmap,
+ pheno_name2id=pheno_name2id,
+ samples=samples,
+ control_data=_control_data,
+ filesdir=_outdir)
+ logger.info("Saved %s new phenotype data rows.", _num_data_rows)
+
+ # 4. If standard errors and N exist, save them too
+ # (use IDs returned in `3. b.` above).
+ if _control_data.get("phenose"):
+ logger.info("Saving new phenotypes standard errors.")
+ _num_se_rows = save_phenotypes_se(conn=conn,
+ dataidmap=dataidmap,
+ pheno_name2id=pheno_name2id,
+ samples=samples,
+ control_data=_control_data,
+ filesdir=_outdir)
+ logger.info("Saved %s new phenotype standard error rows.", _num_se_rows)
+
+ if _control_data.get("phenonum"):
+ logger.info("Saving new phenotypes sample counts.")
+ _num_n_rows = save_phenotypes_n(conn=conn,
+ dataidmap=dataidmap,
+ pheno_name2id=pheno_name2id,
+ samples=samples,
+ control_data=_control_data,
+ filesdir=_outdir)
+ logger.info("Saved %s new phenotype sample counts rows.", _num_n_rows)
+
+ return (_species, _population, _dataset, _xrefs)
+
+
+def update_means(
+ conn: mysqldb.Connection,
+ population_id: int,
+ xref_ids: tuple[int, ...]
+):
+ """Compute the means from the data and update them in the database."""
+ logger.info("Computing means for %02d phenotypes.", len(xref_ids))
+ query = (
+ "UPDATE PublishXRef SET mean = "
+ "(SELECT AVG(value) FROM PublishData"
+ " WHERE PublishData.Id=PublishXRef.DataId) "
+ "WHERE PublishXRef.Id=%(xref_id)s "
+ "AND PublishXRef.InbredSetId=%(population_id)s")
+ _xref_iterator = (_xref_id for _xref_id in xref_ids)
+ with conn.cursor(cursorclass=DictCursor) as cursor:
+ while True:
+ batch = take(_xref_iterator, 10000)
+ if len(batch) == 0:
+ break
+ logger.info("\tComputing means for batch of %02d phenotypes.", len(batch))
+ cursor.executemany(
+ query,
+ tuple({
+ "population_id": population_id,
+ "xref_id": _xref_id
+ } for _xref_id in batch))
+
+
+if __name__ == "__main__":
+ def parse_args():
+ """Setup command-line arguments."""
+ parser = argparse.ArgumentParser(
+ prog="load_phenotypes_to_db",
+ description="Process the phenotypes' data and load it into the database.")
+ parser.add_argument("db_uri", type=str, help="MariaDB/MySQL connection URL")
+ parser.add_argument(
+ "jobs_db_path", type=Path, help="Path to jobs' SQLite database.")
+ parser.add_argument("job_id", type=uuid.UUID, help="ID of the running job")
+ parser.add_argument(
+ "--log-level",
+ type=str,
+ help="Determines what is logged out.",
+ choices=("debug", "info", "warning", "error", "critical"),
+ default="info")
+ return parser.parse_args()
+
+ def setup_logging(log_level: str):
+ """Setup logging for the script."""
+ logger.setLevel(log_level)
+ logging.getLogger("uploader.phenotypes.models").setLevel(log_level)
+
+
+ def __parse_resource_details__(meta) -> dict:
+ """Parse out details regarding the wrapper resource from the metadata."""
+ _key_mappings_ = {
+ # allow both 'data_*' and 'data*' for the metadata.
+ "data_description": "description",
+ "datadescription": "description"
+ }
+ return {
+ "resource_name": meta.get(
+ "dataname",
+ meta.get("data_name",
+ "Unnamed phenotypes - " + datetime.now().isoformat())),
+ "resource_metadata": {
+ rkey: meta[mkey]
+ for mkey, rkey in _key_mappings_.items() if mkey in meta
+ }
+ }
+
+
+ def main():
+ """Entry-point for this script."""
+ args = parse_args()
+ setup_logging(args.log_level.upper())
+
+ with (mysqldb.database_connection(args.db_uri) as conn,
+ conn.cursor(cursorclass=DictCursor) as cursor,
+ sqlite3.connection(args.jobs_db_path) as jobs_conn):
+ job = jobs.job(jobs_conn, args.job_id)
+
+ # Lock the PublishXRef/PublishData/PublishSE/NStrain here: Why?
+ # The `DataId` values are sequential, but not auto-increment
+ # Can't convert `PublishXRef`.`DataId` to AUTO_INCREMENT.
+ # `SELECT MAX(DataId) FROM PublishXRef;`
+ # How do you check for a table lock?
+ # https://oracle-base.com/articles/mysql/mysql-identify-locked-tables
+ # `SHOW OPEN TABLES LIKE 'Publish%';`
+ _db_tables_ = (
+ "Species",
+ "InbredSet",
+ "Strain",
+ "StrainXRef",
+ "Publication",
+ "Phenotype",
+ "PublishXRef",
+ "PublishFreeze",
+ "PublishData",
+ "PublishSE",
+ "NStrain")
+
+ logger.debug(
+ ("Locking database tables for the connection:" +
+ "".join("\n\t- %s" for _ in _db_tables_) + "\n"),
+ *_db_tables_)
+ cursor.execute(# Lock the tables to avoid race conditions
+ "LOCK TABLES " + ", ".join(
+ f"{_table} WRITE" for _table in _db_tables_))
+
+ db_results = load_data(conn, job)
+ _xref_ids = tuple(xref["xref_id"] for xref in db_results[3])
+ jobs.update_metadata(
+ jobs_conn,
+ args.job_id,
+ "xref_ids",
+ json.dumps(_xref_ids))
+
+ logger.info("Unlocking all database tables.")
+ cursor.execute("UNLOCK TABLES")
+
+ logger.info("Updating means.")
+ update_means(conn, db_results[1]["Id"], _xref_ids)
+
+ # Update authorisations (break this down) — maybe loop until it works?
+ logger.info("Updating authorisation.")
+ _job_metadata = job["metadata"]
+
+ return update_auth((_job_metadata["authserver"],
+ _job_metadata["token"]),
+ __parse_resource_details__(_job_metadata),
+ *db_results)
+
+
+ try:
+ sys.exit(main())
+ except Exception as _exc:# pylint: disable=[broad-exception-caught]
+ logger.debug("Data loading failed… Halting!",
+ exc_info=True)
+ sys.exit(1)
diff --git a/scripts/phenotypes/__init__.py b/scripts/phenotypes/__init__.py
new file mode 100644
index 0000000..73ad839
--- /dev/null
+++ b/scripts/phenotypes/__init__.py
@@ -0,0 +1 @@
+"Scripts for dealing with phenotypes."
diff --git a/scripts/phenotypes/delete_phenotypes.py b/scripts/phenotypes/delete_phenotypes.py
new file mode 100644
index 0000000..461f3ec
--- /dev/null
+++ b/scripts/phenotypes/delete_phenotypes.py
@@ -0,0 +1,173 @@
+"""Delete phenotypes."""
+import sys
+import logging
+from pathlib import Path
+from typing import Optional
+from urllib.parse import urljoin
+from argparse import Namespace, ArgumentParser
+
+import requests
+from MySQLdb.cursors import DictCursor, BaseCursor
+
+from gn_libs.mysqldb import database_connection
+
+from uploader.phenotypes.models import delete_phenotypes
+from scripts.cli.logging import setup_logging
+from scripts.cli.options import (add_logging,
+ add_mariadb_uri,
+ add_population_id)
+
+logger = logging.getLogger(__name__)
+
+def read_xref_ids_file(filepath: Optional[Path]) -> tuple[int, ...]:
+ """Read the phenotypes' cross-reference IDS from file."""
+ if filepath is None:
+ return tuple()
+
+ logger.debug("Using file '%s' to retrieve XREF IDs for deletion.",
+ filepath.name)
+ _ids: tuple[int, ...] = tuple()
+ with filepath.open(mode="r") as infile:
+ for line in infile.readlines():
+ try:
+ _ids += (int(line.strip()),)
+ except TypeError:
+ pass
+
+ return _ids
+
+
+def fetch_all_xref_ids(
+ cursor: BaseCursor, population_id: int) -> tuple[int, ...]:
+ """Fetch all cross-reference IDs."""
+ cursor.execute("SELECT Id FROM PublishXRef WHERE InbredSetId=%s",
+ (population_id,))
+ return tuple(int(row["Id"]) for row in cursor.fetchall())
+
+
+def update_auth(
+ auth_details: tuple[str, str],
+ species_id: int,
+ population_id: int,
+ dataset_id: int,
+ xref_ids: tuple[int, ...] = tuple()
+):
+ """Update the authorisation server: remove items to delete."""
+ authserver, token = auth_details
+ resp = requests.post(
+ urljoin(authserver,
+ (f"/auth/data/phenotypes/{species_id}/{population_id}"
+ f"/{dataset_id}/delete")),
+ timeout=(9.13, 20),
+ headers={
+ "Authorization": f"Bearer {token}",
+ "Content-Type": "application/json"
+ },
+ json={"xref_ids": xref_ids})
+ resp.raise_for_status()
+
+
+def delete_the_phenotypes(
+ cursor: BaseCursor,
+ population_id: int,
+ xref_ids: tuple[int, ...] = tuple()) -> int:
+ """Process and delete the phenotypes."""
+ delete_phenotypes(cursor, population_id, xref_ids)
+
+ return 0
+
+if __name__ == "__main__":
+ def parse_args() -> Namespace:
+ """Parse CLI arguments."""
+ parser = add_logging(
+ add_population_id(
+ add_mariadb_uri(
+ ArgumentParser(
+ prog="delete-phenotypes",
+ description=(
+ "Script to delete phenotypes from the database.")))))
+ parser.add_argument(
+ "dataset_id",
+ metavar="DATASET-ID",
+ type=int,
+ help="The dataset identifier for phenotypes to delete.")
+ parser.add_argument(
+ "auth_server_uri",
+ metavar="AUTH-SERVER-URI",
+ type=str,
+ help="URI to the authorisation server.")
+ parser.add_argument(
+ "auth_token",
+ metavar="AUTH-TOKEN",
+ type=str,
+ help=("Token to use to update the authorisation system with the "
+ "deletions done."))
+ parser.add_argument(
+ "--xref_ids_file",
+ metavar="XREF-IDS-FILE",
+ type=Path,
+ help=("Path to a file with phenotypes cross-reference IDs to "
+ "delete."))
+ parser.add_argument(
+ "--delete-all",
+ action="store_true",
+ help=("If no 'XREF-IDS-FILE' is provided, this flag determines "
+ "whether or not all the phenotypes for the given population "
+ "will be deleted."))
+ return parser.parse_args()
+
+
+ def main():
+ """The `delete-phenotypes` script's entry point."""
+ args = parse_args()
+ setup_logging(logger, args.log_level.upper(), tuple())
+ with (database_connection(args.db_uri) as conn,
+ conn.cursor(cursorclass=DictCursor) as cursor):
+ xref_ids = read_xref_ids_file(args.xref_ids_file)
+ try:
+ assert not (len(xref_ids) > 0 and args.delete_all)
+ xref_ids = (fetch_all_xref_ids(cursor, args.population_id)
+ if args.delete_all else xref_ids)
+ logger.debug("Will delete %s phenotypes and related data",
+ len(xref_ids))
+ if len(xref_ids) == 0:
+ print("No cross-reference IDs were provided. Aborting.")
+ return 0
+
+ print("Updating authorisations: ", end="")
+ update_auth((args.auth_server_uri, args.auth_token),
+ args.species_id,
+ args.population_id,
+ args.dataset_id,
+ xref_ids)
+ print("OK.")
+ print("Deleting the data: ", end="")
+ delete_phenotypes(cursor, args.population_id, xref_ids=xref_ids)
+ print("OK.")
+ if args.xref_ids_file is not None:
+ print("Deleting temporary file: ", end="")
+ args.xref_ids_file.unlink()
+ print("OK.")
+
+ return 0
+ except AssertionError:
+ logger.error(
+ "'DELETE-ALL' and 'XREF-IDS' are mutually exclusive. "
+ "If you specify the list of XREF-IDS (in a file) to delete "
+ "and also specify to 'DELETE-ALL' phenotypes in the "
+ "population, we have no way of knowing what it is you want.")
+ return 1
+ except requests.exceptions.HTTPError as _exc:
+ resp = _exc.response
+ resp_data = resp.json()
+ logger.debug("%s: %s",
+ resp_data["error"],
+ resp_data["error_description"],
+ exc_info=True)
+ return 1
+ except Exception as _exc:# pylint: disable=[broad-exception-caught]
+ logger.debug("Failed while attempting to delete phenotypes.",
+ exc_info=True)
+ return 1
+
+ sys.exit(main())
diff --git a/scripts/process_rqtl2_bundle.py b/scripts/process_rqtl2_bundle.py
index 8b7a0fb..e2ce420 100644
--- a/scripts/process_rqtl2_bundle.py
+++ b/scripts/process_rqtl2_bundle.py
@@ -104,7 +104,7 @@ def process_bundle(dbconn: mdb.Connection,
rqtl2bundle=Path(meta["rqtl2-bundle-file"])),
logger)
if genoexit != 0:
- raise Exception("Processing 'geno' file failed.")
+ raise Exception("Processing 'geno' file failed.")# pylint: disable=[broad-exception-raised]
logger.debug(
"geno file processing completed successfully. (ExitCode: %s)",
genoexit)
@@ -122,7 +122,7 @@ def process_bundle(dbconn: mdb.Connection,
rqtl2bundle=Path(meta["rqtl2-bundle-file"])),
logger)
if phenoexit != 0:
- raise Exception("Processing 'pheno' file failed.")
+ raise Exception("Processing 'pheno' file failed.")# pylint: disable=[broad-exception-raised]
logger.debug(
"pheno file processing completed successfully. (ExitCode: %s)",
phenoexit)
diff --git a/scripts/qc_on_rqtl2_bundle.py b/scripts/qc_on_rqtl2_bundle.py
index 9f9248c..4e6ef00 100644
--- a/scripts/qc_on_rqtl2_bundle.py
+++ b/scripts/qc_on_rqtl2_bundle.py
@@ -40,7 +40,7 @@ def add_to_errors(rconn: Redis,
"""Add `errors` to a given list of errors"""
errs = tuple(dict(item) for item in set(
[dict2tuple(old) for old in
- json.loads(rconn.hget(fqjobid, key) or "[]")] +
+ json.loads(rconn.hget(fqjobid, key) or "[]")] +# type: ignore[arg-type]
[dict2tuple({"type": type(error).__name__, **error._asdict()})
for error in errors]))
rconn.hset(fqjobid, key, json.dumps(errs))
@@ -83,7 +83,8 @@ def retrieve_errors_with_progress(rconn: Redis,#pylint: disable=[too-many-locals
count = 0
checked = 0
cdata = rqtl2.control_data(zfile)
- rconn.hset(fqjobid, f"{filetype}-filesize", compute_filesize(zfile, filetype))
+ rconn.hset(
+ fqjobid, f"{filetype}-filesize", str(compute_filesize(zfile, filetype)))
def __update_processed__(value):
nonlocal checked
checked = checked + len(value)
@@ -104,7 +105,7 @@ def retrieve_errors_with_progress(rconn: Redis,#pylint: disable=[too-many-locals
yield error
__update_processed__(value)
- rconn.hset(fqjobid, f"{filetype}-linecount", count)
+ rconn.hset(fqjobid, f"{filetype}-linecount", count)# type: ignore[arg-type]
except rqe.MissingFileException:
fname = cdata.get(filetype)
yield rqfe.MissingFile(filetype, fname, (
@@ -191,7 +192,7 @@ def check_pheno_samples(
return allerrors
-def qc_pheno_errors(# pylint: disable=[too-many-arguments]
+def qc_pheno_errors(# pylint: disable=[too-many-arguments, too-many-positional-arguments]
rconn, fqjobid, dburi, speciesid, zfile, logger) -> bool:
"""Check for errors in `pheno` file(s)."""
cdata = rqtl2.control_data(zfile)
@@ -260,7 +261,7 @@ def run_qc(rconn: Redis,
if qc_missing_files(rconn, fqjobid, zfile, logger):
return 1
- def with_zipfile(# pylint: disable=[too-many-arguments]
+ def with_zipfile(# pylint: disable=[too-many-arguments,too-many-positional-arguments]
rconn, fqjobid, dbconn, speciesid, filename, logger, func
):
with ZipFile(filename, "r") as zfile:
@@ -295,7 +296,7 @@ def run_qc(rconn: Redis,
return 1
def __fetch_errors__(rkey: str) -> tuple:
- return tuple(json.loads(rconn.hget(fqjobid, rkey) or "[]"))
+ return tuple(json.loads(rconn.hget(fqjobid, rkey) or "[]")) # type: ignore[arg-type]
return (1 if any((
bool(__fetch_errors__(key))
diff --git a/scripts/redis_logger.py b/scripts/redis_logger.py
index d3fde5f..a74e5e4 100644
--- a/scripts/redis_logger.py
+++ b/scripts/redis_logger.py
@@ -6,7 +6,7 @@ from redis import Redis
class RedisLogger(logging.Handler):
"""Log out to redis for our worker scripts"""
- def __init__(self,#pylint: disable=[too-many-arguments]
+ def __init__(self,#pylint: disable=[too-many-arguments, too-many-positional-arguments]
rconn: Redis,
fullyqualifiedjobid: str,
messageslistname: str,
diff --git a/scripts/rqtl2/entry.py b/scripts/rqtl2/entry.py
index 327ed2c..7423a4b 100644
--- a/scripts/rqtl2/entry.py
+++ b/scripts/rqtl2/entry.py
@@ -20,27 +20,23 @@ def build_main(
[Redis, Connection, str, Namespace, logging.Logger],
int
],
- loggername: str
+ logger: logging.Logger
) -> Callable[[],int]:
"""Build a function to be used as an entry-point for scripts."""
def main():
- try:
- logging.basicConfig(
- format=(
- "%(asctime)s - %(levelname)s %(name)s: "
- "(%(pathname)s: %(lineno)d) %(message)s"),
- level=args.loglevel)
- logger = logging.getLogger(loggername)
- with (Redis.from_url(args.redisuri, decode_responses=True) as rconn,
- database_connection(args.databaseuri) as dbconn):
- fqjobid = jobs.job_key(args.redisprefix, args.jobid)
+ with (Redis.from_url(args.redisuri, decode_responses=True) as rconn,
+ database_connection(args.databaseuri) as dbconn):
+ logger.setLevel(args.log_level.upper())
+ fqjobid = jobs.job_key(args.redisprefix, args.jobid)
+
+ try:
rconn.hset(fqjobid, "status", "started")
logger.addHandler(setup_redis_logger(
rconn,
fqjobid,
f"{fqjobid}:log-messages",
args.redisexpiry))
- logger.addHandler(StreamHandler(stream=sys.stdout))
+ logger.addHandler(StreamHandler(stream=sys.stderr))
check_db(args.databaseuri)
check_redis(args.redisuri)
@@ -48,15 +44,15 @@ def build_main(
logger.error("File not found: '%s'.", args.rqtl2bundle)
return 2
- returncode = run_fn(rconn, dbconn, fqjobid, args, logger)
+ returncode = run_fn(rconn, dbconn, fqjobid, args)
if returncode == 0:
rconn.hset(fqjobid, "status", "completed:success")
return returncode
rconn.hset(fqjobid, "status", "completed:error")
return returncode
- except Exception as _exc:# pylint: disable=[broad-except]
- logger.error("The process failed!", exc_info=True)
- rconn.hset(fqjobid, "status", "completed:error")
- return 4
+ except Exception as _exc:# pylint: disable=[broad-except]
+ logger.error("The process failed!", exc_info=True)
+ rconn.hset(fqjobid, "status", "completed:error")
+ return 4
return main
diff --git a/scripts/rqtl2/install_genotypes.py b/scripts/rqtl2/install_genotypes.py
index 8762655..5e6abb0 100644
--- a/scripts/rqtl2/install_genotypes.py
+++ b/scripts/rqtl2/install_genotypes.py
@@ -20,7 +20,7 @@ from scripts.rqtl2.entry import build_main
from scripts.rqtl2.cli_parser import add_common_arguments
from scripts.cli_parser import init_cli_parser, add_global_data_arguments
-__MODULE__ = "scripts.rqtl2.install_genotypes"
+logger = getLogger(__name__)
def insert_markers(
dbconn: mdb.Connection,
@@ -191,7 +191,7 @@ def install_genotypes(#pylint: disable=[too-many-locals]
dbconn: mdb.Connection,
fullyqualifiedjobid: str,#pylint: disable=[unused-argument]
args: argparse.Namespace,
- logger: Logger = getLogger(__name__)
+ logger: Logger = logger # pylint: disable=[redefined-outer-name]
) -> int:
"""Load any existing genotypes into the database."""
(speciesid, populationid, datasetid, rqtl2bundle) = (
@@ -257,5 +257,5 @@ if __name__ == "__main__":
return parser.parse_args()
- main = build_main(cli_args(), install_genotypes, __MODULE__)
+ main = build_main(cli_args(), install_genotypes, logger)
sys.exit(main())
diff --git a/scripts/rqtl2/install_phenos.py b/scripts/rqtl2/install_phenos.py
index 9059cd6..11ac8a4 100644
--- a/scripts/rqtl2/install_phenos.py
+++ b/scripts/rqtl2/install_phenos.py
@@ -19,7 +19,7 @@ from r_qtl import r_qtl2_qc as rqc
from functional_tools import take
-__MODULE__ = "scripts.rqtl2.install_phenos"
+logger = getLogger(__name__)
def insert_probesets(dbconn: mdb.Connection,
platformid: int,
@@ -101,7 +101,8 @@ def install_pheno_files(#pylint: disable=[too-many-locals]
dbconn: mdb.Connection,
fullyqualifiedjobid: str,#pylint: disable=[unused-argument]
args: argparse.Namespace,
- logger: Logger = getLogger()) -> int:
+ logger: Logger = logger # pylint: disable=[redefined-outer-name]
+) -> int:
"""Load data in `pheno` files and other related files into the database."""
(speciesid, platformid, datasetid, rqtl2bundle) = (
args.speciesid, args.platformid, args.datasetid, args.rqtl2bundle)
@@ -159,5 +160,5 @@ if __name__ == "__main__":
return parser.parse_args()
- main = build_main(cli_args(), install_pheno_files, __MODULE__)
+ main = build_main(cli_args(), install_pheno_files, logger)
sys.exit(main())
diff --git a/scripts/rqtl2/phenotypes_qc.py b/scripts/rqtl2/phenotypes_qc.py
index 76ecb8d..084c876 100644
--- a/scripts/rqtl2/phenotypes_qc.py
+++ b/scripts/rqtl2/phenotypes_qc.py
@@ -36,8 +36,15 @@ from scripts.cli_parser import init_cli_parser, add_global_data_arguments
from scripts.rqtl2.bundleutils import build_line_joiner, build_line_splitter
__MODULE__ = "scripts.rqtl2.phenotypes_qc"
+logging.basicConfig(
+ format=("%(asctime)s - %(levelname)s %(name)s: "
+ "(%(pathname)s: %(lineno)d) %(message)s"))
+logger = logging.getLogger(__MODULE__)
-def validate(phenobundle: Path, logger: Logger) -> dict:
+def validate(
+ phenobundle: Path,
+ logger: Logger# pylint: disable=[redefined-outer-name]
+) -> dict:
"""Check that the bundle is generally valid"""
try:
rqc.validate_bundle(phenobundle)
@@ -59,7 +66,7 @@ def validate(phenobundle: Path, logger: Logger) -> dict:
def check_for_mandatory_pheno_keys(
phenobundle: Path,
- logger: Logger,
+ logger: Logger,# pylint: disable=[redefined-outer-name]
**kwargs
) -> dict:
"""Check that the mandatory keys exist for phenotypes."""
@@ -86,7 +93,7 @@ def check_for_mandatory_pheno_keys(
def check_for_averages_files(
phenobundle: Path,
- logger: Logger,
+ logger: Logger,# pylint: disable=[redefined-outer-name]
**kwargs
) -> dict:
"""Check that averages files appear together"""
@@ -140,15 +147,15 @@ def redis_logger(
) -> Iterator[logging.Logger]:
"""Build a Redis message-list logger."""
rconn = Redis.from_url(redisuri, decode_responses=True)
- logger = logging.getLogger(loggername)
- logger.propagate = False
+ _logger = logging.getLogger(loggername)
+ _logger.propagate = False
handler = RedisMessageListHandler(
rconn,
fullyqualifiedkey(fqkey, filename))#type: ignore[arg-type]
handler.setFormatter(logging.getLogger().handlers[0].formatter)
- logger.addHandler(handler)
+ _logger.addHandler(handler)
try:
- yield logger
+ yield _logger
finally:
rconn.close()
@@ -175,9 +182,9 @@ def qc_phenocovar_file(
redisuri,
f"{__MODULE__}.qc_phenocovar_file",
filepath.name,
- f"{fqkey}:logs") as logger,
+ f"{fqkey}:logs") as _logger,
Redis.from_url(redisuri, decode_responses=True) as rconn):
- logger.info("Running QC on file: %s", filepath.name)
+ print("Running QC on file: ", filepath.name)
_csvfile = rqtl2.read_csv_file(filepath, separator, comment_char)
_headings = tuple(heading.lower() for heading in next(_csvfile))
_errors: tuple[InvalidValue, ...] = tuple()
@@ -191,11 +198,11 @@ def qc_phenocovar_file(
"-",
"-",
(f"File {filepath.name} is missing the {heading} heading "
- "in the header line."))),)
+ "in the header row/line."))),)
def collect_errors(errors_and_linecount, line):
_errs, _lc = errors_and_linecount
- logger.info("Testing record '%s'", line[0])
+ _logger.info("Testing record '%s'", line[0])
if len(line) != len(_headings):
_errs = _errs + (save_error(InvalidValue(
filepath.name,
@@ -205,12 +212,12 @@ def qc_phenocovar_file(
(f"Record {_lc} in file {filepath.name} has a different "
"number of columns than the number of headings"))),)
_line = dict(zip(_headings, line))
- if not bool(_line["description"]):
+ if not bool(_line.get("description")):
_errs = _errs + (
save_error(InvalidValue(filepath.name,
_line[_headings[0]],
"description",
- _line["description"],
+ _line.get("description"),
"The description is not provided!")),)
rconn.hset(file_fqkey(fqkey, "metadata", filepath),
@@ -236,7 +243,7 @@ def merge_dicts(*dicts):
return reduce(lambda merged, dct: {**merged, **dct}, dicts, {})
-def decimal_points_error(# pylint: disable=[too-many-arguments]
+def decimal_points_error(# pylint: disable=[too-many-arguments,too-many-positional-arguments]
filename: str,
rowtitle: str,
coltitle: str,
@@ -267,7 +274,7 @@ def integer_error(
return InvalidValue(filename, rowtitle, coltitle, cellvalue, message)
-def qc_pheno_file(# pylint: disable=[too-many-locals, too-many-arguments]
+def qc_pheno_file(# pylint: disable=[too-many-locals, too-many-arguments, too-many-positional-arguments]
filepath: Path,
redisuri: str,
fqkey: str,
@@ -283,9 +290,9 @@ def qc_pheno_file(# pylint: disable=[too-many-locals, too-many-arguments]
redisuri,
f"{__MODULE__}.qc_pheno_file",
filepath.name,
- f"{fqkey}:logs") as logger,
+ f"{fqkey}:logs") as _logger,
Redis.from_url(redisuri, decode_responses=True) as rconn):
- logger.info("Running QC on file: %s", filepath.name)
+ print("Running QC on file: ", filepath.name)
save_error = partial(
push_error, rconn, file_fqkey(fqkey, "errors", filepath))
_csvfile = rqtl2.read_csv_file(filepath, separator, comment_char)
@@ -305,12 +312,13 @@ def qc_pheno_file(# pylint: disable=[too-many-locals, too-many-arguments]
"header row",
"-",
", ".join(_absent),
- ("The following phenotype names do not exist in any of the "
- f"provided phenocovar files: ({', '.join(_absent)})"))),)
+ ("The following trait names/identifiers do not exist in any of "
+ "the provided descriptions/covariates files: "
+ f"({', '.join(_absent)})"))),)
def collect_errors(errors_and_linecount, line):
_errs, _lc = errors_and_linecount
- logger.debug("Checking row %s", line[0])
+ _logger.debug("Checking row %s", line[0])
if line[0] not in samples:
_errs = _errs + (save_error(InvalidValue(
filepath.name,
@@ -370,10 +378,10 @@ def run_qc(# pylint: disable=[too-many-locals]
dbconn: mdb.Connection,
fullyqualifiedjobid: str,
args: Namespace,
- logger: Logger
+ logger: Logger = logger # pylint: disable=[redefined-outer-name]
) -> int:
"""Run quality control checks on the bundle."""
- logger.debug("Beginning the quality assurance checks.")
+ print("Beginning the quality assurance checks.")
results = check_for_averages_files(
**check_for_mandatory_pheno_keys(
**validate(args.rqtl2bundle, logger)))
@@ -398,7 +406,7 @@ def run_qc(# pylint: disable=[too-many-locals]
for ftype in ("pheno", "phenocovar", "phenose", "phenonum")))
# - Fetch samples/individuals from database.
- logger.debug("Fetching samples/individuals from the database.")
+ print("Fetching samples/individuals from the database.")
samples = tuple(#type: ignore[var-annotated]
item for item in set(reduce(
lambda acc, item: acc + (
@@ -415,7 +423,7 @@ def run_qc(# pylint: disable=[too-many-locals]
json.dumps(tuple(f"{fullyqualifiedjobid}:phenocovar:{_file}"
for _file in cdata.get("phenocovar", []))))
with mproc.Pool(mproc.cpu_count() - 1) as pool:
- logger.debug("Check for errors in 'phenocovar' file(s).")
+ print("Check for errors in 'phenocovar' file(s).")
_phenocovar_qc_res = merge_dicts(*pool.starmap(qc_phenocovar_file, tuple(
(extractiondir.joinpath(_file),
args.redisuri,
@@ -437,7 +445,7 @@ def run_qc(# pylint: disable=[too-many-locals]
"Expected a non-negative number with at least one decimal "
"place."))
- logger.debug("Check for errors in 'pheno' file(s).")
+ print("Check for errors in 'pheno' file(s).")
_pheno_qc_res = merge_dicts(*pool.starmap(qc_pheno_file, tuple((
extractiondir.joinpath(_file),
args.redisuri,
@@ -456,7 +464,7 @@ def run_qc(# pylint: disable=[too-many-locals]
# - Check the 3 checks above for phenose and phenonum values too
# qc_phenose_files(…)
# qc_phenonum_files(…)
- logger.debug("Check for errors in 'phenose' file(s).")
+ print("Check for errors in 'phenose' file(s).")
_phenose_qc_res = merge_dicts(*pool.starmap(qc_pheno_file, tuple((
extractiondir.joinpath(_file),
args.redisuri,
@@ -472,7 +480,7 @@ def run_qc(# pylint: disable=[too-many-locals]
dec_err_fn
) for _file in cdata.get("phenose", []))))
- logger.debug("Check for errors in 'phenonum' file(s).")
+ print("Check for errors in 'phenonum' file(s).")
_phenonum_qc_res = merge_dicts(*pool.starmap(qc_pheno_file, tuple((
extractiondir.joinpath(_file),
args.redisuri,
@@ -509,5 +517,5 @@ if __name__ == "__main__":
type=Path)
return parser.parse_args()
- main = build_main(cli_args(), run_qc, __MODULE__)
+ main = build_main(cli_args(), run_qc, logger)
sys.exit(main())
diff --git a/scripts/run_qtlreaper.py b/scripts/run_qtlreaper.py
new file mode 100644
index 0000000..a461d9a
--- /dev/null
+++ b/scripts/run_qtlreaper.py
@@ -0,0 +1,238 @@
+"""Script to run rust-qtlreaper and update database with results."""
+import os
+import sys
+import csv
+import time
+import secrets
+import logging
+import subprocess
+import multiprocessing
+from pathlib import Path
+from functools import reduce
+from typing import Union, Iterator
+from argparse import Namespace, ArgumentParser
+
+from gn_libs import mysqldb
+
+from uploader.phenotypes.models import phenotypes_vector_data
+from uploader.population.models import population_by_species_and_id
+from uploader.samples.models import samples_by_species_and_population
+
+from scripts.cli.logging import setup_logging
+from scripts.cli.validators import directory_exists
+from scripts.cli.options import add_logging, add_mariadb_uri, add_population_id
+
+logger = logging.getLogger(__name__)
+
+
+def retrieve_genotype_file(genotypes_dir: Path, population_code: str) -> Path:
+ """Retrieves the genotype file"""
+ _genofile = genotypes_dir.joinpath(f"{population_code}.geno")
+ if _genofile.exists():
+ return _genofile
+ raise FileNotFoundError(f"Could not find the genotype file '{population_code}.geno'")
+
+
+def samples_from_genofile(genofile: Path) -> tuple[str, ...]:
+ """Read samples from the genotype file."""
+ with genofile.open(mode="r", encoding="utf-8") as inptr:
+ while True:
+ line = inptr.readline()
+ if (line.startswith("#") # comment line
+ or line.startswith("@") # allele? heterozygosity?
+ or line.strip() == "" # empty line
+ ):
+ continue
+ return tuple(line.strip().split("\t")[4:])
+
+
+def reconcile_samples(
+ genosamples: tuple[str, ...],
+ dbsamples: tuple[str, ...]
+) -> tuple[tuple[str, ...], tuple[str, ...]]:
+ """merge samples in genosamples and dbsamples and retain order in genosamples."""
+ in_db_not_geno = set(dbsamples).difference(genosamples)
+ return genosamples, tuple(in_db_not_geno)
+
+
+def generate_qtlreaper_traits_file(
+ outdir: Path,
+ samples: tuple[str, ...],
+ traits_data: tuple[dict[str, Union[int, float]], ...],
+ filename_prefix: str = ""
+) -> Path:
+ """Generate a file for use with qtlreaper that contains the traits' data."""
+ _dialect = csv.unix_dialect()
+ _dialect.delimiter="\t"
+ _dialect.quoting=0
+
+ _traitsfile = outdir.joinpath(
+ f"{filename_prefix}_{secrets.token_urlsafe(15)}.tsv")#type: ignore[attr-defined]
+ with _traitsfile.open(mode="w", encoding="utf-8") as outptr:
+ writer = csv.DictWriter(
+ outptr, fieldnames=("Trait",) + samples, dialect=_dialect)
+ writer.writeheader()
+ for row in traits_data:
+ writer.writerow({
+ "Trait": row["xref_id"],
+ **{sample: row.get(sample, "") for sample in samples}
+ })
+
+ return _traitsfile
+
+
+def parse_tsv_file(results_file: Path) -> Iterator[dict]:
+ """Parse the rust-qtlreaper output into usable python objects."""
+ with results_file.open("r", encoding="utf-8") as readptr:
+ _dialect = csv.unix_dialect()
+ _dialect.delimiter = "\t"
+ reader = csv.DictReader(readptr, dialect=_dialect)
+ yield from reader
+
+
+def __qtls_by_trait__(qtls, current):
+ """Organise QTL results by trait"""
+ return {
+ **qtls,
+ current["ID"]: qtls.get(current["ID"], tuple()) + (current,)
+ }
+
+
+def save_qtl_values_to_db(conn, qtls: tuple[dict, ...]):
+ """Save computed QTLs to the database."""
+ with conn.cursor() as cursor:
+ cursor.executemany(
+ "UPDATE PublishXRef SET "
+ "Locus=%(Locus)s, LRS=%(LRS)s, additive=%(Additive)s "
+ "WHERE Id=%(ID)s",
+ qtls)
+
+
+def dispatch(args: Namespace) -> int:# pylint: disable=[too-many-locals]
+ """Dispatch the actual logic."""
+ exitcode = 1
+ with mysqldb.database_connection(args.db_uri) as conn:
+ try:
+ population = population_by_species_and_id(conn, args.species_id, args.population_id)
+ assert population, (f"No population with ID '{args.population_id} for "
+ f"species with ID '{args.species_id}'.")
+ _genofile = retrieve_genotype_file(args.genotypes_dir, population["Name"])
+ logger.debug("Genotype file: %s", _genofile)
+ samples, _samples_not_in_genofile = reconcile_samples(
+ samples_from_genofile(_genofile),
+ tuple(
+ sample["Name"] for sample in
+ samples_by_species_and_population(
+ conn, args.species_id, args.population_id)))
+ if len(_samples_not_in_genofile) > 0:
+ logger.warning(
+ "Ignoring %d samples that are in the database but not in "
+ "the provided genotype file.",
+ len(_samples_not_in_genofile))
+ logger.debug("Ignored the following samples: %s",
+ ", ".join(_samples_not_in_genofile))
+
+ # Fetch traits data: provided list, or all traits in db
+ _traitsdata = tuple(phenotypes_vector_data(
+ conn,
+ args.species_id,
+ args.population_id,
+ xref_ids=tuple(args.xref_ids)).values())
+ logger.debug("Successfully got traits data. Generating the QTLReaper's traits file…")
+ _traitsfile = generate_qtlreaper_traits_file(
+ args.working_dir,
+ samples,
+ _traitsdata,
+ filename_prefix="qtlreaper_input_traits_file")
+ logger.debug("QTLReaper's Traits file: %s", _traitsfile)
+
+ _qtlreaper_main_output = args.working_dir.joinpath(
+ f"main-output-{secrets.token_urlsafe(15)}.tsv")#type: ignore[attr-defined]
+ _qtlreaper_permu_output = args.working_dir.joinpath(
+ f"permu-output-{secrets.token_urlsafe(15)}.tsv")
+ logger.debug("Main output filename: %s", _qtlreaper_main_output)
+ with subprocess.Popen(
+ ("qtlreaper",
+ "--n_permutations", "1000",
+ "--geno", _genofile,
+ "--traits", _traitsfile,
+ "--main_output", _qtlreaper_main_output,
+ "--permu_output", _qtlreaper_permu_output,
+ "--threads", str(int(1+(multiprocessing.cpu_count()/2)))),
+ env=({**os.environ, "RUST_BACKTRACE": "full"}
+ if logger.getEffectiveLevel() == logging.DEBUG
+ else dict(os.environ))) as _qtlreaper:
+ while _qtlreaper.poll() is None:
+ logger.debug("QTLReaper process running…")
+ time.sleep(1)
+ results = (
+ tuple(#type: ignore[var-annotated]
+ max(qtls, key=lambda qtl: qtl["LRS"])
+ for qtls in
+ reduce(__qtls_by_trait__,
+ parse_tsv_file(_qtlreaper_main_output),
+ {}).values())
+ if _qtlreaper_main_output.exists()
+ else tuple())
+ logger.debug("Cleaning up temporary files.")
+
+ # short-circuits to delete file if exists
+ if _traitsfile.exists():
+ _traitsfile.unlink()
+ logger.info("Deleted generated traits' file for QTLReaper.")
+
+ if _qtlreaper_main_output.exists():
+ _qtlreaper_main_output.unlink()
+ logger.info("Deleted QTLReaper's main output file.")
+
+ if _qtlreaper_permu_output.exists():
+ _qtlreaper_permu_output.unlink()
+ logger.info("Deleted QTLReaper's permutations file.")
+
+ if _qtlreaper.returncode != 0:
+ return _qtlreaper.returncode
+
+ save_qtl_values_to_db(conn, results)
+ logger.info("Successfully computed p values for %s traits.", len(_traitsdata))
+ return 0
+ except FileNotFoundError as fnf:
+ logger.error(", ".join(str(arg) for arg in fnf.args), exc_info=False)
+ except AssertionError as aserr:
+ logger.error(", ".join(aserr.args), exc_info=False)
+ except Exception as _exc:# pylint: disable=[broad-exception-caught]
+ logger.debug("Type of exception: %s", type(_exc))
+ logger.error("General exception!", exc_info=True)
+
+ return exitcode
+
+
+if __name__ == "__main__":
+ def main():
+ """run_qtlreaper.py: entry point."""
+ parser = add_logging(add_population_id(add_mariadb_uri(
+ ArgumentParser("run_qtlreaper"))))
+ parser.add_argument(
+ "genotypes_dir",
+ metavar="GENOTYPES-DIRECTORY",
+ type=directory_exists,
+ help="Path to directory with the genotypes.")
+ parser.add_argument(
+ "working_dir",
+ metavar="WORKING-DIRECTORY",
+ type=directory_exists,
+ help="Directory where the script will write temporary files.")
+ parser.add_argument(
+ "xref_ids",
+ metavar="CROSS-REFERENCE-IDS",
+ type=int,
+ nargs="*",
+ help=("Optional list of specific cross-reference IDs to narrow down"
+ " to. If provided, QTLReaper will only run against them. "
+ "If NOT provided, QTLReaper will run against all the traits "
+ "in the population."))
+ args = parser.parse_args()
+ setup_logging(logger, args.log_level)
+
+ return dispatch(args)
+
+ sys.exit(main())
diff --git a/scripts/worker.py b/scripts/worker.py
index 91b0332..3165fe7 100644
--- a/scripts/worker.py
+++ b/scripts/worker.py
@@ -79,7 +79,7 @@ def main():
fqjobid = jobs.job_key(args.redisprefix, args.jobid)
rconn.hset(fqjobid, "stderr", f"No such job. '{args.job_id}'.")
rconn.expire(name=jobs.job_key(args.redisprefix, args.job_id),
- time=timedelta(seconds=(2 * 60 * 60)))
+ time=timedelta(seconds=2 * 60 * 60))
print(f"No such job. '{args.job_id}'.", file=sys.stderr)
return 2
return 3
generated by cgit v1.2.3 (git 2.44.0) at 2026-09-24 10:10:31 +0000