'None' as key in kwargs when trying to merge two tmp BedTools

[jorgsk]

I try to merge two BedTools objects but internally one key in kwargs is 'None' which is not accepted.

import pybedtools

annotated = pybedtools.BedTool('file1.bed')
discovered = pybedtools.BedTool('file2.bed')

hg19genome = pybedtools.chromsizes('hg19')

ann_slop = annotated.slop(g=hg19genome, b=10)
disc_slop = discovered.slop(g=hg19genome, b=10)

mergd = ann_slop.merge(disc_slop, s=True, nms=True)

Gives the error:

/usr/local/lib/python2.7/dist-packages/pybedtools-0.5.5-py2.7-linux-x86_64.egg/pybedtools/bedtool.pyc in wrapped(self, *args, **kwargs)
171             # At runtime, this will parse the kwargs, convert streams to

172             # tempfiles if needed, and return all the goodies

--> 173             cmds, tmp, stdin = self.handle_kwargs(prog=prog, **kwargs)
174 
175             # Do the actual call

TypeError: handle_kwargs() keywords must be strings

And kwargs is

{None: <BedTool(/tmp/pybedtools.Iz9EUg.tmp)>, 's': True, 'nms': True}

The slopped files are fine themselves:

ipdb> disc_slop.head(2)
chr21   40563340        40563360        AATAAA  8       -
chr21   40564162        40564182        AGTAAA  12      -
ipdb> ann_slop.head(2)
chr1    328449  328470  0       0       +
chr1    529943  529964  0       0       +

[daler]

This error is because merge() only operates on a single file (e.g., mergeBed -i infile from the command line). Currently the decorator that performs the wrapping does not have a useful error message. I'll have to think about how best to improve this. Also, for convenience you can pass an assembly name as the "genome" kwarg to save a step:

import pybedtools
annotated = pybedtools.BedTool('file1.bed')
discovered = pybedtools.BedTool('file2.bed')
ann_slop = annotated.slop(genome='hg19', b=10)
disc_slop = discovered.slop(genome='hg19', b=10)

# gives an error, because merge only operates on a single file.
# By default it's the file the BedTool points to (here, file1.bed)
#
# mergd = ann_slop.merge(disc_slop, s=True, nms=True)

# "merge" only merges one file
ann_merged = ann_slop.merge(s=True, nms=True)

# perhaps this is what you want?
in_common = ann_slop.intersect(disc_slop, s=True)

What is the output you'd like to get? You may need to use args like wo=True for the intersection depending on your needs.

[jorgsk]

I see now that I used mergeBed in an unintended way, thanks for the good explanation and the tip about 'hg19'. The output I was after was actually to first concatenate the two files and then merge them:

cat annotated merged | slopBed -i stdin -b 10 -g hg19 | mergeBed -i stdin | sortBed -i stdin > my_output.bed

The solution I ended up with was what you outlined. But since I wanted the entries unique to ann_slop and disc_slop too I added them manually with an intersectBed -v against the in_common file.

What would you would suggest for the pipeline above with only pybedtools? (avoiding the cat method since it truncates at BED3)

[daler]

I just committed a change b64ceb8 that checks to see if the input filenames are 1) the same type and 2) have the same number of fields. If so, then no truncation will happen.

In my use-cases, I typically would cat and then merge, so the BedTool.cat() method does a post-merge unless you tell it not to; additional kwargs are sent to the merge command after cat-ing. Also keep in mind that the latest BEDTools version requires that a file be sorted before merging, but the command in BedTool.cat takes care of that.

However, you're doing a slop call before merging. So something like this should do the same as the command you showed:

import pybedtools
results = pybedtools.BedTool('file1.bed')\
              .cat('file2.bed', post_merge=False)\
              .slop(genome='hg19', b=10)\
              .sort()\
              .merge()

[jorgsk]

Wonderful! Thanks. I'm closing the issue with this.