一、背景
- 有的时候，我们需要手动对结果文件(可能是vcf，xls，或者其他形式)进行bed文件的过滤，这时候如何高效准确的对posi是否在目标区域内，成为关键。
二、关键点
① bed是0-based，vcf等格式是1-based，注意转换时，start+1，end不变
def get_bed_dic(input_bed):
'''将不同染色体的region信息存为列表，整体存为字典'''
with open(input_bed) as ff:
dic={}
for line in ff:
tmps=line.strip().split("\t")
if tmps[0] not in dic:
dic[tmps[0]]=[]
dic[tmps[0]].append([int(tmps[1])+1,int(tmps[2])])
return dic
② 大家一定不要想着用循环遍历，太慢了。建议用二分查找判断目标位置是否在target region里，速度奥奥快，嗖嗖的。
def binary_search(target,bed_list):
'''用二分查找的方法，看看目标位置是不是落于bed区域里'''
target=int(target)
first=0
last=len(bed_list)-1
while first<=last:
mid = (first + last)//2
if target>=int(bed_list[mid][0]) and target<=int(bed_list[mid][1]):
return True
elif target < int(bed_list[mid][0]):
last = mid-1
else:
first= mid+1
return False
如果大家对原理感兴趣，给我留言吧，博主太懒了，就先不写原理了，代码见上啊。
