将fasta序列解析为字典

Question

我需要一个最简单的办法，把一个叫做fasta.txt的文件里的多个核苷酸序列转换成一个字典对象，字典里包含名称和对应的值。

>seq1
TAGATTCTGAGTTATCTCTTGCATTAGCAGGTCATCCTGGTCAAACCGCTACTGTTCCGG
CTTTCTGATAATTGATAGCATACGCTGCGAACCCACGGAAGGGGGTCGAGGACAGTGGTG
>seq2
TCCCTCTAGAGGCTCTTTACCGTGATGCTACATCTTACAGGTATTTCTGAGGCTCTTTCA
AACAGGTGCGCGTGAACAACAACCCACGGCAAACGAGTACAGTGTGTACGCCTGAGAGTA
>seq3
GGTTCCGCTCTAAGCCTCTAACTCCCGCACAGGGAAGAGATGTCGATTAACTTGCGCCCA
TAGAGCTCTGCGCGTGCGTCGAAGGCTCTTTTCGCGATATCTGTGTGGTCTCACTTTGGT

这里的名称就是以>开头的标题，而值则是对应的序列。

下面是我尝试用两个列表来实现这个功能的代码，但这个方法对包含多行的长序列不太管用。

f = open('input2.txt', 'r')
list={}
names=[]
seq=[]
for line in f:
 if line.startswith('>'):
  names.append(line[1:-1])
 elif line.startswith('A') or line.startswith('C') or line.startswith('G') or line.startswith('T'):
  seq.append(line)

list = dict(zip(names, seq))

如果你能告诉我怎么修复这个问题，并给我一个单独函数的示例，我会非常感激。

谢谢你的帮助，

Gleb

Answer 1

对你的代码做个简单的修正：

from collections import defaultdict #this will make your life simpler
f = open('input2.txt','r')
list=defaultdict(str)
name = ''
for line in f:
    #if your line starts with a > then it is the name of the following sequence
    if line.startswith('>'):
        name = line[1:-1]
        continue #this means skips to the next line
    #This code is only executed if it is a sequence of bases and not a name.
    list[name]+=line.strip()

更新：

因为我收到通知说这个旧答案被点赞了，所以我决定分享一下我现在认为的正确解决方案，使用的是Python 3.7。要转换成Python 2.7，只需要去掉类型导入那一行和函数注释：

from collections import OrderedDict
from typing import Dict

NAME_SYMBOL = '>'


def parse_sequences(filename: str,
                    ordered: bool=False) -> Dict[str, str]:
    """
    Parses a text file of genome sequences into a dictionary.
    Arguments:
      filename: str - The name of the file containing the genome info.
      ordered: bool - Set this to True if you want the result to be ordered.
    """
    result = OrderedDict() if ordered else {}

    last_name = None
    with open(filename) as sequences:
        for line in sequences:
            if line.startswith(NAME_SYMBOL):
                last_name = line[1:-1]
                result[last_name] = []
            else:
                result[last_name].append(line[:-1])

    for name in result:
        result[name] = ''.join(result[name])

    return result

现在，我明白提问者想要的是“最简单的解决方案”，不过因为他们在处理基因组数据，所以可以合理地假设每个序列可能会非常大。在这种情况下，稍微优化一下是有意义的，可以先把序列行收集到一个列表中，然后在最后使用str.join方法把这些列表合并成最终结果。

Answer 2

使用biopython库会更好

from Bio import SeqIO
input_file = open("input.fasta")
my_dict = SeqIO.to_dict(SeqIO.parse(input_file, "fasta"))