Changed code parsing marker csv file to increase GEMMA speed

Currently commenting out the vector map/table for GEMMA since they cause the page to hang, I'll try to figure out a way to prevent this

3 files changed, 69 insertions, 26 deletions

diff --git a/wqflask/wqflask/marker_regression/marker_regression.py b/wqflask/wqflask/marker_regression/marker_regression.py
index 37f9351d..b75a30f6 100755
--- a/wqflask/wqflask/marker_regression/marker_regression.py
+++ b/wqflask/wqflask/marker_regression/marker_regression.py
@@ -18,6 +18,7 @@ import numpy as np
 from scipy import linalg
 
 import cPickle as pickle
+import itertools
 
 import simplejson as json
 
@@ -90,10 +91,14 @@ class MarkerRegression(object):
         self.dataset.group.get_markers()
         if self.mapping_method == "gemma":
             self.score_type = "LOD"
-            included_markers, p_values = gemma_mapping.run_gemma(self.dataset, self.samples, self.vals)
-            self.dataset.group.get_specified_markers(markers = included_markers)
-            self.dataset.group.markers.add_pvalues(p_values)
-            results = self.dataset.group.markers.markers
+            with Bench("Running GEMMA"):
+                included_markers, p_values = gemma_mapping.run_gemma(self.dataset, self.samples, self.vals)
+            with Bench("Getting markers from csv"):
+                marker_obs = get_markers_from_csv(included_markers, p_values, self.dataset.group.name)
+            results = marker_obs
+            #self.dataset.group.get_specified_markers(markers = included_markers)
+            #self.dataset.group.markers.add_pvalues(p_values)
+            #results = self.dataset.group.markers.markers
         elif self.mapping_method == "rqtl_plink":
             results = self.run_rqtl_plink()
         elif self.mapping_method == "rqtl_geno":
@@ -183,7 +188,7 @@ class MarkerRegression(object):
                 #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")
+                    #print("changing to X")
                     self.json_data['chr'].append("X")
                 else:
                     self.json_data['chr'].append(str(qtl['chr']))
@@ -1006,6 +1011,52 @@ def create_snp_iterator_file(group):
     with gzip.open(snp_file_base, "wb") as fh:
         pickle.dump(data, fh, pickle.HIGHEST_PROTOCOL)
 
+def get_markers_from_csv(included_markers, p_values, group_name):
+    marker_data_fh = open(os.path.join(webqtlConfig.PYLMM_PATH + group_name + '_markers.csv'))
+    markers = []
+    for marker_name, p_value in itertools.izip(included_markers, p_values):
+        if not p_value or len(included_markers) < 1: 
+            continue   
+        for line in marker_data_fh:
+            splat = line.strip().split()
+            if splat[0] == marker_name:
+                marker = {}
+                marker['name'] = splat[0]
+                marker['chr'] = int(splat[1])
+                marker['Mb'] = float(splat[2])
+                marker['p_value'] = 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'])
+                    marker['lrs_value'] = -math.log10(marker['p_value']) * 4.61
+                markers.append(marker)
+                break
+
+#    for line, p_value in itertools.izip(marker_data_fh, p_values):
+#        if not p_value or len(included_markers) < 1: 
+#            continue
+#        splat = line.strip().split()
+#        if splat[0] in included_markers:
+#            marker = {}
+#            marker['name'] = splat[0]
+#            marker['chr'] = int(splat[1])
+#            marker['Mb'] = float(splat[2])
+#            marker['p_value'] = 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'])
+#                marker['lrs_value'] = -math.log10(marker['p_value']) * 4.61
+#            markers.append(marker)
+#        else:
+#            continue
+
+    return markers
+        
+
 #if __name__ == '__main__':
 #    import cPickle as pickle
 #    import gzip
diff --git a/wqflask/wqflask/templates/marker_regression_gn1.html b/wqflask/wqflask/templates/marker_regression_gn1.html
index 4eb5862a..58478758 100644
--- a/wqflask/wqflask/templates/marker_regression_gn1.html
+++ b/wqflask/wqflask/templates/marker_regression_gn1.html
@@ -15,7 +15,6 @@
         <input type="hidden" name="trait_id" value="{{ this_trait.name }}">
         <input type="hidden" name="dataset" value="{{ dataset.name }}">
         <input type="hidden" name="method" value="{{ mapping_method }}">
-        {% if mapping_method != "gemma" %}
         {% for sample in dataset.group.samplelist %}
         <input type="hidden" name="value:{{ sample }}" value="{{ vals[loop.index - 1] }}">
         {% endfor %}
@@ -29,7 +28,6 @@
         <input type="hidden" name="mapmethod_rqtl_geno" value="{{ mapmethod_rqtl_geno }}">
         <input type="hidden" name="mapmodel_rqtl_geno" value="{{ mapmodel_rqtl_geno }}">
         <input type="hidden" name="pair_scan" value="{{ pair_scan }}">
-        {% endif %}
         
         <div class="container">
           <div class="col-xs-5">
@@ -42,7 +40,6 @@
               <b>Location:</b> Chr {{ this_trait.chr }} @ {{ this_trait.mb }} Mb
               {% endif %}
           </div>
-          {% if mapping_method != "gemma" %}
           <div id="gn1_map_options" class="col-xs-5" style="border:2px solid black; padding: 10px; margin: 10px;">
             <div class="col-xs-7" style="padding-left: 0px;">
               <table>
@@ -69,40 +66,37 @@
               </div>
           </div>
         </div>
-        {% endif %}
-        
  
         <div class="tabbable" style="margin: 10px;">
           <ul class="nav nav-tabs">
-            {% if mapping_method != "gemma" %}
             <li id="gn1_map_tab">
               <a href="#gn1_map" data-toggle="tab" aria-expanded="true">GN1 Map</a>
             </li>
-            {% endif %}
+            {% if mapping_method != "gemma" %}
             <li id="vector_map_tab">
               <a href="#vector_map" data-toggle="tab" {% if mapping_method != "gemma" %}aria-expanded="false"{% else %}aria-expanded="true"{% endif %}>Vector Map</a>
             </li>
+            {% endif %}
           </ul>
           <div class="tab-content">
-            {% if mapping_method != "gemma" %}
             <div class="tab-pane active" id="gn1_map">
               <div class="qtlcharts">
                   {{ gifmap|safe }}
                   <img src="/static/output/{{ filename }}.jpeg" usemap="#WebQTLImageMap">
               </div>
             </div>
-            {% endif %}
+            {% if mapping_method != "gemma" %}
             <div class="tab-pane {% if mapping_method == "gemma" %}active{% endif %}" id="vector_map">
               <div id="chart_container">
                 <div class="qtlcharts" id="topchart"></div>
               </div>
             </div>
+            {% endif %}
           </div>
         </div>
 
         </form>
-        
-        {% if mapping_method != "gemma" %} 
+        {% if mapping_method != "gemma" %}
         <div style="width:48%;">
             <h2>
                 Results
@@ -182,11 +176,13 @@
         js_data = {{ js_data | safe }}
     </script>
 
+    {% if mapping_method != "gemma" %}
     <script language="javascript" type="text/javascript" src="/static/new/javascript/panelutil.js"></script>
     <script language="javascript" type="text/javascript" src="/static/new/javascript/chr_lod_chart.js"></script> 
     <script language="javascript" type="text/javascript" src="/static/new/javascript/lod_chart.js"></script>
     <script language="javascript" type="text/javascript" src="/static/new/javascript/create_lodchart.js"></script>
-    
+    {% endif %}
+
     <script type="text/javascript" charset="utf-8">
         $(document).ready( function () {
             console.time("Creating table");
diff --git a/wqflask/wqflask/views.py b/wqflask/wqflask/views.py
index 17845826..2af03f96 100755
--- a/wqflask/wqflask/views.py
+++ b/wqflask/wqflask/views.py
@@ -380,25 +380,24 @@ def marker_regression_page():
             result = pickle.loads(result)
     else:
         print("Cache miss!!!")
-        template_vars = marker_regression.MarkerRegression(start_vars, temp_uuid)
+        with Bench("Total time in MarkerRegression"):
+            template_vars = marker_regression.MarkerRegression(start_vars, temp_uuid)
 
         template_vars.js_data = json.dumps(template_vars.js_data,
                                            default=json_default_handler,
                                            indent="   ")
 
         result = template_vars.__dict__
-        #print("initial result:", result['qtl_results'])
 
         #for item in template_vars.__dict__.keys():
         #    print("  ---**--- {}: {}".format(type(template_vars.__dict__[item]), item))
         
-        if start_vars['method'] != "gemma":
-            gn1_template_vars = marker_regression_gn1.MarkerRegression(result).__dict__
+        gn1_template_vars = marker_regression_gn1.MarkerRegression(result).__dict__
 
         #qtl_length = len(result['js_data']['qtl_results'])
         #print("qtl_length:", qtl_length)
         pickled_result = pickle.dumps(result, pickle.HIGHEST_PROTOCOL)
-        print("pickled result length:", len(pickled_result))
+        #print("pickled result length:", len(pickled_result))
         Redis.set(key, pickled_result)
         Redis.expire(key, 1*60)
 
@@ -416,10 +415,7 @@ def marker_regression_page():
             rendered_template = render_template("pair_scan_results.html", **result)
         else:
             #rendered_template = render_template("marker_regression.html", **result)
-            if start_vars['method'] != "gemma":
-                rendered_template = render_template("marker_regression_gn1.html", **gn1_template_vars)
-            else:
-                rendered_template = render_template("marker_regression_gn1.html", **result)
+            rendered_template = render_template("marker_regression_gn1.html", **gn1_template_vars)
 
     return rendered_template