about summary refs log tree commit diff
diff options
context:
space:
mode:
-rw-r--r--misc/notes_DA.txt10
-rwxr-xr-xwqflask/base/data_set.py2
-rwxr-xr-xwqflask/wqflask/templates/collections/view.html12
-rw-r--r--wqflask/wqflask/templates/wgcna_results.html23
-rw-r--r--wqflask/wqflask/templates/wgcna_setup.html23
-rwxr-xr-xwqflask/wqflask/views.py21
-rwxr-xr-xwqflask/wqflask/wgcna/__init__.py0
-rw-r--r--wqflask/wqflask/wgcna/wgcna_analysis.py134
8 files changed, 222 insertions, 3 deletions
diff --git a/misc/notes_DA.txt b/misc/notes_DA.txt
new file mode 100644
index 00000000..410e0182
--- /dev/null
+++ b/misc/notes_DA.txt
@@ -0,0 +1,10 @@
+Danny's notes about the genenetwork source
+
+Location of static files:
+
+Location of HTML templates:       wqflask/wqflask/templates/
+
+Entry point of the wqflask app:   wqflask/wqflask/__init__.py
+
+Application routes:               wqflask/wqflask/views.py
+
diff --git a/wqflask/base/data_set.py b/wqflask/base/data_set.py
index e2db9ad7..427fd991 100755
--- a/wqflask/base/data_set.py
+++ b/wqflask/base/data_set.py
@@ -89,7 +89,7 @@ class Dataset_Types(object):
             for group in data['datasets'][species]:
                 for dataset_type in data['datasets'][species][group]:
                     for dataset in data['datasets'][species][group][dataset_type]:
-                        print("dataset is:", dataset)
+                        #print("dataset is:", dataset)
                         
                         short_dataset_name = dataset[0]
                         if dataset_type == "Phenotypes":
diff --git a/wqflask/wqflask/templates/collections/view.html b/wqflask/wqflask/templates/collections/view.html
index 29c65058..5450f903 100755
--- a/wqflask/wqflask/templates/collections/view.html
+++ b/wqflask/wqflask/templates/collections/view.html
@@ -38,6 +38,18 @@
                     <input type="submit" class="btn btn-primary" value="Correlation Matrix" />
                 </div>
             </form>
+            <form action="/wgcna_setup" method="post">
+                {% if uc %}
+                <input type="hidden" name="uc_id" id="uc_id" value="{{ uc.id }}" />
+                {% endif %}
+                <input type="hidden" name="trait_list" id="trait_list" value= "
+                {% for this_trait in trait_obs %}
+                    {{ this_trait.name }}:{{ this_trait.dataset.name }},
+                {% endfor %}" >
+                <div class="col-xs-2 controls">
+                    <input type="submit" class="btn btn-primary" value="WGCNA Analysis" />
+                </div>
+            </form>
             <form action="/heatmap" method="post">
                 {% if uc %}
                 <input type="hidden" name="uc_id" id="uc_id" value="{{ uc.id }}" />
diff --git a/wqflask/wqflask/templates/wgcna_results.html b/wqflask/wqflask/templates/wgcna_results.html
new file mode 100644
index 00000000..b8a4a4f5
--- /dev/null
+++ b/wqflask/wqflask/templates/wgcna_results.html
@@ -0,0 +1,23 @@
+{% extends "base.html" %}
+{% block title %}WCGNA results{% endblock %}
+
+{% block content %} <!-- Start of body -->
+  <div class="container">
+    <h1>WGCNA Results</h1>
+    Phenotype / Module table:
+    <table>
+      <tr><th>Phenotype</th><th>Module</th><th>Weight</th></tr>
+      {% for k in results['network'] %} 
+        <tr><td>{{k}}</td><td>&nbsp;</td><td>&nbsp;</td></tr>
+      {% endfor %}
+    </table>
+    <h2>WGCNA module plot</h2>
+                <a href="/tmp/{{ results['imgurl'] }}">
+                    <img alt="Embedded Image" src="data:image/png;base64,
+                    {% for elem in results['imgdata'] -%}
+                    {% print("%c"|format(elem)) %}
+                    {%- endfor %}
+                    " /></a>
+  </div>
+{% endblock %}
+
diff --git a/wqflask/wqflask/templates/wgcna_setup.html b/wqflask/wqflask/templates/wgcna_setup.html
new file mode 100644
index 00000000..821e4954
--- /dev/null
+++ b/wqflask/wqflask/templates/wgcna_setup.html
@@ -0,0 +1,23 @@
+{% extends "base.html" %}
+{% block title %}WCGNA analysis{% endblock %}
+
+{% block content %} <!-- Start of body -->
+<form action="/wgcna_results" method="post">
+  <div class="container">
+
+    <h1>WGCNA analysis parameters</h1>
+    Phenotypes used as input: {{request.form}}
+    <input type="hidden" name="trait_list" id="trait_list" value= "{{request.form['trait_list']}}">
+
+    <table>
+      <tr><td>Soft threshold:</td><td><input type="text" class="form-inline" id="SoftThreshold", value="1,2,3,4,5,6,7,8,9"></td></tr>
+      <tr><td>Power:</td><td><input type="text" class="form-inline" id="Power", value="6"></td></tr>
+      <tr><td>Minimum module size:</td><td><input type="text" class="form-inline" id="MinModuleSize", value="30"></td></tr>
+      <tr><td>TOMtype:</td><td><input type="text" class="form-inline" id="TOMtype", value="unsigned"></td></tr>
+      <tr><td>mergeCutHeight:</td><td><input type="text" class="form-inline" id="mergeCutHeight", value="0.25"></td></tr>
+    </table>
+    <input type="submit" class="btn btn-primary" value="Run WGCNA using these settings" />
+  </div>
+</form>
+{% endblock %}
+
diff --git a/wqflask/wqflask/views.py b/wqflask/wqflask/views.py
index 7a8a0f03..dc199ab2 100755
--- a/wqflask/wqflask/views.py
+++ b/wqflask/wqflask/views.py
@@ -46,6 +46,9 @@ from wqflask.interval_mapping import interval_mapping
 from wqflask.correlation import show_corr_results
 from wqflask.correlation_matrix import show_corr_matrix
 from wqflask.correlation import corr_scatter_plot
+
+from wqflask.wgcna import wgcna_analysis
+
 from utility import temp_data
 
 from base import webqtlFormData
@@ -93,7 +96,7 @@ def tmp_page(img_path):
     print("img_path:", img_path)
     initial_start_vars = request.form
     print("initial_start_vars:", initial_start_vars)
-    imgfile = open('/home/zas1024/tmp/' + img_path, 'rb')
+    imgfile = open(webqtlConfig.TMPDIR + img_path, 'rb')
     imgdata = imgfile.read()
     imgB64 = imgdata.encode("base64")
     bytesarray = array.array('B', imgB64)
@@ -174,6 +177,20 @@ def help():
     doc = docs.Docs("help")
     return render_template("docs.html", **doc.__dict__)
 
+@app.route("/wgcna_setup", methods=('POST',))
+def wcgna_setup():
+    print("In wgcna, request.form is:", request.form)             # We are going to get additional user input for the analysis
+    return render_template("wgcna_setup.html", **request.form)          # Display them using the template
+
+@app.route("/wgcna_results", methods=('POST',))
+def wcgna_results():
+    print("In wgcna, request.form is:", request.form)
+    wgcna = wgcna_analysis.WGCNA()                                # Start R, load the package and pointers and create the analysis
+    wgcnaA = wgcna.run_analysis(request.form)                     # Start the analysis, a wgcnaA object should be a separate long running thread
+    result = wgcna.process_results(wgcnaA)                        # After the analysis is finished store the result
+    return render_template("wgcna_results.html", **result)        # Display them using the template
+
+
 @app.route("/news")
 def news_route():
     newsobject = news.News()
@@ -328,7 +345,7 @@ def marker_regression_page():
         'mapmethod_rqtl_geno',
         'mapmodel_rqtl_geno'
     )
-    print("Random Print too see if it is running:", initial_start_vars)
+    print("Marker regression called with initial_start_vars:", initial_start_vars)
     start_vars = {}
     for key, value in initial_start_vars.iteritems():
         if key in wanted or key.startswith(('value:')):
diff --git a/wqflask/wqflask/wgcna/__init__.py b/wqflask/wqflask/wgcna/__init__.py
new file mode 100755
index 00000000..e69de29b
--- /dev/null
+++ b/wqflask/wqflask/wgcna/__init__.py
diff --git a/wqflask/wqflask/wgcna/wgcna_analysis.py b/wqflask/wqflask/wgcna/wgcna_analysis.py
new file mode 100644
index 00000000..9ab7950b
--- /dev/null
+++ b/wqflask/wqflask/wgcna/wgcna_analysis.py
@@ -0,0 +1,134 @@
+# WGCNA analysis for GN2
+# Author / Maintainer: Danny Arends <Danny.Arends@gmail.com>
+import sys
+from numpy import *
+import scipy as sp                            # SciPy
+import rpy2.robjects as ro                    # R Objects
+import rpy2.rinterface as ri
+
+from base import webqtlConfig                 # For paths and stuff
+from utility import webqtlUtil                # Random number for the image
+
+import base64
+import array
+
+from utility import helper_functions
+
+from rpy2.robjects.packages import importr
+utils = importr("utils")
+
+## Get pointers to some common R functions
+r_library       = ro.r["library"]             # Map the library function
+r_options       = ro.r["options"]             # Map the options function
+r_read_csv      = ro.r["read.csv"]            # Map the read.csv function
+r_dim           = ro.r["dim"]                 # Map the dim function
+r_c             = ro.r["c"]                   # Map the c function
+r_cat           = ro.r["cat"]                 # Map the cat function
+r_paste         = ro.r["paste"]               # Map the paste function
+r_unlist        = ro.r["unlist"]              # Map the unlist function
+r_unique        = ro.r["unique"]              # Map the unique function
+r_length        = ro.r["length"]              # Map the length function
+r_list          = ro.r.list                   # Map the list function
+r_matrix        = ro.r.matrix                 # Map the matrix function
+r_seq           = ro.r["seq"]                 # Map the seq function
+r_table         = ro.r["table"]               # Map the table function
+r_names         = ro.r["names"]               # Map the names function
+r_sink          = ro.r["sink"]                # Map the sink function
+r_file          = ro.r["file"]                # Map the file function
+r_png           = ro.r["png"]                 # Map the png function for plotting
+r_dev_off       = ro.r["dev.off"]             # Map the dev.off function
+
+class WGCNA(object):
+    def __init__(self):
+        print("Initialization of WGCNA")
+        log = r_file("/tmp/genenetwork_wcgna.log", open = "wt")
+        #r_sink(log)                                  # Uncomment the r_sink() commands to log output from stdout/stderr to a file
+        #r_sink(log, type = "message")
+        r_library("WGCNA")                            # Load WGCNA - Should only be done once, since it is quite expensive
+        r_options(stringsAsFactors = False)
+        print("Initialization of WGCNA done, package loaded in R session")
+        self.r_enableWGCNAThreads    = ro.r["enableWGCNAThreads"]        # Map the enableWGCNAThreads function
+        self.r_pickSoftThreshold     = ro.r["pickSoftThreshold"]         # Map the pickSoftThreshold function
+        self.r_blockwiseModules      = ro.r["blockwiseModules"]          # Map the blockwiseModules function
+        self.r_labels2colors         = ro.r["labels2colors"]             # Map the labels2colors function
+        self.r_plotDendroAndColors   = ro.r["plotDendroAndColors"]       # Map the plotDendroAndColors function
+        print("Obtained pointers to WGCNA functions")
+
+    def run_analysis(self, requestform):
+        print("Starting WGCNA analysis on dataset")
+        self.r_enableWGCNAThreads()                                      # Enable multi threading
+        self.trait_db_list = [trait.strip() for trait in requestform['trait_list'].split(',')]
+        print("Retrieved phenotype data from database", requestform['trait_list'])
+        helper_functions.get_trait_db_obs(self, self.trait_db_list)
+
+        self.input = {}           # self.input contains the phenotype values we need to send to R
+        strains = []              # All the strains we have data for (contains duplicates)
+        traits  = []              # All the traits we have data for (should not contain duplicates)
+        for trait in self.trait_list:
+            traits.append(trait[0].name)
+            self.input[trait[0].name] = {}
+            for strain in trait[0].data:
+                strains.append(strain)
+                self.input[trait[0].name][strain]  = trait[0].data[strain].value
+
+        # Transfer the load data from python to R
+        uStrainsR = r_unique(ro.Vector(strains))    # Unique strains in R vector
+        uTraitsR = r_unique(ro.Vector(traits))      # Unique traits in R vector
+
+        r_cat("The number of unique strains:", r_length(uStrainsR), "\n")
+        r_cat("The number of unique traits:", r_length(uTraitsR), "\n")
+
+        # rM is the datamatrix holding all the data in R /rows = strains columns = traits
+        rM = ro.r.matrix(ri.NA_Real, nrow=r_length(uStrainsR), ncol=r_length(uTraitsR), dimnames = r_list(uStrainsR, uTraitsR))
+        for t in uTraitsR:
+            trait = t[0]                  # R uses vectors every single element is a vector
+            for s in uStrainsR:
+                strain = s[0]             # R uses vectors every single element is a vector
+                rM.rx[strain, trait] = self.input[trait].get(strain)  # Update the matrix location
+                print(trait, strain, " in python: ", self.input[trait].get(strain), "in R:", rM.rx(strain,trait)[0])
+                sys.stdout.flush()
+
+        # TODO: Get the user specified parameters
+
+        self.results = {}
+        # Calculate a good soft threshold
+        powers = r_c(r_c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10), r_seq(12, 20, 2))
+        sft    = self.r_pickSoftThreshold(rM, powerVector = powers, verbose = 5)
+
+        # Create block wise modules using WGCNA
+        network = self.r_blockwiseModules(rM, power = 6, verbose = 3)
+
+        # Save the network for the GUI
+        self.results['network'] = network
+
+        # How many modules and how many gene per module ?
+        print(r_table(network[1]))
+
+        # The iconic WCGNA plot of the modules in the hanging tree
+        self.results['imgurl'] = webqtlUtil.genRandStr("WGCNAoutput_") + ".png"
+        self.results['imgloc'] = webqtlConfig.TMPDIR + self.results['imgurl']
+        r_png(self.results['imgloc'])
+        mergedColors = self.r_labels2colors(network[1])
+        self.r_plotDendroAndColors(network[5][0], mergedColors, "Module colors", dendroLabels = False, hang = 0.03, addGuide = True, guideHang = 0.05)
+        r_dev_off()
+        sys.stdout.flush()
+
+    def render_image(self, results):
+        print("pre-loading imgage results:", self.results['imgloc'])
+        imgfile = open(self.results['imgloc'], 'rb')
+        imgdata = imgfile.read()
+        imgB64 = imgdata.encode("base64")
+        bytesarray = array.array('B', imgB64)
+        self.results['imgdata'] = bytesarray
+
+    def process_results(self, results):
+        print("Processing WGCNA output")
+        template_vars = {}
+        template_vars["input"] = self.input
+        template_vars["results"] = self.results
+        self.render_image(results)
+        sys.stdout.flush()
+        #r_sink(type = "message")                                   # This restores R output to the stdout/stderr
+        #r_sink()                                                   # We should end the Rpy session more or less
+        return(dict(template_vars))
+