From 102523493e2f8a7660c63f117f1d8dfd009eff02 Mon Sep 17 00:00:00 2001 From: Pjotr Prins Date: Mon, 13 Apr 2015 08:14:43 +0000 Subject: Improved assertion message --- wqflask/wqflask/my_pylmm/pyLMM/tsvreader.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/wqflask/wqflask/my_pylmm/pyLMM/tsvreader.py b/wqflask/wqflask/my_pylmm/pyLMM/tsvreader.py index b4027fa3..b24ffe8f 100644 --- a/wqflask/wqflask/my_pylmm/pyLMM/tsvreader.py +++ b/wqflask/wqflask/my_pylmm/pyLMM/tsvreader.py @@ -56,7 +56,8 @@ def geno(fn): print fn with open(fn,'r') as tsvin: - assert(tsvin.readline().strip() == "# Genotype format version 1.0") + line = tsvin.readline().strip() + assert line == "# Genotype format version 1.0", line tsvin.readline() tsvin.readline() tsvin.readline() -- cgit v1.2.3 From a1d8f68d5428a4ceec9a2d9a771b000ecabec5e6 Mon Sep 17 00:00:00 2001 From: Pjotr Prins Date: Sat, 18 Apr 2015 09:36:00 +0000 Subject: pylmm: fix integration problems --- wqflask/runserver.py | 6 ++--- .../wqflask/marker_regression/marker_regression.py | 14 +++++----- wqflask/wqflask/my_pylmm/pyLMM/gn2.py | 30 +++++++++++++++------- wqflask/wqflask/my_pylmm/pyLMM/lmm.py | 19 ++++++++------ wqflask/wqflask/my_pylmm/pyLMM/lmm2.py | 18 ++++++++++--- 5 files changed, 56 insertions(+), 31 deletions(-) diff --git a/wqflask/runserver.py b/wqflask/runserver.py index 9d5686a9..fadae6bf 100755 --- a/wqflask/runserver.py +++ b/wqflask/runserver.py @@ -20,9 +20,9 @@ from wqflask import app import logging #from themodule import TheHandlerYouWant -file_handler = logging.FileHandler("/tmp/flask_gn_log_danny_unsecure") -file_handler.setLevel(logging.DEBUG) -app.logger.addHandler(file_handler) +# file_handler = logging.FileHandler("/tmp/flask_gn_log_danny_unsecure") +# file_handler.setLevel(logging.DEBUG) +# app.logger.addHandler(file_handler) import logging_tree logging_tree.printout() diff --git a/wqflask/wqflask/marker_regression/marker_regression.py b/wqflask/wqflask/marker_regression/marker_regression.py index 7708356b..ae3e062f 100755 --- a/wqflask/wqflask/marker_regression/marker_regression.py +++ b/wqflask/wqflask/marker_regression/marker_regression.py @@ -40,6 +40,7 @@ from utility import temp_data from utility.benchmark import Bench +PYLMM_COMMAND= 'python /home/pjotr/izip/git/opensource/python/gn2/wqflask/wqflask/my_pylmm/pyLMM/lmm.py' class MarkerRegression(object): @@ -272,7 +273,7 @@ class MarkerRegression(object): """) def run_rqtl_geno(self): - print("Calling R/qtl from python") + print("Calling R/qtl") self.geno_to_rqtl_function() @@ -655,8 +656,7 @@ class MarkerRegression(object): Redis.set(key, json_params) Redis.expire(key, 60*60) - command = 'python /home/zas1024/gene/wqflask/wqflask/my_pylmm/pyLMM/lmm.py --key {} --species {}'.format(key, - "other") + command = PYLMM_COMMAND+' --key {} --species {}'.format(key,"other") os.system(command) @@ -713,8 +713,8 @@ class MarkerRegression(object): # "refit": False, # "temp_data": tempdata} - print("genotype_matrix:", str(genotype_matrix.tolist())) - print("pheno_vector:", str(pheno_vector.tolist())) + # print("genotype_matrix:", str(genotype_matrix.tolist())) + # print("pheno_vector:", str(pheno_vector.tolist())) params = dict(pheno_vector = pheno_vector.tolist(), genotype_matrix = genotype_matrix.tolist(), @@ -732,7 +732,7 @@ class MarkerRegression(object): Redis.expire(key, 60*60) print("before printing command") - command = 'python /home/zas1024/gene/wqflask/wqflask/my_pylmm/pyLMM/lmm.py --key {} --species {}'.format(key, + command = PYLMM_COMMAND + ' --key {} --species {}'.format(key, "other") print("command is:", command) print("after printing command") @@ -806,7 +806,7 @@ class MarkerRegression(object): print("Before creating the command") - command = 'python /home/zas1024/gene/wqflask/wqflask/my_pylmm/pyLMM/lmm.py --key {} --species {}'.format(key, + command = PYLMM_COMMAND+' --key {} --species {}'.format(key, "human") print("command is:", command) diff --git a/wqflask/wqflask/my_pylmm/pyLMM/gn2.py b/wqflask/wqflask/my_pylmm/pyLMM/gn2.py index 7bceb089..b128bfab 100644 --- a/wqflask/wqflask/my_pylmm/pyLMM/gn2.py +++ b/wqflask/wqflask/my_pylmm/pyLMM/gn2.py @@ -1,7 +1,10 @@ -# Genenetwork2 specific methods and callback handler +# Standalone specific methods and callback handler # # Copyright (C) 2015 Pjotr Prins (pjotr.prins@thebird.nl) # +# Set the log level with +# +# logging.basicConfig(level=logging.DEBUG) from __future__ import absolute_import, print_function, division @@ -9,10 +12,12 @@ import numpy as np import sys import logging -# logging.basicConfig(level=logging.DEBUG) -# np.set_printoptions() +# logger = logging.getLogger(__name__) +logger = logging.getLogger('lmm2') +logging.basicConfig(level=logging.DEBUG) +np.set_printoptions(precision=3,suppress=True) -progress_location = None +progress_location = None progress_current = None progress_prev_perc = None @@ -20,30 +25,37 @@ def progress_default_func(location,count,total): global progress_current value = round(count*100.0/total) progress_current = value - + progress_func = progress_default_func def progress_set_func(func): global progress_func progress_func = func - + def progress(location, count, total): global progress_location global progress_prev_perc - + perc = round(count*100.0/total) if perc != progress_prev_perc and (location != progress_location or perc > 98 or perc > progress_prev_perc + 5): progress_func(location, count, total) logger.info("Progress: %s %d%%" % (location,perc)) progress_location = location progress_prev_perc = perc - + def mprint(msg,data): """ Array/matrix print function """ m = np.array(data) - print(msg,m.shape,"=\n",m) + if m.ndim == 1: + print(msg,m.shape,"=\n",m[0:3]," ... ",m[-3:]) + if m.ndim == 2: + print(msg,m.shape,"=\n[", + m[0][0:3]," ... ",m[0][-3:],"\n ", + m[1][0:3]," ... ",m[1][-3:],"\n ...\n ", + m[-2][0:3]," ... ",m[-2][-3:],"\n ", + m[-1][0:3]," ... ",m[-1][-3:],"]") def fatal(msg): logger.critical(msg) diff --git a/wqflask/wqflask/my_pylmm/pyLMM/lmm.py b/wqflask/wqflask/my_pylmm/pyLMM/lmm.py index b2067b27..6fff5f1d 100644 --- a/wqflask/wqflask/my_pylmm/pyLMM/lmm.py +++ b/wqflask/wqflask/my_pylmm/pyLMM/lmm.py @@ -42,25 +42,27 @@ from redis import Redis Redis = Redis() import sys -sys.path.append("/home/zas1024/gene/wqflask/") - -has_gn2=True from utility.benchmark import Bench from utility import temp_data -sys.path.append("/home/zas1024/gene/wqflask/wqflask/my_pylmm/pyLMM/") - from kinship import kinship, kinship_full, kvakve import genotype import phenotype import gwas +has_gn2=True +sys.stderr.write("INFO: pylmm system path is "+":".join(sys.path)+"\n") +sys.stderr.write("INFO: pylmm file is "+__file__+"\n") + # ---- A trick to decide on the environment: try: - from wqflask.my_pylmm.pyLMM import chunks + sys.stderr.write("INFO: trying loading module\n") + import utility.formatting # this is never used, just to check the environment + sys.stderr.write("INFO: This is a genenetwork2 environment\n") from gn2 import uses, progress_set_func except ImportError: + # Failed to load gn2 has_gn2=False import standalone as handlers from standalone import uses, progress_set_func @@ -856,7 +858,8 @@ def gwas_with_redis(key,species,new_code=True): print(key) v = params[key] if v is not None: - v = np.array(v) + v = np.array(v).astype(np.float) + print(v) return v def narrayT(key): @@ -969,6 +972,6 @@ if __name__ == '__main__': if has_gn2: gn2_main() else: - print("Run from runlmm.py instead") + fatal("Run from runlmm.py instead") diff --git a/wqflask/wqflask/my_pylmm/pyLMM/lmm2.py b/wqflask/wqflask/my_pylmm/pyLMM/lmm2.py index 358bf27e..c65843ec 100644 --- a/wqflask/wqflask/my_pylmm/pyLMM/lmm2.py +++ b/wqflask/wqflask/my_pylmm/pyLMM/lmm2.py @@ -24,14 +24,24 @@ from scipy import optimize from optmatrix import matrixMult import kinship +sys.stderr.write("INFO: pylmm (lmm2) system path is "+":".join(sys.path)+"\n") +sys.stderr.write("INFO: pylmm (lmm2) file is "+__file__+"\n") + # ---- A trick to decide on the environment: try: - from wqflask.my_pylmm.pyLMM import chunks - from gn2 import uses + sys.stderr.write("INFO: trying loading module\n") + import utility.formatting # this is never used, just to check the environment + sys.stderr.write("INFO: This is a genenetwork2 environment (lmm2)\n") + from gn2 import uses, progress_set_func except ImportError: - sys.stderr.write("WARNING: LMM2 standalone version missing the Genenetwork2 environment\n") + # Failed to load gn2 has_gn2=False - from standalone import uses + import standalone as handlers + from standalone import uses, progress_set_func + sys.stderr.write("WARNING: LMM2 standalone version missing the Genenetwork2 environment\n") + +progress,mprint,debug,info,fatal = uses('progress','mprint','debug','info','fatal') + def calculateKinship(W,center=False): """ -- cgit v1.2.3 From 02660b9406a97943d4c33946250fc3f08b80c556 Mon Sep 17 00:00:00 2001 From: Pjotr Prins Date: Sat, 18 Apr 2015 09:40:57 +0000 Subject: pylmm: fix integration problems --- wqflask/wqflask/my_pylmm/pyLMM/lmm.py | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/wqflask/wqflask/my_pylmm/pyLMM/lmm.py b/wqflask/wqflask/my_pylmm/pyLMM/lmm.py index 6fff5f1d..618f8332 100644 --- a/wqflask/wqflask/my_pylmm/pyLMM/lmm.py +++ b/wqflask/wqflask/my_pylmm/pyLMM/lmm.py @@ -854,22 +854,23 @@ def gwas_with_redis(key,species,new_code=True): debug("Updating REDIS percent_complete=%d" % (round(i*100.0/total))) progress_set_func(update_tempdata) - def narray(key): - print(key) - v = params[key] + def narray(t): + info("Type is "+t) + v = params[t] if v is not None: v = np.array(v).astype(np.float) print(v) return v - def narrayT(key): - m = narray(key) + def narrayT(t): + m = narray(t) if m is not None: return m.T return m # We are transposing before we enter run_gwas - this should happen on the webserver # side (or when reading data from file) + print(params) k = narray('kinship_matrix') g = narrayT('genotype_matrix') y = narray('pheno_vector') -- cgit v1.2.3 From 8706319923b3830a4d8cd63fd9a3f6b9a2b04563 Mon Sep 17 00:00:00 2001 From: Pjotr Prins Date: Sat, 18 Apr 2015 10:58:26 +0000 Subject: Fix NA to float tranforms --- wqflask/wqflask/my_pylmm/pyLMM/lmm.py | 37 +++++++++++++++++++++++++---------- 1 file changed, 27 insertions(+), 10 deletions(-) diff --git a/wqflask/wqflask/my_pylmm/pyLMM/lmm.py b/wqflask/wqflask/my_pylmm/pyLMM/lmm.py index 618f8332..5b06c9ae 100644 --- a/wqflask/wqflask/my_pylmm/pyLMM/lmm.py +++ b/wqflask/wqflask/my_pylmm/pyLMM/lmm.py @@ -856,30 +856,47 @@ def gwas_with_redis(key,species,new_code=True): def narray(t): info("Type is "+t) - v = params[t] + v = params.get(t) if v is not None: - v = np.array(v).astype(np.float) - print(v) + # Note input values can be array of string or float + v1 = [x if x != 'NA' else 'nan' for x in v] + v = np.array(v1).astype(np.float) return v - def narrayT(t): - m = narray(t) + def marray(t): + info("Type is "+t) + v = params.get(t) + if v is not None: + m = [] + for r in v: + # Note input values can be array of string or float + r1 = [x if x != 'NA' else 'nan' for x in r] + m.append(np.array(r1).astype(np.float)) + return np.array(m) + return np.array(v) + + def marrayT(t): + m = marray(t) if m is not None: return m.T return m # We are transposing before we enter run_gwas - this should happen on the webserver # side (or when reading data from file) - print(params) - k = narray('kinship_matrix') - g = narrayT('genotype_matrix') + k = marray('kinship_matrix') + g = marrayT('genotype_matrix') + mprint("geno",g) y = narray('pheno_vector') n = len(y) - m = params['num_genotypes'] - ps,ts = run_gwas(species,n,m,k,y,g,narray('covariate_matrix'),params['restricted_max_likelihood'],params['refit'],params['input_file_name'],new_code) + m = params.get('num_genotypes') + if m is None: + m = g.shape[0] + info("m=%d,n=%d" % (m,n)) + ps,ts = run_gwas(species,n,m,k,y,g,narray('covariate_matrix'),params['restricted_max_likelihood'],params['refit'],params.get('input_file_name'),new_code) results_key = "pylmm:results:" + params['temp_uuid'] + # fatal(results_key) json_results = json.dumps(dict(p_values = ps, t_stats = ts)) -- cgit v1.2.3 From ced6f0c49c155a2ab47adfe93578d4718504566b Mon Sep 17 00:00:00 2001 From: Pjotr Prins Date: Sat, 18 Apr 2015 11:09:20 +0000 Subject: Disable some print statements - will introduce debug levels soon --- wqflask/wqflask/marker_regression/marker_regression.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/wqflask/wqflask/marker_regression/marker_regression.py b/wqflask/wqflask/marker_regression/marker_regression.py index ae3e062f..c80bba8e 100755 --- a/wqflask/wqflask/marker_regression/marker_regression.py +++ b/wqflask/wqflask/marker_regression/marker_regression.py @@ -128,7 +128,7 @@ class MarkerRegression(object): #Need to convert the QTL objects that qtl reaper returns into a json serializable dictionary self.qtl_results = [] for qtl in self.filtered_markers: - print("lod score is:", qtl['lod_score']) + # print("lod score is:", qtl['lod_score']) if qtl['chr'] == highest_chr and highest_chr != "X" and highest_chr != "X/Y": print("changing to X") self.json_data['chr'].append("X") @@ -145,7 +145,7 @@ class MarkerRegression(object): self.json_data['chrnames'].append([self.species.chromosomes.chromosomes[key].name, self.species.chromosomes.chromosomes[key].mb_length]) chromosome_mb_lengths[key] = self.species.chromosomes.chromosomes[key].mb_length - print("json_data:", self.json_data) + # print("json_data:", self.json_data) self.js_data = dict( @@ -745,7 +745,7 @@ class MarkerRegression(object): json_results = Redis.blpop("pylmm:results:" + temp_uuid, 45*60) results = json.loads(json_results[1]) p_values = [float(result) for result in results['p_values']] - print("p_values:", p_values) + print("p_values:", p_values[:10]) #p_values = self.trim_results(p_values) t_stats = results['t_stats'] -- cgit v1.2.3 From 1f6386cbddfd02d8abbd4e9bcb502c06be6864d1 Mon Sep 17 00:00:00 2001 From: Pjotr Prins Date: Sat, 18 Apr 2015 12:41:07 +0000 Subject: Show first 40 LOD scores --- wqflask/wqflask/marker_regression/marker_regression.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/wqflask/wqflask/marker_regression/marker_regression.py b/wqflask/wqflask/marker_regression/marker_regression.py index c80bba8e..fba34b99 100755 --- a/wqflask/wqflask/marker_regression/marker_regression.py +++ b/wqflask/wqflask/marker_regression/marker_regression.py @@ -127,8 +127,9 @@ class MarkerRegression(object): #Need to convert the QTL objects that qtl reaper returns into a json serializable dictionary self.qtl_results = [] - for qtl in self.filtered_markers: - # print("lod score is:", qtl['lod_score']) + for index,qtl in enumerate(self.filtered_markers): + if index<40: + print("lod score is:", qtl['lod_score']) if qtl['chr'] == highest_chr and highest_chr != "X" and highest_chr != "X/Y": print("changing to X") self.json_data['chr'].append("X") -- cgit v1.2.3