Splitting a BAM based on # of reads?

[kga1978]

Hi All,

I'm surprised I can't find anything on this topic, but I would like to take a BAM file of, say 100million reads, and split it into multiple files containing 1million reads each. Is there an easy way of doing that? (I thought BamTools, samtools, GATK, Picard, but none of them seem to be able to do this - these tools can only split based on RG, chromosomes, etc. - I wan't to do it based on number of reads).

Thanks very much in advance.

[adaptivegenome]

Yes. Do a view with samtools and use the unix command "split" to split by number of lines and then pipe this output to a new file. You will get a set of output files that are in SAM format. Next run them back through samtools (view -b) and pipe them to a file with a BAM extension. Now you have BAM files split by # of reads.

[maubp]

And reinsert the SAM header into each BAM file at the end if required using 'samtools reheader'

[adaptivegenome]

And reinsert the SAM header into each BAM file at the end if required using 'samtools reheader'

Yes, sorry for omitting this!

[kga1978]

Ah, yes - that's elegant. Thanks very much for the response - so easy nobody bothered making a script .

[CHRYSES]

And reinsert the SAM header into each BAM file at the end if required using 'samtools reheader'

I had to do:
samtools view -Sb -T reference.fasta out.sam > out.bam

to each split, to get the header...

I don't think that the "samtools reheader" would have worked as you describe it in this thread, as the splits themselves do not have any headers.

[carmeyeii]

Is there any way to split by number of reads, and not lines?

Such that each read, even if it has 10 alignments, is contained in a single new BAM file, along with other reads?