Zooming into chromosome or remapping with GEMMA is now very fast, since it stores the output file

Added option to search for terms in quotes, searching for all the contents within the quotes

Added link to R/qtl documentation when user is in the R/qtl mapping tab

8 files changed, 147 insertions, 128 deletions

diff --git a/wqflask/wqflask/marker_regression/display_mapping_results.py b/wqflask/wqflask/marker_regression/display_mapping_results.py
index bc8f8fd9..41d5c9d0 100644
--- a/wqflask/wqflask/marker_regression/display_mapping_results.py
+++ b/wqflask/wqflask/marker_regression/display_mapping_results.py
@@ -154,6 +154,10 @@ class DisplayMappingResults(object):
 
         self.geno_db_exists = start_vars['geno_db_exists']
 
+        self.first_run = True
+        if 'first_run' in start_vars:
+            self.first_run = start_vars['first_run']
+
         #Needing for form submission when doing single chr mapping or remapping after changing options
         self.samples = start_vars['samples']
         self.vals = start_vars['vals']
@@ -226,6 +230,8 @@ class DisplayMappingResults(object):
             self.maf = start_vars['maf']
         if 'use_loco' in start_vars.keys():
             self.use_loco = start_vars['use_loco']
+            if self.use_loco == "True":
+                self.gwa_filename = start_vars['gwa_filename']
 
         self.selectedChr = int(start_vars['selected_chr'])
 
diff --git a/wqflask/wqflask/marker_regression/gemma_mapping.py b/wqflask/wqflask/marker_regression/gemma_mapping.py
index 9e92838e..0f37e711 100644
--- a/wqflask/wqflask/marker_regression/gemma_mapping.py
+++ b/wqflask/wqflask/marker_regression/gemma_mapping.py
@@ -12,7 +12,7 @@ GEMMAOPTS = "-debug"
 if WEBSERVER_MODE == 'PROD':
   GEMMAOPTS = "-no-check"
 
-def run_gemma(this_trait, this_dataset, samples, vals, covariates, use_loco, maf=0.01):
+def run_gemma(this_trait, this_dataset, samples, vals, covariates, use_loco, maf=0.01, first_run=True, gwa_output_filename=None):
     """Generates p-values for each marker using GEMMA"""
 
     if this_dataset.group.genofile != None:
@@ -20,127 +20,114 @@ def run_gemma(this_trait, this_dataset, samples, vals, covariates, use_loco, maf
     else:
         genofile_name = this_dataset.group.name
 
-    trait_filename = str(this_trait.name) + "_" + str(this_dataset.name) + "_pheno"
-    gen_pheno_txt_file(this_dataset, genofile_name, vals, trait_filename)
-
-    if not os.path.isfile("{}{}_output.assoc.txt".format(webqtlConfig.GENERATED_IMAGE_DIR, genofile_name)):
-        open("{}{}_output.assoc.txt".format(webqtlConfig.GENERATED_IMAGE_DIR, genofile_name), "w+")
-
-    this_chromosomes = this_dataset.species.chromosomes.chromosomes
-    chr_list_string = ""
-    for i in range(len(this_chromosomes)):
-        if i < (len(this_chromosomes) - 1):
-            chr_list_string += this_chromosomes[i+1].name + ","
-        else:
-            chr_list_string += this_chromosomes[i+1].name
-
-    if covariates != "":
-        gen_covariates_file(this_dataset, covariates)
-
-    k_output_filename = this_dataset.group.name + "_K_" + ''.join(random.choice(string.ascii_uppercase + string.digits) for _ in range(6))
-    gwa_output_filename = this_dataset.group.name + "_GWA_" + ''.join(random.choice(string.ascii_uppercase + string.digits) for _ in range(6))
-    if use_loco == "True":
-        generate_k_command = GEMMA_WRAPPER_COMMAND + ' --json --loco ' + chr_list_string + ' -- ' + GEMMAOPTS + ' -g %s/%s_geno.txt -p %s/gn2/%s.txt -a %s/%s_snps.txt -gk > %s/gn2/%s.json' % (flat_files('genotype/bimbam'),
-                                                                                        genofile_name,
-                                                                                        TEMPDIR,
-                                                                                        trait_filename,
-                                                                                        flat_files('genotype/bimbam'),
-                                                                                        genofile_name,
-                                                                                        TEMPDIR,
-                                                                                        k_output_filename)
-        logger.debug("k_command:" + generate_k_command)
-        os.system(generate_k_command)
-
-        gemma_command = GEMMA_WRAPPER_COMMAND + ' --json --loco --input %s/gn2/%s.json -- ' % (TEMPDIR, k_output_filename) + GEMMAOPTS + ' -g %s/%s_geno.txt -p %s/gn2/%s.txt' % (flat_files('genotype/bimbam'),
-                                                                                        genofile_name,
-                                                                                        TEMPDIR,
-                                                                                        trait_filename)
-        if covariates != "":
-            gemma_command += ' -c %s/%s_covariates.txt -a %s/%s_snps.txt -lmm 2 -maf %s > %s/gn2/%s.json' % (flat_files('mapping'),
-                                                                                                              this_dataset.group.name,
-                                                                                                              flat_files('genotype/bimbam'),
-                                                                                                              genofile_name,
-                                                                                                              maf,
-                                                                                                              TEMPDIR,
-                                                                                                              gwa_output_filename)
-        else:
-            gemma_command += ' -a %s/%s_snps.txt -lmm 2 -maf %s > %s/gn2/%s.json' % (flat_files('genotype/bimbam'),
-                                                                                                             genofile_name,
-                                                                                                             maf,
-                                                                                                             TEMPDIR,
-                                                                                                             gwa_output_filename)
-
-    else:
-        generate_k_command = GEMMA_COMMAND + ' ' + GEMMAOPTS + ' -g %s/%s_geno.txt -p %s/gn2/%s.txt -a %s/%s_snps.txt -gk -outdir %s/gn2/ -o %s' % (flat_files('genotype/bimbam'),
-                                                                                        genofile_name,
-                                                                                        TEMPDIR,
-                                                                                        trait_filename,
-                                                                                        flat_files('genotype/bimbam'),
-                                                                                        genofile_name,
-                                                                                        TEMPDIR,
-                                                                                        k_output_filename)
-        #generate_k_command = GEMMA_WRAPPER_COMMAND + ' --json -- ' + GEMMAOPTS + ' -g %s/%s_geno.txt -p %s/%s.txt -a %s/%s_snps.txt -gk > %s/gn2/%s.json' % (flat_files('genotype/bimbam'),
-        #                                                                                genofile_name,
-        #                                                                                flat_files('genotype/bimbam'),
-        #                                                                                trait_filename,
-        #                                                                                flat_files('genotype/bimbam'),
-        #                                                                                genofile_name,
-        #                                                                                TEMPDIR,
-        #                                                                                k_output_filename)
-
-        logger.debug("k_command:" + generate_k_command)
-        os.system(generate_k_command)
-
-        gemma_command = GEMMA_COMMAND + ' ' + GEMMAOPTS + ' -g %s/%s_geno.txt -p %s/gn2/%s.txt -a %s/%s_snps.txt -k %s/gn2/%s.cXX.txt -lmm 2 -maf %s' % (flat_files('genotype/bimbam'),
-                                                                                        genofile_name,
-                                                                                        TEMPDIR,
-                                                                                        trait_filename,
-                                                                                        flat_files('genotype/bimbam'),
-                                                                                        genofile_name,
-                                                                                        TEMPDIR,
-                                                                                        k_output_filename,
-                                                                                        maf)
-
-        #gemma_command = GEMMA_WRAPPER_COMMAND + ' --json --input %s/gn2/%s.json -- ' % (TEMPDIR, k_output_filename) + GEMMAOPTS + ' -g %s/%s_geno.txt -p %s/%s_pheno.txt' % (flat_files('genotype/bimbam'),
-        #                                                                                genofile_name,
-        #                                                                                flat_files('genotype/bimbam'),
-        #                                                                                genofile_name)
-
-
-        if covariates != "":
-            gemma_command += ' -c %s/%s_covariates.txt -outdir %s -o %s_output' % (flat_files('mapping'),
-                                                                                                         this_dataset.group.name,
-                                                                                                         webqtlConfig.GENERATED_IMAGE_DIR,
-                                                                                                         genofile_name)
-        else:
-            gemma_command += ' -outdir %s -o %s_output' % (webqtlConfig.GENERATED_IMAGE_DIR,
-                                                                  genofile_name)
-
-        #if covariates != "":
-        #    gemma_command += ' -c %s/%s_covariates.txt -a %s/%s_snps.txt -lmm 2 -maf %s > %s/gn2/%s.json' % (flat_files('mapping'),
-        #                                                                                                     this_dataset.group.name,
-        #                                                                                                     flat_files('genotype/bimbam'),
-        #                                                                                                     genofile_name,
-        #                                                                                                     maf,
-        #                                                                                                     TEMPDIR,
-        #                                                                                                     gwa_output_filename)
-        #else:
-        #    gemma_command += ' -a %s/%s_snps.txt -lmm 2 -maf %s > %s/gn2/%s.json' % (flat_files('genotype/bimbam'),
-        #                                                                             genofile_name,
-        #                                                                             maf,
-        #                                                                             TEMPDIR,
-        #                                                                             gwa_output_filename)
-
-    logger.debug("gemma_command:" + gemma_command)
-    os.system(gemma_command)
+    if first_run:
+      trait_filename = str(this_trait.name) + "_" + str(this_dataset.name) + "_pheno"
+      gen_pheno_txt_file(this_dataset, genofile_name, vals, trait_filename)
+
+      if not os.path.isfile("{}{}_output.assoc.txt".format(webqtlConfig.GENERATED_IMAGE_DIR, genofile_name)):
+          open("{}{}_output.assoc.txt".format(webqtlConfig.GENERATED_IMAGE_DIR, genofile_name), "w+")
+
+      k_output_filename = this_dataset.group.name + "_K_" + ''.join(random.choice(string.ascii_uppercase + string.digits) for _ in range(6))
+      gwa_output_filename = this_dataset.group.name + "_GWA_" + ''.join(random.choice(string.ascii_uppercase + string.digits) for _ in range(6))
+
+      this_chromosomes = this_dataset.species.chromosomes.chromosomes
+      chr_list_string = ""
+      for i in range(len(this_chromosomes)):
+          if i < (len(this_chromosomes) - 1):
+              chr_list_string += this_chromosomes[i+1].name + ","
+          else:
+              chr_list_string += this_chromosomes[i+1].name
+
+      if covariates != "":
+          gen_covariates_file(this_dataset, covariates)
+
+      if use_loco == "True":
+          generate_k_command = GEMMA_WRAPPER_COMMAND + ' --json --loco ' + chr_list_string + ' -- ' + GEMMAOPTS + ' -g %s/%s_geno.txt -p %s/gn2/%s.txt -a %s/%s_snps.txt -gk > %s/gn2/%s.json' % (flat_files('genotype/bimbam'),
+                                                                                          genofile_name,
+                                                                                          TEMPDIR,
+                                                                                          trait_filename,
+                                                                                          flat_files('genotype/bimbam'),
+                                                                                          genofile_name,
+                                                                                          TEMPDIR,
+                                                                                          k_output_filename)
+          logger.debug("k_command:" + generate_k_command)
+          os.system(generate_k_command)
+
+          gemma_command = GEMMA_WRAPPER_COMMAND + ' --json --loco --input %s/gn2/%s.json -- ' % (TEMPDIR, k_output_filename) + GEMMAOPTS + ' -g %s/%s_geno.txt -p %s/gn2/%s.txt' % (flat_files('genotype/bimbam'),
+                                                                                          genofile_name,
+                                                                                          TEMPDIR,
+                                                                                          trait_filename)
+          if covariates != "":
+              gemma_command += ' -c %s/%s_covariates.txt -a %s/%s_snps.txt -lmm 2 -maf %s > %s/gn2/%s.json' % (flat_files('mapping'),
+                                                                                                                this_dataset.group.name,
+                                                                                                                flat_files('genotype/bimbam'),
+                                                                                                                genofile_name,
+                                                                                                                maf,
+                                                                                                                TEMPDIR,
+                                                                                                                gwa_output_filename)
+          else:
+              gemma_command += ' -a %s/%s_snps.txt -lmm 2 -maf %s > %s/gn2/%s.json' % (flat_files('genotype/bimbam'),
+                                                                                                               genofile_name,
+                                                                                                               maf,
+                                                                                                               TEMPDIR,
+                                                                                                               gwa_output_filename)
+
+      else:
+          generate_k_command = GEMMA_COMMAND + ' ' + GEMMAOPTS + ' -g %s/%s_geno.txt -p %s/gn2/%s.txt -a %s/%s_snps.txt -gk -outdir %s/gn2/ -o %s' % (flat_files('genotype/bimbam'),
+                                                                                          genofile_name,
+                                                                                          TEMPDIR,
+                                                                                          trait_filename,
+                                                                                          flat_files('genotype/bimbam'),
+                                                                                          genofile_name,
+                                                                                          TEMPDIR,
+                                                                                          k_output_filename)
+          #generate_k_command = GEMMA_WRAPPER_COMMAND + ' --json -- ' + GEMMAOPTS + ' -g %s/%s_geno.txt -p %s/%s.txt -a %s/%s_snps.txt -gk > %s/gn2/%s.json' % (flat_files('genotype/bimbam'),
+          #                                                                                genofile_name,
+          #                                                                                flat_files('genotype/bimbam'),
+          #                                                                                trait_filename,
+          #                                                                                flat_files('genotype/bimbam'),
+          #                                                                                genofile_name,
+          #                                                                                TEMPDIR,
+          #                                                                                k_output_filename)
+
+          logger.debug("k_command:" + generate_k_command)
+          os.system(generate_k_command)
+
+          gemma_command = GEMMA_COMMAND + ' ' + GEMMAOPTS + ' -g %s/%s_geno.txt -p %s/gn2/%s.txt -a %s/%s_snps.txt -k %s/gn2/%s.cXX.txt -lmm 2 -maf %s' % (flat_files('genotype/bimbam'),
+                                                                                          genofile_name,
+                                                                                          TEMPDIR,
+                                                                                          trait_filename,
+                                                                                          flat_files('genotype/bimbam'),
+                                                                                          genofile_name,
+                                                                                          TEMPDIR,
+                                                                                          k_output_filename,
+                                                                                          maf)
+
+          #gemma_command = GEMMA_WRAPPER_COMMAND + ' --json --input %s/gn2/%s.json -- ' % (TEMPDIR, k_output_filename) + GEMMAOPTS + ' -g %s/%s_geno.txt -p %s/%s_pheno.txt' % (flat_files('genotype/bimbam'),
+          #                                                                                genofile_name,
+          #                                                                                flat_files('genotype/bimbam'),
+          #                                                                                genofile_name)
+
+
+          if covariates != "":
+              gemma_command += ' -c %s/%s_covariates.txt -outdir %s -o %s_output' % (flat_files('mapping'),
+                                                                                                           this_dataset.group.name,
+                                                                                                           webqtlConfig.GENERATED_IMAGE_DIR,
+                                                                                                           genofile_name)
+          else:
+              gemma_command += ' -outdir %s -o %s_output' % (webqtlConfig.GENERATED_IMAGE_DIR,
+                                                                    genofile_name)
+
+
+      logger.debug("gemma_command:" + gemma_command)
+      os.system(gemma_command)
 
     if use_loco == "True":
         marker_obs = parse_loco_output(this_dataset, gwa_output_filename)
+        return marker_obs, gwa_output_filename
     else:
-        #marker_obs = parse_loco_output(this_dataset, gwa_output_filename)
         marker_obs = parse_gemma_output(genofile_name)
-
-    return marker_obs
+        return marker_obs
 
 def gen_pheno_txt_file(this_dataset, genofile_name, vals, trait_filename):
     """Generates phenotype file for GEMMA"""
diff --git a/wqflask/wqflask/marker_regression/run_mapping.py b/wqflask/wqflask/marker_regression/run_mapping.py
index c0d6dc0b..a0466c0f 100644
--- a/wqflask/wqflask/marker_regression/run_mapping.py
+++ b/wqflask/wqflask/marker_regression/run_mapping.py
@@ -175,10 +175,19 @@ class RunMapping(object):
             self.dataset.group.genofile = self.genofile_string.split(":")[0]
         self.dataset.group.get_markers()
         if self.mapping_method == "gemma":
+            self.first_run = True
+            self.gwa_filename = None
+            if 'first_run' in start_vars: #ZS: check if first run so existing result files can be used if it isn't (for example zooming on a chromosome, etc)
+                self.first_run = False
+                if 'gwa_filename' in start_vars:
+                    self.gwa_filename = start_vars['gwa_filename']
             self.score_type = "-log(p)"
             self.manhattan_plot = True
             with Bench("Running GEMMA"):
-                marker_obs = gemma_mapping.run_gemma(self.this_trait, self.dataset, self.samples, self.vals, self.covariates, self.use_loco, self.maf)
+                if self.use_loco == "True":
+                    marker_obs, self.gwa_filename = gemma_mapping.run_gemma(self.this_trait, self.dataset, self.samples, self.vals, self.covariates, self.use_loco, self.maf, self.first_run, self.gwa_filename)
+                else:
+                    marker_obs = gemma_mapping.run_gemma(self.this_trait, self.dataset, self.samples, self.vals, self.covariates, self.use_loco, self.maf, self.first_run)
             results = marker_obs
         elif self.mapping_method == "rqtl_plink":
             results = self.run_rqtl_plink()
diff --git a/wqflask/wqflask/parser.py b/wqflask/wqflask/parser.py
index ea7ab1b9..1ca5ecff 100644
--- a/wqflask/wqflask/parser.py
+++ b/wqflask/wqflask/parser.py
@@ -33,20 +33,18 @@ def parse(pstring):
     """
     pstring = re.split(r"""(?:(\w+\s*=\s*[\('"\[][^)'"]*[\)\]'"])  |  # LRS=(1 2 3), cisLRS=[4 5 6], etc
                        (\w+\s*[=:\>\<][\w\*]+)  |  # wiki=bar, GO:foobar, etc
+                       (".*?") | ('.*?') | # terms in quotes, i.e. "brain weight"
                        ([\w\*]+))  # shh, brain, etc """, pstring,
                                                     flags=re.VERBOSE)
+
     pstring = [item.strip() for item in pstring if item and item.strip()]
-    print(pstring)
+    logger.debug("pstring:", pstring)
 
     items = []
 
     separators = [re.escape(x) for x in ("<=", ">=", ":", "=", "<", ">")]
     separators = '(%s)' % ("|".join(separators))
 
-    logger.debug("separators:", separators)
-
-
-
     for item in pstring:
         splat = re.split(separators, item)
         logger.debug("splat is:", splat)
@@ -70,6 +68,8 @@ def parse(pstring):
                         separator=separator,
                         search_term=value)
         else:
+            if (item[0] == "\"" and item[-1] == "\"") or (item[0] == "'" and item[-1] == "'"):
+                item = item[1:-1]
             term = dict(key=None,
                         separator=None,
                         search_term=[item])
diff --git a/wqflask/wqflask/static/new/javascript/show_trait_mapping_tools.js b/wqflask/wqflask/static/new/javascript/show_trait_mapping_tools.js
index 28af232d..95d3a537 100644
--- a/wqflask/wqflask/static/new/javascript/show_trait_mapping_tools.js
+++ b/wqflask/wqflask/static/new/javascript/show_trait_mapping_tools.js
@@ -267,6 +267,14 @@
     }
   });
 
+  $("li.mapping_tab").click(function() {
+    if ($(this).hasClass("rqtl")){
+      $(".rqtl_description").css("display", "block");
+    } else {
+      $(".rqtl_description").css("display", "none");
+    }
+  });
+
   toggle_enable_disable = function(elem) {
     return $(elem).prop("disabled", !$(elem).prop("disabled"));
   };
diff --git a/wqflask/wqflask/templates/mapping_results.html b/wqflask/wqflask/templates/mapping_results.html
index 72d9d0dd..317747c1 100644
--- a/wqflask/wqflask/templates/mapping_results.html
+++ b/wqflask/wqflask/templates/mapping_results.html
@@ -16,6 +16,10 @@
         <input type="hidden" name="dataset" value="{{ dataset.name }}">
         <input type="hidden" name="genofile" value="{{ genofile_string }}">
         <input type="hidden" name="geno_db_exists" value="{{ geno_db_exists }}">
+        <input type="hidden" name="first_run" value="{{ first_run }}">
+        {% if gwa_filename is defined %}
+        <input type="hidden" name="gwa_filename" value="{{ gwa_filename }}">
+        {% endif %}
         <input type="hidden" name="results_path" value="{{ mapping_results_path }}">
         <input type="hidden" name="method" value="{{ mapping_method }}">
         {% for sample in samples %}
diff --git a/wqflask/wqflask/templates/show_trait_mapping_tools.html b/wqflask/wqflask/templates/show_trait_mapping_tools.html
index b5ed1c8f..b5c37b8b 100644
--- a/wqflask/wqflask/templates/show_trait_mapping_tools.html
+++ b/wqflask/wqflask/templates/show_trait_mapping_tools.html
@@ -5,23 +5,23 @@
 
             <ul class="nav nav-pills">
                 {% if dataset.group.mapping_id == "1" %}
-                <li class="active">
+                <li class="mapping_tab active">
                     <a href="#gemma" data-toggle="tab">GEMMA</a>
                 </li>
-                <li>
+                <li class="mapping_tab">
                     <a href="#interval_mapping" data-toggle="tab">Interval Mapping</a>
                 </li>
-                <li>
+                <li class="mapping_tab rqtl">
                     <a href="#rqtl_geno" data-toggle="tab">R/qtl</a>
                 </li>
                 {% endif %}
                 {% for mapping_method in dataset.group.mapping_names %}
                 {% if mapping_method == "GEMMA" %}
-                <li class="active">
+                <li class="mapping_tab active">
                     <a href="#gemma_plink" data-toggle="tab">GEMMA</a>
                 </li>
                 {% elif mapping_method == "PLINK" %}
-                <li>
+                <li class="mapping_tab">
                     <a href="#plink" data-toggle="tab">PLINK</a>
                 </li>
                 {% endif %}
@@ -353,6 +353,9 @@
              <dd>R/qtl is an extensible, interactive environment for mapping quantitative trait loci (QTL) in experimental crosses.</dd>
              {% endif %}
         </dl>
+        <div class="rqtl_description" style="padding-top: 40px; display: none;">
+        More information on R/qtl mapping models and methods can be found <a href="http://www.rqtl.org/tutorials/rqtltour.pdf">here</a>.
+        </div>
     </div>
     <div id="mapping_result_holder_wrapper" style="display:none;">
         <div id="mapping_result_holder"></div>
diff --git a/wqflask/wqflask/views.py b/wqflask/wqflask/views.py
index d00f067e..bff3aa2a 100644
--- a/wqflask/wqflask/views.py
+++ b/wqflask/wqflask/views.py
@@ -575,6 +575,8 @@ def mapping_results_page():
         'group',
         'species',
         'vals',
+        'first_run',
+        'gwa_filename',
         'geno_db_exists',
         'method',
         'mapping_results_path',