diff options
| -rw-r--r-- | wqflask/wqflask/marker_regression/plink_mapping.py | 16 | ||||
| -rw-r--r-- | wqflask/wqflask/marker_regression/qtlreaper_mapping.py | 2 | ||||
| -rw-r--r-- | wqflask/wqflask/marker_regression/rqtl_mapping.py | 38 | ||||
| -rw-r--r-- | wqflask/wqflask/model.py | 33 | ||||
| -rw-r--r-- | wqflask/wqflask/network_graph/network_graph.py | 32 | ||||
| -rw-r--r-- | wqflask/wqflask/show_trait/export_trait_data.py | 2 | ||||
| -rw-r--r-- | wqflask/wqflask/show_trait/show_trait.py | 39 | ||||
| -rw-r--r-- | wqflask/wqflask/templates/show_trait_mapping_tools.html | 17 | ||||
| -rw-r--r-- | wqflask/wqflask/tracer.py | 41 | ||||
| -rw-r--r-- | wqflask/wqflask/views.py | 5 |
10 files changed, 29 insertions, 196 deletions
diff --git a/wqflask/wqflask/marker_regression/plink_mapping.py b/wqflask/wqflask/marker_regression/plink_mapping.py index 4de88f00..2f327faf 100644 --- a/wqflask/wqflask/marker_regression/plink_mapping.py +++ b/wqflask/wqflask/marker_regression/plink_mapping.py @@ -10,9 +10,7 @@ logger = utility.logger.getLogger(__name__ ) def run_plink(this_trait, dataset, species, vals, maf): plink_output_filename = webqtlUtil.genRandStr("%s_%s_"%(dataset.group.name, this_trait.name)) - gen_pheno_txt_file(dataset, vals) - #gen_pheno_txt_file_plink(this_trait, dataset, vals, pheno_filename = plink_output_filename) plink_command = PLINK_COMMAND + ' --noweb --bfile %s/%s --no-pheno --no-fid --no-parents --no-sex --maf %s --out %s%s --assoc ' % ( flat_files('mapping'), dataset.group.name, maf, TMPDIR, plink_output_filename) @@ -22,12 +20,6 @@ def run_plink(this_trait, dataset, species, vals, maf): count, p_values = parse_plink_output(plink_output_filename, species) - #for marker in self.dataset.group.markers.markers: - # if marker['name'] not in included_markers: - # logger.debug("marker:", marker) - # self.dataset.group.markers.markers.remove(marker) - # #del self.dataset.group.markers.markers[marker] - logger.debug("p_values:", p_values) dataset.group.markers.add_pvalues(p_values) @@ -108,7 +100,6 @@ def parse_plink_output(output_filename, species): result_fp = open("%s%s.qassoc"% (TMPDIR, output_filename), "rb") - header_line = result_fp.readline()# read header line line = result_fp.readline() value_list = [] # initialize value list, this list will include snp, bp and pvalue info @@ -156,11 +147,6 @@ def parse_plink_output(output_filename, species): else: line = result_fp.readline() - #if p_value_list: - # min_p_value = min(p_value_list) - #else: - # min_p_value = 0 - return count, p_value_dict ###################################################### @@ -173,4 +159,4 @@ def build_line_list(line=None): line_list = [item for item in line_list if item <>''] line_list = map(string.strip, line_list) - return line_list + return line_list \ No newline at end of file diff --git a/wqflask/wqflask/marker_regression/qtlreaper_mapping.py b/wqflask/wqflask/marker_regression/qtlreaper_mapping.py index 6b58190f..ffbfb5c5 100644 --- a/wqflask/wqflask/marker_regression/qtlreaper_mapping.py +++ b/wqflask/wqflask/marker_regression/qtlreaper_mapping.py @@ -26,7 +26,7 @@ def gen_reaper_results(this_trait, dataset, samples_before, trait_vals, json_dat perm_output = genotype.permutation(strains = trimmed_samples, trait = trimmed_values, nperm=num_perm) suggestive = perm_output[int(num_perm*0.37-1)] significant = perm_output[int(num_perm*0.95-1)] - highly_significant = perm_output[int(num_perm*0.99-1)] + #highly_significant = perm_output[int(num_perm*0.99-1)] #ZS: Currently not used, but leaving it here just in case json_data['suggestive'] = suggestive json_data['significant'] = significant diff --git a/wqflask/wqflask/marker_regression/rqtl_mapping.py b/wqflask/wqflask/marker_regression/rqtl_mapping.py index f3694f0b..41d67012 100644 --- a/wqflask/wqflask/marker_regression/rqtl_mapping.py +++ b/wqflask/wqflask/marker_regression/rqtl_mapping.py @@ -5,15 +5,16 @@ from base.webqtlConfig import TMPDIR from utility import webqtlUtil from utility.tools import locate, TEMPDIR +import utility.logger +logger = utility.logger.getLogger(__name__ ) + def run_rqtl_geno(vals, dataset, method, model, permCheck, num_perm, do_control, control_marker, manhattan_plot, pair_scan): geno_to_rqtl_function(dataset) ## Get pointers to some common R functions r_library = ro.r["library"] # Map the library function r_c = ro.r["c"] # Map the c function - r_sum = ro.r["sum"] # Map the sum function plot = ro.r["plot"] # Map the plot function - postscript = ro.r["postscript"] # Map the postscript function png = ro.r["png"] # Map the png function dev_off = ro.r["dev.off"] # Map the device off function @@ -23,17 +24,13 @@ def run_rqtl_geno(vals, dataset, method, model, permCheck, num_perm, do_control, scanone = ro.r["scanone"] # Map the scanone function scantwo = ro.r["scantwo"] # Map the scantwo function calc_genoprob = ro.r["calc.genoprob"] # Map the calc.genoprob function - read_cross = ro.r["read.cross"] # Map the read.cross function - write_cross = ro.r["write.cross"] # Map the write.cross function GENOtoCSVR = ro.r["GENOtoCSVR"] # Map the local GENOtoCSVR function crossname = dataset.group.name genofilelocation = locate(crossname + ".geno", "genotype") crossfilelocation = TMPDIR + crossname + ".cross" - #print("Conversion of geno to cross at location:", genofilelocation, " to ", crossfilelocation) - - cross_object = GENOtoCSVR(genofilelocation, crossfilelocation) # TODO: Add the SEX if that is available + cross_object = GENOtoCSVR(genofilelocation, crossfilelocation) # TODO: Add the SEX if that is available if manhattan_plot: cross_object = calc_genoprob(cross_object) @@ -42,18 +39,14 @@ def run_rqtl_geno(vals, dataset, method, model, permCheck, num_perm, do_control, cross_object = add_phenotype(cross_object, sanitize_rqtl_phenotype(vals)) # Add the phenotype - # for debug: write_cross(cross_object, "csvr", "test.csvr") - # Scan for QTLs - covar = create_covariates(control_marker, cross_object) # Create the additive covariate matrix + covar = create_covariates(control_marker, cross_object) # Create the additive covariate matrix if pair_scan: - if do_control == "true": # If sum(covar) > 0 we have a covariate matrix - print("Using covariate"); result_data_frame = scantwo(cross_object, pheno = "the_pheno", addcovar = covar, model=model, method=method, n_cluster = 16) + if do_control == "true": + logger.info("Using covariate"); result_data_frame = scantwo(cross_object, pheno = "the_pheno", addcovar = covar, model=model, method=method, n_cluster = 16) else: - print("No covariates"); result_data_frame = scantwo(cross_object, pheno = "the_pheno", model=model, method=method, n_cluster = 16) - - #print("Pair scan results:", result_data_frame) + logger.info("No covariates"); result_data_frame = scantwo(cross_object, pheno = "the_pheno", model=model, method=method, n_cluster = 16) pair_scan_filename = webqtlUtil.genRandStr("scantwo_") + ".png" png(file=TEMPDIR+pair_scan_filename) @@ -63,9 +56,9 @@ def run_rqtl_geno(vals, dataset, method, model, permCheck, num_perm, do_control, return process_pair_scan_results(result_data_frame) else: if do_control == "true": - print("Using covariate"); result_data_frame = scanone(cross_object, pheno = "the_pheno", addcovar = covar, model=model, method=method) + logger.info("Using covariate"); result_data_frame = scanone(cross_object, pheno = "the_pheno", addcovar = covar, model=model, method=method) else: - print("No covariates"); result_data_frame = scanone(cross_object, pheno = "the_pheno", model=model, method=method) + logger.info("No covariates"); result_data_frame = scanone(cross_object, pheno = "the_pheno", model=model, method=method) if num_perm > 0 and permCheck == "ON": # Do permutation (if requested by user) if do_control == "true": @@ -79,7 +72,6 @@ def run_rqtl_geno(vals, dataset, method, model, permCheck, num_perm, do_control, return process_rqtl_results(result_data_frame) def geno_to_rqtl_function(dataset): # TODO: Need to figure out why some genofiles have the wrong format and don't convert properly - ro.r(""" trim <- function( x ) { gsub("(^[[:space:]]+|[[:space:]]+$)", "", x) } @@ -117,15 +109,13 @@ def add_phenotype(cross, pheno_as_string): def create_covariates(control_marker, cross): ro.globalenv["the_cross"] = cross ro.r('genotypes <- pull.geno(the_cross)') # Get the genotype matrix - userinputS = control_marker.replace(" ", "").split(",") # TODO: sanitize user input, Never Ever trust a user + userinputS = control_marker.replace(" ", "").split(",") # TODO: sanitize user input, Never Ever trust a user covariate_names = ', '.join('"{0}"'.format(w) for w in userinputS) - #print("Marker names of selected covariates:", covariate_names) ro.r('covnames <- c(' + covariate_names + ')') ro.r('covInGeno <- which(covnames %in% colnames(genotypes))') ro.r('covnames <- covnames[covInGeno]') ro.r("cat('covnames (purged): ', covnames,'\n')") ro.r('covariates <- genotypes[,covnames]') # Get the covariate matrix by using the marker name as index to the genotype file - #print("R/qtl matrix of covariates:", ro.r["covariates"]) return ro.r["covariates"] def sanitize_rqtl_phenotype(vals): @@ -149,7 +139,6 @@ def process_pair_scan_results(result): result = result[1] output = [tuple([result[j][i] for j in range(result.ncol)]) for i in range(result.nrow)] - #print("R/qtl scantwo output:", output) for i, line in enumerate(result.iter_row()): marker = {} @@ -175,9 +164,7 @@ def process_rqtl_perm_results(num_perm, results): def process_rqtl_results(result): # TODO: how to make this a one liner and not copy the stuff in a loop qtl_results = [] - output = [tuple([result[j][i] for j in range(result.ncol)]) for i in range(result.nrow)] - #print("R/qtl scanone output:", output) for i, line in enumerate(result.iter_row()): marker = {} @@ -187,5 +174,4 @@ def process_rqtl_results(result): # TODO: how to make this a one liner an marker['lod_score'] = output[i][2] qtl_results.append(marker) - return qtl_results - + return qtl_results \ No newline at end of file diff --git a/wqflask/wqflask/model.py b/wqflask/wqflask/model.py index 5321e420..38117a8e 100644 --- a/wqflask/wqflask/model.py +++ b/wqflask/wqflask/model.py @@ -6,31 +6,16 @@ import datetime import simplejson as json from flask import request -from flask.ext.sqlalchemy import SQLAlchemy from wqflask import app import sqlalchemy - -from sqlalchemy import (Column, Integer, String, Table, ForeignKey, Unicode, Boolean, DateTime, +from sqlalchemy import (Column, ForeignKey, Unicode, Boolean, DateTime, Text, Index) -from sqlalchemy.orm import relationship, backref +from sqlalchemy.orm import relationship from wqflask.database import Base, init_db - - -# Define models -#roles_users = Table('roles_users', -# Column('user_id', Integer(), ForeignKey('user.the_id')), -# Column('role_id', Integer(), ForeignKey('role.the_id'))) - -#class Role(Base): -# __tablename__ = "role" -# id = Column(Unicode(36), primary_key=True, default=lambda: unicode(uuid.uuid4())) -# name = Column(Unicode(80), unique=True, nullable=False) -# description = Column(Unicode(255)) - class User(Base): __tablename__ = "user" id = Column(Unicode(36), primary_key=True, default=lambda: unicode(uuid.uuid4())) @@ -133,11 +118,6 @@ class User(Base): except IndexError: return None - - #roles = relationship('Role', secondary=roles_users, - # backref=backref('users', lazy='dynamic')) - - class Login(Base): __tablename__ = "login" id = Column(Unicode(36), primary_key=True, default=lambda: unicode(uuid.uuid4())) @@ -177,22 +157,15 @@ class UserCollection(Base): except: return 0 - #@property - #def display_num_members(self): - # return display_collapsible(self.num_members) - def members_as_set(self): return set(json.loads(self.members)) - def display_collapsible(number): if number: return number else: return "" - def user_uuid(): """Unique cookie for a user""" - user_uuid = request.cookies.get('user_uuid') - + user_uuid = request.cookies.get('user_uuid') \ No newline at end of file diff --git a/wqflask/wqflask/network_graph/network_graph.py b/wqflask/wqflask/network_graph/network_graph.py index b42904a4..63273a29 100644 --- a/wqflask/wqflask/network_graph/network_graph.py +++ b/wqflask/wqflask/network_graph/network_graph.py @@ -21,9 +21,7 @@ from __future__ import absolute_import, print_function, division import sys -# sys.path.append(".") Never do this in a webserver! -import gc import string import cPickle import os @@ -95,7 +93,6 @@ class NetworkGraph(object): self.lowest_overlap = 8 #ZS: Variable set to the lowest overlapping samples in order to notify user, or 8, whichever is lower (since 8 is when we want to display warning) - self.network_data = {} self.nodes_list = [] self.edges_list = [] for trait_db in self.trait_list: @@ -107,9 +104,9 @@ class NetworkGraph(object): corr_result_row = [] is_spearman = False #ZS: To determine if it's above or below the diagonal - + max_corr = 0 #ZS: Used to determine whether node should be hidden when correlation coefficient slider is used - + for target in self.trait_list: target_trait = target[0] target_db = target[1] @@ -141,7 +138,7 @@ class NetworkGraph(object): continue else: sample_r, sample_p = scipy.stats.spearmanr(this_trait_vals, target_vals) - + if -1 <= sample_r < -0.7: color = "#0000ff" width = 3 @@ -163,10 +160,10 @@ class NetworkGraph(object): else: color = "#000000" width = 0 - + if abs(sample_r) > max_corr: max_corr = abs(sample_r) - + edge_data = {'id' : str(this_trait.name) + '_to_' + str(target_trait.name), 'source' : str(this_trait.name) + ":" + str(this_trait.dataset.name), 'target' : str(target_trait.name) + ":" + str(target_trait.dataset.name), @@ -176,11 +173,11 @@ class NetworkGraph(object): 'overlap' : num_overlap, 'color' : color, 'width' : width } - + edge_dict = { 'data' : edge_data } - + self.edges_list.append(edge_dict) - + if trait_db[1].type == "ProbeSet": node_dict = { 'data' : {'id' : str(this_trait.name) + ":" + str(this_trait.dataset.name), 'label' : this_trait.symbol, @@ -197,19 +194,13 @@ class NetworkGraph(object): 'label' : this_trait.name, 'max_corr' : max_corr } } self.nodes_list.append(node_dict) - - #self.network_data['dataSchema'] = {'nodes' : [{'name' : "label" , 'type' : "string"}], - # 'edges' : [{'name' : "label" , 'type' : "string"}] } - - #self.network_data['data'] = {'nodes' : self.nodes_list, - # 'edges' : self.edges_list } self.elements = json.dumps(self.nodes_list + self.edges_list) - + groups = [] for sample in self.all_sample_list: groups.append(1) - + self.js_data = dict(traits = [trait.name for trait in self.traits], groups = groups, cols = range(len(self.traits)), @@ -217,7 +208,6 @@ class NetworkGraph(object): samples = self.all_sample_list, sample_data = self.sample_data, elements = self.elements,) - # corr_results = [result[1] for result in result_row for result_row in self.corr_results]) def get_trait_db_obs(self, trait_db_list): self.trait_list = [] @@ -229,4 +219,4 @@ class NetworkGraph(object): trait_ob = GeneralTrait(dataset=dataset_ob, name=trait_name, cellid=None) - self.trait_list.append((trait_ob, dataset_ob)) + self.trait_list.append((trait_ob, dataset_ob)) \ No newline at end of file diff --git a/wqflask/wqflask/show_trait/export_trait_data.py b/wqflask/wqflask/show_trait/export_trait_data.py index ac3cd366..558372bb 100644 --- a/wqflask/wqflask/show_trait/export_trait_data.py +++ b/wqflask/wqflask/show_trait/export_trait_data.py @@ -1,7 +1,5 @@ from __future__ import print_function, division -import operator - import simplejson as json from pprint import pformat as pf diff --git a/wqflask/wqflask/show_trait/show_trait.py b/wqflask/wqflask/show_trait/show_trait.py index 1f000564..8b801396 100644 --- a/wqflask/wqflask/show_trait/show_trait.py +++ b/wqflask/wqflask/show_trait/show_trait.py @@ -6,7 +6,6 @@ import datetime import cPickle import uuid import json as json -#import pyXLWriter as xl from collections import OrderedDict @@ -73,14 +72,6 @@ class ShowTrait(object): cellid=None) self.trait_vals = Redis.get(self.trait_id).split() - #self.dataset.group.read_genotype_file() - - #if this_trait: - # if this_trait.dataset and this_trait.dataset.type and this_trait.dataset.type == 'ProbeSet': - # self.cursor.execute("SELECT h2 from ProbeSetXRef WHERE DataId = %d" % - # this_trait.mysqlid) - # heritability = self.cursor.fetchone() - #ZS: Get verify/rna-seq link URLs try: blatsequence = self.this_trait.blatseq @@ -192,7 +183,7 @@ class ShowTrait(object): self.sample_group_types['samples_primary'] = self.dataset.group.name sample_lists = [group.sample_list for group in self.sample_groups] - self.get_mapping_methods() + self.genofiles = get_genofiles(self.dataset) self.stats_table_width, self.trait_table_width = get_table_widths(self.sample_groups) @@ -211,28 +202,6 @@ class ShowTrait(object): temp_uuid = self.temp_uuid) self.js_data = js_data - def get_mapping_methods(self): - '''Only display mapping methods when the dataset group's genotype file exists''' - def check_plink_gemma(): - if flat_file_exists("mapping"): - MAPPING_PATH = flat_files("mapping")+"/" - if (os.path.isfile(MAPPING_PATH+self.dataset.group.name+".bed") and - (os.path.isfile(MAPPING_PATH+self.dataset.group.name+".map") or - os.path.isfile(MAPPING_PATH+self.dataset.group.name+".bim"))): - return True - return False - - def check_pylmm_rqtl(): - if os.path.isfile(webqtlConfig.GENODIR+self.dataset.group.name+".geno") and (os.path.getsize(webqtlConfig.JSON_GENODIR+self.dataset.group.name+".json") > 0): - return True - else: - return False - - self.genofiles = get_genofiles(self.dataset) - self.use_plink_gemma = check_plink_gemma() - self.use_pylmm_rqtl = check_pylmm_rqtl() - - def build_correlation_tools(self): if self.temp_trait == True: this_group = self.temp_group @@ -245,7 +214,6 @@ class ShowTrait(object): this_group = 'BXD' if this_group: - #dataset_menu = self.dataset.group.datasets() if self.temp_trait == True: dataset_menu = data_set.datasets(this_group) else: @@ -263,7 +231,6 @@ class ShowTrait(object): return_results_menu = return_results_menu, return_results_menu_selected = return_results_menu_selected,) - def make_sample_lists(self): all_samples_ordered = self.dataset.group.all_samples_ordered() @@ -315,10 +282,6 @@ class ShowTrait(object): sample_group_type='primary', header="%s Only" % (self.dataset.group.name)) self.sample_groups = (primary_samples,) - #TODO: Figure out why this if statement is written this way - Zach - #if (other_sample_names or (fd.f1list and this_trait.data.has_key(fd.f1list[0])) - # or (fd.f1list and this_trait.data.has_key(fd.f1list[1]))): - # logger.debug("hjs") self.dataset.group.allsamples = all_samples_ordered def get_nearest_marker(this_trait, this_db): diff --git a/wqflask/wqflask/templates/show_trait_mapping_tools.html b/wqflask/wqflask/templates/show_trait_mapping_tools.html index 03590c2c..0ecf1eb9 100644 --- a/wqflask/wqflask/templates/show_trait_mapping_tools.html +++ b/wqflask/wqflask/templates/show_trait_mapping_tools.html @@ -35,7 +35,6 @@ </ul> <div class="tab-content"> - {# if use_pylmm_rqtl and not use_plink_gemma and dataset.group.species != "human" #} {% if dataset.group.mapping_id == "1" %} <div class="tab-pane active" id="gemma"> <div style="padding-top: 20px;" class="form-horizontal"> @@ -70,10 +69,6 @@ </label> </div> </div> - <!-- - </div> - <div style="padding-top: 5px; padding-bottom: 5px; padding-left: 20px;" class="form-horizontal"> - --> <div class="mapping_method_fields form-group"> <label style="text-align: right;" class="col-xs-3 control-label">Covariates</label> <div style="margin-left:20px;" class="col-xs-7"> @@ -99,18 +94,6 @@ </div> </div> </div> -<!-- - <div class="form-group"> - <div class="col-xs-4 controls"> - <label class="col-xs-2 control-label"></label> - <div class="col-xs-4"> - <button id="gemma_bimbam_compute" class="btn submit_special btn-success" data-url="/marker_regression" title="Compute Marker Regression"> - Compute - </button> - </div> - </div> - </div> ---> </div> <div class="tab-pane" id="interval_mapping"> <div style="margin-top: 20px" class="form-horizontal"> diff --git a/wqflask/wqflask/tracer.py b/wqflask/wqflask/tracer.py deleted file mode 100644 index a1043d28..00000000 --- a/wqflask/wqflask/tracer.py +++ /dev/null @@ -1,41 +0,0 @@ -from __future__ import absolute_import, division, print_function - -print("At top of tracer") - -import sys - -#################################################################################### - -# Originally based on http://stackoverflow.com/a/8315566 -def tracefunc(frame, event, arg, indent=[0]): - - func = dict(funcname = frame.f_code.co_name, - filename = frame.f_code.co_filename, - lineno = frame.f_lineno) - - #These are too common to bother printing... - too_common = ( - '/home/sam/ve27/local/lib/python2.7/site-packages/werkzeug/', - '/home/sam/ve27/local/lib/python2.7/site-packages/jinja2/', - ) - - - if func['filename'].startswith(too_common): - return tracefunc - - info = "{funcname} [{filename}: {lineno}]".format(**func) - - if event == "call": - indent[0] += 2 - #print("-" * indent[0] + "> call function", frame.f_code.co_name) - print("-" * indent[0] + "> call function:", info) - elif event == "return": - print("<" + "-" * indent[0], "exit function:", info) - indent[0] -= 2 - return tracefunc - -def turn_on(): - sys.settrace(tracefunc) - print("Tracing turned on!!!!") -#################################################################################### - diff --git a/wqflask/wqflask/views.py b/wqflask/wqflask/views.py index 8ff359a7..49b47123 100644 --- a/wqflask/wqflask/views.py +++ b/wqflask/wqflask/views.py @@ -89,11 +89,6 @@ def shutdown_session(exception=None): db_session.remove() g.db = None -#@app.before_request -#def trace_it(): -# from wqflask import tracer -# tracer.turn_on() - @app.errorhandler(Exception) def handle_bad_request(e): err_msg = str(e) |