使用BLASTp计算蛋白序列相似性时匹配率始终为100%的问题

核心问题 #

你的代码完全偏离了需求：你要计算的是列表中两两序列之间的相似性，但当前逻辑是把其中一个序列拿去BLAST nr数据库，得到的是该序列与数据库中同源序列的比对结果，根本没用到配对的另一个序列query_sequence2。返回的100%匹配率，其实是query和数据库里几乎完全相同的序列的比对结果，和你要的两两序列相似性毫无关系。

正确解决方案：用双序列比对直接计算 #

用Biopython的pairwise2模块可以直接完成两个序列的比对，高效得到你要的匹配率，代码如下：

1. 导入依赖模块 #

import itertools
from Bio import pairwise2
from Bio.Align import substitution_matrices

2. 修正后的完整代码 #

# filtered_list是你的蛋白质序列列表，元素可以是字符串或Seq对象
pairs = itertools.combinations(filtered_list, 2)

for seq1, seq2 in pairs:
    # 用BLOSUM62矩阵做蛋白质全局比对，参数符合常规蛋白比对标准
    matrix = substitution_matrices.load("BLOSUM62")
    # gap_open=-10：打开缺口的罚分；gap_extend=-0.5：延伸缺口的罚分
    alignments = pairwise2.align.globalds(seq1, seq2, matrix, -10, -0.5)
    
    # 取得分最高的最优比对结果
    best_align = alignments[0]
    aligned_seq1, aligned_seq2 = best_align[0], best_align[1]
    
    # 统计有效匹配数（排除缺口'-'的位置）
    match_count = 0
    for a, b in zip(aligned_seq1, aligned_seq2):
        if a == b and a != '-':
            match_count += 1
    
    # 计算匹配率：匹配数除以两个原始序列的最小长度
    min_len = min(len(seq1), len(seq2))
    match_ratio = match_count / min_len
    
    # 输出结果
    print(f"****序列对比对结果****")
    print(f"序列1长度: {len(seq1)}, 序列2长度: {len(seq2)}")
    print(f"匹配残基数: {match_count}, 最小序列长度: {min_len}")
    print(f"匹配率: {match_ratio:.2f}\n")

补充说明 #

• 如果你不需要精准的蛋白矩阵，也可以用简单的全局比对方法pairwise2.align.globalxx(seq1, seq2)，但对于蛋白质序列，用BLOSUM矩阵的结果更可靠。
• 原代码里的BLAST调用逻辑只适合搜索数据库找同源序列，不适合做两两序列的相似性计算，完全没必要调用NCBI的远程BLAST服务，既慢又达不到你的需求。

内容的提问来源于stack exchange，提问作者y.a