diff options
Diffstat (limited to 'wqflask/base')
| -rw-r--r-- | wqflask/base/data_set.py | 1301 | ||||
| -rw-r--r-- | wqflask/base/data_set/__init__.py | 124 | ||||
| -rw-r--r-- | wqflask/base/data_set/dataset.py | 292 | ||||
| -rw-r--r-- | wqflask/base/data_set/datasetgroup.py | 193 | ||||
| -rw-r--r-- | wqflask/base/data_set/datasettype.py | 122 | ||||
| -rw-r--r-- | wqflask/base/data_set/genotypedataset.py | 75 | ||||
| -rw-r--r-- | wqflask/base/data_set/markers.py | 96 | ||||
| -rw-r--r-- | wqflask/base/data_set/mrnaassaydataset.py | 178 | ||||
| -rw-r--r-- | wqflask/base/data_set/phenotypedataset.py | 133 | ||||
| -rw-r--r-- | wqflask/base/data_set/probably_unused.py | 34 | ||||
| -rw-r--r-- | wqflask/base/data_set/tempdataset.py | 23 | ||||
| -rw-r--r-- | wqflask/base/data_set/utils.py | 77 |
12 files changed, 1347 insertions, 1301 deletions
diff --git a/wqflask/base/data_set.py b/wqflask/base/data_set.py deleted file mode 100644 index 9685e33e..00000000 --- a/wqflask/base/data_set.py +++ /dev/null @@ -1,1301 +0,0 @@ -# Copyright (C) University of Tennessee Health Science Center, Memphis, TN. -# -# This program is free software: you can redistribute it and/or modify it -# under the terms of the GNU Affero General Public License -# as published by the Free Software Foundation, either version 3 of the -# License, or (at your option) any later version. -# -# This program is distributed in the hope that it will be useful, -# but WITHOUT ANY WARRANTY; without even the implied warranty of -# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. -# See the GNU Affero General Public License for more details. -# -# This program is available from Source Forge: at GeneNetwork Project -# (sourceforge.net/projects/genenetwork/). -# -# Contact Drs. Robert W. Williams and Xiaodong Zhou (2010) -# at rwilliams@uthsc.edu and xzhou15@uthsc.edu -# -# This module is used by GeneNetwork project (www.genenetwork.org) -from dataclasses import dataclass -from dataclasses import field -from dataclasses import InitVar -from typing import Optional, Dict, List -from utility.tools import USE_REDIS, flat_file_exists, GN2_BASE_URL -from utility.db_tools import escape -from utility.db_tools import mescape -from utility.db_tools import create_in_clause -from maintenance import get_group_samplelists -from utility.tools import locate, locate_ignore_error, flat_files -from utility import gen_geno_ob -from utility import chunks -from utility import webqtlUtil -from db import webqtlDatabaseFunction -from base import species -from base import webqtlConfig -from base.webqtlConfig import TMPDIR -from urllib.parse import urlparse -from utility.tools import SQL_URI -from wqflask.database import database_connection -import os -import math -import collections -import codecs - -import json -import requests -import pickle as pickle -import hashlib -from redis import Redis - - -r = Redis() - -# Used by create_database to instantiate objects -# Each subclass will add to this -DS_NAME_MAP = {} - - -def create_dataset(dataset_name, dataset_type=None, - get_samplelist=True, group_name=None): - if dataset_name == "Temp": - dataset_type = "Temp" - - if not dataset_type: - dataset_type = Dataset_Getter(dataset_name) - - dataset_ob = DS_NAME_MAP[dataset_type] - dataset_class = globals()[dataset_ob] - if dataset_type == "Temp": - return dataset_class(dataset_name, get_samplelist, group_name) - else: - return dataset_class(dataset_name, get_samplelist) - - -@dataclass -class DatasetType: - """Create a dictionary of samples where the value is set to Geno, - Publish or ProbeSet. E.g. - - {'AD-cases-controls-MyersGeno': 'Geno', - 'AD-cases-controls-MyersPublish': 'Publish', - 'AKXDGeno': 'Geno', - 'AXBXAGeno': 'Geno', - 'AXBXAPublish': 'Publish', - 'Aging-Brain-UCIPublish': 'Publish', - 'All Phenotypes': 'Publish', - 'B139_K_1206_M': 'ProbeSet', - 'B139_K_1206_R': 'ProbeSet' ... - } - """ - redis_instance: InitVar[Redis] - datasets: Optional[Dict] = field(init=False, default_factory=dict) - data: Optional[Dict] = field(init=False) - - def __post_init__(self, redis_instance): - self.redis_instance = redis_instance - data = redis_instance.get("dataset_structure") - if data: - self.datasets = json.loads(data) - else: - # ZS: I don't think this should ever run unless Redis is - # emptied - try: - data = json.loads(requests.get( - GN2_BASE_URL + "/api/v_pre1/gen_dropdown", - timeout=5).content) - for _species in data['datasets']: - for group in data['datasets'][_species]: - for dataset_type in data['datasets'][_species][group]: - for dataset in data['datasets'][_species][group][dataset_type]: - short_dataset_name = dataset[1] - if dataset_type == "Phenotypes": - new_type = "Publish" - elif dataset_type == "Genotypes": - new_type = "Geno" - else: - new_type = "ProbeSet" - self.datasets[short_dataset_name] = new_type - except Exception: # Do nothing - pass - - self.redis_instance.set("dataset_structure", - json.dumps(self.datasets)) - self.data = data - - def set_dataset_key(self, t, name): - """If name is not in the object's dataset dictionary, set it, and - update dataset_structure in Redis - args: - t: Type of dataset structure which can be: 'mrna_expr', 'pheno', - 'other_pheno', 'geno' - name: The name of the key to inserted in the datasets dictionary - - """ - sql_query_mapping = { - 'mrna_expr': ("SELECT ProbeSetFreeze.Id FROM " - "ProbeSetFreeze WHERE " - "ProbeSetFreeze.Name = %s "), - 'pheno': ("SELECT InfoFiles.GN_AccesionId " - "FROM InfoFiles, PublishFreeze, InbredSet " - "WHERE InbredSet.Name = %s AND " - "PublishFreeze.InbredSetId = InbredSet.Id AND " - "InfoFiles.InfoPageName = PublishFreeze.Name"), - 'other_pheno': ("SELECT PublishFreeze.Name " - "FROM PublishFreeze, InbredSet " - "WHERE InbredSet.Name = %s AND " - "PublishFreeze.InbredSetId = InbredSet.Id"), - 'geno': ("SELECT GenoFreeze.Id FROM GenoFreeze WHERE " - "GenoFreeze.Name = %s ") - } - - dataset_name_mapping = { - "mrna_expr": "ProbeSet", - "pheno": "Publish", - "other_pheno": "Publish", - "geno": "Geno", - } - - group_name = name - if t in ['pheno', 'other_pheno']: - group_name = name.replace("Publish", "") - - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute(sql_query_mapping[t], (group_name,)) - if cursor.fetchone(): - self.datasets[name] = dataset_name_mapping[t] - self.redis_instance.set( - "dataset_structure", json.dumps(self.datasets)) - return True - - - def __call__(self, name): - if name not in self.datasets: - for t in ["mrna_expr", "pheno", "other_pheno", "geno"]: - # This has side-effects, with the end result being a - # truth-y value - if(self.set_dataset_key(t, name)): - break - # Return None if name has not been set - return self.datasets.get(name, None) - - -# Do the intensive work at startup one time only -Dataset_Getter = DatasetType(r) - - -def create_datasets_list(): - if USE_REDIS: - key = "all_datasets" - result = r.get(key) - - if result: - datasets = pickle.loads(result) - - if result is None: - datasets = list() - type_dict = {'Publish': 'PublishFreeze', - 'ProbeSet': 'ProbeSetFreeze', - 'Geno': 'GenoFreeze'} - - for dataset_type in type_dict: - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute("SELECT Name FROM %s", - (type_dict[dataset_type],)) - results = cursor.fetchall(query) - if results: - for result in results: - datasets.append( - create_dataset(result.Name, dataset_type)) - if USE_REDIS: - r.set(key, pickle.dumps(datasets, pickle.HIGHEST_PROTOCOL)) - r.expire(key, 60 * 60) - - return datasets - - -class Markers: - """Todo: Build in cacheing so it saves us reading the same file more than once""" - - def __init__(self, name): - json_data_fh = open(locate(name + ".json", 'genotype/json')) - - markers = [] - with open("%s/%s_snps.txt" % (flat_files('genotype/bimbam'), name), 'r') as bimbam_fh: - if len(bimbam_fh.readline().split(", ")) > 2: - delimiter = ", " - elif len(bimbam_fh.readline().split(",")) > 2: - delimiter = "," - elif len(bimbam_fh.readline().split("\t")) > 2: - delimiter = "\t" - else: - delimiter = " " - for line in bimbam_fh: - marker = {} - marker['name'] = line.split(delimiter)[0].rstrip() - marker['Mb'] = float(line.split(delimiter)[ - 1].rstrip()) / 1000000 - marker['chr'] = line.split(delimiter)[2].rstrip() - markers.append(marker) - - for marker in markers: - if (marker['chr'] != "X") and (marker['chr'] != "Y") and (marker['chr'] != "M"): - marker['chr'] = int(marker['chr']) - marker['Mb'] = float(marker['Mb']) - - self.markers = markers - - def add_pvalues(self, p_values): - if isinstance(p_values, list): - # THIS IS only needed for the case when we are limiting the number of p-values calculated - # if len(self.markers) > len(p_values): - # self.markers = self.markers[:len(p_values)] - - for marker, p_value in zip(self.markers, p_values): - if not p_value: - continue - marker['p_value'] = float(p_value) - if math.isnan(marker['p_value']) or marker['p_value'] <= 0: - marker['lod_score'] = 0 - marker['lrs_value'] = 0 - else: - marker['lod_score'] = -math.log10(marker['p_value']) - # Using -log(p) for the LRS; need to ask Rob how he wants to get LRS from p-values - marker['lrs_value'] = -math.log10(marker['p_value']) * 4.61 - elif isinstance(p_values, dict): - filtered_markers = [] - for marker in self.markers: - if marker['name'] in p_values: - marker['p_value'] = p_values[marker['name']] - if math.isnan(marker['p_value']) or (marker['p_value'] <= 0): - marker['lod_score'] = 0 - marker['lrs_value'] = 0 - else: - marker['lod_score'] = -math.log10(marker['p_value']) - # Using -log(p) for the LRS; need to ask Rob how he wants to get LRS from p-values - marker['lrs_value'] = - \ - math.log10(marker['p_value']) * 4.61 - filtered_markers.append(marker) - self.markers = filtered_markers - - -class HumanMarkers(Markers): - - def __init__(self, name, specified_markers=[]): - marker_data_fh = open(flat_files('mapping') + '/' + name + '.bim') - self.markers = [] - for line in marker_data_fh: - splat = line.strip().split() - if len(specified_markers) > 0: - if splat[1] in specified_markers: - marker = {} - marker['chr'] = int(splat[0]) - marker['name'] = splat[1] - marker['Mb'] = float(splat[3]) / 1000000 - else: - continue - else: - marker = {} - marker['chr'] = int(splat[0]) - marker['name'] = splat[1] - marker['Mb'] = float(splat[3]) / 1000000 - self.markers.append(marker) - - def add_pvalues(self, p_values): - super(HumanMarkers, self).add_pvalues(p_values) - - -class DatasetGroup: - """ - Each group has multiple datasets; each species has multiple groups. - - For example, Mouse has multiple groups (BXD, BXA, etc), and each group - has multiple datasets associated with it. - - """ - - def __init__(self, dataset, name=None): - """This sets self.group and self.group_id""" - with database_connection() as conn, conn.cursor() as cursor: - if not name: - cursor.execute(dataset.query_for_group, - (dataset.name,)) - else: - cursor.execute( - "SELECT InbredSet.Name, " - "InbredSet.Id, " - "InbredSet.GeneticType, " - "InbredSet.InbredSetCode " - "FROM InbredSet WHERE Name = %s", - (dataset.name,)) - results = cursor.fetchone() - if results: - (self.name, self.id, self.genetic_type, self.code) = results - else: - self.name = name or dataset.name - if self.name == 'BXD300': - self.name = "BXD" - - self.f1list = None - self.parlist = None - self.get_f1_parent_strains() - - self.mapping_id, self.mapping_names = self.get_mapping_methods() - - self.species = webqtlDatabaseFunction.retrieve_species(self.name) - - self.incparentsf1 = False - self.allsamples = None - self._datasets = None - self.genofile = None - - def get_mapping_methods(self): - mapping_id = () - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute( - "SELECT MappingMethodId FROM " - "InbredSet WHERE Name= %s", - (self.name,)) - results = cursor.fetchone() - if results and results[0]: - mapping_id = results[0] - if mapping_id == "1": - mapping_names = ["GEMMA", "QTLReaper", "R/qtl"] - elif mapping_id == "2": - mapping_names = ["GEMMA"] - elif mapping_id == "3": - mapping_names = ["R/qtl"] - elif mapping_id == "4": - mapping_names = ["GEMMA", "PLINK"] - else: - mapping_names = [] - - return mapping_id, mapping_names - - def get_markers(self): - def check_plink_gemma(): - if flat_file_exists("mapping"): - MAPPING_PATH = flat_files("mapping") + "/" - if os.path.isfile(MAPPING_PATH + self.name + ".bed"): - return True - return False - - if check_plink_gemma(): - marker_class = HumanMarkers - else: - marker_class = Markers - - if self.genofile: - self.markers = marker_class(self.genofile[:-5]) - else: - self.markers = marker_class(self.name) - - def get_f1_parent_strains(self): - try: - # NL, 07/27/2010. ParInfo has been moved from webqtlForm.py to webqtlUtil.py; - f1, f12, maternal, paternal = webqtlUtil.ParInfo[self.name] - except KeyError: - f1 = f12 = maternal = paternal = None - - if f1 and f12: - self.f1list = [f1, f12] - if maternal and paternal: - self.parlist = [maternal, paternal] - - def get_study_samplelists(self): - study_sample_file = locate_ignore_error( - self.name + ".json", 'study_sample_lists') - try: - f = open(study_sample_file) - except: - return [] - study_samples = json.load(f) - return study_samples - - def get_genofiles(self): - jsonfile = "%s/%s.json" % (webqtlConfig.GENODIR, self.name) - try: - f = open(jsonfile) - except: - return None - jsondata = json.load(f) - return jsondata['genofile'] - - def get_samplelist(self): - result = None - key = "samplelist:v3:" + self.name - if USE_REDIS: - result = r.get(key) - - if result is not None: - self.samplelist = json.loads(result) - else: - genotype_fn = locate_ignore_error(self.name + ".geno", 'genotype') - if genotype_fn: - self.samplelist = get_group_samplelists.get_samplelist( - "geno", genotype_fn) - else: - self.samplelist = None - - if USE_REDIS: - r.set(key, json.dumps(self.samplelist)) - r.expire(key, 60 * 5) - - def all_samples_ordered(self): - result = [] - lists = (self.parlist, self.f1list, self.samplelist) - [result.extend(l) for l in lists if l] - return result - - def read_genotype_file(self, use_reaper=False): - '''Read genotype from .geno file instead of database''' - # genotype_1 is Dataset Object without parents and f1 - # genotype_2 is Dataset Object with parents and f1 (not for intercross) - - # reaper barfs on unicode filenames, so here we ensure it's a string - if self.genofile: - if "RData" in self.genofile: # ZS: This is a temporary fix; I need to change the way the JSON files that point to multiple genotype files are structured to point to other file types like RData - full_filename = str( - locate(self.genofile.split(".")[0] + ".geno", 'genotype')) - else: - full_filename = str(locate(self.genofile, 'genotype')) - else: - full_filename = str(locate(self.name + '.geno', 'genotype')) - genotype_1 = gen_geno_ob.genotype(full_filename) - - if genotype_1.type == "group" and self.parlist: - genotype_2 = genotype_1.add( - Mat=self.parlist[0], Pat=self.parlist[1]) # , F1=_f1) - else: - genotype_2 = genotype_1 - - # determine default genotype object - if self.incparentsf1 and genotype_1.type != "intercross": - genotype = genotype_2 - else: - self.incparentsf1 = 0 - genotype = genotype_1 - - self.samplelist = list(genotype.prgy) - - return genotype - - -def datasets(group_name, this_group=None): - key = "group_dataset_menu:v2:" + group_name - dataset_menu = [] - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute(''' - (SELECT '#PublishFreeze',PublishFreeze.FullName,PublishFreeze.Name - FROM PublishFreeze,InbredSet - WHERE PublishFreeze.InbredSetId = InbredSet.Id - and InbredSet.Name = '%s' - ORDER BY PublishFreeze.Id ASC) - UNION - (SELECT '#GenoFreeze',GenoFreeze.FullName,GenoFreeze.Name - FROM GenoFreeze, InbredSet - WHERE GenoFreeze.InbredSetId = InbredSet.Id - and InbredSet.Name = '%s') - UNION - (SELECT Tissue.Name, ProbeSetFreeze.FullName,ProbeSetFreeze.Name - FROM ProbeSetFreeze, ProbeFreeze, InbredSet, Tissue - WHERE ProbeSetFreeze.ProbeFreezeId = ProbeFreeze.Id - and ProbeFreeze.TissueId = Tissue.Id - and ProbeFreeze.InbredSetId = InbredSet.Id - and InbredSet.Name like %s - ORDER BY Tissue.Name, ProbeSetFreeze.OrderList DESC) - ''' % (group_name, - group_name, - "'" + group_name + "'")) - the_results = cursor.fetchall() - - sorted_results = sorted(the_results, key=lambda kv: kv[0]) - - # ZS: This is kind of awkward, but need to ensure Phenotypes show up before Genotypes in dropdown - pheno_inserted = False - geno_inserted = False - for dataset_item in sorted_results: - tissue_name = dataset_item[0] - dataset = dataset_item[1] - dataset_short = dataset_item[2] - if tissue_name in ['#PublishFreeze', '#GenoFreeze']: - if tissue_name == '#PublishFreeze' and (dataset_short == group_name + 'Publish'): - dataset_menu.insert( - 0, dict(tissue=None, datasets=[(dataset, dataset_short)])) - pheno_inserted = True - elif pheno_inserted and tissue_name == '#GenoFreeze': - dataset_menu.insert( - 1, dict(tissue=None, datasets=[(dataset, dataset_short)])) - geno_inserted = True - else: - dataset_menu.append( - dict(tissue=None, datasets=[(dataset, dataset_short)])) - else: - tissue_already_exists = False - for i, tissue_dict in enumerate(dataset_menu): - if tissue_dict['tissue'] == tissue_name: - tissue_already_exists = True - break - - if tissue_already_exists: - dataset_menu[i]['datasets'].append((dataset, dataset_short)) - else: - dataset_menu.append(dict(tissue=tissue_name, - datasets=[(dataset, dataset_short)])) - - if USE_REDIS: - r.set(key, pickle.dumps(dataset_menu, pickle.HIGHEST_PROTOCOL)) - r.expire(key, 60 * 5) - - if this_group != None: - this_group._datasets = dataset_menu - return this_group._datasets - else: - return dataset_menu - - -class DataSet: - """ - DataSet class defines a dataset in webqtl, can be either Microarray, - Published phenotype, genotype, or user input dataset(temp) - - """ - - def __init__(self, name, get_samplelist=True, group_name=None): - - assert name, "Need a name" - self.name = name - self.id = None - self.shortname = None - self.fullname = None - self.type = None - self.data_scale = None # ZS: For example log2 - self.accession_id = None - - self.setup() - - if self.type == "Temp": # Need to supply group name as input if temp trait - # sets self.group and self.group_id and gets genotype - self.group = DatasetGroup(self, name=group_name) - else: - self.check_confidentiality() - self.retrieve_other_names() - # sets self.group and self.group_id and gets genotype - self.group = DatasetGroup(self) - self.accession_id = self.get_accession_id() - if get_samplelist == True: - self.group.get_samplelist() - self.species = species.TheSpecies(self) - - def as_dict(self): - return { - 'name': self.name, - 'shortname': self.shortname, - 'fullname': self.fullname, - 'type': self.type, - 'data_scale': self.data_scale, - 'group': self.group.name, - 'accession_id': self.accession_id - } - - def get_accession_id(self): - results = None - with database_connection() as conn, conn.cursor() as cursor: - if self.type == "Publish": - cursor.execute( - "SELECT InfoFiles.GN_AccesionId FROM " - "InfoFiles, PublishFreeze, InbredSet " - "WHERE InbredSet.Name = %s AND " - "PublishFreeze.InbredSetId = InbredSet.Id " - "AND InfoFiles.InfoPageName = PublishFreeze.Name " - "AND PublishFreeze.public > 0 AND " - "PublishFreeze.confidentiality < 1 " - "ORDER BY PublishFreeze.CreateTime DESC", - (self.group.name,) - ) - results = cursor.fetchone() - elif self.type == "Geno": - cursor.execute( - "SELECT InfoFiles.GN_AccesionId FROM " - "InfoFiles, GenoFreeze, InbredSet " - "WHERE InbredSet.Name = %s AND " - "GenoFreeze.InbredSetId = InbredSet.Id " - "AND InfoFiles.InfoPageName = GenoFreeze.ShortName " - "AND GenoFreeze.public > 0 AND " - "GenoFreeze.confidentiality < 1 " - "ORDER BY GenoFreeze.CreateTime DESC", - (self.group.name,) - ) - results = cursor.fetchone() - - if results: - return str(results[0]) - return "None" - - def retrieve_other_names(self): - """This method fetches the the dataset names in search_result. - - If the data set name parameter is not found in the 'Name' field of - the data set table, check if it is actually the FullName or - ShortName instead. - - This is not meant to retrieve the data set info if no name at - all is passed. - - """ - with database_connection() as conn, conn.cursor() as cursor: - try: - if self.type == "ProbeSet": - cursor.execute( - "SELECT ProbeSetFreeze.Id, ProbeSetFreeze.Name, " - "ProbeSetFreeze.FullName, ProbeSetFreeze.ShortName, " - "ProbeSetFreeze.DataScale, Tissue.Name " - "FROM ProbeSetFreeze, ProbeFreeze, Tissue " - "WHERE ProbeSetFreeze.ProbeFreezeId = ProbeFreeze.Id " - "AND ProbeFreeze.TissueId = Tissue.Id " - "AND (ProbeSetFreeze.Name = %s OR " - "ProbeSetFreeze.FullName = %s " - "OR ProbeSetFreeze.ShortName = %s)", - (self.name,)*3) - (self.id, self.name, self.fullname, self.shortname, - self.data_scale, self.tissue) = cursor.fetchone() - else: - self.tissue = "N/A" - cursor.execute( - "SELECT Id, Name, FullName, ShortName " - f"FROM {self.type}Freeze " - "WHERE (Name = %s OR FullName = " - "%s OR ShortName = %s)", - (self.name,)*3) - (self.id, self.name, self.fullname, - self.shortname) = cursor.fetchone() - except TypeError: - pass - - def chunk_dataset(self, dataset, n): - - results = {} - traits_name_dict = () - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute( - "SELECT ProbeSetXRef.DataId,ProbeSet.Name " - "FROM ProbeSet, ProbeSetXRef, ProbeSetFreeze " - "WHERE ProbeSetFreeze.Name = %s AND " - "ProbeSetXRef.ProbeSetFreezeId = ProbeSetFreeze.Id " - "AND ProbeSetXRef.ProbeSetId = ProbeSet.Id", - (self.name,)) - # should cache this - traits_name_dict = dict(cursor.fetchall()) - - for i in range(0, len(dataset), n): - matrix = list(dataset[i:i + n]) - trait_name = traits_name_dict[matrix[0][0]] - - my_values = [value for (trait_name, strain, value) in matrix] - results[trait_name] = my_values - return results - - def get_probeset_data(self, sample_list=None, trait_ids=None): - - # improvement of get trait data--->>> - if sample_list: - self.samplelist = sample_list - - else: - self.samplelist = self.group.samplelist - - if self.group.parlist != None and self.group.f1list != None: - if (self.group.parlist + self.group.f1list) in self.samplelist: - self.samplelist += self.group.parlist + self.group.f1list - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute( - "SELECT Strain.Name, Strain.Id FROM " - "Strain, Species WHERE Strain.Name IN " - f"{create_in_clause(self.samplelist)} " - "AND Strain.SpeciesId=Species.Id AND " - "Species.name = %s", (self.group.species,) - ) - results = dict(cursor.fetchall()) - sample_ids = [results[item] for item in self.samplelist] - - sorted_samplelist = [strain_name for strain_name, strain_id in sorted( - results.items(), key=lambda item: item[1])] - - cursor.execute( - "SELECT * from ProbeSetData WHERE StrainID IN " - f"{create_in_clause(sample_ids)} AND id IN " - "(SELECT ProbeSetXRef.DataId FROM " - "(ProbeSet, ProbeSetXRef, ProbeSetFreeze) " - "WHERE ProbeSetXRef.ProbeSetFreezeId = ProbeSetFreeze.Id " - "AND ProbeSetFreeze.Name = %s AND " - "ProbeSet.Id = ProbeSetXRef.ProbeSetId)", - (self.name,) - ) - - query_results = list(cursor.fetchall()) - data_results = self.chunk_dataset(query_results, len(sample_ids)) - self.samplelist = sorted_samplelist - self.trait_data = data_results - - def get_trait_data(self, sample_list=None): - if sample_list: - self.samplelist = sample_list - else: - self.samplelist = self.group.samplelist - - if self.group.parlist != None and self.group.f1list != None: - if (self.group.parlist + self.group.f1list) in self.samplelist: - self.samplelist += self.group.parlist + self.group.f1list - - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute( - "SELECT Strain.Name, Strain.Id FROM Strain, Species " - f"WHERE Strain.Name IN {create_in_clause(self.samplelist)} " - "AND Strain.SpeciesId=Species.Id " - "AND Species.name = %s", - (self.group.species,) - ) - results = dict(cursor.fetchall()) - sample_ids = [ - sample_id for sample_id in - (results.get(item) for item in self.samplelist - if item is not None) - if sample_id is not None - ] - - # MySQL limits the number of tables that can be used in a join to 61, - # so we break the sample ids into smaller chunks - # Postgres doesn't have that limit, so we can get rid of this after we transition - chunk_size = 50 - number_chunks = int(math.ceil(len(sample_ids) / chunk_size)) - - cached_results = fetch_cached_results(self.name, self.type, self.samplelist) - - if cached_results is None: - trait_sample_data = [] - for sample_ids_step in chunks.divide_into_chunks(sample_ids, number_chunks): - if self.type == "Publish": - dataset_type = "Phenotype" - else: - dataset_type = self.type - temp = ['T%s.value' % item for item in sample_ids_step] - if self.type == "Publish": - query = "SELECT {}XRef.Id".format(escape(self.type)) - else: - query = "SELECT {}.Name".format(escape(dataset_type)) - data_start_pos = 1 - if len(temp) > 0: - query = query + ", " + ', '.join(temp) - query += ' FROM ({}, {}XRef, {}Freeze) '.format(*mescape(dataset_type, - self.type, - self.type)) - - for item in sample_ids_step: - query += """ - left join {}Data as T{} on T{}.Id = {}XRef.DataId - and T{}.StrainId={}\n - """.format(*mescape(self.type, item, item, self.type, item, item)) - - if self.type == "Publish": - query += """ - WHERE {}XRef.InbredSetId = {}Freeze.InbredSetId - and {}Freeze.Name = '{}' - and {}.Id = {}XRef.{}Id - order by {}.Id - """.format(*mescape(self.type, self.type, self.type, self.name, - dataset_type, self.type, dataset_type, dataset_type)) - else: - query += """ - WHERE {}XRef.{}FreezeId = {}Freeze.Id - and {}Freeze.Name = '{}' - and {}.Id = {}XRef.{}Id - order by {}.Id - """.format(*mescape(self.type, self.type, self.type, self.type, - self.name, dataset_type, self.type, self.type, dataset_type)) - cursor.execute(query) - results = cursor.fetchall() - trait_sample_data.append([list(result) for result in results]) - - trait_count = len(trait_sample_data[0]) - self.trait_data = collections.defaultdict(list) - - data_start_pos = 1 - for trait_counter in range(trait_count): - trait_name = trait_sample_data[0][trait_counter][0] - for chunk_counter in range(int(number_chunks)): - self.trait_data[trait_name] += ( - trait_sample_data[chunk_counter][trait_counter][data_start_pos:]) - - cache_dataset_results( - self.name, self.type, self.samplelist, self.trait_data) - else: - self.trait_data = cached_results - - -class PhenotypeDataSet(DataSet): - DS_NAME_MAP['Publish'] = 'PhenotypeDataSet' - - def setup(self): - # Fields in the database table - self.search_fields = ['Phenotype.Post_publication_description', - 'Phenotype.Pre_publication_description', - 'Phenotype.Pre_publication_abbreviation', - 'Phenotype.Post_publication_abbreviation', - 'PublishXRef.mean', - 'Phenotype.Lab_code', - 'Publication.PubMed_ID', - 'Publication.Abstract', - 'Publication.Title', - 'Publication.Authors', - 'PublishXRef.Id'] - - # Figure out what display_fields is - self.display_fields = ['name', 'group_code', - 'pubmed_id', - 'pre_publication_description', - 'post_publication_description', - 'original_description', - 'pre_publication_abbreviation', - 'post_publication_abbreviation', - 'mean', - 'lab_code', - 'submitter', 'owner', - 'authorized_users', - 'authors', 'title', - 'abstract', 'journal', - 'volume', 'pages', - 'month', 'year', - 'sequence', 'units', 'comments'] - - # Fields displayed in the search results table header - self.header_fields = ['Index', - 'Record', - 'Description', - 'Authors', - 'Year', - 'Max LRS', - 'Max LRS Location', - 'Additive Effect'] - - self.type = 'Publish' - self.query_for_group = """ -SELECT InbredSet.Name, InbredSet.Id, InbredSet.GeneticType, InbredSet.InbredSetCode FROM InbredSet, PublishFreeze WHERE PublishFreeze.InbredSetId = InbredSet.Id AND PublishFreeze.Name = %s""" - - def check_confidentiality(self): - # (Urgently?) Need to write this - pass - - def get_trait_info(self, trait_list, species=''): - for this_trait in trait_list: - - if not this_trait.haveinfo: - this_trait.retrieve_info(get_qtl_info=True) - - description = this_trait.post_publication_description - - # If the dataset is confidential and the user has access to confidential - # phenotype traits, then display the pre-publication description instead - # of the post-publication description - if this_trait.confidential: - this_trait.description_display = "" - continue # for now, because no authorization features - - if not webqtlUtil.hasAccessToConfidentialPhenotypeTrait( - privilege=self.privilege, - userName=self.userName, - authorized_users=this_trait.authorized_users): - - description = this_trait.pre_publication_description - - if len(description) > 0: - this_trait.description_display = description.strip() - else: - this_trait.description_display = "" - - if not this_trait.year.isdigit(): - this_trait.pubmed_text = "N/A" - else: - this_trait.pubmed_text = this_trait.year - - if this_trait.pubmed_id: - this_trait.pubmed_link = webqtlConfig.PUBMEDLINK_URL % this_trait.pubmed_id - - # LRS and its location - this_trait.LRS_score_repr = "N/A" - this_trait.LRS_location_repr = "N/A" - - if this_trait.lrs: - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute( - "SELECT Geno.Chr, Geno.Mb FROM " - "Geno, Species WHERE " - "Species.Name = %s AND " - "Geno.Name = %s AND " - "Geno.SpeciesId = Species.Id", - (species, this_trait.locus,) - ) - if result := cursor.fetchone(): - if result[0] and result[1]: - LRS_Chr, LRS_Mb = result[0], result[1] - this_trait.LRS_score_repr = LRS_score_repr = '%3.1f' % this_trait.lrs - this_trait.LRS_location_repr = LRS_location_repr = 'Chr%s: %.6f' % ( - LRS_Chr, float(LRS_Mb)) - - def retrieve_sample_data(self, trait): - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute( - "SELECT Strain.Name, PublishData.value, " - "PublishSE.error, NStrain.count, " - "Strain.Name2 FROM (PublishData, Strain, " - "PublishXRef, PublishFreeze) LEFT JOIN " - "PublishSE ON " - "(PublishSE.DataId = PublishData.Id " - "AND PublishSE.StrainId = PublishData.StrainId) " - "LEFT JOIN NStrain ON " - "(NStrain.DataId = PublishData.Id AND " - "NStrain.StrainId = PublishData.StrainId) " - "WHERE PublishXRef.InbredSetId = PublishFreeze.InbredSetId " - "AND PublishData.Id = PublishXRef.DataId AND " - "PublishXRef.Id = %s AND PublishFreeze.Id = %s " - "AND PublishData.StrainId = Strain.Id " - "ORDER BY Strain.Name", (trait, self.id)) - return cursor.fetchall() - - -class GenotypeDataSet(DataSet): - DS_NAME_MAP['Geno'] = 'GenotypeDataSet' - - def setup(self): - # Fields in the database table - self.search_fields = ['Name', - 'Chr'] - - # Find out what display_fields is - self.display_fields = ['name', - 'chr', - 'mb', - 'source2', - 'sequence'] - - # Fields displayed in the search results table header - self.header_fields = ['Index', - 'ID', - 'Location'] - - # Todo: Obsolete or rename this field - self.type = 'Geno' - self.query_for_group = """ -SELECT InbredSet.Name, InbredSet.Id, InbredSet.GeneticType, InbredSet.InbredSetCode -FROM InbredSet, GenoFreeze WHERE GenoFreeze.InbredSetId = InbredSet.Id AND -GenoFreeze.Name = %s""" - - def check_confidentiality(self): - return geno_mrna_confidentiality(self) - - def get_trait_info(self, trait_list, species=None): - for this_trait in trait_list: - if not this_trait.haveinfo: - this_trait.retrieveInfo() - - if this_trait.chr and this_trait.mb: - this_trait.location_repr = 'Chr%s: %.6f' % ( - this_trait.chr, float(this_trait.mb)) - - def retrieve_sample_data(self, trait): - results = [] - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute( - "SELECT Strain.Name, GenoData.value, " - "GenoSE.error, 'N/A', Strain.Name2 " - "FROM (GenoData, GenoFreeze, Strain, Geno, " - "GenoXRef) LEFT JOIN GenoSE ON " - "(GenoSE.DataId = GenoData.Id AND " - "GenoSE.StrainId = GenoData.StrainId) " - "WHERE Geno.SpeciesId = %s AND " - "Geno.Name = %s AND GenoXRef.GenoId = Geno.Id " - "AND GenoXRef.GenoFreezeId = GenoFreeze.Id " - "AND GenoFreeze.Name = %s AND " - "GenoXRef.DataId = GenoData.Id " - "AND GenoData.StrainId = Strain.Id " - "ORDER BY Strain.Name", - (webqtlDatabaseFunction.retrieve_species_id(self.group.name), - trait, self.name,)) - results = list(cursor.fetchall()) - - if self.group.name in webqtlUtil.ParInfo: - f1_1, f1_2, ref, nonref = webqtlUtil.ParInfo[self.group.name] - results.append([f1_1, 0, None, "N/A", f1_1]) - results.append([f1_2, 0, None, "N/A", f1_2]) - results.append([ref, -1, None, "N/A", ref]) - results.append([nonref, 1, None, "N/A", nonref]) - - return results - - -class MrnaAssayDataSet(DataSet): - ''' - An mRNA Assay is a quantitative assessment (assay) associated with an mRNA trait - - This used to be called ProbeSet, but that term only refers specifically to the Affymetrix - platform and is far too specific. - - ''' - DS_NAME_MAP['ProbeSet'] = 'MrnaAssayDataSet' - - def setup(self): - # Fields in the database table - self.search_fields = ['Name', - 'Description', - 'Probe_Target_Description', - 'Symbol', - 'Alias', - 'GenbankId', - 'UniGeneId', - 'RefSeq_TranscriptId'] - - # Find out what display_fields is - self.display_fields = ['name', 'symbol', - 'description', 'probe_target_description', - 'chr', 'mb', - 'alias', 'geneid', - 'genbankid', 'unigeneid', - 'omim', 'refseq_transcriptid', - 'blatseq', 'targetseq', - 'chipid', 'comments', - 'strand_probe', 'strand_gene', - 'proteinid', 'uniprotid', - 'probe_set_target_region', - 'probe_set_specificity', - 'probe_set_blat_score', - 'probe_set_blat_mb_start', - 'probe_set_blat_mb_end', - 'probe_set_strand', - 'probe_set_note_by_rw', - 'flag'] - - # Fields displayed in the search results table header - self.header_fields = ['Index', - 'Record', - 'Symbol', - 'Description', - 'Location', - 'Mean', - 'Max LRS', - 'Max LRS Location', - 'Additive Effect'] - - # Todo: Obsolete or rename this field - self.type = 'ProbeSet' - self.query_for_group = """ -SELECT InbredSet.Name, InbredSet.Id, InbredSet.GeneticType, InbredSet.InbredSetCode -FROM InbredSet, ProbeSetFreeze, ProbeFreeze WHERE ProbeFreeze.InbredSetId = InbredSet.Id AND -ProbeFreeze.Id = ProbeSetFreeze.ProbeFreezeId AND ProbeSetFreeze.Name = %s""" - - def check_confidentiality(self): - return geno_mrna_confidentiality(self) - - def get_trait_info(self, trait_list=None, species=''): - - # Note: setting trait_list to [] is probably not a great idea. - if not trait_list: - trait_list = [] - with database_connection() as conn, conn.cursor() as cursor: - for this_trait in trait_list: - - if not this_trait.haveinfo: - this_trait.retrieveInfo(QTL=1) - - if not this_trait.symbol: - this_trait.symbol = "N/A" - - # XZ, 12/08/2008: description - # XZ, 06/05/2009: Rob asked to add probe target description - description_string = str( - str(this_trait.description).strip(codecs.BOM_UTF8), 'utf-8') - target_string = str( - str(this_trait.probe_target_description).strip(codecs.BOM_UTF8), 'utf-8') - - if len(description_string) > 1 and description_string != 'None': - description_display = description_string - else: - description_display = this_trait.symbol - - if (len(description_display) > 1 and description_display != 'N/A' - and len(target_string) > 1 and target_string != 'None'): - description_display = description_display + '; ' + target_string.strip() - - # Save it for the jinja2 template - this_trait.description_display = description_display - - if this_trait.chr and this_trait.mb: - this_trait.location_repr = 'Chr%s: %.6f' % ( - this_trait.chr, float(this_trait.mb)) - - # Get mean expression value - cursor.execute( - "SELECT ProbeSetXRef.mean FROM " - "ProbeSetXRef, ProbeSet WHERE " - "ProbeSetXRef.ProbeSetFreezeId = %s " - "AND ProbeSet.Id = ProbeSetXRef.ProbeSetId " - "AND ProbeSet.Name = %s", - (str(this_trait.dataset.id), this_trait.name,) - ) - result = cursor.fetchone() - - mean = result[0] if result else 0 - - if mean: - this_trait.mean = "%2.3f" % mean - - # LRS and its location - this_trait.LRS_score_repr = 'N/A' - this_trait.LRS_location_repr = 'N/A' - - # Max LRS and its Locus location - if this_trait.lrs and this_trait.locus: - cursor.execute( - "SELECT Geno.Chr, Geno.Mb FROM " - "Geno, Species WHERE " - "Species.Name = %s AND " - "Geno.Name = %s AND " - "Geno.SpeciesId = Species.Id", - (species, this_trait.locus,) - ) - if result := cursor.fetchone(): - lrs_chr, lrs_mb = result - this_trait.LRS_score_repr = '%3.1f' % this_trait.lrs - this_trait.LRS_location_repr = 'Chr%s: %.6f' % ( - lrs_chr, float(lrs_mb)) - - return trait_list - - def retrieve_sample_data(self, trait): - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute( - "SELECT Strain.Name, ProbeSetData.value, " - "ProbeSetSE.error, NStrain.count, " - "Strain.Name2 FROM (ProbeSetData, " - "ProbeSetFreeze, Strain, ProbeSet, " - "ProbeSetXRef) LEFT JOIN ProbeSetSE ON " - "(ProbeSetSE.DataId = ProbeSetData.Id AND " - "ProbeSetSE.StrainId = ProbeSetData.StrainId) " - "LEFT JOIN NStrain ON " - "(NStrain.DataId = ProbeSetData.Id AND " - "NStrain.StrainId = ProbeSetData.StrainId) " - "WHERE ProbeSet.Name = %s AND " - "ProbeSetXRef.ProbeSetId = ProbeSet.Id " - "AND ProbeSetXRef.ProbeSetFreezeId = ProbeSetFreeze.Id " - "AND ProbeSetFreeze.Name = %s AND " - "ProbeSetXRef.DataId = ProbeSetData.Id " - "AND ProbeSetData.StrainId = Strain.Id " - "ORDER BY Strain.Name", - (trait, self.name,) - ) - return cursor.fetchall() - - def retrieve_genes(self, column_name): - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute( - f"SELECT ProbeSet.Name, ProbeSet.{column_name} " - "FROM ProbeSet,ProbeSetXRef WHERE " - "ProbeSetXRef.ProbeSetFreezeId = %s " - "AND ProbeSetXRef.ProbeSetId=ProbeSet.Id", - (str(self.id),)) - return dict(cursor.fetchall()) - - -class TempDataSet(DataSet): - """Temporary user-generated data set""" - DS_NAME_MAP['Temp'] = 'TempDataSet' - - def setup(self): - self.search_fields = ['name', - 'description'] - - self.display_fields = ['name', - 'description'] - - self.header_fields = ['Name', - 'Description'] - - self.type = 'Temp' - - # Need to double check later how these are used - self.id = 1 - self.fullname = 'Temporary Storage' - self.shortname = 'Temp' - - -def geno_mrna_confidentiality(ob): - with database_connection() as conn, conn.cursor() as cursor: - cursor.execute( - "SELECT confidentiality, " - f"AuthorisedUsers FROM {ob.type}Freeze WHERE Name = %s", - (ob.name,) - ) - result = cursor.fetchall() - if len(result) > 0 and result[0]: - return True - - -def parse_db_url(): - parsed_db = urlparse(SQL_URI) - - return (parsed_db.hostname, parsed_db.username, - parsed_db.password, parsed_db.path[1:]) - - -def query_table_timestamp(dataset_type: str): - """function to query the update timestamp of a given dataset_type""" - - # computation data and actions - with database_connection() as conn, conn.cursor() as cursor: - fetch_db_name = parse_db_url() - cursor.execute( - "SELECT UPDATE_TIME FROM " - "information_schema.tables " - f"WHERE TABLE_SCHEMA = '{fetch_db_name[-1]}' " - f"AND TABLE_NAME = '{dataset_type}Data'") - date_time_obj = cursor.fetchone()[0] - return date_time_obj.strftime("%Y-%m-%d %H:%M:%S") - - -def generate_hash_file(dataset_name: str, dataset_type: str, dataset_timestamp: str, samplelist: str): - """given the trait_name generate a unique name for this""" - string_unicode = f"{dataset_name}{dataset_timestamp}{samplelist}".encode() - md5hash = hashlib.md5(string_unicode) - return md5hash.hexdigest() - - -def cache_dataset_results(dataset_name: str, dataset_type: str, samplelist: List, query_results: List): - """function to cache dataset query results to file - input dataset_name and type query_results(already processed in default dict format) - """ - # data computations actions - # store the file path on redis - - table_timestamp = query_table_timestamp(dataset_type) - samplelist_as_str = ",".join(samplelist) - - file_name = generate_hash_file(dataset_name, dataset_type, table_timestamp, samplelist_as_str) - file_path = os.path.join(TMPDIR, f"{file_name}.json") - - with open(file_path, "w") as file_handler: - json.dump(query_results, file_handler) - - -def fetch_cached_results(dataset_name: str, dataset_type: str, samplelist: List): - """function to fetch the cached results""" - - table_timestamp = query_table_timestamp(dataset_type) - samplelist_as_str = ",".join(samplelist) - - file_name = generate_hash_file(dataset_name, dataset_type, table_timestamp, samplelist_as_str) - file_path = os.path.join(TMPDIR, f"{file_name}.json") - try: - with open(file_path, "r") as file_handler: - - return json.load(file_handler) - - except Exception: - pass diff --git a/wqflask/base/data_set/__init__.py b/wqflask/base/data_set/__init__.py new file mode 100644 index 00000000..eaf80b19 --- /dev/null +++ b/wqflask/base/data_set/__init__.py @@ -0,0 +1,124 @@ +"The data_set package ..." + +# builtins imports +import json +import pickle as pickle + +# 3rd-party imports +from redis import Redis + +# local imports +from .dataset import DataSet +from base import webqtlConfig +from utility.tools import USE_REDIS +from .datasettype import DatasetType +from .tempdataset import TempDataSet +from .datasetgroup import DatasetGroup +from .utils import query_table_timestamp +from .markers import Markers, HumanMarkers +from .genotypedataset import GenotypeDataSet +from .phenotypedataset import PhenotypeDataSet +from .mrnaassaydataset import MrnaAssayDataSet +from wqflask.database import database_connection + +# Used by create_database to instantiate objects +# Each subclass will add to this + +DS_NAME_MAP = { + "Temp": "TempDataSet", + "Geno": "GenotypeDataSet", + "Publish": "PhenotypeDataSet", + "ProbeSet": "MrnaAssayDataSet" +} + +# Do the intensive work at startup one time only +# TODO: Pass in the Redis conniction from elsewhere to allow fo flexible +# configuration +Dataset_Getter = DatasetType(Redis()) + +def create_dataset(dataset_name, dataset_type=None, + get_samplelist=True, group_name=None): + if dataset_name == "Temp": + dataset_type = "Temp" + + if not dataset_type: + dataset_type = Dataset_Getter(dataset_name) + + dataset_ob = DS_NAME_MAP[dataset_type] + dataset_class = globals()[dataset_ob] + if dataset_type == "Temp": + return dataset_class(dataset_name, get_samplelist, group_name) + else: + return dataset_class(dataset_name, get_samplelist) + +def datasets(group_name, this_group=None, redis_conn=Redis()): + key = "group_dataset_menu:v2:" + group_name + dataset_menu = [] + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute(''' + (SELECT '#PublishFreeze',PublishFreeze.FullName,PublishFreeze.Name + FROM PublishFreeze,InbredSet + WHERE PublishFreeze.InbredSetId = InbredSet.Id + and InbredSet.Name = '%s' + ORDER BY PublishFreeze.Id ASC) + UNION + (SELECT '#GenoFreeze',GenoFreeze.FullName,GenoFreeze.Name + FROM GenoFreeze, InbredSet + WHERE GenoFreeze.InbredSetId = InbredSet.Id + and InbredSet.Name = '%s') + UNION + (SELECT Tissue.Name, ProbeSetFreeze.FullName,ProbeSetFreeze.Name + FROM ProbeSetFreeze, ProbeFreeze, InbredSet, Tissue + WHERE ProbeSetFreeze.ProbeFreezeId = ProbeFreeze.Id + and ProbeFreeze.TissueId = Tissue.Id + and ProbeFreeze.InbredSetId = InbredSet.Id + and InbredSet.Name like %s + ORDER BY Tissue.Name, ProbeSetFreeze.OrderList DESC) + ''' % (group_name, + group_name, + "'" + group_name + "'")) + the_results = cursor.fetchall() + + sorted_results = sorted(the_results, key=lambda kv: kv[0]) + + # ZS: This is kind of awkward, but need to ensure Phenotypes show up before Genotypes in dropdown + pheno_inserted = False + geno_inserted = False + for dataset_item in sorted_results: + tissue_name = dataset_item[0] + dataset = dataset_item[1] + dataset_short = dataset_item[2] + if tissue_name in ['#PublishFreeze', '#GenoFreeze']: + if tissue_name == '#PublishFreeze' and (dataset_short == group_name + 'Publish'): + dataset_menu.insert( + 0, dict(tissue=None, datasets=[(dataset, dataset_short)])) + pheno_inserted = True + elif pheno_inserted and tissue_name == '#GenoFreeze': + dataset_menu.insert( + 1, dict(tissue=None, datasets=[(dataset, dataset_short)])) + geno_inserted = True + else: + dataset_menu.append( + dict(tissue=None, datasets=[(dataset, dataset_short)])) + else: + tissue_already_exists = False + for i, tissue_dict in enumerate(dataset_menu): + if tissue_dict['tissue'] == tissue_name: + tissue_already_exists = True + break + + if tissue_already_exists: + dataset_menu[i]['datasets'].append((dataset, dataset_short)) + else: + dataset_menu.append(dict(tissue=tissue_name, + datasets=[(dataset, dataset_short)])) + + if USE_REDIS: + redis_conn.set(key, pickle.dumps(dataset_menu, pickle.HIGHEST_PROTOCOL)) + redis_conn.expire(key, 60 * 5) + + if this_group != None: + this_group._datasets = dataset_menu + return this_group._datasets + else: + return dataset_menu diff --git a/wqflask/base/data_set/dataset.py b/wqflask/base/data_set/dataset.py new file mode 100644 index 00000000..f035e028 --- /dev/null +++ b/wqflask/base/data_set/dataset.py @@ -0,0 +1,292 @@ +"Base Dataset class ..." + +import math +import collections + + +from redis import Redis + + +from base import species +from utility import chunks +from .datasetgroup import DatasetGroup +from wqflask.database import database_connection +from utility.db_tools import escape, mescape, create_in_clause +from .utils import fetch_cached_results, cache_dataset_results + +class DataSet: + """ + DataSet class defines a dataset in webqtl, can be either Microarray, + Published phenotype, genotype, or user input dataset(temp) + + """ + + def __init__(self, name, get_samplelist=True, group_name=None, redis_conn=Redis()): + + assert name, "Need a name" + self.name = name + self.id = None + self.shortname = None + self.fullname = None + self.type = None + self.data_scale = None # ZS: For example log2 + self.accession_id = None + + self.setup() + + if self.type == "Temp": # Need to supply group name as input if temp trait + # sets self.group and self.group_id and gets genotype + self.group = DatasetGroup(self, name=group_name) + else: + self.check_confidentiality() + self.retrieve_other_names() + # sets self.group and self.group_id and gets genotype + self.group = DatasetGroup(self) + self.accession_id = self.get_accession_id() + if get_samplelist == True: + self.group.get_samplelist(redis_conn) + self.species = species.TheSpecies(self) + + def as_dict(self): + return { + 'name': self.name, + 'shortname': self.shortname, + 'fullname': self.fullname, + 'type': self.type, + 'data_scale': self.data_scale, + 'group': self.group.name, + 'accession_id': self.accession_id + } + + def get_accession_id(self): + results = None + with database_connection() as conn, conn.cursor() as cursor: + if self.type == "Publish": + cursor.execute( + "SELECT InfoFiles.GN_AccesionId FROM " + "InfoFiles, PublishFreeze, InbredSet " + "WHERE InbredSet.Name = %s AND " + "PublishFreeze.InbredSetId = InbredSet.Id " + "AND InfoFiles.InfoPageName = PublishFreeze.Name " + "AND PublishFreeze.public > 0 AND " + "PublishFreeze.confidentiality < 1 " + "ORDER BY PublishFreeze.CreateTime DESC", + (self.group.name,) + ) + results = cursor.fetchone() + elif self.type == "Geno": + cursor.execute( + "SELECT InfoFiles.GN_AccesionId FROM " + "InfoFiles, GenoFreeze, InbredSet " + "WHERE InbredSet.Name = %s AND " + "GenoFreeze.InbredSetId = InbredSet.Id " + "AND InfoFiles.InfoPageName = GenoFreeze.ShortName " + "AND GenoFreeze.public > 0 AND " + "GenoFreeze.confidentiality < 1 " + "ORDER BY GenoFreeze.CreateTime DESC", + (self.group.name,) + ) + results = cursor.fetchone() + + if results: + return str(results[0]) + return "None" + + def retrieve_other_names(self): + """This method fetches the the dataset names in search_result. + + If the data set name parameter is not found in the 'Name' field of + the data set table, check if it is actually the FullName or + ShortName instead. + + This is not meant to retrieve the data set info if no name at + all is passed. + + """ + with database_connection() as conn, conn.cursor() as cursor: + try: + if self.type == "ProbeSet": + cursor.execute( + "SELECT ProbeSetFreeze.Id, ProbeSetFreeze.Name, " + "ProbeSetFreeze.FullName, ProbeSetFreeze.ShortName, " + "ProbeSetFreeze.DataScale, Tissue.Name " + "FROM ProbeSetFreeze, ProbeFreeze, Tissue " + "WHERE ProbeSetFreeze.ProbeFreezeId = ProbeFreeze.Id " + "AND ProbeFreeze.TissueId = Tissue.Id " + "AND (ProbeSetFreeze.Name = %s OR " + "ProbeSetFreeze.FullName = %s " + "OR ProbeSetFreeze.ShortName = %s)", + (self.name,)*3) + (self.id, self.name, self.fullname, self.shortname, + self.data_scale, self.tissue) = cursor.fetchone() + else: + self.tissue = "N/A" + cursor.execute( + "SELECT Id, Name, FullName, ShortName " + f"FROM {self.type}Freeze " + "WHERE (Name = %s OR FullName = " + "%s OR ShortName = %s)", + (self.name,)*3) + (self.id, self.name, self.fullname, + self.shortname) = cursor.fetchone() + except TypeError: + pass + + def chunk_dataset(self, dataset, n): + + results = {} + traits_name_dict = () + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute( + "SELECT ProbeSetXRef.DataId,ProbeSet.Name " + "FROM ProbeSet, ProbeSetXRef, ProbeSetFreeze " + "WHERE ProbeSetFreeze.Name = %s AND " + "ProbeSetXRef.ProbeSetFreezeId = ProbeSetFreeze.Id " + "AND ProbeSetXRef.ProbeSetId = ProbeSet.Id", + (self.name,)) + # should cache this + traits_name_dict = dict(cursor.fetchall()) + + for i in range(0, len(dataset), n): + matrix = list(dataset[i:i + n]) + trait_name = traits_name_dict[matrix[0][0]] + + my_values = [value for (trait_name, strain, value) in matrix] + results[trait_name] = my_values + return results + + def get_probeset_data(self, sample_list=None, trait_ids=None): + + # improvement of get trait data--->>> + if sample_list: + self.samplelist = sample_list + + else: + self.samplelist = self.group.samplelist + + if self.group.parlist != None and self.group.f1list != None: + if (self.group.parlist + self.group.f1list) in self.samplelist: + self.samplelist += self.group.parlist + self.group.f1list + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute( + "SELECT Strain.Name, Strain.Id FROM " + "Strain, Species WHERE Strain.Name IN " + f"{create_in_clause(self.samplelist)} " + "AND Strain.SpeciesId=Species.Id AND " + "Species.name = %s", (self.group.species,) + ) + results = dict(cursor.fetchall()) + sample_ids = [results[item] for item in self.samplelist] + + sorted_samplelist = [strain_name for strain_name, strain_id in sorted( + results.items(), key=lambda item: item[1])] + + cursor.execute( + "SELECT * from ProbeSetData WHERE StrainID IN " + f"{create_in_clause(sample_ids)} AND id IN " + "(SELECT ProbeSetXRef.DataId FROM " + "(ProbeSet, ProbeSetXRef, ProbeSetFreeze) " + "WHERE ProbeSetXRef.ProbeSetFreezeId = ProbeSetFreeze.Id " + "AND ProbeSetFreeze.Name = %s AND " + "ProbeSet.Id = ProbeSetXRef.ProbeSetId)", + (self.name,) + ) + + query_results = list(cursor.fetchall()) + data_results = self.chunk_dataset(query_results, len(sample_ids)) + self.samplelist = sorted_samplelist + self.trait_data = data_results + + def get_trait_data(self, sample_list=None): + if sample_list: + self.samplelist = sample_list + else: + self.samplelist = self.group.samplelist + + if self.group.parlist != None and self.group.f1list != None: + if (self.group.parlist + self.group.f1list) in self.samplelist: + self.samplelist += self.group.parlist + self.group.f1list + + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute( + "SELECT Strain.Name, Strain.Id FROM Strain, Species " + f"WHERE Strain.Name IN {create_in_clause(self.samplelist)} " + "AND Strain.SpeciesId=Species.Id " + "AND Species.name = %s", + (self.group.species,) + ) + results = dict(cursor.fetchall()) + sample_ids = [ + sample_id for sample_id in + (results.get(item) for item in self.samplelist + if item is not None) + if sample_id is not None + ] + + # MySQL limits the number of tables that can be used in a join to 61, + # so we break the sample ids into smaller chunks + # Postgres doesn't have that limit, so we can get rid of this after we transition + chunk_size = 50 + number_chunks = int(math.ceil(len(sample_ids) / chunk_size)) + + cached_results = fetch_cached_results(self.name, self.type, self.samplelist) + + if cached_results is None: + trait_sample_data = [] + for sample_ids_step in chunks.divide_into_chunks(sample_ids, number_chunks): + if self.type == "Publish": + dataset_type = "Phenotype" + else: + dataset_type = self.type + temp = ['T%s.value' % item for item in sample_ids_step] + if self.type == "Publish": + query = "SELECT {}XRef.Id".format(escape(self.type)) + else: + query = "SELECT {}.Name".format(escape(dataset_type)) + data_start_pos = 1 + if len(temp) > 0: + query = query + ", " + ', '.join(temp) + query += ' FROM ({}, {}XRef, {}Freeze) '.format(*mescape(dataset_type, + self.type, + self.type)) + + for item in sample_ids_step: + query += """ + left join {}Data as T{} on T{}.Id = {}XRef.DataId + and T{}.StrainId={}\n + """.format(*mescape(self.type, item, item, self.type, item, item)) + + if self.type == "Publish": + query += """ + WHERE {}XRef.InbredSetId = {}Freeze.InbredSetId + and {}Freeze.Name = '{}' + and {}.Id = {}XRef.{}Id + order by {}.Id + """.format(*mescape(self.type, self.type, self.type, self.name, + dataset_type, self.type, dataset_type, dataset_type)) + else: + query += """ + WHERE {}XRef.{}FreezeId = {}Freeze.Id + and {}Freeze.Name = '{}' + and {}.Id = {}XRef.{}Id + order by {}.Id + """.format(*mescape(self.type, self.type, self.type, self.type, + self.name, dataset_type, self.type, self.type, dataset_type)) + cursor.execute(query) + results = cursor.fetchall() + trait_sample_data.append([list(result) for result in results]) + + trait_count = len(trait_sample_data[0]) + self.trait_data = collections.defaultdict(list) + + data_start_pos = 1 + for trait_counter in range(trait_count): + trait_name = trait_sample_data[0][trait_counter][0] + for chunk_counter in range(int(number_chunks)): + self.trait_data[trait_name] += ( + trait_sample_data[chunk_counter][trait_counter][data_start_pos:]) + + cache_dataset_results( + self.name, self.type, self.samplelist, self.trait_data) + else: + self.trait_data = cached_results diff --git a/wqflask/base/data_set/datasetgroup.py b/wqflask/base/data_set/datasetgroup.py new file mode 100644 index 00000000..a13d2f93 --- /dev/null +++ b/wqflask/base/data_set/datasetgroup.py @@ -0,0 +1,193 @@ +"Dataset Group class ..." + +import os +import json + + +from base import webqtlConfig +from utility import webqtlUtil +from utility import gen_geno_ob +from db import webqtlDatabaseFunction +from maintenance import get_group_samplelists +from wqflask.database import database_connection +from utility.tools import ( + locate, + USE_REDIS, + flat_files, + flat_file_exists, + locate_ignore_error) + +class DatasetGroup: + """ + Each group has multiple datasets; each species has multiple groups. + + For example, Mouse has multiple groups (BXD, BXA, etc), and each group + has multiple datasets associated with it. + + """ + + def __init__(self, dataset, name=None): + """This sets self.group and self.group_id""" + with database_connection() as conn, conn.cursor() as cursor: + if not name: + cursor.execute(dataset.query_for_group, + (dataset.name,)) + else: + cursor.execute( + "SELECT InbredSet.Name, " + "InbredSet.Id, " + "InbredSet.GeneticType, " + "InbredSet.InbredSetCode " + "FROM InbredSet WHERE Name = %s", + (dataset.name,)) + results = cursor.fetchone() + if results: + (self.name, self.id, self.genetic_type, self.code) = results + else: + self.name = name or dataset.name + if self.name == 'BXD300': + self.name = "BXD" + + self.f1list = None + self.parlist = None + self.get_f1_parent_strains() + + self.mapping_id, self.mapping_names = self.get_mapping_methods() + + self.species = webqtlDatabaseFunction.retrieve_species(self.name) + + self.incparentsf1 = False + self.allsamples = None + self._datasets = None + self.genofile = None + + def get_mapping_methods(self): + mapping_id = () + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute( + "SELECT MappingMethodId FROM " + "InbredSet WHERE Name= %s", + (self.name,)) + results = cursor.fetchone() + if results and results[0]: + mapping_id = results[0] + if mapping_id == "1": + mapping_names = ["GEMMA", "QTLReaper", "R/qtl"] + elif mapping_id == "2": + mapping_names = ["GEMMA"] + elif mapping_id == "3": + mapping_names = ["R/qtl"] + elif mapping_id == "4": + mapping_names = ["GEMMA", "PLINK"] + else: + mapping_names = [] + + return mapping_id, mapping_names + + def get_markers(self): + def check_plink_gemma(): + if flat_file_exists("mapping"): + MAPPING_PATH = flat_files("mapping") + "/" + if os.path.isfile(MAPPING_PATH + self.name + ".bed"): + return True + return False + + if check_plink_gemma(): + marker_class = HumanMarkers + else: + marker_class = Markers + + if self.genofile: + self.markers = marker_class(self.genofile[:-5]) + else: + self.markers = marker_class(self.name) + + def get_f1_parent_strains(self): + try: + # NL, 07/27/2010. ParInfo has been moved from webqtlForm.py to webqtlUtil.py; + f1, f12, maternal, paternal = webqtlUtil.ParInfo[self.name] + except KeyError: + f1 = f12 = maternal = paternal = None + + if f1 and f12: + self.f1list = [f1, f12] + if maternal and paternal: + self.parlist = [maternal, paternal] + + def get_study_samplelists(self): + study_sample_file = locate_ignore_error( + self.name + ".json", 'study_sample_lists') + try: + f = open(study_sample_file) + except: + return [] + study_samples = json.load(f) + return study_samples + + def get_genofiles(self): + jsonfile = "%s/%s.json" % (webqtlConfig.GENODIR, self.name) + try: + f = open(jsonfile) + except: + return None + jsondata = json.load(f) + return jsondata['genofile'] + + def get_samplelist(self, redis_conn): + result = None + key = "samplelist:v3:" + self.name + if USE_REDIS: + result = redis_conn.get(key) + + if result is not None: + self.samplelist = json.loads(result) + else: + genotype_fn = locate_ignore_error(self.name + ".geno", 'genotype') + if genotype_fn: + self.samplelist = get_group_samplelists.get_samplelist( + "geno", genotype_fn) + else: + self.samplelist = None + + if USE_REDIS: + redis_conn.set(key, json.dumps(self.samplelist)) + redis_conn.expire(key, 60 * 5) + + def all_samples_ordered(self): + result = [] + lists = (self.parlist, self.f1list, self.samplelist) + [result.extend(l) for l in lists if l] + return result + + def read_genotype_file(self, use_reaper=False): + '''Read genotype from .geno file instead of database''' + # genotype_1 is Dataset Object without parents and f1 + # genotype_2 is Dataset Object with parents and f1 (not for intercross) + + # reaper barfs on unicode filenames, so here we ensure it's a string + if self.genofile: + if "RData" in self.genofile: # ZS: This is a temporary fix; I need to change the way the JSON files that point to multiple genotype files are structured to point to other file types like RData + full_filename = str( + locate(self.genofile.split(".")[0] + ".geno", 'genotype')) + else: + full_filename = str(locate(self.genofile, 'genotype')) + else: + full_filename = str(locate(self.name + '.geno', 'genotype')) + genotype_1 = gen_geno_ob.genotype(full_filename) + + if genotype_1.type == "group" and self.parlist: + genotype_2 = genotype_1.add( + Mat=self.parlist[0], Pat=self.parlist[1]) # , F1=_f1) + else: + genotype_2 = genotype_1 + + # determine default genotype object + if self.incparentsf1 and genotype_1.type != "intercross": + genotype = genotype_2 + else: + self.incparentsf1 = 0 + genotype = genotype_1 + + self.samplelist = list(genotype.prgy) + + return genotype diff --git a/wqflask/base/data_set/datasettype.py b/wqflask/base/data_set/datasettype.py new file mode 100644 index 00000000..ca6515b6 --- /dev/null +++ b/wqflask/base/data_set/datasettype.py @@ -0,0 +1,122 @@ +# builtins imports + +import json +import requests +from dataclasses import field +from dataclasses import InitVar +from typing import Optional, Dict +from dataclasses import dataclass + + +from redis import Redis + + +from utility.tools import GN2_BASE_URL +from wqflask.database import database_connection + +@dataclass +class DatasetType: + """Create a dictionary of samples where the value is set to Geno, + Publish or ProbeSet. E.g. + + {'AD-cases-controls-MyersGeno': 'Geno', + 'AD-cases-controls-MyersPublish': 'Publish', + 'AKXDGeno': 'Geno', + 'AXBXAGeno': 'Geno', + 'AXBXAPublish': 'Publish', + 'Aging-Brain-UCIPublish': 'Publish', + 'All Phenotypes': 'Publish', + 'B139_K_1206_M': 'ProbeSet', + 'B139_K_1206_R': 'ProbeSet' ... + } + """ + redis_instance: InitVar[Redis] + datasets: Optional[Dict] = field(init=False, default_factory=dict) + data: Optional[Dict] = field(init=False) + + def __post_init__(self, redis_instance): + self.redis_instance = redis_instance + data = redis_instance.get("dataset_structure") + if data: + self.datasets = json.loads(data) + else: + # ZS: I don't think this should ever run unless Redis is + # emptied + try: + data = json.loads(requests.get( + GN2_BASE_URL + "/api/v_pre1/gen_dropdown", + timeout=5).content) + for _species in data['datasets']: + for group in data['datasets'][_species]: + for dataset_type in data['datasets'][_species][group]: + for dataset in data['datasets'][_species][group][dataset_type]: + short_dataset_name = dataset[1] + if dataset_type == "Phenotypes": + new_type = "Publish" + elif dataset_type == "Genotypes": + new_type = "Geno" + else: + new_type = "ProbeSet" + self.datasets[short_dataset_name] = new_type + except Exception: # Do nothing + pass + + self.redis_instance.set("dataset_structure", + json.dumps(self.datasets)) + self.data = data + + def set_dataset_key(self, t, name): + """If name is not in the object's dataset dictionary, set it, and + update dataset_structure in Redis + args: + t: Type of dataset structure which can be: 'mrna_expr', 'pheno', + 'other_pheno', 'geno' + name: The name of the key to inserted in the datasets dictionary + + """ + sql_query_mapping = { + 'mrna_expr': ("SELECT ProbeSetFreeze.Id FROM " + "ProbeSetFreeze WHERE " + "ProbeSetFreeze.Name = %s "), + 'pheno': ("SELECT InfoFiles.GN_AccesionId " + "FROM InfoFiles, PublishFreeze, InbredSet " + "WHERE InbredSet.Name = %s AND " + "PublishFreeze.InbredSetId = InbredSet.Id AND " + "InfoFiles.InfoPageName = PublishFreeze.Name"), + 'other_pheno': ("SELECT PublishFreeze.Name " + "FROM PublishFreeze, InbredSet " + "WHERE InbredSet.Name = %s AND " + "PublishFreeze.InbredSetId = InbredSet.Id"), + 'geno': ("SELECT GenoFreeze.Id FROM GenoFreeze WHERE " + "GenoFreeze.Name = %s ") + } + + dataset_name_mapping = { + "mrna_expr": "ProbeSet", + "pheno": "Publish", + "other_pheno": "Publish", + "geno": "Geno", + } + + group_name = name + if t in ['pheno', 'other_pheno']: + group_name = name.replace("Publish", "") + + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute(sql_query_mapping[t], (group_name,)) + if cursor.fetchone(): + self.datasets[name] = dataset_name_mapping[t] + self.redis_instance.set( + "dataset_structure", json.dumps(self.datasets)) + return True + + + def __call__(self, name): + if name not in self.datasets: + for t in ["mrna_expr", "pheno", "other_pheno", "geno"]: + # This has side-effects, with the end result being a + # truth-y value + if(self.set_dataset_key(t, name)): + break + # Return None if name has not been set + return self.datasets.get(name, None) diff --git a/wqflask/base/data_set/genotypedataset.py b/wqflask/base/data_set/genotypedataset.py new file mode 100644 index 00000000..2381a42a --- /dev/null +++ b/wqflask/base/data_set/genotypedataset.py @@ -0,0 +1,75 @@ +"GenotypeDataSet class ..." + +from .dataset import DataSet +from utility import webqtlUtil +from db import webqtlDatabaseFunction +from .utils import geno_mrna_confidentiality +from wqflask.database import database_connection + +class GenotypeDataSet(DataSet): + + def setup(self): + # Fields in the database table + self.search_fields = ['Name', + 'Chr'] + + # Find out what display_fields is + self.display_fields = ['name', + 'chr', + 'mb', + 'source2', + 'sequence'] + + # Fields displayed in the search results table header + self.header_fields = ['Index', + 'ID', + 'Location'] + + # Todo: Obsolete or rename this field + self.type = 'Geno' + self.query_for_group = """ +SELECT InbredSet.Name, InbredSet.Id, InbredSet.GeneticType, InbredSet.InbredSetCode +FROM InbredSet, GenoFreeze WHERE GenoFreeze.InbredSetId = InbredSet.Id AND +GenoFreeze.Name = %s""" + + def check_confidentiality(self): + return geno_mrna_confidentiality(self) + + def get_trait_info(self, trait_list, species=None): + for this_trait in trait_list: + if not this_trait.haveinfo: + this_trait.retrieveInfo() + + if this_trait.chr and this_trait.mb: + this_trait.location_repr = 'Chr%s: %.6f' % ( + this_trait.chr, float(this_trait.mb)) + + def retrieve_sample_data(self, trait): + results = [] + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute( + "SELECT Strain.Name, GenoData.value, " + "GenoSE.error, 'N/A', Strain.Name2 " + "FROM (GenoData, GenoFreeze, Strain, Geno, " + "GenoXRef) LEFT JOIN GenoSE ON " + "(GenoSE.DataId = GenoData.Id AND " + "GenoSE.StrainId = GenoData.StrainId) " + "WHERE Geno.SpeciesId = %s AND " + "Geno.Name = %s AND GenoXRef.GenoId = Geno.Id " + "AND GenoXRef.GenoFreezeId = GenoFreeze.Id " + "AND GenoFreeze.Name = %s AND " + "GenoXRef.DataId = GenoData.Id " + "AND GenoData.StrainId = Strain.Id " + "ORDER BY Strain.Name", + (webqtlDatabaseFunction.retrieve_species_id(self.group.name), + trait, self.name,)) + results = list(cursor.fetchall()) + + if self.group.name in webqtlUtil.ParInfo: + f1_1, f1_2, ref, nonref = webqtlUtil.ParInfo[self.group.name] + results.append([f1_1, 0, None, "N/A", f1_1]) + results.append([f1_2, 0, None, "N/A", f1_2]) + results.append([ref, -1, None, "N/A", ref]) + results.append([nonref, 1, None, "N/A", nonref]) + + return results diff --git a/wqflask/base/data_set/markers.py b/wqflask/base/data_set/markers.py new file mode 100644 index 00000000..6f56445e --- /dev/null +++ b/wqflask/base/data_set/markers.py @@ -0,0 +1,96 @@ +"Base Class: Markers - " + +import math + +from utility.tools import locate, flat_files + +class Markers: + """Todo: Build in cacheing so it saves us reading the same file more than once""" + + def __init__(self, name): + json_data_fh = open(locate(name + ".json", 'genotype/json')) + + markers = [] + with open("%s/%s_snps.txt" % (flat_files('genotype/bimbam'), name), 'r') as bimbam_fh: + if len(bimbam_fh.readline().split(", ")) > 2: + delimiter = ", " + elif len(bimbam_fh.readline().split(",")) > 2: + delimiter = "," + elif len(bimbam_fh.readline().split("\t")) > 2: + delimiter = "\t" + else: + delimiter = " " + for line in bimbam_fh: + marker = {} + marker['name'] = line.split(delimiter)[0].rstrip() + marker['Mb'] = float(line.split(delimiter)[ + 1].rstrip()) / 1000000 + marker['chr'] = line.split(delimiter)[2].rstrip() + markers.append(marker) + + for marker in markers: + if (marker['chr'] != "X") and (marker['chr'] != "Y") and (marker['chr'] != "M"): + marker['chr'] = int(marker['chr']) + marker['Mb'] = float(marker['Mb']) + + self.markers = markers + + def add_pvalues(self, p_values): + if isinstance(p_values, list): + # THIS IS only needed for the case when we are limiting the number of p-values calculated + # if len(self.markers) > len(p_values): + # self.markers = self.markers[:len(p_values)] + + for marker, p_value in zip(self.markers, p_values): + if not p_value: + continue + marker['p_value'] = float(p_value) + if math.isnan(marker['p_value']) or marker['p_value'] <= 0: + marker['lod_score'] = 0 + marker['lrs_value'] = 0 + else: + marker['lod_score'] = -math.log10(marker['p_value']) + # Using -log(p) for the LRS; need to ask Rob how he wants to get LRS from p-values + marker['lrs_value'] = -math.log10(marker['p_value']) * 4.61 + elif isinstance(p_values, dict): + filtered_markers = [] + for marker in self.markers: + if marker['name'] in p_values: + marker['p_value'] = p_values[marker['name']] + if math.isnan(marker['p_value']) or (marker['p_value'] <= 0): + marker['lod_score'] = 0 + marker['lrs_value'] = 0 + else: + marker['lod_score'] = -math.log10(marker['p_value']) + # Using -log(p) for the LRS; need to ask Rob how he wants to get LRS from p-values + marker['lrs_value'] = - \ + math.log10(marker['p_value']) * 4.61 + filtered_markers.append(marker) + self.markers = filtered_markers + + +class HumanMarkers(Markers): + "Markers for humans ..." + + def __init__(self, name, specified_markers=[]): + marker_data_fh = open(flat_files('mapping') + '/' + name + '.bim') + self.markers = [] + for line in marker_data_fh: + splat = line.strip().split() + if len(specified_markers) > 0: + if splat[1] in specified_markers: + marker = {} + marker['chr'] = int(splat[0]) + marker['name'] = splat[1] + marker['Mb'] = float(splat[3]) / 1000000 + else: + continue + else: + marker = {} + marker['chr'] = int(splat[0]) + marker['name'] = splat[1] + marker['Mb'] = float(splat[3]) / 1000000 + self.markers.append(marker) + + def add_pvalues(self, p_values): + super(HumanMarkers, self).add_pvalues(p_values) diff --git a/wqflask/base/data_set/mrnaassaydataset.py b/wqflask/base/data_set/mrnaassaydataset.py new file mode 100644 index 00000000..b76c5a74 --- /dev/null +++ b/wqflask/base/data_set/mrnaassaydataset.py @@ -0,0 +1,178 @@ +"MrnaAssayDataSet class ..." + +import codecs + + +from .dataset import DataSet +from .utils import geno_mrna_confidentiality +from wqflask.database import database_connection + +class MrnaAssayDataSet(DataSet): + ''' + An mRNA Assay is a quantitative assessment (assay) associated with an mRNA trait + + This used to be called ProbeSet, but that term only refers specifically to the Affymetrix + platform and is far too specific. + + ''' + + def setup(self): + # Fields in the database table + self.search_fields = ['Name', + 'Description', + 'Probe_Target_Description', + 'Symbol', + 'Alias', + 'GenbankId', + 'UniGeneId', + 'RefSeq_TranscriptId'] + + # Find out what display_fields is + self.display_fields = ['name', 'symbol', + 'description', 'probe_target_description', + 'chr', 'mb', + 'alias', 'geneid', + 'genbankid', 'unigeneid', + 'omim', 'refseq_transcriptid', + 'blatseq', 'targetseq', + 'chipid', 'comments', + 'strand_probe', 'strand_gene', + 'proteinid', 'uniprotid', + 'probe_set_target_region', + 'probe_set_specificity', + 'probe_set_blat_score', + 'probe_set_blat_mb_start', + 'probe_set_blat_mb_end', + 'probe_set_strand', + 'probe_set_note_by_rw', + 'flag'] + + # Fields displayed in the search results table header + self.header_fields = ['Index', + 'Record', + 'Symbol', + 'Description', + 'Location', + 'Mean', + 'Max LRS', + 'Max LRS Location', + 'Additive Effect'] + + # Todo: Obsolete or rename this field + self.type = 'ProbeSet' + self.query_for_group = """ +SELECT InbredSet.Name, InbredSet.Id, InbredSet.GeneticType, InbredSet.InbredSetCode +FROM InbredSet, ProbeSetFreeze, ProbeFreeze WHERE ProbeFreeze.InbredSetId = InbredSet.Id AND +ProbeFreeze.Id = ProbeSetFreeze.ProbeFreezeId AND ProbeSetFreeze.Name = %s""" + + def check_confidentiality(self): + return geno_mrna_confidentiality(self) + + def get_trait_info(self, trait_list=None, species=''): + + # Note: setting trait_list to [] is probably not a great idea. + if not trait_list: + trait_list = [] + with database_connection() as conn, conn.cursor() as cursor: + for this_trait in trait_list: + + if not this_trait.haveinfo: + this_trait.retrieveInfo(QTL=1) + + if not this_trait.symbol: + this_trait.symbol = "N/A" + + # XZ, 12/08/2008: description + # XZ, 06/05/2009: Rob asked to add probe target description + description_string = str( + str(this_trait.description).strip(codecs.BOM_UTF8), 'utf-8') + target_string = str( + str(this_trait.probe_target_description).strip(codecs.BOM_UTF8), 'utf-8') + + if len(description_string) > 1 and description_string != 'None': + description_display = description_string + else: + description_display = this_trait.symbol + + if (len(description_display) > 1 and description_display != 'N/A' + and len(target_string) > 1 and target_string != 'None'): + description_display = description_display + '; ' + target_string.strip() + + # Save it for the jinja2 template + this_trait.description_display = description_display + + if this_trait.chr and this_trait.mb: + this_trait.location_repr = 'Chr%s: %.6f' % ( + this_trait.chr, float(this_trait.mb)) + + # Get mean expression value + cursor.execute( + "SELECT ProbeSetXRef.mean FROM " + "ProbeSetXRef, ProbeSet WHERE " + "ProbeSetXRef.ProbeSetFreezeId = %s " + "AND ProbeSet.Id = ProbeSetXRef.ProbeSetId " + "AND ProbeSet.Name = %s", + (str(this_trait.dataset.id), this_trait.name,) + ) + result = cursor.fetchone() + + mean = result[0] if result else 0 + + if mean: + this_trait.mean = "%2.3f" % mean + + # LRS and its location + this_trait.LRS_score_repr = 'N/A' + this_trait.LRS_location_repr = 'N/A' + + # Max LRS and its Locus location + if this_trait.lrs and this_trait.locus: + cursor.execute( + "SELECT Geno.Chr, Geno.Mb FROM " + "Geno, Species WHERE " + "Species.Name = %s AND " + "Geno.Name = %s AND " + "Geno.SpeciesId = Species.Id", + (species, this_trait.locus,) + ) + if result := cursor.fetchone(): + lrs_chr, lrs_mb = result + this_trait.LRS_score_repr = '%3.1f' % this_trait.lrs + this_trait.LRS_location_repr = 'Chr%s: %.6f' % ( + lrs_chr, float(lrs_mb)) + + return trait_list + + def retrieve_sample_data(self, trait): + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute( + "SELECT Strain.Name, ProbeSetData.value, " + "ProbeSetSE.error, NStrain.count, " + "Strain.Name2 FROM (ProbeSetData, " + "ProbeSetFreeze, Strain, ProbeSet, " + "ProbeSetXRef) LEFT JOIN ProbeSetSE ON " + "(ProbeSetSE.DataId = ProbeSetData.Id AND " + "ProbeSetSE.StrainId = ProbeSetData.StrainId) " + "LEFT JOIN NStrain ON " + "(NStrain.DataId = ProbeSetData.Id AND " + "NStrain.StrainId = ProbeSetData.StrainId) " + "WHERE ProbeSet.Name = %s AND " + "ProbeSetXRef.ProbeSetId = ProbeSet.Id " + "AND ProbeSetXRef.ProbeSetFreezeId = ProbeSetFreeze.Id " + "AND ProbeSetFreeze.Name = %s AND " + "ProbeSetXRef.DataId = ProbeSetData.Id " + "AND ProbeSetData.StrainId = Strain.Id " + "ORDER BY Strain.Name", + (trait, self.name,) + ) + return cursor.fetchall() + + def retrieve_genes(self, column_name): + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute( + f"SELECT ProbeSet.Name, ProbeSet.{column_name} " + "FROM ProbeSet,ProbeSetXRef WHERE " + "ProbeSetXRef.ProbeSetFreezeId = %s " + "AND ProbeSetXRef.ProbeSetId=ProbeSet.Id", + (str(self.id),)) + return dict(cursor.fetchall()) diff --git a/wqflask/base/data_set/phenotypedataset.py b/wqflask/base/data_set/phenotypedataset.py new file mode 100644 index 00000000..bdf28a7a --- /dev/null +++ b/wqflask/base/data_set/phenotypedataset.py @@ -0,0 +1,133 @@ +"PhenotypeDataSet class ..." + +from .dataset import DataSet +from base import webqtlConfig +from wqflask.database import database_connection + +class PhenotypeDataSet(DataSet): + + def setup(self): + # Fields in the database table + self.search_fields = ['Phenotype.Post_publication_description', + 'Phenotype.Pre_publication_description', + 'Phenotype.Pre_publication_abbreviation', + 'Phenotype.Post_publication_abbreviation', + 'PublishXRef.mean', + 'Phenotype.Lab_code', + 'Publication.PubMed_ID', + 'Publication.Abstract', + 'Publication.Title', + 'Publication.Authors', + 'PublishXRef.Id'] + + # Figure out what display_fields is + self.display_fields = ['name', 'group_code', + 'pubmed_id', + 'pre_publication_description', + 'post_publication_description', + 'original_description', + 'pre_publication_abbreviation', + 'post_publication_abbreviation', + 'mean', + 'lab_code', + 'submitter', 'owner', + 'authorized_users', + 'authors', 'title', + 'abstract', 'journal', + 'volume', 'pages', + 'month', 'year', + 'sequence', 'units', 'comments'] + + # Fields displayed in the search results table header + self.header_fields = ['Index', + 'Record', + 'Description', + 'Authors', + 'Year', + 'Max LRS', + 'Max LRS Location', + 'Additive Effect'] + + self.type = 'Publish' + self.query_for_group = """ +SELECT InbredSet.Name, InbredSet.Id, InbredSet.GeneticType, InbredSet.InbredSetCode FROM InbredSet, PublishFreeze WHERE PublishFreeze.InbredSetId = InbredSet.Id AND PublishFreeze.Name = %s""" + + def check_confidentiality(self): + # (Urgently?) Need to write this + pass + + def get_trait_info(self, trait_list, species=''): + for this_trait in trait_list: + + if not this_trait.haveinfo: + this_trait.retrieve_info(get_qtl_info=True) + + description = this_trait.post_publication_description + + # If the dataset is confidential and the user has access to confidential + # phenotype traits, then display the pre-publication description instead + # of the post-publication description + if this_trait.confidential: + this_trait.description_display = "" + continue # for now, because no authorization features + + if not webqtlUtil.hasAccessToConfidentialPhenotypeTrait( + privilege=self.privilege, + userName=self.userName, + authorized_users=this_trait.authorized_users): + + description = this_trait.pre_publication_description + + if len(description) > 0: + this_trait.description_display = description.strip() + else: + this_trait.description_display = "" + + if not this_trait.year.isdigit(): + this_trait.pubmed_text = "N/A" + else: + this_trait.pubmed_text = this_trait.year + + if this_trait.pubmed_id: + this_trait.pubmed_link = webqtlConfig.PUBMEDLINK_URL % this_trait.pubmed_id + + # LRS and its location + this_trait.LRS_score_repr = "N/A" + this_trait.LRS_location_repr = "N/A" + + if this_trait.lrs: + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute( + "SELECT Geno.Chr, Geno.Mb FROM " + "Geno, Species WHERE " + "Species.Name = %s AND " + "Geno.Name = %s AND " + "Geno.SpeciesId = Species.Id", + (species, this_trait.locus,) + ) + if result := cursor.fetchone(): + if result[0] and result[1]: + LRS_Chr, LRS_Mb = result[0], result[1] + this_trait.LRS_score_repr = LRS_score_repr = '%3.1f' % this_trait.lrs + this_trait.LRS_location_repr = LRS_location_repr = 'Chr%s: %.6f' % ( + LRS_Chr, float(LRS_Mb)) + + def retrieve_sample_data(self, trait): + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute( + "SELECT Strain.Name, PublishData.value, " + "PublishSE.error, NStrain.count, " + "Strain.Name2 FROM (PublishData, Strain, " + "PublishXRef, PublishFreeze) LEFT JOIN " + "PublishSE ON " + "(PublishSE.DataId = PublishData.Id " + "AND PublishSE.StrainId = PublishData.StrainId) " + "LEFT JOIN NStrain ON " + "(NStrain.DataId = PublishData.Id AND " + "NStrain.StrainId = PublishData.StrainId) " + "WHERE PublishXRef.InbredSetId = PublishFreeze.InbredSetId " + "AND PublishData.Id = PublishXRef.DataId AND " + "PublishXRef.Id = %s AND PublishFreeze.Id = %s " + "AND PublishData.StrainId = Strain.Id " + "ORDER BY Strain.Name", (trait, self.id)) + return cursor.fetchall() diff --git a/wqflask/base/data_set/probably_unused.py b/wqflask/base/data_set/probably_unused.py new file mode 100644 index 00000000..4e54bcff --- /dev/null +++ b/wqflask/base/data_set/probably_unused.py @@ -0,0 +1,34 @@ +"Functions that are probably unused in the code" + +import pickle as pickle + +from wqflask.database import database_connection + +def create_datasets_list(): + if USE_REDIS: + key = "all_datasets" + result = redis_conn.get(key) + + if result: + datasets = pickle.loads(result) + + if result is None: + datasets = list() + type_dict = {'Publish': 'PublishFreeze', + 'ProbeSet': 'ProbeSetFreeze', + 'Geno': 'GenoFreeze'} + + for dataset_type in type_dict: + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute("SELECT Name FROM %s", + (type_dict[dataset_type],)) + results = cursor.fetchall(query) + if results: + for result in results: + datasets.append( + create_dataset(result.Name, dataset_type)) + if USE_REDIS: + redis_conn.set(key, pickle.dumps(datasets, pickle.HIGHEST_PROTOCOL)) + redis_conn.expire(key, 60 * 60) + + return datasets diff --git a/wqflask/base/data_set/tempdataset.py b/wqflask/base/data_set/tempdataset.py new file mode 100644 index 00000000..b1c26a3b --- /dev/null +++ b/wqflask/base/data_set/tempdataset.py @@ -0,0 +1,23 @@ +"TempDataSet class ..." + +from .dataset import DataSet + +class TempDataSet(DataSet): + """Temporary user-generated data set""" + + def setup(self): + self.search_fields = ['name', + 'description'] + + self.display_fields = ['name', + 'description'] + + self.header_fields = ['Name', + 'Description'] + + self.type = 'Temp' + + # Need to double check later how these are used + self.id = 1 + self.fullname = 'Temporary Storage' + self.shortname = 'Temp' diff --git a/wqflask/base/data_set/utils.py b/wqflask/base/data_set/utils.py new file mode 100644 index 00000000..0077c292 --- /dev/null +++ b/wqflask/base/data_set/utils.py @@ -0,0 +1,77 @@ +"data_set package utilities" + +import os +import json +import hashlib +from typing import List + + +from utility.tools import SQL_URI +from base.webqtlConfig import TMPDIR +from wqflask.database import parse_db_url, database_connection + +def geno_mrna_confidentiality(ob): + with database_connection() as conn, conn.cursor() as cursor: + cursor.execute( + "SELECT confidentiality, " + f"AuthorisedUsers FROM {ob.type}Freeze WHERE Name = %s", + (ob.name,) + ) + result = cursor.fetchall() + if len(result) > 0 and result[0]: + return True + +def query_table_timestamp(dataset_type: str): + """function to query the update timestamp of a given dataset_type""" + + # computation data and actions + with database_connection() as conn, conn.cursor() as cursor: + fetch_db_name = parse_db_url(SQL_URI) + cursor.execute( + "SELECT UPDATE_TIME FROM " + "information_schema.tables " + f"WHERE TABLE_SCHEMA = '{fetch_db_name[3]}' " + f"AND TABLE_NAME = '{dataset_type}Data'") + date_time_obj = cursor.fetchone()[0] + return date_time_obj.strftime("%Y-%m-%d %H:%M:%S") + + +def generate_hash_file(dataset_name: str, dataset_type: str, dataset_timestamp: str, samplelist: str): + """given the trait_name generate a unique name for this""" + string_unicode = f"{dataset_name}{dataset_timestamp}{samplelist}".encode() + md5hash = hashlib.md5(string_unicode) + return md5hash.hexdigest() + + +def cache_dataset_results(dataset_name: str, dataset_type: str, samplelist: List, query_results: List): + """function to cache dataset query results to file + input dataset_name and type query_results(already processed in default dict format) + """ + # data computations actions + # store the file path on redis + + table_timestamp = query_table_timestamp(dataset_type) + samplelist_as_str = ",".join(samplelist) + + file_name = generate_hash_file(dataset_name, dataset_type, table_timestamp, samplelist_as_str) + file_path = os.path.join(TMPDIR, f"{file_name}.json") + + with open(file_path, "w") as file_handler: + json.dump(query_results, file_handler) + + +def fetch_cached_results(dataset_name: str, dataset_type: str, samplelist: List): + """function to fetch the cached results""" + + table_timestamp = query_table_timestamp(dataset_type) + samplelist_as_str = ",".join(samplelist) + + file_name = generate_hash_file(dataset_name, dataset_type, table_timestamp, samplelist_as_str) + file_path = os.path.join(TMPDIR, f"{file_name}.json") + try: + with open(file_path, "r") as file_handler: + + return json.load(file_handler) + + except Exception: + pass |