Skip to content

Commit f83f4f4

Browse files
committed
added min-readcount option for bamrc2af
1 parent 1b46003 commit f83f4f4

File tree

1 file changed

+6
-2
lines changed

1 file changed

+6
-2
lines changed

hometools/hometools.py

+6-2
Original file line numberDiff line numberDiff line change
@@ -1699,6 +1699,7 @@ def run_bam_readcount(tmp_com, outfile):
16991699
o = p.communicate()
17001700
# if o[1] != b'Minimum mapping quality is set to 40':
17011701
# sys.exit("Error in running bam-readcount:\n{}".format(o[1].decode()))
1702+
return
17021703
# END
17031704

17041705

@@ -1783,6 +1784,7 @@ def bamrc2af(args):
17831784
rcfin = args.bamrc.name
17841785
vcffin = args.vcf.name
17851786
outfin = 'bamrc_af.txt' if args.out is None else args.out.name
1787+
minrc = args.min_rc
17861788

17871789
logger.info('Reading VCF')
17881790
posdict = dict()
@@ -1802,7 +1804,8 @@ def bamrc2af(args):
18021804
with open(rcfin, 'r') as rc, open(outfin, 'w') as out:
18031805
for line in rc:
18041806
line = line.strip().split()
1805-
if line[3] == '0': continue
1807+
# if line[3] == '0': continue
1808+
if int(line[3]) < min_rc: continue
18061809
try:
18071810
ref, alt = posdict[(line[0], line[1])]
18081811
except KeyError:
@@ -1811,7 +1814,7 @@ def bamrc2af(args):
18111814
refi = basedict[ref]
18121815
alti = basedict[alt]
18131816
out.write(f'{line[0]}\t{line[1]}\t{ref}\t{alt}\t{round(int(line[refi])/int(line[3]) , 2)}\t{round(int(line[alti])/int(line[3]), 2)}\n')
1814-
logger.info('Finishe')
1817+
logger.info('Finished')
18151818
# END
18161819

18171820

@@ -2661,6 +2664,7 @@ def main(cmd):
26612664
parser_bamrc2af.add_argument("bamrc", help="BAM readcount file generated using bamrc", type=argparse.FileType('r'))
26622665
parser_bamrc2af.add_argument("vcf", help="VCF file", type=argparse.FileType('r'))
26632666
parser_bamrc2af.add_argument("out", help="Output file", type=argparse.FileType('w'))
2667+
parser_bamrc2af.add_argument("--min_rc", help="Minimum required read count. Position with lower number of reads would be filtered out", type=int, default=1)
26642668

26652669
# xls2csv
26662670
parser_xls2tsv.set_defaults(func=xls2csv)

0 commit comments

Comments
 (0)