来自 FASTA 序列的字典

Question

我正在尝试从 FASTA 序列创建字典。但是我一直收到同样的错误，说变量没有定义，即使它被定义了。有人可以帮我理解这段代码是否正确吗？

#!/usr/bin/pyhton
# Build a dictionary containing all sequences from a FASTA file from typing import List



try:
    f= open("FASTA.rtf",'r')
except IOError:
    print('The file "FASTA" does not exist')

seqs = {}
for line in f:
    # Let's discard the newline at the end (if any)
    line = line.rstrip()
    if line[0]==">": #or line.startswith('>')
        words=line.split()
        name= words[0][1:]
        seqs[name]=""
    else: #sequence, not header
        seqs[name]=seqs[name]+line

for name,seq in seqs.items():
    print(name,seq)

Answer 1

首先，您为什么使用 .rtf（富文本格式）的 FASTA 而不是 .fa/.fasta 或 .fq/.fastq？这可能会导致你的一些问题。下面是我用来编辑 FASTA 文件 in-place.

基因标签的代码片段

def remove_trailing_inplace(filename):

    with fileinput.input(files=str(filename), inplace=True, backup=".bak") as f
        for line in f:
            if line[0] == ">": # only engages with lines that have genes
                a = line.split(" ") # series of splits to properly modify each
                idx = a[0].split(".")[0]
                gene = a[3].split(".")[0]
                rest = " ".join(a[4:])
                print(idx + " " + a[1] + " " + a[2] + " " + gene + " " + rest.str)
            else:
                print(line, end="")

这段代码工作得很好，捕获了我需要的所有基因。最大的区别是它直接采用 FASTA 文件，而不是 fasta 格式的 rtf。另外，我不使用任何类型的剥离：据我所知，FASTA 格式不需要它。

此外，正如@Frank Yellin 评论的那样，您的 if/else 语句在逻辑上是不正确的，因为在 else 部分中，您正试图使用 name，它仅在if 部分。 if 和 else 语句是互斥的，因此您不能依赖其中一个声明的变量。

幸运的是修复不是太困难。我建议在循环开始时声明 name（在 if/else 之外）。

seqs = {}
names = []
for line in f:
    
    if line[0]==">":
        words=line.split()
        name= words[0][1:]
        names.append(name)
        seqs[name]=""
    else: #sequence, not header
        seqs[names[-1]]=seqs[name[-1]]+line

在乱码的情况下，应该避免名字错误。此外，因为 FASTA 就是它们的样子，我们可以假设我们总是以“>”开头。

编辑：为了进一步说明那些不熟悉 FASTA 格式的人，这里有一个维基百科的例子：

>SEQUENCE_1
SATVSEINSETDFVAKNDQFIALTKDTTAHIQSNSLQSVEELHSSTINGVKFEEYLKSQI
ATIGENLVVRRFATLKAGANGVVNGYIHTNGRVGVVIAAACDSAEVASKSRDLLRQICMH

其中带有 > 的行包含序列信息，而之后的所有行（直到下一个“>”）代表序列本身。因此，我们总是可以假设在完成第一个序列之前我们永远不会命中另一个序列。这就是上面代码有效的原因，也是使用 names[-1] 没问题的原因。

来自 FASTA 序列的字典

Dictionary from a FASTA sequence

python

bioinformatics