-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathfilterSitesInBed.py
53 lines (43 loc) · 1.64 KB
/
filterSitesInBed.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
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
# -*- coding: utf-8 -*-
"""
Created on Tue Oct 16 15:12:22 2018
@author: YudongCai
@Email: [email protected]
"""
import os
import gzip
import click
import numpy as np
import pandas as pd
def count_header(gzvcf):
with gzip.open(gzvcf, 'rb') as f:
for nline, line in enumerate(f):
if line.decode()[0] != '#':
return nline
@click.command()
@click.option('--sitesfile', help='Chrom\\tPos', default=None)
@click.option('--gzvcf', help='bgzip vcffile', default=None)
@click.option('--bedfile', help='Chrom\\tStart\\tEnd')
@click.option('--outsites', help='output file')
def main(sitesfile, gzvcf, bedfile, outsites):
"""
output sites (in sitesfile or gzvcffile) within the bedfile
"""
try:
os.remove(outsites)
except Exception:
pass
if gzvcf:
nheader = count_header(gzvcf)
sdf = pd.read_csv(gzvcf, skiprows=nheader, header=None, names=['chrom','pos'],
usecols=[0,1], sep='\t', dtype={'chrom': str, 'pos': int})
elif sitesfile:
sdf = pd.read_csv(sitesfile, sep='\t', header=None, names=['chrom', 'pos'],
dtype={'chrom': str, 'pos': int})
bdf = pd.read_csv(bedfile, sep='\t', header=None, names=['chrom', 'start', 'end'],
dtype={'chrom': str, 'start': int, 'end': int})
for chrom, start, end in bdf[['chrom', 'start', 'end']].values:
sdf.loc[(sdf['chrom']==chrom) & (sdf['pos']>=start) & (sdf['pos']<=end), :]\
.to_csv(outsites, mode='a', header=False, sep='\t', index=False)
if __name__ == '__main__':
main()