参考资料,了解如何使用Bash和/或Biopython管理表格数据(来自BLAST+6格式)

2024-06-25 05:18:51 发布

您现在位置:Python中文网/ 问答频道 /正文

我正在运行BLAST,希望使用BLAST+6格式处理输出

例如,我想获取每个命中率的E值、查询覆盖率和标识,然后将它们插入一个等式,将所有三者加权为一个“分数”。然后我想把所有的分数放在一个表中,这样我就可以通过减少“分数”来对每个命中进行排序

我还想为数据库中的每个爆炸命中生成一个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  

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

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


Tags: 数据格式覆盖率资源分数fasta表格blast
1条回答
网友
1楼 · 发布于 2024-06-25 05:18:51

对于使用表格数据,我真的建议使用pandas

首先,您需要将输出转换为^{},这是一种非常适合存储以表格形式出现的数据的数据结构

对于这个示例,我使用了tblastn和示例文件^{}^{}

>>> 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

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

现在可以轻松地在此数据帧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列轻松地对数据帧进行排序

>>> 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

相关问题 更多 >