python - 学习如何使用 Bash 和/或 Biopython 管理表格数据(来自 BLAST+6 格式)的资源

Question

我正在运行 BLAST 并希望使用 BLAST+6 格式操作输出。

例如，我想获取每个命中的 E 值、查询覆盖率和身份，然后将它们插入一个等式中，该等式将所有三个加权成一个“分数”。然后我想获取所有分数并将它们放入表格中，以便我可以通过减少“分数”来对每个命中进行排序。

我还想为数据库上的每个 BLAST 命中生成一个 ORF，并将它们添加到表中的相应位置。

有人可以指出我可以搜索的任何资源/关键字以了解如何操作表格数据吗？

例子:

blastn  -query genes.fasta  -subject genome.fasta  -outfmt "6 qseqid pident qcovs  evalue"

输出:

    qseqid pident qcovs evalue   
0   moaC     100.00    0.0       161.0      
1   moaC     99.38     1.0       161.0  

我想从每一列中获取值并将它们用作等式中的变量，然后在相应的行中打印该值。我将使用 bash 脚本或 BioPython 进行 BLAST，所以我想将数据操作作为其中的一部分

而不是这个例子的解决方案，我想看看是否有一个资源可以让我了解这个主题(到目前为止，我会使用电子表格程序来操作表格数据)

Answer 1

对于处理表格数据，我真的推荐 pandas .

首先，您要将输出转换为 pandas DataFrame ，这是一种数据结构，非常适合存储表格形式的数据。

对于这个例子，我使用 tblastn和示例文件 four_human_proteins.fasta 和 rhodopsin_nucs.fasta .

>>> import pandas as pd
>>> from Bio.Blast.Applications import NcbiblastnCommandline
>>> cline = NcbiblastnCommandline(cmd='/path/to/BLAST+/2.8.1/bin/tblastn', 
                          query='four_human_proteins.fasta', 
                          subject='rhodopsin_nucs.fasta', 
                          evalue='1e-10',
                          outfmt='"6 qseqid pident qcovs evalue"')
>>> print(cline)
/path/to/BLAST+/2.8.1/bin/tblastn -outfmt "6 qseqid pident qcovs evalue" -query four_human_proteins.fasta -evalue 1e-10 -subject rhodopsin_nucs.fasta

>>> blast_output = cline()[0].strip()
>>> print(blast_output) 
sp|P08100|OPSD_HUMAN    96.552  100 0.0
sp|P08100|OPSD_HUMAN    93.391  100 0.0
sp|P08100|OPSD_HUMAN    95.092  94  0.0
sp|P08100|OPSD_HUMAN    84.795  98  0.0
sp|P08100|OPSD_HUMAN    82.164  98  0.0
sp|P08100|OPSD_HUMAN    96.396  89  2.65e-67
sp|P08100|OPSD_HUMAN    92.308  89  7.50e-36
sp|P08100|OPSD_HUMAN    93.220  89  1.81e-32
sp|P08100|OPSD_HUMAN    96.296  89  6.37e-32
sp|P08100|OPSD_HUMAN    88.462  89  4.64e-12


>>> headers = ['qseqid', 'pident', 'qcovs', 'evalue']
>>> rows = [line.split() for line in blast_output.splitlines()]    
>>> df = pd.DataFrame(rows, columns=headers)
>>> print(df)
    qseqid  pident  qcovs   evalue
0   sp|P08100|OPSD_HUMAN    96.552  100     0.0
1   sp|P08100|OPSD_HUMAN    93.391  100     0.0
2   sp|P08100|OPSD_HUMAN    95.092  94  0.0
3   sp|P08100|OPSD_HUMAN    84.795  98  0.0
4   sp|P08100|OPSD_HUMAN    82.164  98  0.0
5   sp|P08100|OPSD_HUMAN    96.396  89  2.65e-67
6   sp|P08100|OPSD_HUMAN    92.308  89  7.50e-36
7   sp|P08100|OPSD_HUMAN    93.220  89  1.81e-32
8   sp|P08100|OPSD_HUMAN    96.296  89  6.37e-32
9   sp|P08100|OPSD_HUMAN    88.462  89  4.64e-12

首先我们需要告诉 pandas哪些列包含 float s。

>>> convert = {'pident': float, 
              'qcovs': float, 
              'evalue': float,
              'qseqid': str}
>>> df = df.astype(convert) 

现在可以轻松地对这个 DataFrame 执行列操作 df .
定义您的评分函数并将结果添加为额外的列:

>>> df['score'] = df['qcovs'] / df['pident']  # adapt to your own needs
>>> print(df)
    qseqid  pident  qcovs   evalue  score
0   sp|P08100|OPSD_HUMAN    96.552  100.0   0.000000e+00    1.035711
1   sp|P08100|OPSD_HUMAN    93.391  100.0   0.000000e+00    1.070767
2   sp|P08100|OPSD_HUMAN    95.092  94.0    0.000000e+00    0.988516
3   sp|P08100|OPSD_HUMAN    84.795  98.0    0.000000e+00    1.155729
4   sp|P08100|OPSD_HUMAN    82.164  98.0    0.000000e+00    1.192736
5   sp|P08100|OPSD_HUMAN    96.396  89.0    2.650000e-67    0.923275
6   sp|P08100|OPSD_HUMAN    92.308  89.0    7.500000e-36    0.964163
7   sp|P08100|OPSD_HUMAN    93.220  89.0    1.810000e-32    0.954731
8   sp|P08100|OPSD_HUMAN    96.296  89.0    6.370000e-32    0.924234
9   sp|P08100|OPSD_HUMAN    88.462  89.0    4.640000e-12    1.006082

您可以通过此 score 轻松地对 DataFrame 进行排序柱子

>>> df.sort_values(['score'], inplace=True)
>>> print(df)
    qseqid  pident  qcovs   evalue  score
5   sp|P08100|OPSD_HUMAN    96.396  89.0    2.650000e-67    0.923275
8   sp|P08100|OPSD_HUMAN    96.296  89.0    6.370000e-32    0.924234
7   sp|P08100|OPSD_HUMAN    93.220  89.0    1.810000e-32    0.954731
6   sp|P08100|OPSD_HUMAN    92.308  89.0    7.500000e-36    0.964163
2   sp|P08100|OPSD_HUMAN    95.092  94.0    0.000000e+00    0.988516
9   sp|P08100|OPSD_HUMAN    88.462  89.0    4.640000e-12    1.006082
0   sp|P08100|OPSD_HUMAN    96.552  100.0   0.000000e+00    1.035711
1   sp|P08100|OPSD_HUMAN    93.391  100.0   0.000000e+00    1.070767
3   sp|P08100|OPSD_HUMAN    84.795  98.0    0.000000e+00    1.155729
4   sp|P08100|OPSD_HUMAN    82.164  98.0    0.000000e+00    1.192736