parallel-python and pysam

[MBender]

Hi next-gen sequencing community,

as you all know, the analysis of next-generation sequencing data can be very time consuming. Therefore, i'd like to integrate parallelization into my python scripts. What i use is the python-pp package which makes it very convenient to put independent tasks into different processes which can be run in parallel. On the other hand i use the pysam package which gives me the desired functionality to access, and manipulate sequencing data. Unfortunately, it turned out that python-pp and pysam are not the best friends, and thus trouble was bound to occur.

A small example what can go wrong:
If you want to sort a set bam files within a python script in parallel, the following code should be sufficient:

import pp
import pysam

def sortBam(bam_filepath, out_prefix):
  pysam.sort(bam_filepath, out_prefix)

job_server = pp.Server(ncpus=2, ppservers = (), secret='secret')

job1 = job_server.submit(sortBam, ('<in1.bam>', '<out_prefix1>'), modules=('pysam',))
job2 = job_server.submit(sortBam, ('<in2.bam>', '<out_prefix2>'), modules=('pysam',))

result1 = job1()
result1 = job2()

Unfortunately, i got the following error message:

An error has occured during the function execution
Traceback (most recent call last):
File "/home/mbieg/.local/lib/python2.6/site-packages/ppworker.py", line 90, in run
__result = __f(*__args)
File "<string>", line 2, in sortBam
File "/home/mbieg/project/tools/python_packages/lib/python2.6/site-packages/pysam-0.6-py2.6-linux-x86_64.egg/pysam/__init__.py", line 50, in __call__
retval, stderr, stdout = csamtools._samtools_dispatch( self.dispatch, args )
File "csamtools.pyx", line 2860, in csamtools._samtools_dispatch (pysam/csamtools.c:26160)
AttributeError: StringIO instance has no attribute 'fileno'

Does anybody have an idea how to fix this?

Best,

Matthias

[rflrob]

I'm not sure exactly what the issue is there (I'm assuming it works when you take the parallelization out), but you might try checking out GNU parallel: http://www.gnu.org/software/parallel/
From the command line, you could do:

$ parallel 'samtools sort {} {.}_sorted' ::: in1.bam in2.bam in3.bam ... inN.bam

It should also be possible to do this call from within a python script using, e.g., subprocess.Popen().

[MBender]

Hi rflrob,

thanks for your reply. First of all, the example i posted is just an example. I want to integrate the features of pysam into a larger software project with much more complicated analyses. Therefore, i need native parallelization in python, combined with the features of pysam. Second, the python packages thread, or threading are not suitable for me, since they are limited by the global interpreter lock of python, which means, even if i fork some threads there will be no real parallelization. The python-pp package circumvents the GIL and thus seems to be a good choice. The problem is, that it seems that python-pp has some problems with pysam (especially the samtools function calls like pysam.sort, or pysam.index, etc.).

I would be really glad, if somebody could tell me about his/her experience with pysam and python-pp.

Best,

Matthias

[rflrob]

Hi Matthias,

If I had to guess, I'd say what pp is doing is turning all file buffers into StringIO buffers, which is making pysam unhappy. Are you locked into pp? Using the multiprocessing module (which also gets around the GIL) seems to work fine for me:

import multiprocessing as mp
import pysam

def sortBam(args):
    pysam.sort(args[0], args[1])

pool = mp.Pool()
pool.map(sortBam, [('<in1.bam>', '<sorted_prefix1>'),
                   ('<in2.bam>', '<sorted_prefix2>')])

EDIT: Looking more at the code of pysam and pp, it actually looks like they're both competing for access to stdout or stderr. in ppworker.py, pp sets stdout and stderr to StringIOs, whereas I'd bet that pysam is expecting them to be bona fide filedescriptors.