用python匹配fasta文件中的头文件

Question

我有两个文件：第一个是带有标题和序列的 fasta 文件，第二个仅由标题组成。

文件_1：

>DF94KKQ1|265|D0M1LACXX|3|2103|4637|10742|1|N|0|TGACCA
TTCCAAAGAAACATGGAAGACCCAGGACTTGGAGGCACCAGGCACCAGCACACAGGGGTA
GGCACATGGCATGGTGTTGGTTGAAGTCTACTTTTCCCACC
>DF94KKQ1|265|D0M1LACXX|3|2103|4565|10742|1|N|0|TGACCA
TTAATTTTTTCAGGCAAGTTTTGTGGATTTCAGTGTGTAAGTCTTTCACCTCTTTGGTTA
AATTTATTCCTATGTATTTTATTCCTTTAGATGCTATTATG
>DF94KKQ1|265|D0M1LACXX|3|2103|4565|10742|2|N|0|TGACCA
TTAATTTTTTCAGGCAAGTTTTGTGGATTTCAGTGTGTAAGTCTTTCACCTCTTTGGTTA
AATTTATTCCTATGTATTTTATTCCTTTAGATGCTATTATG

文件_2：

>DF94KKQ1|265|D0M1LACXX|3|2103|4565|10742|2|N|0|TGACCA
>DF94KKQ1|265|D0M1LACXX|3|2103|4668|10746|1|N|0|TGACCA
>DF94KKQ1|265|D0M1LACXX|3|2103|4668|10746|2|N|0|TGACCA
>DF94KKQ1|265|D0M1LACXX|1|2207|10852|3331|2|N|0|TGACCA

我想将 File_2 中的标头与 File_1 中直到第 7 个“|”之前具有相同确切字符的任何内容相匹配。

我拆分了 File_1 中的项目（标题的每个部分都被索引到一个列表中）。任何以 '>' 开头的行都被放入一个变量中：

#!/usr/bin/env python
import sys
from Bio import SeqIO

#Function, split header line into a list
def getHeaderInfo(blastLine):
   myFields = blastLine.strip("\n").split("|") 
   HeaderInfo = myFields[:6]
   return HeaderInfo

input_file = sys.argv[1]

#Get input file from the command line
inFileName = sys.argv[1]

#open the input file
inFileHandle = open(inFileName)

#loop over the input file line by line
for thisLine in inFileHandle.readlines():
    if thisLine [0] == '>': 
       print getHeaderInfo(thisLine)
       HeaderInfo = getHeaderInfo(thisLine)

我一直在尝试找到一种方法，在该方法中我可以比较 File_2 中的这些相同索引以返回以下输出：

>DF94KKQ1|265|D0M1LACXX|3|2103|4565|10742|1|N|0|TGACCA
TTAATTTTTTCAGGCAAGTTTTGTGGATTTCAGTGTGTAAGTCTTTCACCTCTTTGGTTA
AATTTATTCCTATGTATTTTATTCCTTTAGATGCTATTATG
>DF94KKQ1|265|D0M1LACXX|3|2103|4565|10742|2|N|0|TGACCA
TTAATTTTTTCAGGCAAGTTTTGTGGATTTCAGTGTGTAAGTCTTTCACCTCTTTGGTTA
AATTTATTCCTATGTATTTTATTCCTTTAGATGCTATTATG

我尝试过的几种方法都使用索引，但是，我的键不是唯一的。如何获取前六个元素并将它们作为我的关键，或者有没有比我正在尝试的当前方法更好的方法？谢谢。

Answer 1

这是你想要的吗？

def make_key(line):
    return "|".join(line.split("|", 7)[ : 7]) + "|"

header_set = set()
with open("file_2.txt") as in_f:
    for line in in_f:
        header_set.add(make_key(line))

with open("file_1.txt") as in_f, open("file_3.txt", "w") as out_f:
    accept = False
    for line in in_f:
        if line.startswith(">"):
            key = make_key(line)
            accept = key in header_set

        if accept:
            out_f.write(line)

Answer 2

这种方法涉及内存映射 file_1.txt，通过它运行 re.finditer 以将行加载到由 header 键入的 defaultdict 中

import re
import mmap
from collections import defaultdict
pat = re.compile('>.+?(?=(>|$))', re.DOTALL)
with open('file_1.txt') as f:
    map = mmap.mmap(f.fileno(), 0, access=mmap.ACCESS_READ)

line_1 = defaultdict(list)
for line in pat.finditer(map):
    fields = line.group().split('|')
    key = '|'.join(fields[:6])
    line_1[key].append(line.group())

map.close()
line_2 = []
with open('file_2.txt') as f:
    for line in f:
      fields = line.split('|')
      key = '|'.join(fields[:6])
      line_2.append(key)

line_2 = set(line_2)
for key in line_2.intersection(line_1.keys()):
    print "".join(line_1[key])