#!/usr/bin/env python # Copyright (C) 2007 Monash University # # This program is free software; you can redistribute it and/or modify # it under the terms of the GNU General Public License as published by # the Free Software Foundation; either version 2 of the License, or # (at your option) any later version. # # This program is distributed in the hope that it will be useful, # but WITHOUT ANY WARRANTY; without even the implied warranty of # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the # GNU General Public License for more details. # # You should have received a copy of the GNU General Public License # along with this program; if not, write to the Free Software # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA # TODO: base plots on genes only, reverse complement where appropriate # TODO: Should offer independant components rather than PCA in scatter plot # (PCA doesn't make sense except with first two components) # TODO: More scatter plot variables import matplotlib matplotlib.use('TkAgg') import sys, numpy, pylab, Tkinter, tkFileDialog, tkMessageBox, gzip, sets, traceback, code, optparse from matplotlib import patches, collections from numpy import random from Bio import SeqIO SHRINK = 20 AMINOS = 'ACDEFGHIKLMNPQRSTVWY' def cached(func): """ Caching decorator """ cache = { } def outer_func(*args): if args not in cache: cache[args] = func(*args) return cache[args] return outer_func # Status windows ======================================================== # TODO: use these more def make_status_window(message): window = Tkinter.Tk() window.wm_title('Gene Spaghetti is noodling...') Tkinter.Label(window, text=message).pack() window.update() return window # Plot window base-class ======================================================== class Controller: """ Base class """ figsize = None description = 'plot' def __init__(self, *args, **kwargs): status_window = make_status_window('Drawing your ' + self.description) try: self.figure = pylab.figure(figsize=self.figsize) pylab.ioff() self.axes = pylab.gca() self.manager = pylab.get_current_fig_manager() self.canvas = self.figure.canvas Tkinter.Label(master=self.manager.toolbar, text='<- DPI:').pack(side=Tkinter.LEFT,padx=5) self.dpi_var = Tkinter.DoubleVar(self.manager.toolbar) self.dpi_var.set(80.0) entry = Tkinter.Entry(master=self.manager.toolbar, width=6, textvariable=self.dpi_var) entry.pack(side=Tkinter.LEFT) self.manager.toolbar.bsave['command'] = self.handle_save_figure self._setup(*args, **kwargs) finally: status_window.withdraw() pylab.ion() self.manager.show() def set_title(self, title): self.manager.toolbar.window.wm_title(title) def handle_save_figure(self): from tkFileDialog import asksaveasfilename fname = asksaveasfilename( master=self.manager.toolbar.window, title='Save the figure at %.1f dpi' % self.dpi_var.get(), filetypes=[ ('Portable Network Graphics','*.png'), ('Encapsulated Postscript File','*.eps'), ('Scalable Vector Graphics','*.svg'), ]) if fname == "" : return else: self.canvas.print_figure(fname, dpi=self.dpi_var.get()) class Scatter_controller(Controller): description = 'scatter plot' def _setup_scatter(self): self.selected = None self.markers = [] self.handlers = [ ] self.dragged = False self.canvas.mpl_connect('button_press_event', self.handle_press) self.canvas.mpl_connect('motion_notify_event', self.handle_motion) self.canvas.mpl_connect('button_release_event', self.handle_release) def handle_press(self, event): self.dragged = False def handle_motion(self, event): self.dragged = True def handle_release(self, event): if event.xdata is None or \ self.manager.toolbar.mode or \ self.dragged: return xlim = self.axes.get_xlim() ylim = self.axes.get_ylim() dx = (self.x - event.xdata) / (xlim[1]-xlim[0]) dy = (self.y - event.ydata) / (ylim[1]-ylim[0]) d = (dx*dx)+(dy*dy) #i = self.position[numpy.argmin(d)] n = numpy.argmin(d) for handler in self.handlers: handler(self.start[n],self.end[n]) #i-self.window_size*2,i+self.window_size*2) def refresh_marker(self): while self.markers: self.canvas.get_tk_widget().delete(self.markers[-1]) del self.markers[-1] if self.selected is None: return height = self.figure.bbox.height() x, y = self.axes.transData.xy_tup((self.x[self.selected],self.y[self.selected])) y = height-y self.markers.append( self.canvas.get_tk_widget().create_line(x-20,y,x+20,y) ) self.markers.append( self.canvas.get_tk_widget().create_line(x,y-20,x,y+20) ) def goto(self, start, end): if start is None or end is None: self.selected = None else: self.selected = numpy.argmin(numpy.maximum(0.0, self.start-start) + numpy.maximum(0.0, end-self.end) + 1e-30*numpy.absolute(self.start+self.end-(start+end))) self.refresh_marker() # Biopython stuff ======================================================== def gzip_open(filename): if filename.lower().endswith('.gz'): file = gzip.open(filename) filename = filename[:-3] else: file = open(filename, 'rb') return file def load(filename): window = make_status_window('Loading '+filename) try: for format in ['embl', 'genbank', 'fasta']: result = list(SeqIO.parse(gzip_open(filename), format)) if result: print filename, 'is', format, 'with', len(result), len(result) == 1 and 'record' or 'records' break else: raise Exception("Couldn't load file") #for record in result: # for letter in sets.Set(record.seq.data): # if letter not in 'ATCGNatcgn': # raise Exception("Not DNA:" + repr(list(sets.Set(record.seq.data)))) return result finally: window.withdraw() # Feature viewer stuff =================================================== class Gene(patches.Rectangle): def get_verts(self): bbox1 = self._transform.get_bbox1() bbox2 = self._transform.get_bbox2() aspect = (bbox2.width() / bbox2.height()) / (bbox1.width() / bbox1.height()) x, y = self.xy retreat = min(abs(self.height/aspect*0.5), abs(self.width*0.5)) if self.width < 0.0: retreat = -retreat return ( (x, y), (x, y+self.height), (x+self.width-retreat, y+self.height), (x+self.width, y+self.height*0.5), (x+self.width-retreat, y), ) def patch_draw(axes, pairs): old_draw = axes.draw def new_draw(renderer): axes.set_ylim((-0.5,3)) left = axes.transData.get_bbox1().xmin() right = axes.transData.get_bbox1().xmax() scale = axes.transData.get_bbox1().width() / axes.transData.get_bbox2().width() visible = scale < 50.0 for a,b in pairs: this_left = min(a.xy[0],a.xy[0]+a.width) this_right = max(a.xy[0],a.xy[0]+a.width) this_visible = (this_right >= left and this_left <= right) a.set_visible( this_visible ) b.set_visible( (visible or abs(a.width)/scale > 10.0) and this_visible ) #extent1.width() > 10 ) old_draw(renderer) axes.draw = new_draw class Feature_controller(Controller): figsize = (12,3) description = 'gene browser' def _setup(self, record): self.record = record self.axes.set_yticks(()) self.axes.set_position((0.005,0.2,0.99,0.79)) pairs = [ ] for feature in record.features: if feature.type != 'CDS': continue x = feature.location.nofuzzy_start width = feature.location.nofuzzy_end - x name = ' '.join((feature.qualifiers.get('gene',[]) or feature.qualifiers.get('locus_tag',[])) + feature.qualifiers.get('product',[])) color = feature.qualifiers.get('colour',['0 0 0'])[0] color = tuple([ int(i)/255.0 for i in color.split() ]) if feature.strand == -1: x += width width = -width patch = Gene((x,0),width,1, linewidth=0, facecolor=color) text = self.axes.text(min(x+width*0.25,x+width*0.75),1,name, clip_on=True, ha='left', rotation=45.0) self.axes.add_patch( patch ) pairs.append((patch, text)) x += width self.axes.autoscale_view() patch_draw(self.axes, pairs) Tkinter.Label(master=self.manager.toolbar, text='Search:').pack(side=Tkinter.LEFT,padx=5) self.find_var = Tkinter.StringVar(self.manager.toolbar) self.find_entry = Tkinter.Entry(master=self.manager.toolbar, textvariable=self.find_var) self.find_entry.pack(side=Tkinter.LEFT) self.find_entry.bind('', self.handle_find) Tkinter.Button(master=self.manager.toolbar, text='Find', command=self.handle_find).pack(side=Tkinter.LEFT,padx=5) self.find_status = Tkinter.Label(master=self.manager.toolbar, text='', fg='red') self.find_status.pack(side=Tkinter.LEFT,padx=5) self.cursor = patches.Rectangle((0,-1), 200, 10, linewidth=1, edgecolor=[0.75,0.75,0.75], facecolor=[0.75,0.75,0.75], zorder=-1, visible=False) self.axes.add_patch(self.cursor) self.handlers = [ ] def my_handler(event): if event.xdata is None or \ self.manager.toolbar.mode: return for handler in self.handlers: handler(event.xdata,event.xdata) self.figure.canvas.mpl_connect('button_press_event', my_handler) self.manager.toolbar.push_current() def goto(self, start, end): if start is None or end is None: self.cursor.set_visible(False) return self.axes.set_xlim((start-10000,end+10000)) self.cursor.set_x(start) self.cursor.set_width(end-start) self.cursor.set_visible(True) self.figure.canvas.draw() def handle_find(self, *args): text = self.find_var.get().strip().lower() self.find_status['text'] = '' self.find_entry.select_range(0,'end') for field in ['locus_tag','gene','product']: for feature in self.record.features: if feature.type != 'CDS': continue for field_text in feature.qualifiers.get(field,()): if text in field_text.lower(): break else: continue break else: continue break else: self.find_status['text'] = 'Not found' return for handler in self.handlers: handler(feature.location.nofuzzy_start,feature.location.nofuzzy_end) # Base usage stuff ======================================================= def extend(sequence, amount=1): return numpy.concatenate([sequence[:amount][::-1], sequence, sequence[-amount:][::-1]]) def bump(sequence, stddev, n=5): variance = stddev ** 2 factor = 1.0 while n: radius = 0 r_sum = 0 while r_sum*n < variance*(radius*2+1): radius += 1 r_sum += 2*radius*radius variance -= r_sum/float(radius*2+1) n -= 1 if not radius: continue sequence = extend(sequence, radius+1)[:-1] cumseq = numpy.cumsum(sequence) sequence = (cumseq[radius*2+1:] - cumseq[:-radius*2-1]) factor /= (radius*2+1) return sequence*factor def shrink(sequence, amount): sequence = sequence[:len(sequence)//amount*amount] sequence.shape = (len(sequence)//amount, amount) return numpy.average(sequence, 1) def show_base_plot_axes(mode, sc=1.0): pylab.axis('equal') pylab.axis('off') if sc == 1.0: prefix = '' else: prefix = '%.0f%% '%(sc*100.0) if mode in ('at-g-c','gc-a-t'): pylab.gca().set_position((0.0,0.04,1.0,0.92)) pylab.fill( [0.87*sc,0.87*sc,-0.87*sc], [1.0*sc,-1.0*sc,0.0*sc], facecolor='w', zorder=-2, clip_on=False ) at_g_c = (mode == 'at-g-c') pylab.text(-0.87*sc, 0.0*sc, prefix+(at_g_c and 'A+T' or 'C+G'), ha='right',va='center') pylab.text(0.87*sc, 1.0*sc, prefix+(at_g_c and 'G' or 'A'), ha='left',va='bottom') pylab.text(0.87*sc, -1.0*sc, prefix+(at_g_c and 'C' or 'T'), ha='left',va='top') pylab.text(-1.0*sc, -1.0*sc, { 'at-g-c' : 'Colour is AT skew (T-A)/(T+A)\nRed = T-rich\nBlue = A-rich', 'gc-a-t' : 'Colour is GC skew (C-G)/(C+G)\nRed = C-rich\nBlue = G-rich' }[mode]) else: pylab.fill( [sc,sc,-sc,-sc], [sc,-sc,-sc,sc], facecolor='w', zorder=-2, lw=2, clip_on=False ) pylab.text( -sc, -sc, prefix+'A', ha='right', va='top') pylab.text( -sc, sc, prefix+'T', ha='right', va='bottom') pylab.text( sc, sc, prefix+'G', ha='left', va='bottom') pylab.text( sc, -sc, prefix+'C', ha='left', va='top') pylab.text( 0, -sc*1.05, 'Colour is pyramidine-purine skew\nRed = CT-rich Blue = AG-rich', ha='center', va='top' ) class Base_plot_controller(Scatter_controller): description = 'base usage plot' #Modes available: at-g-c, gc-a-t, a-t-g-c def _setup(self, sequence_starts, sequences, mode='at-g-c', window_size=200.0, scale=1.0): self.set_title('Base usage plot') show_base_plot_axes(mode, sc=scale) self.window_size = window_size self.x = [] self.y = [] self.position = [] for sequence_start, sequence in zip(sequence_starts, sequences): sequence = numpy.fromstring(sequence.upper(), dtype=numpy.uint8) def get(letter): return bump(shrink(sequence == ord(letter), SHRINK), window_size/SHRINK) a = get('A') c = get('C') g = get('G') t = get('T') total = a+c+g+t a /= total c /= total g /= total t /= total del total if mode == 'at-g-c': # # G # AT # C # skew = (t-a)/(t+a) y = g-c x = 0.87-(a+t)*1.74 elif mode == 'gc-a-t': # # A # GC # T # skew = (c-g)/(c+g) y = a-t x = 0.87-(c+g)*1.74 elif mode == 'a-t-g-c': # # T G # # A C # x = g+c-a-t y = t+g-a-c skew = c+t-a-g norm_x = extend(y[:-2]-y[2:])*0.5 norm_y = extend(x[2:]-x[:-2])*0.5 lengths = numpy.sqrt(norm_x*norm_x+norm_y*norm_y) + 1e-300 #Avoid divide by zero angle = numpy.arctan2(norm_x,norm_y) angle_change = extend(numpy.absolute((angle[2:]-angle[:-2]+numpy.pi)%(2*numpy.pi)-numpy.pi)) c = 0.00000005 * SHRINK max_width = numpy.sqrt(c/(numpy.pi/2)) widths = numpy.minimum(c/lengths, max_width) adjustment = widths/lengths norm_x *= adjustment norm_y *= adjustment X = numpy.array((x-norm_x,x+norm_x)) Y = numpy.array((y-norm_y,y+norm_y)) mid_skew = (skew[:-1]+skew[1:])*0.5 lw = 0.1 pylab.plot(X[0],Y[0],'k',linewidth=lw,zorder=-2) pylab.plot(X[1],Y[1],'k',linewidth=lw,zorder=-2) polys = numpy.array(( (X[0,:-1],Y[0,:-1]),(X[1,:-1],Y[1,:-1]),(X[1,1:],Y[1,1:]),(X[0,1:],Y[0,1:]) )) polys = numpy.transpose(polys, axes=(2,0,1)) collection = collections.PolyCollection(polys, linewidths=(0,)) collection.set_array(mid_skew) collection.set_cmap(None) collection.set_norm(None) collection.set_clim(-0.6, 0.6) collection.update_scalarmappable() collection.update_scalarmappable = lambda: None #Don't do it on each redraw pylab.gca().add_collection(collection) self.x.append(x) self.y.append(y) self.position.append(numpy.arange(len(x))*SHRINK) self.x = numpy.concatenate(self.x) self.y = numpy.concatenate(self.y) position = numpy.concatenate(self.position) self.start = position - window_size*2 self.end = position + window_size*2 pylab.colorbar(collection) self._setup_scatter() # Gene scatter plots ===================================================== @cached def amino_usage(records): hists = [ ] i = 0 for record in records: a = i*2.0*numpy.pi/len(records) color = [ numpy.sin(a)*0.5+0.5, numpy.sin(a+numpy.pi*2.0/3)+0.5*0.5, numpy.sin(a-numpy.pi*2.0/3)+0.5*0.5 ] i += 1 for feature in record.features: if feature.type != 'CDS': continue if 'translation' not in feature.qualifiers: continue sequence = feature.qualifiers['translation'][0] if sequence[-1:] == '*': sequence = sequence[:-1] #Remove end marker sequence = numpy.fromstring(sequence.upper(), dtype=numpy.uint8) def get(letter): return numpy.sum(sequence == ord(letter)) hist = numpy.array([ float(get(amino)) for amino in AMINOS ]) hist /= numpy.sum(hist) hists.append(hist) return numpy.array(hists) @cached def codon_usage(records): hists = [ ] i = 0 for record in records: a = i*2.0*numpy.pi/len(records) color = [ numpy.sin(a)*0.5+0.5, numpy.sin(a+numpy.pi*2.0/3)+0.5*0.5, numpy.sin(a-numpy.pi*2.0/3)+0.5*0.5 ] i += 1 for feature in record.features: if feature.type != 'CDS': continue if 'translation' not in feature.qualifiers: continue sequence = record.seq[feature.location.nofuzzy_start:feature.location.nofuzzy_end] if feature.strand == -1: sequence = sequence.reverse_complement() sequence = sequence.data sequence = numpy.fromstring(sequence.upper(), dtype=numpy.uint8) s = (sequence==ord('C'))*1+(sequence==ord('G'))*2+(sequence==ord('T'))*3 if len(s)%3: print 'Weird CDS, length not a multiple of 3', feature.location s = s[:len(s)//3*3] codons = s[::3]*16+s[1::3]*4+s[2::3] hist = numpy.zeros(64) for codon in codons: hist[codon] += 1.0 hist /= numpy.sum(hist) hists.append(hist) return numpy.array(hists) @cached def base_usage(records): hists = [ ] i = 0 for record in records: a = i*2.0*numpy.pi/len(records) color = [ numpy.sin(a)*0.5+0.5, numpy.sin(a+numpy.pi*2.0/3)+0.5*0.5, numpy.sin(a-numpy.pi*2.0/3)+0.5*0.5 ] i += 1 for feature in record.features: if feature.type != 'CDS': continue if 'translation' not in feature.qualifiers: continue sequence = record.seq[feature.location.nofuzzy_start:feature.location.nofuzzy_end] if feature.strand == -1: sequence = sequence.reverse_complement() sequence = sequence.data sequence = numpy.fromstring(sequence.upper(), dtype=numpy.uint8) hist = numpy.array([ numpy.sum(sequence==ord('A')), numpy.sum(sequence==ord('C')), numpy.sum(sequence==ord('G')), numpy.sum(sequence==ord('T')) ], dtype=numpy.double) hist /= numpy.sum(hist) hists.append(hist) return numpy.array(hists) @cached def svd(records, func): hists = func(records) #import multinomial_eb # Pre-process a bit #hists = multinomial_eb.multinomial_eb(hists) #hists /= numpy.sum(hists,1)[:,None] mean = numpy.average(hists,0) hists = hists - mean[None,:] from numpy import linalg return linalg.svd(hists, full_matrices=0) + (mean,) def gene_stat(func): @cached def inner_func(records): result = [ ] for record in records: for feature in record.features: if feature.type != 'CDS': continue if 'translation' not in feature.qualifiers: continue result.append(func(feature)) return numpy.array(result) return inner_func @gene_stat def gene_stat_position(feature): return (feature.location.nofuzzy_start+feature.location.nofuzzy_end) * 0.5 @gene_stat def gene_stat_number_of_codons(feature): return len(feature.qualifiers['translation'][0]) @gene_stat def gene_stat_isoelectric_point(feature): from Bio.SeqUtils import ProtParam a = ProtParam.ProteinAnalysis(feature.qualifiers['translation'][0].replace('*','')) return a.isoelectric_point() @gene_stat def gene_stat_molecular_weight(feature): from Bio.SeqUtils import ProtParam a = ProtParam.ProteinAnalysis(feature.qualifiers['translation'][0].replace('*','')) return a.molecular_weight() def gene_stat_gc_content(records): b = base_usage(records) # ACGT return b[:,1]+b[:,2] def gene_stat_gc_skew(records): b = base_usage(records) # ACGT return (b[:,1]-b[:,2])/(b[:,1]+b[:,2]) def gene_stat_at_skew(records): b = base_usage(records) # ACGT return (b[:,3]-b[:,0])/(b[:,3]+b[:,0]) gene_stats = [ ('Position', gene_stat_position), ('Number of codons', gene_stat_number_of_codons), ('Molecular weight', gene_stat_molecular_weight), ('Isoelectric point', gene_stat_isoelectric_point), ('GC content', gene_stat_gc_content), ('GC skew (C-G)/(C+G)', gene_stat_gc_skew), ('AT skew (T-A)/(T+A)', gene_stat_at_skew), ('Amino principle component 1', lambda records: svd(records,amino_usage)[0][:,0]), ('Amino principle component 2', lambda records: svd(records,amino_usage)[0][:,1]), ('Amino principle component 3', lambda records: svd(records,amino_usage)[0][:,2]), ('Amino principle component 4', lambda records: svd(records,amino_usage)[0][:,3]), ('Amino principle component 5', lambda records: svd(records,amino_usage)[0][:,4]), ('Codon principle component 1', lambda records: svd(records,codon_usage)[0][:,0]), ('Codon principle component 2', lambda records: svd(records,codon_usage)[0][:,1]), ('Codon principle component 3', lambda records: svd(records,codon_usage)[0][:,2]), ('Codon principle component 4', lambda records: svd(records,codon_usage)[0][:,3]), ('Codon principle component 5', lambda records: svd(records,codon_usage)[0][:,4]), ] gene_stat_names = [ item[0] for item in gene_stats ] gene_stats = dict(gene_stats) class Gene_plot_controller(Scatter_controller): description = 'scatter plot' def _setup(self, records, x_axis=None, y_axis=None): self.set_title('Scatter plot') start = [ ] end = [ ] colors = [ ] sizes = [ ] i = 0 for record in records: a = i*2.0*numpy.pi/len(records) color = [ numpy.sin(a)*0.5+0.5, numpy.sin(a+numpy.pi*2.0/3)+0.5*0.5, numpy.sin(a-numpy.pi*2.0/3)+0.5*0.5 ] i += 1 for feature in record.features: if feature.type != 'CDS': continue if 'translation' not in feature.qualifiers: continue start.append(feature.location.nofuzzy_start) end.append(feature.location.nofuzzy_end) if len(records) <= 1: color = feature.qualifiers.get('colour',['0 0 0'])[0] color = tuple([ int(i)/255.0 *0.75 for i in color.split() ]) colors.append(color) sizes.append(end[-1]-start[-1]) self.start = numpy.array(start) self.end = numpy.array(end) self.sizes = numpy.array(sizes)*0.03 self.colors = colors self.x = gene_stats[x_axis](records) self.y = gene_stats[y_axis](records) self.labels = [] pylab.xlabel(x_axis) pylab.ylabel(y_axis) #TODO: describe principle components... self.collection = pylab.scatter(self.x,self.y,self.sizes, linewidth=0.1) self.collection.set_facecolor(self.colors) self._setup_scatter() def component_report(records): window = Tkinter.Toplevel() window.wm_title('Principal components') scrollbar = Tkinter.Scrollbar(window) scrollbar.pack(side='right', fill='y') text = Tkinter.Text(window) text.pack(fill='both', expand=True) scrollbar.config(command=text.yview) text.config(yscrollcommand=scrollbar.set) u, s, vh, mean = svd(records,amino_usage) #Alternatively... #from numpy import linalg #usage = amino_usage(records) #covar = numpy.zeros((usage.shape[1],usage.shape[1])) #for item in usage: # item = item - mean # covar += item[:,None] * item[None,:] #W, V = linalg.eig(covar) #print W[:4] ** 0.5 #print V[:,:4] report = '' report += 'Amino acid principal components\n\n' report += ' Mean 1 2 3 4 5 ...\n\n' report += 'Std dev ' for j in xrange(5): report += ' %6.3f' % s[j] report += '\n\n' for i, amino in enumerate(AMINOS): report += ' ' + amino + ' %6.3f ' % mean[i] for j in xrange(5): report += ' %6.3f' % vh[j,i] report += '\n' u, s, vh, mean = svd(records,codon_usage) report += '\n\n\nCodon principal components\n\n' report += ' Mean 1 2 3 4 5 ...\n\n' report += 'Std dev ' for j in xrange(5): report += ' %6.3f' % s[j] report += '\n\n' for i1, n1 in enumerate('ACGT'): for i2, n2 in enumerate('ACGT'): for i3, n3 in enumerate('ACGT'): i = i1*16+i2*4+i3 report += ' ' + n1+n2+n3 + ' %6.3f ' % mean[i] for j in xrange(5): report += ' %6.3f' % vh[j,i] report += '\n' text.insert('end', report) # Main =================================================================== class Main_controller: def __init__(self, records, feature_controller=None): self.feature_controller = feature_controller self.records = records if feature_controller: window = feature_controller.figure.canvas.get_tk_widget().master else: window = Tkinter.Tk() window.wm_title('Gene Spaghetti') all_frame = Tkinter.Frame(window) all_frame.pack(fill=Tkinter.X) all_frame.columnconfigure(0, weight=1) frame = Tkinter.LabelFrame(all_frame,text='Base usage plots',relief=Tkinter.SUNKEN) frame.grid(row=0,column=0,padx=5,pady=5, sticky='NESW') Tkinter.Label(frame,text='Gaussian window std dev:')\ .grid(row=0,column=0,padx=5,pady=5) self.window_var = Tkinter.DoubleVar(frame) self.window_var.set(200.0) Tkinter.Spinbox(frame,from_=10.0,to=10000.0,increment=10.0,textvariable=self.window_var)\ .grid(row=0,column=1,padx=5,pady=5,sticky='WE') Tkinter.Label(frame,text='Zoom to:')\ .grid(row=0,column=2,padx=5,pady=5) self.scale_var = Tkinter.DoubleVar(frame) self.scale_var.set(1.0) Tkinter.Spinbox(frame,from_=0.1,to=1.0,increment=0.01,textvariable=self.scale_var)\ .grid(row=0,column=3,padx=5,pady=5,sticky='WE') Tkinter.Button(frame,text='AT-G-C plot',command=lambda: self.base_plot('at-g-c'))\ .grid(row=0,column=4,padx=5,pady=5) Tkinter.Button(frame,text='GC-A-T plot',command=lambda: self.base_plot('gc-a-t'))\ .grid(row=0,column=5,padx=5,pady=5) Tkinter.Button(frame,text='A-T-G-C plot',command=lambda: self.base_plot('a-t-g-c'))\ .grid(row=0,column=6,padx=5,pady=5) frame = Tkinter.LabelFrame(all_frame,text='Gene plots',relief=Tkinter.SUNKEN) frame.grid(row=1,column=0,padx=5,pady=5, sticky='NESW') frame.columnconfigure(0, weight=1) frame.columnconfigure(2, weight=1) self.gene_x = Tkinter.StringVar(frame) self.gene_x.set(gene_stat_names[0]) Tkinter.OptionMenu(frame,self.gene_x, *gene_stat_names)\ .grid(row=0,column=0,padx=5,pady=5, sticky='WE') Tkinter.Label(frame,text='vs')\ .grid(row=0,column=1,padx=5,pady=5) self.gene_y = Tkinter.StringVar(frame) self.gene_y.set(gene_stat_names[1]) Tkinter.OptionMenu(frame,self.gene_y, *gene_stat_names)\ .grid(row=0,column=2,padx=5,pady=5, sticky='WE') Tkinter.Button(frame,text='Scatter plot', command=lambda: self.gene_plot(self.gene_x.get(),self.gene_y.get()))\ .grid(row=0,column=3,padx=5,pady=5, sticky='WE') Tkinter.Button(frame,text='Principal component report', command=lambda: component_report(self.records))\ .grid(row=1,column=3,padx=5,pady=5, sticky='WE') self.controllers = [ ] self.start = None self.end = None if feature_controller: self.add_controller(feature_controller) def handle_destroy(*args): window.quit() all_frame.bind('', handle_destroy) def add_controller(self, controller): self.controllers.append(controller) controller.handlers.append(self.goto) controller.goto(self.start, self.end) def handle_destroy(*args): self.controllers = [ item for item in self.controllers if item is not controller ] controller.figure.canvas.get_tk_widget().bind('', handle_destroy) def base_plot(self, mode): plot_controller = Base_plot_controller( [0]*len(self.records), [ record.seq.data for record in self.records ], mode, window_size=self.window_var.get(), scale=self.scale_var.get() ) self.add_controller(plot_controller) return plot_controller def gene_plot(self, *args): plot_controller = Gene_plot_controller(self.records, *args) self.add_controller(plot_controller) return plot_controller def goto(self, start, end): self.start = start self.end = end self.refresh() def refresh(self): for controller in self.controllers: controller.goto(self.start, self.end) def main(): root = Tkinter.Tk() root.withdraw() if len(sys.argv) > 1: filenames = sys.argv[1:] else: filetypes = [ ('All supported', '*.gbk *.gbk.gz *.embl *.embl.gz *.fna *.fna.gz'), ('GenBank', '*.gbk *.gbk.gz'), ('EMBL', '*.embl *.embl.gz'), ('Fasta', '*.fna *.fna.gz'), ('All', '*'), ] filenames = tkFileDialog.askopenfilenames(title='Gene Spaghetti - select files', filetypes=filetypes) if not filenames: return records = tuple(sum([ load(filename) for filename in filenames ], [])) if len(records) == 1 and records[0].features: feature_controller = Feature_controller(records[0]) else: feature_controller = None controller = Main_controller(records, feature_controller) Tkinter.mainloop() return if __name__ == '__main__': try: main() except KeyboardInterrupt: pass except Exception, exception: tkMessageBox.showerror('Problem', ''.join(traceback.format_exception(sys.exc_type, sys.exc_value, sys.exc_traceback)))