Segons llegeixo a Depth-First, després de la mort del desenvolupador de Pymol, el manteniment del programa passarà a mans de la coneguda empresa de programari quimioinformàtic Schrödinger.

On January 8, 2010 Schrödinger reached an agreement with the estate of the late Dr. Warren L. DeLano to acquire PyMOL. Schrödinger will take over continued development and maintenance, as well as support and sales of PyMOL, including all current subscriptions. Schrödinger will also continue to actively support the open-source community of PyMOL.

Prior to Warren’s tragic and unexpected passing, he had been working closely with Schrödinger to progressively integrate PyMOL with Schrödinger’s graphical interface, Maestro. With a great sense of humility, we will work hard to pick up as best we can where Warren left off and will strive to honor his memory by continuing in the spirit and tradition of PyMOL.

Com he dit algunes vegades, empro el pybel a diari.Aquest conjunt de classes que permeten cridar l’openbabel des de qualsevol script python, van molt bé per fer diverses tasques en quimioinformàtica (manejar molècules, càlculs de descriptors, dibuixar molècules,…).

Va molt bé però, cal vigilar amb la memòria quan es treballa amb un gran nombre de molècules (parlo de l’ordre de milions de molècules). Si es fa servir a “petita escala” no hi ha cap problema, ara bé, si s’empra a gran escala cal tenir cura de l’ús de la memòria.

He estat dos dies lluitant amb això. Tenia un programet per analitzar uns quants milions de molècules, totes les proves m’anaven bé fins que en feia ús a gran escala. Aleshores la memòria s’anava carregant gradualment fins a saturar la màquina. M’he mirat el codi de dalt a baix, esborrat totes les variables després de fer-ne ús, escriure tots els fitxers de sortida al moment,… però res, memòria a “tope”. Finalment he vist que l’error venia del pybel (no sé si del propi pybel o de l’openbabel).

Pybel és un conjunt de classes escrites en Python per treballar amb l’OpenBabel. És una eina molt potent que et permet fer ús de totes les utilitats de l’OpenBabel des de Python. Jo l’empro molt i a diari!

Per, per exemple, llegir un smile i escriure’n un sdf només cal fer el següent:

1
2
3
4
5
6
7
8
9
10
11
12

#Importar el pybel
from pybel import *
 
#Ara llegirem l'smile
mol = readstring('smi', 'CCCC(Cl)CCF')
 
#Obrim el fitxer de sortida
out =Outputfile('sdf','NomDelFitxer.sdf')
#Hi escrivim i tanquem
 
out.write(mol)
out.close()

Això va molt bé però m’he trobat que escriu sdfs sense posar cap títol (nom) a les molècules. Pel que ho estic emprant ara però, necessito posar els InchiK com a títol de cada molècula. Hi he estat hores i hores, intentant entendre la classe “Outputfile”, però he desistit, no crec que es pugui fer, així, directament.

Com a alternativa he decidit emprar els mètodes natius de Python per escriure fitxers, barrejat amb el Pybel (l’opció molecula.write). Així doncs:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

#Importar el pybel
from pybel import *
 
#Ara llegirem l'smile
mol = readstring('smi', 'CCCC(Cl)CCF')
 
#Obrim el fitxer de sortida (tal com es fa amb Python). Li diem que és per escriure
out = open('NomDelFitxer.sdf','w')
 
#Hi escrivim el títol
out.write('TítolDeLaMolecula')
 
#I ara, la molècula
out.write(mol.write('sdf'))
 
#Finalment, tanquem
out.close()

De moment, és l’única solució que he trobat i va prou bé.

pDynamo és una llibreria en python per a fer simulacions moleculars mitjançant mecànica quàntica, mecànica molecular i mètodes híbrids mecànica quàntica/mecànica molecular. Està basat en Dynamo, unes llibreries en Fortran 90 de les que per desgràcia, no en tinc massa bones referències. Tot i això i, encara que fa més o menys un any que no s’actualitza (desconec si el projecte segueix viu), caldrà fer-li una ullada

Els “descriptors de Lipinski”, de la norma del cinc, inclouen del pes molecular, el coeficient octanol/aigua (cLogP), els punts donadors d’enllaç d’hidrogen (HBD) i els punts acceptors d’enllaç d’hidrogen (HBA).

Aquests dies he hagut de treballar amb diferents molècules d’una base de dades gegant i em calia calcular diferents descriptors. Per als de Lipinski, he fet aquest petit script (en python). Requereix pybel i es basa en un que hi ha publicat a l’article que presentapybel.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36

"""
Script to compute the Lipibski descriptors
Alfons Nonell-Canals - October 2009
Chemogenomics Lab
 
USAGE:
python Lipinski.py -i sdfFile
 
"""
 
import optparse
from pybel import *
from openbabel import * 
 
p = optparse.OptionParser()
p.add_option('--infile','-i',default='default')
options, arguments = p.parse_args()
#files
infile = options.infile
 
#function to compute the Lipinski descriptors
HBD = Smarts("[#7,#8;!H0]")
HBA = Smarts("[#7,#8]")
 
def lipinski(mol):
   desc = []
   desc.append(mol.molwt)
   desc.append(mol.calcdesc(['LogP']) ['LogP'])
   desc.append(len(HBA.findall(mol)))
   desc.append(len(HBD.findall(mol)))
   return desc
 
for mol in readfile("sdf",infile):
    molName = mol.title
    desc = lipinski(mol)
    print molName,desc

Necessito comparar diverses molècules. Encara no hi entenc massa, però, crec que una de les millors maneres és fer servir el coeficient de Tanimoto [en]. Això, el càlcul d’aquest coeficient l’he de fer de forma reiterada, necessito anar calculant el coeficient per cada molècula que generi i comparar-la així amb una molècula final. A més a més, vull integrar-ho dins un programa que tinc fet amb Python.

La millor manera, doncs, és fent ús de Pybel, un mòdul Python per fer servir eines de l’Openbabel des de codi Python. La manera de fer-ho és ben senzilla. Des de Python:

>>> import pybel #importarem el mòdul
>>> smiles = ['CCCC', 'CCCN'] #Array de molècules
>>> mols = [pybel.readstring("smi", x) for x in smiles] #Llegim les molècules
>>> fps = [x.calcfp() for x in mols] #calculem els fingerprints (descriptors)
>>> print fps[0] | fps[1] #finalment, calculem el coeficient
0.333333333333

I… llestos! suposo que hi deuen haver maneres més acurades i que, com tota cosa estadística, com més descriptors es tinguin per cada molècula, millor, però, per a començar, no va gens malament i el mètode és molt fàcil d’implementar.

Interessant utitlitat que fins ahir desconeixia totalment: els “sets”. Bàsicament, és un conjunt d’eines per a comparar, unir,… conjunts de dades. Amb les funcions de la classe podem, per exemple, obtenir de maneta senzilla la diferència entre dos conjunts de dades, obtenir les diferències, actualitzar-ne un a partir d’un segon.

Un exemple ben senzill és quan es té una llista A:

a = ['carboni', 'nitrogen', 'hidrogen']

I ho volem comparar amb una llista B:

b = ['oxigen', 'nitrogen', 'molibdè', 'calci', 'sofre']

Ara volem comparar les dues llistes, primer generem els dos sets

as = set(a)
bs = set(b)

Ara, per exemple, obtindrem una llista consens, és a dir, amb els elements comuns d’ambdues llistes:

as.intersection(b)

Això és un petit exemple, però, a partir d’aquí, es poden fer meravelles. Aquí hi ha un petit tutorial, molt didàctic.

Fins ahir tenia un problema prou complicat: estic treballant en un programa en Python que requereix múltiples i constants connexions a una base de dades mySQL. Quan treballava en local, no tenia cap problema però, quan executava el programa en màquines diferents apuntant cap a un servidor de bases de dades (com passaria, per exemple, en un cluster), de tant en tant i de manera aleatòria, el servidor mySQL es quedava “tonto” i no responia. Aquesta pèrdua de connexió feia que el programa “petés”.

Amb l’administrador de sistemes del grup hem treballat per mirar què passa amb el mySQL ja que sembla que el problema és allà, però no hem trobat solució. No hi ha ni problemes ni amb el límit màxim de connexions, ni amb la memòria, cpu,… en fi, per algun motiu que desconeixem, el servidor es satura i té petites penjades. Així dons, la solució passa per fer alguna cosa al programa perquè no peti quan la base de dades no respongui.

La solució està en les “excepcions” (de l’anglès exception, desconec si aquesta traducció literal és vàlida). És a dir, per defecte, quan python intenta accedir a la base de dades, si aquesta no respon, es crea una excepció que, si no se li diu res, provoca la sortida del programa. La idea està en modificar aquesta excepció perquè quan python vegi que la base de dades no respon, la solució no sigui sortir del programa sinó reintentar-ho. He jugat amb el try i except.

db = DB()
cursor = db.cursor()
sql = ’select count(*) from molecules’
number = cursor.execute(sql)

class DB:
def cursor(self):
try:
bd = MySQLdb.connect(host=HOSTNAME, user=”USERNAMEl”, passwd=”PASSWD”,db=”DB_NAME)
cursor = bd.cursor()
return cursor
except MySQLdb.OperationalError, message:
return self.cursor()

La primera part del codi inicia la classe DB i n’obté el “cursor”, operador amb el que realitzarem les crides sql. La part més interessant però és la classe DB, que conté la funció “cursor”. Aquesta funció connecta amb la base de dades i retorna el “cursor”, comanda dins l’ordre try. Però, el fet de jugar amb el try i except fa que el programa primer provi el que hi ha a try però, si falla i troba una excepció que coincideix amb la indicada a except, en l’exemple MySQLdb.OperationalError, en comptes d’aturar l’execució del programa (el que faria per defecte) fa el que se li diu, en aquest cas que torni a executar la funció.

Obprop és una aplicació inclosa al paquet Openbabel i que calcula alguns descriptors per una molècula donada. Cal passar-li la molècula en un fitxer sdf (hi pot haver diverses molècules dins el fitxer) i et retorna la seva fórmula, l’SMILE canònic, el pes molecular, el logP,…

Per fer-la servir simplement cal executar l’ordre:

obprop fitxer_amb_molècules.sdf

I, com a resultat, s’obté quelcom com ara:

name [Name]
formula [Formula]
mol_weight [Molecular Weight]
exact_mass [Isotopic Mass]
canonical_SMILES [String]
num_atoms [Number]
num_bonds [Number]
num_residues [Number]
sequence [Residue Sequence]
num_rings [Number of Rings (by SSSR)]
logP [Number (octanol-water partition)]
PSA [Number (topological polar surface area)]
MR [Number (molar refractivity)]
$$$

Eina útil, senzilla i ràpida per obtenir les quatre dades bàsiques d’una molècula.

Fa un parell de dies parlava d’una opció que tenim per “silenciar” part d’un residu en una dinàmica. El que comentava era que, ajustàvem la càrrega d’aquella part a zero i llestos. Per fer això, cada programa té les seves maneres.

La cosa està però en que allò que vaig dir no és ben bé cert i, si només ajustem la càrrega a zero, seguirem tenint tot de boles (àtoms) sense càrrega voltant per la proteïna. La dinàmica té en compte la càrrega i el diàmetre de cada àtom de manera que si només ajustem la càrrega a zero, seguim tenint les “boles” fent nosa.

Doncs, com fer-ho bé (ara espero que sí!)? Ajustant el radi de Van der Waals també a zero. La manera més senzilla per fer-ho és modificant el tipus atòmic dels àtoms que volem “silenciar”, fent-ne un a l’atzar, per exemple Oz, per silenciar un oxigen. Després, a la llibreria de tipus atòmic hi afegim aquest nou “tipus” assignant-li un radi de Van der Waals de zero. I llestos!