0

我将深入解释我的编程问题：我有两个文件；文件 #1 是基因注释文件，文件 #2 是按碱基位置文件计数的文件（只是试图为问题提供一些上下文）。

我想在第 6 列中有“+”的行中提取“start_codon”位置，然后转到文件#2 中的那个位置。例如，我想从文件 #1 中的第 3 列中提取 954，然后转到文件 #2 中的第 954 行。然后，我想计算文件 #2 中第 954 行上方产生 70 或更大计数值的行数。

File#1

    Chromosome  exon    337 774 0.0 -   .   gene_id "A";    
    Chromosome  start_codon 954 956 0.0 +   0   gene_id "B"; 
    Chromosome  stop_codon  2502    2504    0.0 +   0   gene_id "B"; 

 File#2
 .      .
 .      . 
 942    71
 943    63
 944    88
 945    80
 946    80
 947    85
 948    86
 949    97
 950    97
 951    97
 952    104
 953    105
 954    104
 955    108

我的最终输出文件将是gene_id 的制表符分隔文件，后跟产生70 或更大计数值的行数。对于我给出的示例文件，输出如下：

 Gene_id  Count_before_start_codon
 B     10

我想遍历大文件以生成一个长输出文件。

谢谢，我希望这很清楚。我很感激任何指导！

4

2 回答 2

0

如果您想保持简单并假设文件 2 不是很大，您可以将文件 2 加载到一个 numpy 数组中，这样您就可以非常快速地访问任何位置。例如，如果您将其加载到数组 ar 中并希望在位置 p 之前搜索，您可以这样做：

numpy.sum(ar[:p]>=70)

这会给你你正在寻找的号码

然后，您可以通过文件 #1 并即时进行计算。这样你只需要读取每个文件一次，它应该很快。

于 2013-04-09T01:01:19.577 回答
0

这应该可以工作... 第一部分获取文件 1 中的基因信息并填充字典 第二部分打开文件 2，检查字典并生成输出。

D={}
with open("file1.txt","rU") as f1:
    for line in f1:
        line=line.rstrip().rsplit("\t")
        if line[6]=="+" and line[2]=="start_codon":
            D[line[3]] = line[8].rstrip('"')[9]
            keys = D.keys()

count=[]
results=[]
number = 12
with open("file2.txt","rU") as f2:
    for line in f2:
        line=line.rstrip().rsplit("\t")
        if int(line[1]) >= 70:
            count.append(line[1])
            if line[0] in D:
                results.append(D[line[0]])
                if len(count) > number:
                    results.append(str(number))
                else:
                    results.append(str(len(count)-1))
                print "\t".join(results)
                count=[]

        else:
            count=[]

附言。我复制粘贴了你的例子。我将文件编辑为制表符分隔。所以你可能需要玩“切片”

于 2013-04-09T16:08:06.820 回答