fast python with shedskin

There's a new release of the shedskin compiler. It is able to generate fast shared libraries that can be run from CPython. It can also create binaries so I thought I'd see how it did on some code from this BMC bioinformatcs article compared to psyco and CPython.
I took this iterative, brute-force (Needleman-Wunsch?) alignment code and modified it slightly. That's pasted here. (Notice the first line! that's how it appears in the original code). The modifications allow shedskin to infer the function and variable types. Plus, there's a couple changes I made that improve the run-time for all cases. The max() function is also in the original, but unnecessary because of python's builtin max(), however, pysco does run much faster using their hand-coded max(). For the shedskin run, I removed that extra code and used shedskin's builtin 'cause it made me feel better.

The python code was run as
$ time python -c "import alignment; alignment.imain()"

To run with psyco. The only change from the pasted script is to add
'import psyco; psyco.full()' at the top.

To run with shedskin built executable (after removing the max()):
$ shedskin alignment.py ### generates Makefile
$ make
$ time ./alignment

To run with shedskin as a shared (.so) library
$ shedskin -e -n alignment.py # generates Makefile
$ make
$ time python -c "import alignment; print alignment.imain()"
# python finds and imports the shared library (.so) before the .py file.

Timing:
time reported is real.
Python 2.5.1: 19.030s
Psyco: 1.1336s
Shedskin shared: 0.921
Shedskin binary: 0.818

CPython is much slower at 19seconds compared to ~1 for pysco and shedskin. The output format for the shedskin binary is slightly different because it calls the stuff in __main__, but they all generate the same alignments. It's interesting to look at the alignment.ss.py that shedskin generates, as you can see all the inferred types. The alignment.cpp contains the generated cpp code, which is also quite readable. It's also nice to be able to get a binary executable without jumping through any extra hoops.

Shedskin was very easy to use and faster than psyco for this case, I just pulled from SVN and it worked out of the box. It now has support for sets and regular expressions, and seems quite active. I can see using shedskin for purely brute force stuff where numpy is no help and I might otherwise have to resort to cython. I kinda like cython, but it's nice just to be able to get fast code with little to no modifications.


EDIT:
jython 2.3a0 runs this unaltered in 42 seconds.

Comments

Anonymous said…
More comments, and more benchmarking in D:

http://leonardo-m.livejournal.com/57439.html
7:51 AM
Post a Comment

Popular posts from this blog

filtering paired end reads (high throughput sequencing)

NOTE: I don't recommend using this code. It is not supported and currently does not work for some sets of reads. If you use it, be prepared to fix it. I wrote last time about a pipeline for high-throughput sequence data. In it, I mentioned that the fastx toolkit works well for filtering but does not handle paired end reads. The problem is that you can filter each end (file) of reads independently, but most aligners expect that the nth record in file 1 will be the pair of the nth record in file 2. That may not be the case if one end of the pair is completely removed while the other remains. At the end of this post is the code for a simple python script that clips an adaptor sequences and trims low-quality bases from paired end reads. It simply calls fastx toolkit (which is assumed to be on your path). It uses fastx_clipper if an adaptor sequence is specified and then pipes the output to fastq_quality_trimmer for each file then loops through the filtered output and keeps onl
Read more

python interval tree

EDIT: added a couple points inline. I'm obsessed with trees lately -- of the CS variety, not the plant variety. Although we are studying poplar, so I'll be using trees to study trees . I'd tried a couple times to implement an interval tree from scratch following the Wikipedia entry . Today I more or less did that in python. It's the simplest possible form. There's no insertion (though that's trivial to add), it just takes a list of 'things' with start and stop attributes and creates a tree with a .find() method. The wikipedia entry that baffled me was about storing 2 copies of each node's intervals --one sorted by start and the other by stop. I didn't do that as I think in most cases it won't improve lookup time. I figure if you have 1 million elements and a tree of depth 16, then you have on average 15 intervals per node (actually fewer since there are the non-leaf nodes). So I just brute force each of those nodes and move to the next. I th
Read more

needleman-wunsch global sequence alignment -- updates and optimizations to nwalign

I've written previously about nwalign , a python package I wrote to do fast (thanks to cython) global sequence alignment. Over break, I spent some time fixing some bugs and improving performance. Ack It's actually nice to get bug reports for a relatively obscure bit of software like this as it shows it's getting used. Thanks especially to R. Christen for his patience in (repeatedly) showing me places where nwalign was not doing the right thing. Bugs Placement of Gaps Some of the "Bugs" were actually just "unspecified behavior". For example, given input text of "AGEBAMAM" and "AGEBAM", the alignments: AGEBAMAM AGEBAM-- and AGEBAMAM AGEB--AM Have the same score. However, the first version is generally more, um, polite. Previously, the behavior was deterministic, but depended on the length and order of the input sequence. As of version 0.3, nwalign will always put matches together when given a choice between 2 (or more) paths through t
Read more