Did some work towards getting Nick's code working with our plink

human genotype files

Tried to get the automatic searches working when points are highlighted
on the marker regression manhattan plot

6 files changed, 114 insertions, 64 deletions

diff --git a/wqflask/wqflask/do_search.py b/wqflask/wqflask/do_search.py
index a2eddfc6..b390bdc4 100644
--- a/wqflask/wqflask/do_search.py
+++ b/wqflask/wqflask/do_search.py
@@ -20,17 +20,18 @@ class DoSearch(object):
     # Used to translate search phrases into classes
     search_types = dict()
 
-    def __init__(self, search_term, search_operator, dataset):
+    def __init__(self, search_term, search_operator=None, dataset=None):
         self.search_term = search_term
         # Make sure search_operator is something we expect
         assert search_operator in (None, "=", "<", ">", "<=", ">="), "Bad search operator"
         self.search_operator = search_operator
         self.dataset = dataset
-        print("self.dataset is boo: ", type(self.dataset), pf(self.dataset))
-        print("self.dataset.group is: ", pf(self.dataset.group))
-
-        #Get group information for dataset and the species id
-        self.species_id = webqtlDatabaseFunction.retrieve_species_id(self.dataset.group.name)           
+        
+        if self.dataset:
+            print("self.dataset is boo: ", type(self.dataset), pf(self.dataset))
+            print("self.dataset.group is: ", pf(self.dataset.group))
+            #Get group information for dataset and the species id
+            self.species_id = webqtlDatabaseFunction.retrieve_species_id(self.dataset.group.name)           
 
     def execute(self, query):
         """Executes query and returns results"""
@@ -59,6 +60,47 @@ class DoSearch(object):
     def get_search(cls, search_type):
         return cls.search_types[search_type]
 
+class QuickMrnaAssaySearch(DoSearch):
+    """A general search for mRNA assays"""
+    
+    DoSearch.search_types['quick_mrna_assay'] = "QuickMrnaAssaySearch"
+    
+    base_query = """SELECT ProbeSet.Name as ProbeSet_Name,
+                ProbeSet.Chr_num as ProbeSet_Chr_Num,
+                ProbeSet.Mb as ProbeSet_Mb,
+                ProbeSet.Symbol as ProbeSet_Symbol,
+                ProbeSet.name_num as ProbeSet_name_num
+                FROM ProbeSet,
+                ProbeSetXRef,
+                ProbeSetFreeze,
+                ProbeFreeze,
+                InbredSet,
+                Species """
+                
+    header_fields = ['',
+                     'Record ID',
+                     'Symbol',
+                     'Location']
+    
+    def run(self):
+        """Generates and runs a search for assays across all mRNA expression datasets"""
+
+        print("Running ProbeSetSearch")
+        query = self.base_query + """WHERE (MATCH (ProbeSet.Name,
+                    ProbeSet.description,
+                    ProbeSet.symbol)
+                    AGAINST ('%s' IN BOOLEAN MODE)) and
+                    ProbeSet.Id = ProbeSetXRef.ProbeSetId and
+                    ProbeSetXRef.ProbeSetFreezeId = ProbeSetFreeze.Id and
+                    ProbeSetFreeze.ProbeFreezeId = ProbeFreeze.Id and
+                    ProbeFreeze.InbredSetId = InbredSet.Id and
+                    InbredSet.SpeciesId = Species.Id
+                            """ % (escape(self.search_term[0]))
+
+        print("final query is:", pf(query))
+
+        return self.execute(query)
+
 
 class MrnaAssaySearch(DoSearch):
     """A search within an mRNA expression dataset"""
@@ -125,6 +167,8 @@ class MrnaAssaySearch(DoSearch):
         return self.execute(query)
 
 
+#class QuickPhenotypeSearch(DoSearch):
+
 class PhenotypeSearch(DoSearch):
     """A search within a phenotype dataset"""
 
diff --git a/wqflask/wqflask/my_pylmm/pyLMM/lmm.py b/wqflask/wqflask/my_pylmm/pyLMM/lmm.py
index 66830662..1e689e49 100644
--- a/wqflask/wqflask/my_pylmm/pyLMM/lmm.py
+++ b/wqflask/wqflask/my_pylmm/pyLMM/lmm.py
@@ -26,42 +26,42 @@ from scipy import stats
 
 from pprint import pformat as pf
 
-#from utility.benchmark import Bench
-#
-##np.seterr('raise')
-#
-#def run(pheno_vector,
-#        genotype_matrix,
-#        restricted_max_likelihood=True,
-#        refit=False,
-#        temp_data=None):
-#    """Takes the phenotype vector and genotype matrix and returns a set of p-values and t-statistics
-#    
-#    restricted_max_likelihood -- whether to use restricted max likelihood; True or False
-#    refit -- whether to refit the variance component for each marker
-#    temp_data -- TempData object that stores the progress for each major step of the
-#    calculations ("calculate_kinship" and "GWAS" take the majority of time)
-#    
-#    """
-#    
-#    with Bench("Calculate Kinship"):
-#        kinship_matrix = calculate_kinship(genotype_matrix, temp_data)
-#    
-#    with Bench("Create LMM object"):
-#        lmm_ob = LMM(pheno_vector, kinship_matrix)
-#    
-#    with Bench("LMM_ob fitting"):
-#        lmm_ob.fit()
-#
-#    with Bench("Doing GWAS"):
-#        t_stats, p_values = GWAS(pheno_vector,
-#                                genotype_matrix,
-#                                kinship_matrix,
-#                                restricted_max_likelihood=True,
-#                                refit=False,
-#                                temp_data=temp_data)
-#    Bench().report()
-#    return t_stats, p_values
+from utility.benchmark import Bench
+
+#np.seterr('raise')
+
+def run(pheno_vector,
+        genotype_matrix,
+        restricted_max_likelihood=True,
+        refit=False,
+        temp_data=None):
+    """Takes the phenotype vector and genotype matrix and returns a set of p-values and t-statistics
+    
+    restricted_max_likelihood -- whether to use restricted max likelihood; True or False
+    refit -- whether to refit the variance component for each marker
+    temp_data -- TempData object that stores the progress for each major step of the
+    calculations ("calculate_kinship" and "GWAS" take the majority of time)
+    
+    """
+    
+    with Bench("Calculate Kinship"):
+        kinship_matrix = calculate_kinship(genotype_matrix, temp_data)
+    
+    with Bench("Create LMM object"):
+        lmm_ob = LMM(pheno_vector, kinship_matrix)
+    
+    with Bench("LMM_ob fitting"):
+        lmm_ob.fit()
+
+    with Bench("Doing GWAS"):
+        t_stats, p_values = GWAS(pheno_vector,
+                                genotype_matrix,
+                                kinship_matrix,
+                                restricted_max_likelihood=True,
+                                refit=False,
+                                temp_data=temp_data)
+    Bench().report()
+    return t_stats, p_values
 
 
 def matrixMult(A,B):
@@ -99,8 +99,8 @@ def calculate_kinship(genotype_matrix, temp_data):
     """
     n = genotype_matrix.shape[0]
     m = genotype_matrix.shape[1]
-    print("n is:", n)
-    print("m is:", m)
+    #print("n is:", n)
+    #print("m is:", m)
     keep = []
     for counter in range(m):
         print("type of genotype_matrix[:,counter]:", pf(genotype_matrix[:,counter]))
@@ -123,7 +123,6 @@ def calculate_kinship(genotype_matrix, temp_data):
         genotype_matrix[:,counter] = (genotype_matrix[:,counter] - values_mean) / np.sqrt(vr)
         
         percent_complete = int(round((counter/m)*45))
-        print("Percent complete: ", percent_complete)
         temp_data.store("percent_complete", percent_complete)
         
     genotype_matrix = genotype_matrix[:,keep]
@@ -220,7 +219,6 @@ def GWAS(pheno_vector,
             ts, ps = lmm_ob.association(x, REML=restricted_max_likelihood)
             
         percent_complete = 45 + int(round((counter/m)*55))
-        print("Percent complete: ", percent_complete)
         temp_data.store("percent_complete", percent_complete)
 
         p_values.append(ps)
@@ -279,6 +277,8 @@ class LMM:
       self.K = K
       self.Kva = Kva
       self.Kve = Kve
+      print("self.Kva is: ", pf(self.Kva))
+      print("self.Kve is: ", pf(self.Kve))
       self.Y = Y
       self.X0 = X0
       self.N = self.K.shape[0]
diff --git a/wqflask/wqflask/search_results.py b/wqflask/wqflask/search_results.py
index 080568f0..2d39a711 100644
--- a/wqflask/wqflask/search_results.py
+++ b/wqflask/wqflask/search_results.py
@@ -45,7 +45,7 @@ class SearchResultPage():
         ###########################################
         #   Names and IDs of group / F2 set
         ###########################################
-        
+
         # All Phenotypes is a special case we'll deal with later
         #if kw['dataset'] == "All Phenotypes":
         #    self.cursor.execute("""
@@ -57,15 +57,16 @@ class SearchResultPage():
         #    self.dataset_groups = map(lambda x: x[1], results)
         #    self.dataset_group_ids = map(lambda x: x[2], results)
         #else:
-        
+
         self.results = []
-        
+
         if 'q' in kw:
-            self.quick_search = True
+            #self.quick_search = True
             self.search_terms = kw['q']
+            print("self.search_terms is: ", self.search_terms)
             self.quick_search()
         else:
-            self.quick_search = False
+            #self.quick_search = False
             self.search_terms = kw['search_terms']
             self.dataset = create_dataset(kw['dataset'])
             self.search()
@@ -98,10 +99,23 @@ class SearchResultPage():
             self.trait_list.append(this_trait)
 
         self.dataset.get_trait_info(self.trait_list, species)
-        
+
     def quick_search(self):
+        self.search_terms = parser.parse(self.search_terms)
+        print("After parsing:", self.search_terms)
 
-        return True
+        for a_search in self.search_terms:
+            search_term = a_search['search_term']
+            #Do mRNA assay search
+            search_ob = do_search.DoSearch.get_search("quick_mrna_assay")
+            search_class = getattr(do_search, search_ob)
+            the_search = search_class(search_term)
+            
+            self.results.extend(the_search.run())
+            print("in the search results are:", self.results)
+
+
+        #return True
 
         #search_gene
         #search_geno
@@ -114,8 +128,6 @@ class SearchResultPage():
         self.search_terms = parser.parse(self.search_terms)
         print("After parsing:", self.search_terms)
 
-        
-        
         for a_search in self.search_terms:
             print("[kodak] item is:", pf(a_search))
             search_term = a_search['search_term']
diff --git a/wqflask/wqflask/static/new/javascript/marker_regression.coffee b/wqflask/wqflask/static/new/javascript/marker_regression.coffee
index 4d678bde..6e605fa7 100644
--- a/wqflask/wqflask/static/new/javascript/marker_regression.coffee
+++ b/wqflask/wqflask/static/new/javascript/marker_regression.coffee
@@ -85,7 +85,7 @@ $ ->
                 marker_name = marker_info[2]
             else
                 marker_name = ""
-            $("#qtl_results_filter").find("input:first").val(marker_name).keyup()
+            $("#qtl_results_filter").find("input:first").val(marker_name).keypress()
 
         create_svg: () ->
             svg = d3.select("#manhattan_plots")
@@ -279,8 +279,6 @@ $ ->
                 )
                 .classed("circle", true)
                 .on("mouseover", (d) =>
-                    console.log("foodie")
-                    console.log("this:", this)
                     console.log("d3.event is:", d3.event)
                     console.log("d is:", d)
                     this_id = "point_" + String(d[2])
diff --git a/wqflask/wqflask/static/new/javascript/marker_regression.js b/wqflask/wqflask/static/new/javascript/marker_regression.js
index f5bb0c41..cb3c09cb 100644
--- a/wqflask/wqflask/static/new/javascript/marker_regression.js
+++ b/wqflask/wqflask/static/new/javascript/marker_regression.js
@@ -103,7 +103,7 @@
         } else {
           marker_name = "";
         }
-        return $("#qtl_results_filter").find("input:first").val(marker_name).keyup();
+        return $("#qtl_results_filter").find("input:first").val(marker_name).keypress();
       };
 
       Manhattan_Plot.prototype.create_svg = function() {
@@ -252,8 +252,6 @@
           return "point_" + String(d[2]);
         }).classed("circle", true).on("mouseover", function(d) {
           var this_id;
-          console.log("foodie");
-          console.log("this:", _this);
           console.log("d3.event is:", d3.event);
           console.log("d is:", d);
           this_id = "point_" + String(d[2]);
diff --git a/wqflask/wqflask/templates/index_page.html b/wqflask/wqflask/templates/index_page.html
index 46a307b7..09172705 100644
--- a/wqflask/wqflask/templates/index_page.html
+++ b/wqflask/wqflask/templates/index_page.html
@@ -8,12 +8,10 @@
     <header class="jumbotron subhead" id="overview">
         <div class="container">
             <h1>GeneNetwork</h1>
-            <p class="lead">Open source bioinformatics for systems genetics<br />
-            Brought to you by the University of Tennessee</p>
+            <p class="lead">Open source bioinformatics for systems genetics</p>
         </div>
     </header>
 
-
     <div class="container">
         <div class="row">
             <div class="span3 bs-docs-sidebar">