8

我有一个 FASTA 文件，可以很容易地被SeqIO.parse.

我对提取序列 ID 和序列长度感兴趣。我用这些线来做，但我觉得它太重了（两次迭代，转换等）

from Bio import SeqIO
import pandas as pd


# parse sequence fasta file
identifiers = [seq_record.id for seq_record in SeqIO.parse("sequence.fasta",
                                                           "fasta")]
lengths = [len(seq_record.seq) for seq_record in SeqIO.parse("sequence.fasta",
                                                             "fasta")]
#converting lists to pandas Series    
s1 = Series(identifiers, name='ID')
s2 = Series(lengths, name='length')
#Gathering Series into a pandas DataFrame and rename index as ID column
Qfasta = DataFrame(dict(ID=s1, length=s2)).set_index(['ID'])

我只需要一次迭代就可以做到，但我得到了一个 dict ：

records = SeqIO.parse(fastaFile, 'fasta')

我不知何故无法DataFrame.from_dict上班...

我的目标是迭代 FASTA 文件，并DataFrame通过每次迭代将 id 和序列长度放入 a 中。

这是一个简短的 FASTA 文件，供那些想要提供帮助的人使用。

4

2 回答 2

7

你是正确的 - 你绝对不应该对文件进行两次解析，并且将数据存储在字典中会浪费计算资源，因为稍后你只是将其转换为numpy数组。

SeqIO.parse()返回一个生成器，因此您可以逐个记录迭代，构建一个列表，如下所示：

with open('sequences.fasta') as fasta_file:  # Will close handle cleanly
    identifiers = []
    lengths = []
    for seq_record in SeqIO.parse(fasta_file, 'fasta'):  # (generator)
        identifiers.append(seq_record.id)
        lengths.append(len(seq_record.seq))

有关从 FASTA 文件中仅解析 ID 和序列的更有效方法，请参阅Peter Cock 的答案。

您的其余代码对我来说看起来不错。但是，如果您真的想优化与 一起使用pandas，您可以阅读以下内容：

关于最小化内存使用

查阅的来源panda.Series，我们可以看到它data内部存储为numpy ndarray：

class Series(np.ndarray, Picklable, Groupable):
    """Generic indexed series (time series or otherwise) object.

    Parameters
    ----------
    data:  array-like
        Underlying values of Series, preferably as numpy ndarray

如果你创建identifiers一个ndarray，它可以直接使用Series而无需构造一个新数组（参数copy，默认值False）将防止在不需要时创建一个新数组ndarray。通过将序列存储在列表中，您将强制 Series 将所述列表强制为ndarray.

避免初始化列表

如果你事先知道你有多少序列（以及最长的 ID 有多长），你可以初始化一个空ndarray来保存标识符，如下所示：

num_seqs = 50
max_id_len = 60
numpy.empty((num_seqs, 1), dtype='S{:d}'.format(max_id_len))

当然，很难确切地知道您将拥有多少个序列，或者最大的 ID 是多少，因此最简单的方法是让numpy从现有列表进行转换。但是，从技术上讲，这是存储您的数据以供在pandas.

于 2013-10-17T21:15:05.253 回答
6

大卫在pandas侧面给了你一个很好的答案，在 Biopython 方面，如果你想要的只是记录标识符和它们的序列长度，你不需要使用SeqRecord对象Bio.SeqIO- 这应该更快：

from Bio.SeqIO.FastaIO import SimpleFastaParser
with open('sequences.fasta') as fasta_file:  # Will close handle cleanly
    identifiers = []
    lengths = []
    for title, sequence in SimpleFastaParser(fasta_file):
        identifiers.append(title.split(None, 1)[0])  # First word is ID
        lengths.append(len(sequence))
于 2013-10-18T15:08:31.317 回答