PageRenderTime 43ms CodeModel.GetById 20ms app.highlight 19ms RepoModel.GetById 1ms app.codeStats 0ms

/tools/mutation/visualize.py

https://bitbucket.org/cistrome/cistrome-harvard/
Python | 391 lines | 381 code | 6 blank | 4 comment | 0 complexity | d82a8d907b004b44f87141ac629babfe MD5 | raw file
  1#!/usr/bin/env python
  2
  3'''
  4Mutation Visualizer tool
  5'''
  6
  7from __future__ import division
  8
  9import sys, csv, os, math
 10import optparse
 11
 12from galaxy import eggs
 13import pkg_resources
 14pkg_resources.require( "SVGFig" )
 15import svgfig as svg
 16
 17
 18SVGPan = """
 19/**
 20 *  SVGPan library 1.2
 21 * ====================
 22 *
 23 * Given an unique existing element with id "viewport", including the
 24 * the library into any SVG adds the following capabilities:
 25 *
 26 *  - Mouse panning
 27 *  - Mouse zooming (using the wheel)
 28 *  - Object dargging
 29 *
 30 * Known issues:
 31 *
 32 *  - Zooming (while panning) on Safari has still some issues
 33 *
 34 * Releases:
 35 *
 36 * 1.2, Sat Mar 20 08:42:50 GMT 2010, Zeng Xiaohui
 37 *      Fixed a bug with browser mouse handler interaction
 38 *
 39 * 1.1, Wed Feb  3 17:39:33 GMT 2010, Zeng Xiaohui
 40 *      Updated the zoom code to support the mouse wheel on Safari/Chrome
 41 *
 42 * 1.0, Andrea Leofreddi
 43 *      First release
 44 *
 45 * This code is licensed under the following BSD license:
 46 *
 47 * Copyright 2009-2010 Andrea Leofreddi (a.leofreddi@itcharm.com). All rights reserved.
 48 *
 49 * Redistribution and use in source and binary forms, with or without modification, are
 50 * permitted provided that the following conditions are met:
 51 *
 52 *    1. Redistributions of source code must retain the above copyright notice, this list of
 53 *       conditions and the following disclaimer.
 54 *
 55 *    2. Redistributions in binary form must reproduce the above copyright notice, this list
 56 *       of conditions and the following disclaimer in the documentation and/or other materials
 57 *       provided with the distribution.
 58 *
 59 * THIS SOFTWARE IS PROVIDED BY Andrea Leofreddi ``AS IS'' AND ANY EXPRESS OR IMPLIED
 60 * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
 61 * FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL Andrea Leofreddi OR
 62 * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
 63 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
 64 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
 65 * ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
 66 * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
 67 * ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 68 *
 69 * The views and conclusions contained in the software and documentation are those of the
 70 * authors and should not be interpreted as representing official policies, either expressed
 71 * or implied, of Andrea Leofreddi.
 72 */
 73
 74var root = document.documentElement;
 75
 76var state = 'none', stateTarget, stateOrigin, stateTf;
 77
 78setupHandlers(root);
 79
 80/**
 81 * Register handlers
 82 */
 83function setupHandlers(root){
 84        setAttributes(root, {
 85                "onmouseup" : "add(evt)",
 86                "onmousedown" : "handleMouseDown(evt)",
 87                "onmousemove" : "handleMouseMove(evt)",
 88                "onmouseup" : "handleMouseUp(evt)",
 89                //"onmouseout" : "handleMouseUp(evt)", // Decomment this to stop the pan functionality when dragging out of the SVG element
 90        });
 91
 92        if(navigator.userAgent.toLowerCase().indexOf('webkit') >= 0)
 93                window.addEventListener('mousewheel', handleMouseWheel, false); // Chrome/Safari
 94        else
 95                window.addEventListener('DOMMouseScroll', handleMouseWheel, false); // Others
 96}
 97
 98/**
 99 * Instance an SVGPoint object with given event coordinates.
100 */
101function getEventPoint(evt) {
102        var p = root.createSVGPoint();
103
104        p.x = evt.clientX;
105        p.y = evt.clientY;
106
107        return p;
108}
109
110/**
111 * Sets the current transform matrix of an element.
112 */
113function setCTM(element, matrix) {
114        var s = "matrix(" + matrix.a + "," + matrix.b + "," + matrix.c + "," + matrix.d + "," + matrix.e + "," + matrix.f + ")";
115
116        element.setAttribute("transform", s);
117}
118
119/**
120 * Dumps a matrix to a string (useful for debug).
121 */
122function dumpMatrix(matrix) {
123        var s = "[ " + matrix.a + ", " + matrix.c + ", " + matrix.e + "\\n  " + matrix.b + ", " + matrix.d + ", " + matrix.f + "\\n  0, 0, 1 ]";
124
125        return s;
126}
127
128/**
129 * Sets attributes of an element.
130 */
131function setAttributes(element, attributes){
132        for (i in attributes)
133                element.setAttributeNS(null, i, attributes[i]);
134}
135
136/**
137 * Handle mouse move event.
138 */
139function handleMouseWheel(evt) {
140        if(evt.preventDefault)
141                evt.preventDefault();
142
143        evt.returnValue = false;
144
145        var svgDoc = evt.target.ownerDocument;
146
147        var delta;
148
149        if(evt.wheelDelta)
150                delta = evt.wheelDelta / 3600; // Chrome/Safari
151        else
152                delta = evt.detail / -90; // Mozilla
153
154        var z = 1 + delta; // Zoom factor: 0.9/1.1
155
156        var g = svgDoc.getElementById("viewport");
157       
158        var p = getEventPoint(evt);
159
160        p = p.matrixTransform(g.getCTM().inverse());
161
162        // Compute new scale matrix in current mouse position
163        var k = root.createSVGMatrix().translate(p.x, p.y).scale(z).translate(-p.x, -p.y);
164
165        setCTM(g, g.getCTM().multiply(k));
166
167        stateTf = stateTf.multiply(k.inverse());
168}
169
170/**
171 * Handle mouse move event.
172 */
173function handleMouseMove(evt) {
174        if(evt.preventDefault)
175                evt.preventDefault();
176
177        evt.returnValue = false;
178
179        var svgDoc = evt.target.ownerDocument;
180
181        var g = svgDoc.getElementById("viewport");
182
183        if(state == 'pan') {
184                // Pan mode
185                var p = getEventPoint(evt).matrixTransform(stateTf);
186
187                setCTM(g, stateTf.inverse().translate(p.x - stateOrigin.x, p.y - stateOrigin.y));
188        } else if(state == 'move') {
189                // Move mode
190                var p = getEventPoint(evt).matrixTransform(g.getCTM().inverse());
191
192                setCTM(stateTarget, root.createSVGMatrix().translate(p.x - stateOrigin.x, p.y - stateOrigin.y).multiply(g.getCTM().inverse()).multiply(stateTarget.getCTM()));
193
194                stateOrigin = p;
195        }
196}
197
198/**
199 * Handle click event.
200 */
201function handleMouseDown(evt) {
202        if(evt.preventDefault)
203                evt.preventDefault();
204
205        evt.returnValue = false;
206
207        var svgDoc = evt.target.ownerDocument;
208
209        var g = svgDoc.getElementById("viewport");
210
211        if(evt.target.tagName == "svg") {
212                // Pan mode
213                state = 'pan';
214
215                stateTf = g.getCTM().inverse();
216
217                stateOrigin = getEventPoint(evt).matrixTransform(stateTf);
218        }
219        /*else {
220                // Move mode
221                state = 'move';
222
223                stateTarget = evt.target;
224
225                stateTf = g.getCTM().inverse();
226
227                stateOrigin = getEventPoint(evt).matrixTransform(stateTf);
228        }*/
229}
230/**
231 * Handle mouse button release event.
232 */
233function handleMouseUp(evt) {
234        if(evt.preventDefault)
235                evt.preventDefault();
236
237        evt.returnValue = false;
238
239        var svgDoc = evt.target.ownerDocument;
240
241        if(state == 'pan' || state == 'move') {
242                // Quit pan mode
243                state = '';
244        }
245}
246"""
247
248COLS_PER_SAMPLE = 7
249HEADER_COLS = 4
250
251HEIGHT = 6
252WIDTH = 12
253BAR_WIDTH = 1.5
254GAP = 2
255
256
257colors = {'A':'blue', 'C':'green', 'G':'orange', 'T':'red'}
258bases = ['A', 'C', 'G', 'T' ]
259
260def stop_error(message):
261    print >> sys.stderr, message
262    sys.exit(1)
263
264def validate_bases(n_a, n_c, n_g, n_t, total):
265    if n_a > total:
266        return 'A'
267    elif n_c > total:
268        return 'C'
269    elif n_g > total:
270        return 'G'
271    elif n_t > total:
272        return 'T'
273    return None
274
275def main(opts, args):
276    s = svg.SVG('g', id='viewport')
277    
278    # display legend
279    for i, b in enumerate( bases ):
280        bt = svg.SVG("tspan", b, style="font-family:Verdana;font-size:20%")
281        s.append(svg.SVG("text", bt, x=12+(i*10), y=3, stroke="none", fill="black"))
282        s.append(svg.SVG("rect", x=14+(i*10), y=0, width=4, height=3, 
283                         stroke="none", fill=colors[b], fill_opacity=0.5))
284
285    reader = open(opts.input_file, 'U')
286
287    samples = []
288    for i in range(int(len(args)/3)):
289        index = i*3
290        samples.append(dict(name=args[index],
291                            a_col=args[index+1],
292                            totals_col=args[index+2]))
293
294    if opts.zoom == 'interactive':
295        y = 35
296    else:
297        y = 25
298    for i, sample in enumerate(samples):
299        x = 23+(i*(WIDTH+GAP))
300        t = svg.SVG("text", svg.SVG("tspan", sample['name'], style="font-family:Verdana;font-size:25%"), 
301                    x=x, y=y, transform="rotate(-90 %i,%i)" % (x, y), stroke="none", fill="black")
302        s.append(t)
303    
304    count=1
305    for line in reader:
306        row = line.split('\t')
307        highlighted_position = False
308        show_pos = True
309        position = row[int(opts.position_col)-1]
310        ref = row[int(opts.ref_col)-1].strip().upper()
311        # validate
312        if ref not in bases: 
313            stop_error( "The reference column (col%s) contains invalid character '%s' at row %i of the dataset." % ( opts.ref_col, ref, count ) )
314        # display positions
315        if opts.zoom == 'interactive':
316            textx = 0
317        else:
318            textx = 7
319        bt = svg.SVG("tspan", str(position), style="font-family:Verdana;font-size:25%")
320        s.append(svg.SVG("text", bt, x=textx, y=34+(count*(HEIGHT+GAP)), stroke="none", fill="black"))
321        s.append(svg.SVG("rect", x=0, y=30+(count*(HEIGHT+GAP)), width=14, height=HEIGHT, 
322                         stroke='none', fill=colors[ref.upper()], fill_opacity=0.2))
323        
324        for sample_index, sample in enumerate(samples):
325            n_a = int(row[int(sample['a_col'])-1])
326            n_c = int(row[int(sample['a_col'])+1-1])
327            n_g = int(row[int(sample['a_col'])+2-1])
328            n_t = int(row[int(sample['a_col'])+3-1])
329            total = int(row[int(sample['totals_col'])-1])
330            # validate
331            base_error = validate_bases(n_a, n_c, n_g, n_t, total)
332            if base_error:
333                stop_error("For sample %i (%s), the number of base %s reads is more than the coverage on row %i." % (sample_index+1, 
334                                                                                                                     sample['name'], 
335                                                                                                                     base_error, 
336                                                                                                                     count))
337 
338            if total:
339                x = 16+(sample_index*(WIDTH+GAP))
340                y = 30+(count*(HEIGHT+GAP))
341                width = WIDTH
342                height = HEIGHT
343                if count%2:
344                    s.append(svg.SVG("rect", x=x, y=y, width=width, height=height, 
345                                     stroke='none', fill='grey', fill_opacity=0.25))
346                else:
347                    s.append(svg.SVG("rect", x=x, y=y, width=width, height=height, 
348                                     stroke='none', fill='grey', fill_opacity=0.25))
349                
350                for base, value in enumerate([n_a, n_c, n_g, n_t]):
351                    width = int(math.ceil(value / total * WIDTH))
352                    s.append(svg.SVG("rect", x=x, y=y, width=width, height=BAR_WIDTH, 
353                                     stroke='none', fill=colors[bases[base]], fill_opacity=0.6))
354                    y = y + BAR_WIDTH
355
356        count=count+1
357        
358    if opts.zoom == 'interactive':
359        canv = svg.canvas(s)
360        canv.save(opts.output_file)
361        import fileinput
362        flag = False
363        for line in fileinput.input(opts.output_file, inplace=1):
364            if line.startswith('<svg'):
365                print '<svg xmlns="http://www.w3.org/2000/svg" xmlns:xlink="http://www.w3.org/1999/xlink" version="1.1">'
366                flag = True
367                continue
368            else:
369                if flag:
370                    print '<script type="text/javascript">%s</script>' % SVGPan
371                flag = False
372            print line,
373    else:
374        zoom = int(opts.zoom)
375        w = "%ipx" % (x*(10+zoom))
376        h = "%ipx" % (y*(2+zoom))
377        canv = svg.canvas(s, width=w, height=h, viewBox="0 0 %i %i" %(x+100, y+100)) 
378        canv.save(opts.output_file)
379
380if __name__ == '__main__':
381    parser = optparse.OptionParser()
382    parser.add_option('-i', '--input-file', dest='input_file', action='store')
383    parser.add_option('-o', '--output-file', dest='output_file', action='store')
384    parser.add_option('-z', '--zoom', dest='zoom', action='store', default='1')
385    parser.add_option('-p', '--position_col', dest='position_col', action='store', default='c0')
386    parser.add_option('-r', '--ref_col', dest='ref_col', action='store', default='c1')
387    (opts, args) = parser.parse_args()
388    main(opts, args)
389    sys.exit(1)
390
391