about summary refs log tree commit diff
diff options
context:
space:
mode:
authorZachary Sloan2013-10-08 12:13:56 -0500
committerZachary Sloan2013-10-08 12:13:56 -0500
commit3325184b1dd310619626dd31852ab84cae6dc7fc (patch)
treebf0d3893e25c14094066611f5e7d075d8b0fd092
parent4c37e3962d6f34a935ea23511e5b760be4208474 (diff)
downloadgenenetwork2-3325184b1dd310619626dd31852ab84cae6dc7fc.tar.gz
Did some work with interval mapping page; will use qtl reaper to do
the mapping, since it is reliable/fast and avoids us having to rewrite
from scratch using something like r/qtl
-rwxr-xr-xwqflask/base/data_set.py20
-rw-r--r--wqflask/maintenance/quick_search_table.py4
-rw-r--r--wqflask/wqflask/interval_mapping/interval_mapping.py53
3 files changed, 67 insertions, 10 deletions
diff --git a/wqflask/base/data_set.py b/wqflask/base/data_set.py
index 96e04df0..befbd518 100755
--- a/wqflask/base/data_set.py
+++ b/wqflask/base/data_set.py
@@ -168,13 +168,13 @@ class Markers(object):
         
         for marker, p_value in itertools.izip(self.markers, p_values):
             marker['p_value'] = p_value
-            if marker['p_value'] == 0:
-                marker['lod_score'] = 0
-                marker['lrs_value'] = 0
-            else:
-                marker['lod_score'] = -math.log10(marker['p_value'])
-                #Using -log(p) for the LRS; need to ask Rob how he wants to get LRS from p-values
-                marker['lrs_value'] = -math.log10(marker['p_value']) * 4.61
+            if math.isnan(marker['p_value']):
+                print("p_value is:", marker['p_value'])
+            marker['lod_score'] = -math.log10(marker['p_value'])
+            #Using -log(p) for the LRS; need to ask Rob how he wants to get LRS from p-values
+            marker['lrs_value'] = -math.log10(marker['p_value']) * 4.61
+        
+        
 
 
 class HumanMarkers(Markers):
@@ -189,6 +189,8 @@ class HumanMarkers(Markers):
             marker['name'] = splat[1]
             marker['Mb'] = float(splat[3]) / 1000000
             self.markers.append(marker)
+            
+        #print("markers is: ", pf(self.markers))
 
 
     def add_pvalues(self, p_values):
@@ -315,10 +317,10 @@ class DatasetGroup(object):
 
         #determine default genotype object
         if self.incparentsf1 and genotype_1.type != "intercross":
-            genotype = genotype_2
+            self.genotype = genotype_2
         else:
             self.incparentsf1 = 0
-            genotype = genotype_1
+            self.genotype = genotype_1
 
         self.samplelist = list(genotype.prgy)
 
diff --git a/wqflask/maintenance/quick_search_table.py b/wqflask/maintenance/quick_search_table.py
index 9cd792ef..23bd505c 100644
--- a/wqflask/maintenance/quick_search_table.py
+++ b/wqflask/maintenance/quick_search_table.py
@@ -11,8 +11,10 @@ each trait, its dataset, and several columns determined by its trait type (pheno
 
 """
 
-from __future__ import print_function, division, absolute_import
+from __future__ import absolute_import, division, print_function
 
+# We do this here so we can use zach_settings
+# Not to avoid other absoulte_imports
 import sys
 sys.path.append("../../..")
 
diff --git a/wqflask/wqflask/interval_mapping/interval_mapping.py b/wqflask/wqflask/interval_mapping/interval_mapping.py
index 48e8018e..5d660224 100644
--- a/wqflask/wqflask/interval_mapping/interval_mapping.py
+++ b/wqflask/wqflask/interval_mapping/interval_mapping.py
@@ -12,6 +12,7 @@ import collections
 
 import numpy as np
 from scipy import linalg
+import rpy2.robjects
 
 import simplejson as json
 
@@ -83,6 +84,28 @@ class IntervalMapping(object):
         """Generates qtl results for plotting interval map"""
 
         self.dataset.group.get_markers()
+        self.dataset.read_genotype_file()
+
+        samples, values, variances = self.trait.export_informative()
+        if self.control_locus:
+            if self.weighted_regression:
+                qtl_result = self.dataset.genotype.regression(strains = samples,
+                                                              trait = values,
+                                                              variance = variances,
+                                                              control = self.control_locus)
+            else:
+                qtl_result = self.dataset.genotype.regression(strains = samples,
+                                                              trait = values,
+                                                              control = self.control_locus)
+        else:
+            if self.weighted_regression:
+                qtl_result = self.dataset.genotype.regression(strains = samples,
+                                                              trait = values,
+                                                              variance = variances)
+            else:
+                qtl_result = self.dataset.genotype.regression(strains = samples,
+                                                              trait = values)
+
 
         pheno_vector = np.array([val == "x" and np.nan or float(val) for val in self.vals])
 
@@ -108,6 +131,36 @@ class IntervalMapping(object):
         
         self.qtl_results = self.dataset.group.markers.markers
 
+    #def gen_qtl_results_2(self, tempdata):
+    #    """Generates qtl results for plotting interval map"""
+    #
+    #    self.dataset.group.get_markers()
+    #    self.dataset.read_genotype_file()
+    #
+    #    pheno_vector = np.array([val == "x" and np.nan or float(val) for val in self.vals])
+    #
+    #    #if self.dataset.group.species == "human":
+    #    #    p_values, t_stats = self.gen_human_results(pheno_vector, tempdata)
+    #    #else:
+    #    genotype_data = [marker['genotypes'] for marker in self.dataset.group.markers.markers]
+    #    
+    #    no_val_samples = self.identify_empty_samples()
+    #    trimmed_genotype_data = self.trim_genotypes(genotype_data, no_val_samples)
+    #    
+    #    genotype_matrix = np.array(trimmed_genotype_data).T
+    #    
+    #    t_stats, p_values = lmm.run(
+    #        pheno_vector,
+    #        genotype_matrix,
+    #        restricted_max_likelihood=True,
+    #        refit=False,
+    #        temp_data=tempdata
+    #    )
+    #    
+    #    self.dataset.group.markers.add_pvalues(p_values)
+    #    
+    #    self.qtl_results = self.dataset.group.markers.markers
+
 
     def identify_empty_samples(self):
         no_val_samples = []