-
Notifications
You must be signed in to change notification settings - Fork 15
/
reduce_bam.py
executable file
·35 lines (23 loc) · 1022 Bytes
/
reduce_bam.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
#!/usr/bin/python
from subprocess import call
import sys
print "Usage: reduce_bam.py BamFile"
try:
file = sys.argv[1]
except:
file = raw_input("Introduce name of your bam file: ")
name = file.split(".")
name = name[0]
print "\nReducing file %s..." % (file)
#to get all the reads where both mapped.
call("samtools view -b -F 12 %s > temp1.bam" % file, shell=True)
#to get all the reads that did not map, but whose mate mapped
call("samtools view -b -f 4 -F 8 %s > temp2.bam" % file, shell=True)
#to get all the reads that mapped, but whose mates did not.
call("samtools view -b -f 8 -F 4 %s > temp3.bam" % file, shell=True)
#to merge bam files and sort by name
#call("samtools merge -u - temp[123].bam | samtools sort - %s_mapped" % (name), shell=True)
call("samtools merge -u - temp[123].bam | samtools sort -T aln.sorted - -o %s_mapped.bam" % (name), shell=True)
#to index sorted bam file
call("samtools index %s_mapped.bam" % (name), shell=True)
call("rm temp1.bam temp2.bam temp3.bam", shell=True)