Noel O'Blog

Least Publishable Unit: Selective GPCR Agonists

2013-05-09T20:56:00.000+01:00

After a certain length of time, there comes a point when you realise that that idea of yours you always planned to pursue, well, you're never going to get around to it. Such is life and all that.

But maybe if I post some ideas here, it might spark someone else's imagination and lead to something. A mention in your Nobel speech is all I ask (preferably close to the start).

So here is something I worked on last (academic) year...

In collaboration with Dr JJ Keating (University College Cork) I was interested in finding some selective 5-HT_2A receptor agonists. Turns out everyone else wants to find selective 5-HT_2C agonists and so 2A is neglected. As far as I could tell, there are no known 2A-selective agonists. Such a compound would be useful both as a tool compound but also has some therapeutic potential for glaucoma.

To begin with I turned to the literature to see what had been done. In short, I never found any paper describing the search for a selective 2A agonist.

But then I thought of ChEMBL. Despite the fact that no-one had been looking for a selective 2A agonist, maybe there are some accidental examples to be found in papers on selective 2C agonists? So I downloaded the entire database, and searched for all instances where activity (EC50) for the same compound in the same paper was measured against both 2A and 2C. And lo and behold I found a few examples, the best of which was compound 5 in Bioorg. Med. Chem. Lett., 2005, 15, 4555.

So in short, your mission (if you choose to accept it) is to repeat the analysis, take these compounds as a starting point and to develop selective 5-HT_2A agonists through a medicinal chemistry strategy. Note that a similar approach towards receptor subtype selectivity might yield useful results for other receptors.

Image credit: Xavier Béjar on Flickr

How fast is chemfp? We investigate!

2013-04-25T22:24:00.000+01:00

That is, um, I investigate. I described chemfp in an earlier blogpost. Simply put, it comprises software by Python/C guru Andrew Dalke to handle and generate binary fingerprints for molecules. It enables easy comparison of fingerprints from different toolkits, and provides super-fast similarity methods. With a new even faster version out in Feb 2013, I thought I should put it through its paces.

Pre-processing

The standard way to carry out a fast similarity search with Open Babel is the so-called fastsearch method, which precalculates an index file of fingerprints for an sdf file, and then linearly-scans through this for hits against a single query structure. With chemfp, the initial step is similar; you convert everything to fps files, a standard file format for fingerprints developed by Andrew. Previously you needed to use chemfp's ob2fps to do this, but OB 2.3.2 added support for fps files and can generate them directly 2 to 3 times faster.

$ obabel chembl_15mod.sdf -O chembl.fs # for fastsearch
$ obabel chembl_15mod.sdf -O chembl.fps # for chemfp

Single-molecule query

Let's take the first molecule in the database as a query, and find those molecules in ChEMBL that are within 0.8 similarity (Tanimoto). simsearch is the name of chemfp's tool for similarity searching:

$ obabel chembl.fs -osmi -s query.sdf -aa -at 0.8 > threshold.fs.txt
$ simsearch --threshold 0.8 -q query.fps chembl.fps > threshold.fps.txt

They both take about 1.2 seconds. I'm being a bit vague because the exact value doesn't matter; for some query molecules fastsearch takes longer, for some less. It's just to give an idea. But in short, fastsearch appears to be in the same ballpark as simsearch for single molecule queries.

Multi-molecule query

fastsearch doesn't provide an easy way to do multiple queries. If you want to do it yourself, you would just have to do each search one-by-one. So 1000 searches would take 1000s, let's say. In contrast, with simsearch 1000 searches on ChEMBL takes only 8s:

$ simsearch --threshold 0.8 -q largequery.fps chembl.fps > thresholdb.fps.txt

How does it manage this? Well, first of all, it runs in parallel on 4 available CPUs by default (...see Andrew's comment below). These figures come from a server which has 4 hyperthreading CPUs. The rest of the improvement comes from the use of in-memory data structures and algorithmic magic, some of which Andrew has described on his blog over the last year or so.

With speeds like this, it brings all-against-all similarity searches within reach, and simsearch provides an option just for this. The following finds the nearest 8 molecules for each molecule in the dataset:

$ simsearch -k 8 --NxN largequery.fps > NxN.fps.txt

I tried this for the first 10000 molecules in ChEMBL, and it took about 1s; 50K took 11s; 100K took 32s; 250K took 144s; and the whole of ChEMBLdb (1.2 million) took 2900s (48m 20s). To put this in context, an all-against-all search of ChEMBLdb using fastsearch would take something like 14 days.

Single-molecule query revisited

The alert reader may be wondering why, if the multi-molecule query is so darned fast, the single-molecule query is no faster than Open Babel's fastsearch. The answer is rather simple (...but not quite right - see Andrew's comment below): 99% of the time spent on the single-molecule query is setup. Once the search datastructure has been initialised, then the response for a query is within milliseconds. You can see timings for these over at the chemfp website. To achieve these timings in practice, you would need to write a Python script that used the chemfp library and was accessed via some sort of client/server architective, e.g. a local webservice. As Andrew points out, this would allow search results to be returned instantly as a chemist sketches a structure. Trying this out is left as an exercise for the reader. :-)

In conclusion

So, in short, if you have any interest in running multiple queries against a database, comparing two large datasets, or finding pairwise similarity within a dataset, check out chemfp.

Notes:
1. I should point out that Andrew is making chemfp available as open source but with delayed release. Commercial licensees get support and the latest features.
2. For ease of use with Open Babel, the ChEMBLdb SDF file was modified to add the chembl_id to the title. Oh ye gods of ChEMBL please consider doing this for your humble users.

Image credit: Russell J Watkins on Flickr

Talk on Universal SMILES at New Orleans ACS

2013-04-15T14:50:00.000+01:00

Early on Wednesday I presented my recent paper on Universal SMILES at the New Orleans ACS. This is a canonical SMILES string that uses the InChI canonical labels. Usually I tell the audience that the slides will be made available, but this time there was someone in the audience who was standing up every so often and taking photos; I thought this was so awesome I said nothing.

Anyway, here are the slides. They are very wordy, partly because I used up all my visualisation skills fiddling around with the other formal talk I was giving (appearing soon on the NextMove blog), and partly because I was aware that for web readers a picture of a donkey surfing might not spell out how to create canonical SMILES quite as well as traditional bulletpoints.

Universal (and Inchified) SMILES are available right now in Open Babel. Rumour has it that the CDK and RDKit are considering supporting Universal SMILES. If you use these or any other toolkits, and think that having support for Universal SMILES would be nothing short of paradigm-shifting awesome, ask them to add support.

Talk on Open Babel at New Orleans ACS

2013-04-11T21:35:00.000+01:00

Now that Open Babel 2.3.2 has been released 6 months, I thought it might almost be time to talk about what's new, and also mention what's under development.

Here's the talk I gave at Rajarshi's CINF Flash session last Sunday. It was my first time at the flash session but I'm definitely going to make a point of attending and presenting at this in future; it is preceded by a free lunch!

See you at ACSNoelA?

2013-03-29T22:00:00.000Z

I'll be presenting at the Spring ACS National Meeting in New Orleans in just over a week. The last ACS I was at was three years ago so I'm looking forward to catching up with what's been going on, and meeting up with some familiar faces.

I've got three talks lined up, the slides for which I'll post after the event:

1. I'll be talking about "What's new and cooking in Open Babel?", as part of Rajarshi's CINF Flash talks on April 7 (Sunday) sometime between 12.30 and 2.00 in Room 350 (Morial Convention Center). I think Rajarshi is still accepting CINF Flash talks so get them in.

2. "Universal SMILES: Finally, a canonical SMILES string?"
This presents the work described in a recent paper.

DIVISION: CINF: Division of Chemical Information
SESSION: Public Databases Serving the Chemistry Community
DAY & TIME OF PRESENTATION: April 10 (Wed) from 8:35 am to 9:05 am
LOCATION: Morial Convention Center, Room: 350

3. "Roundtripping between small-molecule and biopolymer representations".
This describes some of the issues and challenges that have arisen in the development of the Sugar & Splice software, a toolkit for perceiving biopolymer structures, depicting them and converting them between formats.

DIVISION: CINF: Division of Chemical Information
SESSION: Linking Bioinformatic Data and Cheminformatic Data
DAY & TIME OF PRESENTATION: April 09 (Tues) from 3:10 pm to 3:35 pm
LOCATION: Morial Convention Center, Room: 349

To any first-time ACSers reading this, don't miss the various CINF functions. (I seem to recall that they are not exactly easy to find in the printed programme but details as ever are on the webs.)

Time and the InChI

2013-03-25T21:57:00.002Z

Notes:
1. Created using a Google spreadsheet and Timeline JS.
2. If anyone wants to send me updates (e.g. for commercial software like ChemDraw I found it hard to find dates and versions), feel free.

Police your code with fuzz testing

2013-03-18T20:49:00.000Z

I've just become a convert to Fuzz testing. This is a rather simple idea for testing software that processes some user input: just send in some random junk and see what happens (and keep repeating until something does).

I tried this for some software I've been working on, thinking it was so dumb it couldn't possibly flush out anything useful, but quickly changed my mind. Even more interesting, where the space of all possible user input can be stepped through systematically, fuzz testing can be used to map out the space of allowed input (which may or may not be what you were expecting). For example, you could find all 4-letter words acceptable to OPSIN by using the code below to generate fuzz and use OPSIN through Cinfony (for example) to find out whether an error is raised.

import random

class Fuzzer:
    def __init__(self, allowed, length):
        self.allowed = allowed
        self.length = length

    def systematic(self, text=""):
        if len(text) == self.length:
            yield text
        else:
            for x in self.allowed:
                for y in self.systematic(x+text):
                    yield y

    def random(self):
        return "".join([random.choice(allowed) for x in range(self.length)])

if __name__ == "__main__":
    allowed = "abcdefghijklmnopqrstuvwxyzABCDEFGHIJKLMNOPQRSTUVCWXYZ"
    fuzzer = Fuzzer(allowed, 3)
    for fuzz in fuzzer.systematic():
        print fuzz
    fuzzer = Fuzzer(allowed, 5)
    for i in range(10):
        print fuzzer.random()

Avogadro hydrogen-clicking craziness with Sikuli

2013-03-03T19:08:00.002Z

There's something that's always annoyed me when using Avogadro. On Windows, it often takes several clicks on a hydrogen before I hit the sweet spot and it sprouts a methyl group. This is only mildly annoying when building a structure, but very frustrating for operations where misclicking cause it to lose selection (e.g. rotating around a dihedral). It's awkward to report a bug about this though because it's hard to prove.

Enter Sikuli. It's a Java application running Jython that's used for testing GUIs (among other things). It uses a sort of visual programming style to match regions of the screen based on screenshots and carry out operations based on whether or not a match exists. (At this point you may want to check out the videos on the Sikuli website to figure out what I'm talking about.)

My original idea was to get Sikuli to click all over the offending hydrogen and see which points worked and which didn't. In the interests of time, I reduced this down to all points on a horizontal line through the centre of the hydrogen. First of all, here is the result, with the 'allowed area' indicated by the red line:

I would argue that that red line is neither long enough nor in the expected location. In any case, check out the sort of code needed to run this test:

You mightn't understand everything here, but you should get the gist. On every iteration of the loop it clears the Avogadro window, clicks to create a methane, and then clicks at point x, y on a particular hydrogen (where x, y are relative to the centre of the image of the matched hydrogen) trying to sprout a methyl. The "success test" checks that a methyl sprouted. The "fail test" checks whether a methane was added (i.e. the region to the right of the red bar in the diagram of the atom above). Note that fuzzy matching is used to match screenshots against the actual contents of the screen or a region.

cclib 1.1 and GaussSum 2.2.6 released

2013-02-23T10:14:00.002Z

cclib is a Python library for parsing for analysing comp chem log files from many different QM packages. GaussSum is a GUI that uses cclib to monitor the progress of comp chem calculations and calculate predicted spectra for comparison with experimental results. cclib is the work of Adam Tenderholt, Karol Langner and myself, while GaussSum is just by me.

New releases of both of these are now available for download at their respective websites: http://cclib.sf.net and http://gausssum.sf.net. For help, email cclib-users@lists.sf.net or gausssum-help@lists.sf.net. Here's what's new in these releases:

cclib 1.1
New Features:

Add progress info for all parsers
Support ONIOM calculations in Gaussian (Karen Hemelsoet)
New attribute atomcharges extracts Mulliken and Lowdin atomic charges if present
New attribute atomspins extracts Mulliken and Lowdin atomic spin densities if present
New thermodynamic attributes: freeenergy, temperature, enthalpy (Edward Holland)
Extract PES information: scanenergies, scancoords, scanparm, scannames (Edward Holland)

Bugfixes:

Handle coupled cluster energies in Gaussian 09 (Björn Dahlgren)
Vibrational displacement vectors missing for Gaussian 09 (Björn Dahlgren)
Fix problem parsing vibrational frequencies in some GAMESS-US files
Fix missing final scfenergy in ADF geometry optimisations
Fix missing final scfenergy for ORCA where a specific number of SCF cycles has been specified
ORCA scfenergies not parsed if COSMO solvent effects included
Allow spin unrestricted calculations to use the fragment MO overlaps correctly for the MPA and CDA calculations
Handle Gaussian MO energies that are printed as a row of asterisks (Jerome Kieffer)
Add more explicit license notices, and allow LGPL versions after 2.1
Support Firefly calculations where nmo != nbasis (Pavel Solntsev)
Fix problem parsing vibrational frequency information in recent GAMESS (US) files (Chengju Wang)
Apply patch from Chengju Wang to handle GAMESS calculations with more than 99 atoms
Handle Gaussian files with more than 99 atoms having pseudopotentials (Björn Baumeier)

GaussSum 2.2.6
New Features:

A patch from Thomas Pijper was integrated to enable calculation of Raman intensities (from Raman activity).
Support has been added for calculating charge density changes for unrestricted calculations (requested by Phil Schauer).
Blank lines in Groups.txt are now ignored.
Parser updated to cclib 1.1

A new home for Linux4Chemistry?

2013-02-17T18:50:00.000Z

Are you interested in taking over the stewardship of the Linux4Chemistry website? This website was setup by Nikodem Kuznik in 2001 to promote the usage of Linux for chemistry by listing chemistry software available for Linux, both commercial, free, and open source. In 2005 I took over its maintenance, moved it to its present location and made some changes including the awesome logo. I am currently on the lookout for someone interested and willing to take over the running of this site. Update (25/02/2103): A new home has been found.

Why am I looking for someone new? Well, I think I am no longer the right person to maintain this, as the goals I wished to achieve have already been met. To me the idea of promoting Linux as a platform for chemistry software today seems obsolete when almost all chemistry software vendors target Linux at a minimum. Every docking program is on Linux, every QM package runs on Linux, etc.

In taking over the website I also had the goal of clarifying the distinction between Open Source software and software available for free. To this end, I added license information to all of the software, and made it possible to filter the results by license type (I was very proud of my tongue-in-cheek logo for Shareware software). I think that today more people are aware of this distinction and in any case, I'm not sure that L4C is really playing a role here.

It is also true that while at the time my desktop machine was running Linux (Debian Sarge I think), since then I do most of my work on Windows but use a VM for Linux. And it's harder to maintain something which is outside my day-to-day usage.

So I'm looking for someone to bring a new vision to Linux4Chemistry and shake it up a bit. If you are this person, get in touch.

A compilation of speeds - Compiler face-off

2013-02-06T22:00:00.000Z

Let's get right into this one. I've compiled Open Babel with g++ in various ways, and am going to compare the speed with the MSVC++ release. Specifically I'm going to compare the wallclock time to convert 10000 molecules (the first 10000 in ChEMBL 13) from an SDF file to SMILES.

Our starting point is the time for the MSVC++ compiled release:
29.6s (MSVC++ 2010 Express 32-bit)

I have a Linux Mint 12 VM (VMWare) on the same machine, so let's run the same executable under Wine on Linux:
37.3s (MSVC++ 32-bit under Wine/Linux)
...so it's slower, pretty much as expected. The not-an-emulation layer slows things down a bit.

How about the MinGW compilation described in the previous post?:
24.1s (MinGW g++ 4.6.2 32-bit)
g++ beats MSVC++. To be honest, I was a bit surprised to see this, although I understand from Roger that g++ is surprisingly highly-optimised for cheminformatics toolkits. Maybe we should look into an official MinGW release in future.

What about Open Babel compiled with Cygwin's g++?:
39.5s (Cygwin g++ 4.5.3 32-bit)
As expected it runs like a pig compared to the MinGW version. Cygwin's handy, but when you're in a hurry it's maybe not the best choice.

So far, so not very unexpected. Now we will enter the realm of weirdness. Let's compile it on Linux in the VM and run it there:
14.8s (Linux Mint 12 g++ 4.6.1 64-bit)

So, in short, the fastest way to run Open Babel on Windows is to use a VM to run Linux. Huh? The like-with-like comparison of MinGW's 24.1 versus Linux's 14.8 is the most intriguing. It suggests that the slowdown is either due to rubbish file I/O by Windows, or sub-optimal platform-specific code in Open Babel's I/O handling code.

Either way, it's a pretty interesting result.

Notes:
1. Hardware was a Dell Latitude E6400 bought 3 years ago (Core 2 Duo 2.4 Ghz, 4GB Ram) running Win 7 64-bit. The timing was the best of three after timings had stabilised (the first one or two is usually a second or two slower).
2. After the initial post, I compiled clang on Linux, and then used it to compile Open Babel. Running the conversion took 15.3s.
3. Also, I ran the MinGW compiled version under Linux, and it took 30.7s.

Image credit: Venkatesh Srinivas (Extrudedaluminiu on Flickr)

Compiling Open Babel with MinGW on Windows

2013-02-06T09:31:00.000Z

If you want to compile on Windows using GCC, you have two alternatives: Cygwin's GCC and MinGW's. The one from Cygwin is easier to use (easier installation) but has the disadvantage that the resulting software does not run natively on Windows, various system calls go through Cygwin's emulation layer which slows things down. Here I'll show how to compile Open Babel with MinGW.

Installing MinGW

I've previously found this a bit confusing. This time I did a manual installation by creating a folder C:\MinGW, and then downloading all the relevant dlls on the installation page. To do this quickly just middle click on several links, wait a few seconds, and then hit Save on all the dialog boxes. Once they are all downloaded, move them to C:\MinGW and unzip them there.

Installing MSYS

No need to install MSYS (a kind of build environment for MinGW) for a project such as Open Babel that uses CMake to build. Why do I mention it then? Because the MinGW page talks all about it.

Compiling Open Babel

1. Add C:\MinGW\bin to the PATH
2. Get Cygwin's stuff off the PATH (if it's there). This is most easily accomplished by renaming C:\Cygwin to C:\oldCygwin or so.
3. Configure CMake to create makefiles for MinGW. I had some problems (at runtime) with a shared library version, so I went with the static one:

cmake -G "MinGW Makefiles" ../openbabel-2.3.2 -DWITH_INCHI=FALSE -DBUILD_SHARED=FALSE

4. Build it with MinGW's make.

mingw32-make

Hmmm...I wonder if it's as fast as the MSVC-compiled version we distribute?

You can QSAR that again - Reproducible research with IPython

2013-01-24T21:51:00.003Z

I've mentioned the IPython Notebook before (here and here). It's an interactive Python session that runs in the web browser, and can capture and display the output including plots. It can be saved, loaded and exported to a static HTML page. Entries in the notebook can be edited, and the whole notebook can be run in order to regenerate the output.

In other words, it's the perfect tool for documenting and presenting an analysis of data, thus bringing us one step closer to the goal of reproducible research. There is one area in which it is a particularly good fit for cheminformatics, and that's QSAR.

Greg Landrum and Nikolas Fechner of Novartis have led the way here. Check out this series of IPython notebooks originally presented at the RDKit UGM in 2012, and in particular the one on Using SciKit-Learn and Descriptors to Build Regression Models. Here's an excerpt:

It's pretty much a complete record of how they went about analysing a particular dataset from start to finish. The only thing that I would add is that I would ask the software used (RDKit, ipython, matplotlib and scikits-learn) to print out their version numbers of the top of the notebook (and add some pretty pictures of outliers too of course).

Hopefully others will follow in these footsteps. It would certainly be something to see such a Notebook included as part of the Methods section in a QSAR paper. Almost makes me want to do some QSAR work again...(almost). :-)

Chemistrify your Raspberry Pi Part III

2013-01-19T18:45:00.001Z

Following on from Parts I and II, now for the chemistry bit. It turns out this is the easiest part:

apt-get install python-cinfony python-imaging python-imaging-tk openbabel openbabel-gui indigo-utils python-chemfp python-cclib gausssum pymol jmol rasmol avogadro

Pretty easy huh? A single line install for 13 or so chemistry packages. Note that this install command should work on any other Linux distribution based on Debian (e.g. Linux Mint or Ubuntu).

Specifically, this installs Cinfony 1.1 and all its dependencies (Open Babel, RDKit, CDK, Indigo, OPSIN). Then there's Andrew Dalke's chemfp. Not to mention the 'mols' (Jmol, PyMol, Rasmol) and Avogadro. And let's not forget shameless self-promotion of cclib and GaussSum.

After installation here are some examples of things you could do:

$ obabel -:"CC(=O)Cl" -O testOB.png
$ indigo-depict - "CC(=O)Cl" testIndigo.png
$ gpicview test*.png # Display images
$
$ obgui   # Runs fine
$
$ obabel -:"CC(=O)Cl MyMol" -O tmp.mol --gen2d
$ ob2fps tmp.mol # Run ChemFP
$
$ wget http://www.rcsb.org/pdb/files/1PTQ.pdb
$ pymol 1PTQ.pdb  # Fails to start (no OpenGL GLX extension on RPi)
$ rasmol 1PTQ.pdb # Runs fine
$ jmol 1PTQ.pdb   # Runs fine
$ avogadro        # Fails to start (OpenGL problem)
$
$ cclib-get --list mycompchemfile.log
$ gausssum

To use Cinfony, you need to set some variables first as the Java parts don't work out-of-the-box:

$ export JPYPE_JVM=/usr/lib/jvm/java-6-openjdk-armhf/jre/lib/arm/server/libjvm.so
$ export CLASSPATH=/usr/share/java/cdk-nonotify.jar:/usr/share/java/cdk-io.jar:/usr/share/java/cdk-formula.jar:/usr/share/java/cdk-forcefield.jar:/usr/share/java/cdk-atomtype.jar:/usr/share/java/cdk-pdb.jar:/usr/share/java/cdk-fingerprint.jar:/usr/share/java/cdk-qsar.jar:/usr/share/java/cdk-ionpot.jar:/usr/share/java/cdk-annotation.jar:/usr/share/java/cdk-builder3d.jar:/usr/share/java/cdk-libiocml.jar:/usr/share/java/cdk-libiomd.jar:/usr/share/java/cdk-pcore.jar:/usr/share/java/cdk-ioformats.jar:/usr/share/java/cdk-qsarmolecular.jar:/usr/share/java/cdk-qsaratomic.jar:/usr/share/java/cdk-valencycheck.jar:/usr/share/java/cdk-extra.jar:/usr/share/java/cdk-structgen.jar:/usr/share/java/cdk-dict.jar:/usr/share/java/cdk-smarts.jar:/usr/share/java/cdk-control.jar:/usr/share/java/cdk-render.jar:/usr/share/java/cdk-builder3dtools.jar:/usr/share/java/cdk-qsarprotein.jar:/usr/share/java/cdk-data.jar:/usr/share/java/cdk-charges.jar:/usr/share/java/cdk-qm.jar:/usr/share/java/cdk-qsarionpot.jar:/usr/share/java/cdk-standard.jar:/usr/share/java/cdk-interfaces.jar:/usr/share/java/cdk-core.jar:/usr/share/java/cdk-sdg.jar:/usr/share/java/cdk-isomorphism.jar:/usr/share/java/cdk-qsarbond.jar:/usr/share/java/cdk-reaction.jar:/usr/share/java/cdk-diff.jar:/usr/share/java/cdk-smiles.jar:/usr/share/java/jaxen.jar:/usr/share/java/opsin-1.2.0.jar:/usr/share/java/opsin.jar
$ python
>>> from cinfony import opsin, webel
>>> webel.readstring("name", "aspirin").write("iupac")
'2-acetyloxybenzoic acid'
>>> opsin.readstring("iupac", "2-acetyloxybenzoic acid").write("smi")
'C(C)(=O)OC1=C(C(=O)O)C=CC=C1'

There are many other packages of interest; see under Science category in synaptic (see Notes below). Some examples include autodock, ballview, bkchem, and kalzium. Or to max out on chemistry just install the package science-chemistry.

Notes:
(1) If using apt-get to install software is too hard-core for you, there's also a GUI called synaptic. To install, use "apt-get install synaptic".
(2) After installation, to actually see what has been installed, use "dpkg -L 'package-name'". For example, anything that was installed in /usr/bin is a new command.
(3) The list of CDK jars was created using the following Python script:

import glob

cdk = glob.glob("/usr/share/java/cdk-*.jar")
jar = [jar for jar in cdk if not jar.endswith("1.2.10.jar")]
print ":".join(jar)

(4) Compiling Open Babel oneself works fine but takes 3 or so hours. A similar experience for RDKit has been reported by Jan Holst Jensen.

Compiling RDKit with MSVC 2012

2013-01-15T09:32:00.000Z

Compiling RDKit is a bit like the recipe for Elephant Soup. It's straightforward, but first we have to compile Boost (as there are no binaries provided for MSVC 2012). Unfortunately, the boost build instructions are very poor. The HTML instructions are full of text, none of which will simply tell you how to get the job done.

Preparation
1. Just to be safe, as both Boost and RDKit compile against Python, I deleted all my Python install folders except C:\Python2.7.
2. Make sure that bison and flex are installed in Cygwin, and that they are on the PATH.

Compiling Boost 1.49
1. Choose the right version. Too new, and the API will have changed and RDKit will not compile; too old, and it won't compile with MSVC 2012. I'm using boost 1.49.
2. Unzip into C:\Boost\boost_1_49_0. Do not bother using a different folder as it will install into C:\Boost\lib in any case.
3. Compile bjam as follows:

cd C:\Boost\boost_1_49_0
bootstrap

4. Start the MSVC2012 command prompt (or else bjam won't find 'cl'). Now we're going to compile the bits of Boost that RDKit needs. Some of these are shared libraries and some are dynamically linked libraries.

bjam.exe --with-regex --with-python --with-date_time --with-thread link=shared toolset=msvc-11.0 release install -j4
bjam.exe --with-thread --with-date_time toolset=msvc-11.0 release stage -j4

5. Copy the files from C:\Boost\boost_1_49_0\stage\lib to C:\Boost\lib

Compiling RDKit Q3 2012
1. The setup is all done by one command:

C:\Tools\RDKit\newbuild>cmake -G "Visual Studio 11" ..\RDKit_2012_09_1
-- Check for working C compiler using: Visual Studio 11
-- Check for working C compiler using: Visual Studio 11 -- works
-- Detecting C compiler ABI info
-- Detecting C compiler ABI info - done
-- Check for working CXX compiler using: Visual Studio 11
-- Check for working CXX compiler using: Visual Studio 11 -- works
-- Detecting CXX compiler ABI info
-- Detecting CXX compiler ABI info - done
-- Check if the system is big endian
-- Searching 16 bit integer
-- Looking for sys/types.h
-- Looking for sys/types.h - found
-- Looking for stdint.h
-- Looking for stdint.h - found
-- Looking for stddef.h
-- Looking for stddef.h - found
-- Check size of unsigned short
-- Check size of unsigned short - done
-- Using unsigned short
-- Check if the system is big endian - little endian
-- Found PythonLibs: C:/Python27/libs/python27.lib (found version "2.7.3")
-- Found PythonInterp: C:/Python27/python.exe (found version "2.7.3")
-- Boost version: 1.49.0
-- Found the following Boost libraries:
--   python
-- Found BISON: C:/cygwin/bin/bison.exe
-- Found FLEX: C:/cygwin/bin/flex.exe
-- Looking for include file pthread.h
-- Looking for include file pthread.h - not found.
-- Found Threads: TRUE
-- Boost version: 1.49.0
-- Found the following Boost libraries:
--   regex
-- Configuring done
-- Generating done
-- Build files have been written to: C:/Tools/RDKit/newbuild

2. Type "start RDKit.sln", change to Release build, and build the ALL_BUILD target, followed by the INSTALL target.

Running the RDKit tests
1. Close Visual Studio, and at the command line type:

set RDBASE=C:\Tools\RDKit\RDKit_2012_09_1
set PYTHONPATH=%RDBASE%
set PATH=%RDBASE%\lib;C:\Boost\lib;%PATH%
start RDKit.sln

2. Now you can run the tests by 'building' the RUN_TESTS target

1>  100% tests passed, 0 tests failed out of 76
1>  
1>  Total Test time (real) =  82.36 sec

Notes: For a debug build, you need the debug build of Boost. Just replace release by debug in the bjam command-lines above (to speed things up, use 'stage' for both).

IPython notebook and animated FiPy simulations

2013-01-05T22:29:00.001Z

I was back home this Christmas, and met up with a friend, Johan Hjelm. Naturally the conversation turned to the awesomeness of Python, and whether there was any way to create animated FiPy simulations directly in the IPython Notebook.

FiPy is "an object oriented, partial differential equation (PDE) solver, written in Python, based on a standard finite volume (FV) approach." Fair enough. More usefully, there are a couple of examples on the website that model diffusion, electrodeposition and convection. I focussed on the mesh20x20 diffusion example.

If you run the example at the command-line it pops up a matplotlib window showing the progress of the simulation. However, direct entry of the example into an IPython Notebook just results in a single graph for the simulation. To adapt it, I added a call to clear_output(), and used IPython's display() command to directly display the matplotlib figure associated with the simulation.

In short, here are the results as an IPython notebook, a Python script, and as an HTML page (and another HTML page, created on-the-fly from the notebook URL by http://nbviewer.ipython.org).

Notes: I used IPython 0.13.1 on Windows. To downgrade the notebook to earlier versions, see this discussion. Also, the statement "from IPython.display import clear_output" may need to be changed to "from IPython.core.display import clear_output".

Chemistrify your Raspberry Pi Part II

2013-01-02T21:50:00.000Z

In Part I, we got this very lean mean machine up and running. Now we want to take a look at what's going on in its tiny tiny silicon brain. If you have a monitor/TV and USB keyboard and mouse it's easy - just plug them in. In my case I don't so...

...I'm going to log in remotely using my laptop over the network. The good news is that there's an ssh server running by default on the RPi. The username is pi and password is raspberry. All we need is the RPi's IP address.

Connect the RPi to your router using an ethernet cable. If both are turned on, the router will assign the RPi an IP address. You can find out the value by logging into your router and looking at the details (or you can just guess the IP address by changing the number at the end of your laptop's IP address). Once you have the IP address, you can log in with Putty or Cygwin's ssh (remember, username pi).

The first time you log in, it asks you to run 'sudo raspi-config'. I did so, to set the timezone, expand the root filesystem (otherwise it doesn't use the whole SD card), and reduce the video memory to 32MB from 64MB (under "Memory split"). When you hit Finish it reboots, killing ssh, so you have to wait a minute before logging back in.

While some believe that the Unix command line is the perfect user interface, let's see what the Raspbian GUI looks like. To do so, we are going to use VNC (Virtual Network Computing), and specifically a piece of software called TightVNC. We will set up a server on the RPi, and a viewer on the laptop.

On the RPi:

$ sudo apt-get install tightvncserver
$ tightvncserver :1

If you set a password, make a note of it.

On the Windows laptop, install TightVNC. Please note that when you run the installer, you should untick the box that sets TightVNC running as a Windows service. This would be a BAD idea, as it would mean that your desktop is being broadcast over the network.

Now run "TightVNC Viewer" and connect to the RPi by entering the IP address of the RPi followed by ":5901", e.g. 128.128.0.1:5901. If you set a password, you will need to enter it. Finally, you should see something like this:

I still haven't done any chemistry but I guess that's all in Part III...

Notes: From time to time the router changes the IP address it allocates. If you want to assign a fixed IP address to the RPi, see the information here (untested). If you want the RPi to automatically start a TightVNC server on booting, see the information in the same article.

Chemistrify your Raspberry Pi

2013-01-01T20:47:00.001Z

So you've just gotten a Raspberry Pi Model B for Christmas. Yay - let's install all known cheminformatics software! But you don't have a monitor, TV or a USB keyboard. What are we going to do?

Well, it turns out we don't need them. We do need a couple of things though first:

1. A 5V power supply with a micro USB cable and supplying >700mA. £15 at Maplins - a bit overpriced. You'll probably do better on Amazon, or recycling a phone or Kindle charger.

2. An SDHC card. I bought a 32GB one at ASDA for £15, but you can get by with a 4GB one. Apparently they come in different speeds - I didn't realise this, but your Raspberry Pi might be faster if you get a class 10 SDHC (see list of compatible cards).

Getting the Linux distro image ("Raspbian" - a customised Debian) onto the SD card is the hardest part, at least from Windows, as the recommended method is using Win32DiskImager which (as of version 0.6) has a critical bug that might bork your hard-drive (YMMV). I created a Knoppix Live CD, booted off that, and used dd to copy the image to the SD card*. If you don't know what you're doing you might be better off buying an SD card with Raspbian pre-installed.

Your RPi is now good to go. Stick the SD card into the appropriate place, plug in the micro USB, and turn on the power. Eh walla!

Well, actually, the only thing you'll see are lights blinking on the RPi. But Part II will continue this epic narrative...

(* I used gparted to figure out which device was which.)

My year in review - Reviewing reviewing

2012-12-18T22:12:00.000Z

Like most scientists, from time to time I get asked to review papers, and a few years ago I decided to keep track of the reviews I did for various journals. In 2010, 2011, and 2012 I was asked to review 7, 6 and 7 times respectively.

I've been trying to figure out is this a reasonable level of reviewing? I guess the question is, am I doing more reviewing for the chemistry community than they are doing for me? Not that I mind those freeloading slackers dumping it all on me - I'm just curious.

So what's the other side of the equation? In the same three years I've had my name attached to 8 peer-reviewed publications. If each was reviewed 2.5 times (a reasonable guesstimate), that's also 20 reviews.

And so the balance is maintained.

Note: I know what you're thinking...(eerie isn't it!)...you want to adjust the figures for multiple authors. But both of the values would have to be adjusted in the same way and so it just cancels out. (I think.)

Update (02/01/2013): The previous note is a load of rubbish. As Felix points out in the comments, I should be correcting for multiple authors. Oh well...

Image credit: Laura Ritchie (LMRitchie on Flickr)

Cinfony 1.2 released

2012-12-11T09:53:00.001Z

Cinfony presents a common API to several cheminformatics toolkits. It uses the Python programming language, and builds on top of Open Babel, the RDKit, the CDK, Indigo, OPSIN, JChem and cheminformatics webservices.

Cinfony 1.2 is now available for download.

The two major additions in this release are support for the JChem commercial cheminformatics toolkit, and the ability to specify options (via an 'opt' dictionary) for format conversion and some other operations. There were also some under-the-hood changes to consolidate source files for ease of maintenance.

These additions were principally contributed by Adrià Cereto Massagué (of Universitat Rovira i Virgili and DecoyFinder) who joined the Cinfony team (i.e. me) earlier this year.

As usual, Cinfony has been updated to use the latest stable releases of each toolkit: Open Babel 2.3.2, CDK 1.4.15, RDKit 2012.09, Indigo 1.1, JChem 5.11 and OPSIN 1.3.

The Cinfony website has a somewhat condensed example showing the use of many of these resources in a dozen lines of Python. Here's a smaller example showing part of the new functionality:

>>> from cinfony import pybel, jchem
>>> mol = pybel.readstring("smi", "CC(=O)Cl")
>>> print mol.write("smi", opt={"f": 2, "l": 1}) # Make atom 2 the first and atom 1 the last
C(=O)(Cl)C
>>> fp = jchem.Molecule(mol).calcfp("ECFP")
>>> fp.bits
[39, 47, 55, 246, 397, 429, 700, 908]

To support Cinfony, please cite:

N.M. O'Boyle, G.R. Hutchison, Chem. Cent. J., 2008, 2, 24. [link]

Intro to Open Babel

2012-12-04T11:27:00.000Z

Recently I was asked (for the first time!) to provide introductory training for Open Babel. Here are the slides I put together:
If you're interested in this, you should follow up with the hands-on tutorial in the docs which I've mentioned previously. For further cheminformatics teaching material, see this post.

What's taking so long? - Profiling Open Babel

2012-11-20T16:54:00.000Z

Profiling code shows where all the runtime is spent. In the case of Open Babel, profiling is a bit awkward due to its use of dynamically-loaded libraries (the format plugins, etc.). So here's how you do it on Linux...

[openbabel/build]$ rm CMakeCache.txt 
[openbabel/build]$ CXXFLAGS="-pg" LDFLAGS="-pg" cmake ../trunk -DCMAKE_INSTALL_PREFIX=../tree -DCMAKE_BUILD_TYPE=DEBUG -DBUILD_SHARED=OFF
[openbabel/build]$ make

This should successfully compile the library and plugins with profiling information, but will fail when it comes to linking one of the executables:

[openbabel/build]$ make obabel
[100%] Built target openbabel
Scanning dependencies of target obabel
[100%] Building CXX object tools/CMakeFiles/obabel.dir/obabel.o
Linking CXX executable ../bin/obabel
/usr/bin/ld: dynamic STT_GNU_IFUNC symbol `strcmp' with pointer equality in `/usr/lib/../lib64/libc.a(strcmp.o)' can not be used when making an executable; recompile with -fPIE and relink with -pie
collect2: error: ld returned 1 exit status
make[3]: *** [bin/obabel] Error 1
make[2]: *** [tools/CMakeFiles/obabel.dir/all] Error 2
make[1]: *** [tools/CMakeFiles/obabel.dir/rule] Error 2
make: *** [obabel] Error 2

Yikes! Using VERBOSE=1 we can see the offending command:

[openbabel/build]$ VERBOSE=1 make obabel 
cd /home/noel/Tools/openbabel/profile/build/tools && /usr/local/bin/cmake -E cmake_link_script CMakeFiles/obabel.dir/link.txt --verbose=1
/usr/local/bin/c++   -static -pg  -g -g3 -fno-inline   -pg CMakeFiles/obabel.dir/obabel.o  -o ../bin/obabel -rdynamic ../src/libopenbabel.a -Wl,-Bstatic -lpthread -Wl,-Bdynamic -lm -lz -Wl,-Bstatic 
/usr/bin/ld: dynamic STT_GNU_IFUNC symbol `strcmp' with pointer equality in `/usr/lib/../lib64/libc.a(strcmp.o)' can not be used when making an executable; recompile with -fPIE and relink with -pie
collect2: error: ld returned 1 exit status

With Roger's help, I was able to change this to something simpler which will compile:

[build/tools]$ /usr/local/bin/c++ -pg  -g -g3 -fno-inline   -pg CMakeFiles/obabel.dir/obabel.o  -o ../bin/obabel ../src/libopenbabel.a -lpthread -lm -lz

Success is mine AT LAST!! Now let's profile it:

[build/bin]$ export BABEL_DATADIR=wherever
[build/bin]$./obabel bigfile.smi -onul
[build/bin]$ gprof ./obabel > gprof.out

Now time to read the gprof manual.

We can design molecular wires for you wholesale Part II

2012-11-20T10:01:00.000Z

As described previously, last year I published a paper in J. Phys. Chem. C with Geoff Hutchison on Computational Design and Selection of Optimal Organic Photovoltaic Materials.

A week or two ago someone emailed me for a copy of the paper. Fortunately, after a one year embargo, if you have an open access mandate you can request permission from the editor of an ACS journal to deposit the PDF in an institutional repository.

I went through this process a little while ago, and I am pleased to say that a copy of the PDF can now be found in University College Cork's institutional repository at http://hdl.handle.net/10468/748.

I'm not sure how people are going to find this though - it doesn't seem to be very prominent in Google. In particular, the PDF is not indexed (google search with filetype:pdf).

Tricks with SMILES and SMARTS Part II

2012-11-09T11:42:00.000Z

A much-underused feature of SMILES is the ability to apply 'atom classes' to atoms using a colon (inside square brackets). So, for example, CC and C[CH3:6] both represent ethane, but in the latter case one of the carbons is labelled as being a member of atom class 6.

So what's the meaning of atom class 6? Well, it's whatever you want - it's simply a label that you use to indicate some related information. For example, you might want to record reaction locations, or locations of common substitutions, or mappings between different molecules (reactant/product, or sub/superstructures).

Anyhoo, here's how you access the atom class information in Open Babel from Python:

>>> import pybel
>>> ob = pybel.ob
>>> mol = pybel.readstring("smi", "C[CH3:6]")
>>> print mol.write("smi")
CC

>>> print mol.write("smi", opt={"a":True})
C[CH3:6]

>>> data = ob.toAtomClassData(mol.OBMol.GetData("Atom Class"))
>>> data.HasClass(1)
False
>>> data.GetClassString(1)
''
>>> data.HasClass(2)
True
>>> data.GetClassString(2)
':6'
>>>

Plotting accesses on the axis Part III

2012-11-08T16:44:00.000Z

Following on from Parts I and II, this is the last in a series of posts exploring how journal access statistics provided by Open Access journals such as Journal of Cheminformatics give an insight into relative impact of papers.

It's now over a year since the Blue Obelisk and Open Babel papers were published so let's look at the accesses over that period: ...close to straight lines with a defined slope.

And here's an update of the view over the first month, this time including the Universal SMILES paper: ...it seems like the accesses to this paper mirror almost exactly those of the Blue Obelisk paper.

In conclusion, it would be nice if journals provided these sorts of graphs, or if some third-party website (e.g. one of altmetrics ones) did it. All of the data is on the website; it just needs to be collated as I've done here.