Convert paired-end BAM into a single-end BAM and keep all the reads I want to call peaks using MACS2 on my paired-end ATAC-seq BAM and I am interested in cutting sites. So I make use of the -BAM option (opposed to the -BAMPE) option. However this throws away all the mates, and all the reads of which the mate did not align (that's more than 50%!!!). My idea was to align paired-end, and then remove all the information from the BAM that it is paired-end. So my question is how do I:

Convert paired-end BAM into a single-end BAM and keep all the reads

I made a start based on the SAM-format specification , but when I run the result of this script through markduplicates I get the error:

ERROR: Record 1, Read name SRR6245286.60336908., MRNM should not be set for unpaired read.


I can not find what this MRNM is...

This is my code I used to generate the faulty 'single-end' bam, but I am open to suggestions (e.g. bash, preferably no R). 

import pysam

paired =         1
proper_pair =    2
mate_unmapped =  8
mate_reverse =   32
first_in_pair =  64
second_in_pair = 128

bam_in = pysam.AlignmentFile("/vol/atacpipeline/bwa/GSM2837489-danRer7.samtools-coordinate.bam", "rb")
bam_out = pysam.AlignmentFile("/vol/atacpipeline/bwa/GSM2837489-danRer7_SE.samtools-coordinate.bam", "wb", template=bam_in)
for line in bam_in:
    if line.flag &amp; second_in_pair:
        line.pos += line.template_length
    line.next_reference_id = 0
    line.next_reference_start = 0
    line.flag &amp;= ~(paired + proper_pair + mate_unmapped + mate_reverse + first_in_pair + second_in_pair)

    bam_out.write(line)
Like Devon pointed out, most likely you should sort out whether the files have been marked for duplicates correctly. 

You can also use samtools rmdup too

samtools rmdup &lt;bamfile&gt; &lt;output&gt;


Also, I find it very odd that 50% of the reads are lost when you keep only properly paired reads.

I would suggest removing the duplicates and then converting the bam file to a bed file. Macs2 can take in the bed format. Any extremely crude python script that converts a bam into a bed:

import pysam
samfile = pysam.AlignmentFile("test.bam", "rb")
open("test.bed", 'w').close()
f = open("test.bed", "a")
STRAND = ["+","-"]
for read in samfile.fetch():
    STR = STRAND[int(read.is_reverse)]
    BED =[read.reference_name,read.pos,read.reference_end,".",read.mapq,STR]
    f.write('\t'.join([str(i) for i in BED]))
    f.write('\n')
f.close()


You can also try the conversion with bedtools or definitely refine that script above with pandas etc.