搜索 python 中的 Fasta 文件，return 高效读取

Question

这里是初学者。我想在 python 中编写一个函数，在 fasta 文件中搜索基因名称，然后 return 相应的读取。

FASTA 文件示例：

>name1
AATTCCGG
>name2
ATCGATCG

到目前为止我的代码（非常简陋）：

def findseq(name):
    with open('cel39.fa', 'rb') as csv_file:
        csv_reader = csv.reader(csv_file)
        for i in csv_reader:
            if i == '>' + name:
                return i+1
                break

这实际上不起作用，因为我不能 return 'i+1'。我也可以迭代 len(csv_reader) 因为 'len' 不是一个属性。另外我不确定是否有更有效（但简单）的搜索系统，这样我就不需要每次都遍历整个文件（将是数千行）。

具体来说，有没有更好的读取Fasta文件的方法？有什么方法可以让我 return 阅读吗？

findseq(name1)

应该return'AATTCCGG'

谢谢！！

Answer 1

看看 python 库：Biothon

它包含大量有用的工具和方法。

这是他们解析 fasta 文件的示例：

from Bio import SeqIO
for seq_record in SeqIO.parse("ls_orchid.fasta", "fasta"):
    print(seq_record.id)
    print(repr(seq_record.seq))
    print(len(seq_record))

本例打印出fasta文件中的所有记录。

为了您的目的：

from Bio import SeqIO
for seq_record in SeqIO.parse("ls_orchid.fasta", "fasta"):
    if seq_record.id == name:
        return seq_record.seq

Answer 2

由于 FASTA 文件序列扩展到多行，您将不得不连接行直到找到 > 的下一个实例。下面的代码创建了一个以基因名称为键，以基因序列为值的字典。

with open('cel39.fa', 'rb') as fp:
    lines = fp.read().splitlines()

geneDict = {}

# Just to start populating the dictionary later
geneName = 'dummy'
fastaSeq = ''

for line in lines:
    if line[0] == '>':
        geneDict.update({geneName: fastaSeq})
        geneName = line[1:]
        fastaSeq = ''
    else:
        fastaSeq += line

geneDict.update({geneName: fastaSeq}) # Putting the values from the last loop
geneDict.pop('dummy') # Now removing the dummy

print geneDict['name1']
print geneDict['name2']

并打印：

AATTCCGG
ATCGATCG

搜索 python 中的 Fasta 文件，return 高效读取

Search Fasta file in python, return read efficiently

python

csv

fasta