Sebastian Bassi (aka Bloguero Connor): Gráfico "científico"

sábado, febrero 28, 2009

Gráfico "científico"

Siempre me gustaron los gráficos científicos, me refiero a los que salen en revistas como Nature, Science, PNAS, etc. Por lo general usan SigmaPlot, MathCad y esas cosas, aunque a veces se puede ser un Excel (y suele quedar horrible). En este caso me tocó hacer un gráfico "a mano", esto es, sin usar ningún programa de terceros, sino lo hice yo. La idea era representar uno de estos eventos una alineamiento de 2 secuencias nucleotidicas: Mutacion, Inserción y Deleción. Las mutaciones (cambio de un nucleotido a otro) se tiene que representar por una barra de un color segun el nucleotido. Las inserciones y deleciones son barras grises con un tamaño y orientación en función de la magntud del evento.

Los datos los tenia que sacar de planillas, como por ejemplo esta: Planilla 36

Muestro esa porque el siguiente gráfico tiene (entre otras) a la barra generada en base a esa planilla:

Gráfico de comparación de secuencias nucleotidicas

Aqui está el código, es el que produjo esto, pero no es final porque aun quedan hacer algunas modificaciones con respecto al orden de las figuras, pero no con los graficos en si.


#!/home/user/MMC-SD/partition1/py252/bin/python

import xlrd
import glob
#f_inXLS = '/home/sbassi/bioinfo/INTA/ejemploSNPs.xls'
#outfile = '/tmp/out.svg'
outfile = '/mnt/hda2/all.svg'
f_inXLS = 'ejemplo SNPs.xls'

outfh = open(outfile,'w')

svgini = '''<?xml version="1.0" standalone="no"?>
    <!DOCTYPE svg PUBLIC "-//W3C//DTD SVG 1.1//EN"
    "http://www.w3.org/Graphics/SVG/1.1/DTD/svg11.dtd">
    <svg width="4000" height="1000" version="1.1" 
    xmlns="http://www.w3.org/2000/svg">'''
svgend = '</svg>'    
    
def grafica(xlsfin,outfh,j,tdetodos):
    book = xlrd.open_workbook(xlsfin)
    sh = book.sheet_by_index(0)
    datos = []
    largototal = sh.nrows
    for rx in range(1,largototal):
        n_orden = int(sh.row(rx)[0].value)
        color = ''
        signo = '+'
        valor = 1
        if sh.row(rx)[1].value == 1:
            color = 'red'
        elif sh.row(rx)[2].value == 1:
            color = 'green'
        elif sh.row(rx)[3].value == 1:
            color = 'blue'
        elif sh.row(rx)[4].value == 1:
            color = 'black'
        if sh.row(rx)[5].value != 0:
            valor = abs(int(sh.row(rx)[5].value))
            signo = '+' if sh.row(rx)[5].value>0 else '-'
            color = 'grey'
        elif sh.row(rx)[6].value != 0:
            valor = abs(int(sh.row(rx)[6].value))
            signo = '+' if sh.row(rx)[6].value>0 else '-'
            color = 'grey'
        elif sh.row(rx)[7].value != 0:
            valor = abs(int(sh.row(rx)[7].value))
            signo = '+' if sh.row(rx)[7].value>0 else '-'
            color = 'grey'
        elif sh.row(rx)[8].value != 0:
            valor = abs(int(sh.row(rx)[8].value))
            signo = '+' if sh.row(rx)[8].value>0 else '-'
            color = 'grey'
        if not color:
            color = 'white'
        datos.append((n_orden, color, valor, signo))
    print 'largototal',largototal
    print 'tdetodos',tdetodos
    outfh.write('<g transform="translate(%s,%s)">\n'%
                 (((tdetodos)),j*50))
    tdetodos += largototal+20
    outfh.write('<text style="font-size:12;font-weight:bold;\
                font-family:Consolas;fill:black" x="100" \
                y="175" >%s</text>'%(xlsfin))
    outfh.write('<line x1="100" y1="150" x2="%s" y2="150" \
                style="stroke:black;stroke-width:1"/>\n'%
                (100+largototal))
    #out += htmlini
    multiplic_alto = 5
    for dato in datos:
        #print dato
        x = 100+dato[0]
        #h = 10
        if dato[3]=='+':
            y = 150-(dato[2]*multiplic_alto)
            h = 150-y
        else:
            y = 150
            h = multiplic_alto*dato[2]
        if dato[1]!='white':
            outfh.write('<rect x="%s" y="%s" width="1" \
            height="%s" style="fill:%s;stroke-width:0"/>\n'
            % (x,y,h,dato[1]))
    outfh.write('</g>\n')
    return outfh,tdetodos

outfh.write(svgini)
i = 0
j = 0
tdetodos = 0
for fn in glob.glob('*.xls'):
    if 'ejemplo' not in str(fn) and 'my' not in str(fn):
        #print ' a'
        print fn
        outfh, tdetodos = grafica(fn,outfh,j,tdetodos)
        i+=1
        if i==5:
            i=0
            tdetodos = 0
            j+=1
    

outfh.write(svgend)
outfh.close()