使用熊猫或枚举查找两个数据集之间的重叠/范围
问题描述:
我想对两个文件执行一些区间范围操作,并且出现以下条件 检查chrom是否相等,然后检查我的co0rdinatefile的开始和结束是否在相等或相等开始和结束gene_annotation文件(条件如果链是“+”,开始和结束将是例如10-20,如果它的“ - ”将是20-10),如果匹配打印开始末尾链从坐标和来自geneannotation文件的gene_id,gene_name。 (用于表示目的,我有头annoataion文件)使用熊猫或枚举查找两个数据集之间的重叠/范围
中的行数注释文件〜50000 数协调的文件中的行〜200000
gene_annotationfile
chrom start end gene_id gene_name strand
17 71223692 71274336 ENSMUSG00000085299 Gm16627 -
17 18186448 18211184 ENSMUSG00000067978 Vmn2r-ps113 +
11 84645863 84684319 ENSMUSG00000020530 Ggnbp2 -
7 51097639 51106551 ENSMUSG00000074155 Klk5 +
13 31711037 31712238 ENSMUSG00000087276 Gm11378 +
coordinates_file
chrom start end strand
1 4247322 4247912 -
1 4427449 4432604 +
1 4763414 4764404 -
1 4764597 4767606 -
1 4764597 4766491 -
1 4766882 4767606 -
1 4767729 4772649 -
1 4767729 4768829 -
1 4767729 4775654 -
1 4772382 4772649 -
1 4772814 4774032 -
1 4772814 4774159 -
1 4772814 4775654 -
1 4772814 4774032 +
1 4774186 4775654 -
1 4774186 4775654
1 4774186 4775699 -
想要的输出
chrom, start, end,strand, gene_id, gene_name
1 4427432 4432686 + ENSMUSG0001 abcd
另一个问题是,在某些情况下,如果有匹配它可以映射到在这种情况下,gene_id我还想写
chrom, start, end,strand, gene_id, gene_name
1 4427432 4432686 + ENSMUSG0001,ENSMUSG0002 abcd,efgh
我迄今为止代码:
import csv
with open('coordinates.txt', 'r') as source:
coordinates = list(csv.reader(source, delimiter="\t"))
with open('/gene_annotations.txt', 'rU') as source:
#if i do not use 'rU' i get this error Error: new-line character seen in unquoted field - do you need to open the file in universal-newline mode?
annotations = list(csv.reader(source, delimiter="\t"))
for index,line in enumerate(coordinates):
for index2, line2 in enumerate(annotations):
if coordinates[line][0] == annotations[line2][0] and coordinates[line][1] <= annotations[line2][1] and annotations[line2][2] >= coordinates[line][2] :
print "%s\t%s\t%s\t%s\t%s" % (coordinates[line][0],coordinates[line][1],coordinates[line][2], annotations[line2][3], annotations[line2][4])
break
错误我得到
---> 15 if coordinates[line][0] == annotations[line2][0] and coordinates[line][1] <= annotations[line2][1] and annotations[line2][2] >= coordinates[line][2] :
16 print "%s\t%s\t%s\t%s\t%s" % (coordinates[line][0],coordinates[line][1],coordinates[line][2], annotations[line2][3], annotations[line2][4])
17 break
TypeError: list indices must be integers, not list
大熊猫会是一个很好的方法吗?
答
我假设坐标是列表的列表,如[[1,2],[3,4]]。 线
for index,line in enumerate(coordinates):
迭代坐标,返回坐标线的各条线和索引为索引。
if coordinates[line][0] == annotations[line2][0] and coordinates[line][1] <= annotations[line2][1] and annotations[line2][2] >= coordinates[line][2] :
错误消息表示您在此处使用列表(行)作为索引。你可能想改用线的指标:
if coordinates[index][0] == annotations[index2][0] and coordinates[index][1] <= annotations[index2][1] and annotations[index2][2] >= coordinates[index][2] :
更妙的是只使用线:
if line[0] == line2[0] and line[1] <= line2[1] and line2[2] >= line[2] :
看到https://docs.python.org/2.7/reference/compound_stmts.html?highlight=for_stmt#grammar-token-for_stmt