如何用Python实现重叠肽段高亮匹配及序列覆盖度计算

问题原因 #

原有方案的问题出在两点：

1. 直接使用re.sub做替换时，正则匹配到一段序列后会跳过已匹配的字符，因此重叠的肽段无法被识别
2. 没有预先统计所有匹配的字符位置，无法计算覆盖度

完整实现代码 #

通用逻辑说明 #

采用位掩码标记法解决重叠匹配和覆盖度统计问题：

• 初始化一个和长序列长度一致的布尔数组，默认值为False
• 遍历所有肽段，逐个找出该肽段在长序列中所有的匹配位置（支持重叠），将匹配到的位置对应的布尔数组设为True
• 最终统计布尔数组中True的总数就是高亮字符数，直接计算覆盖度即可
• 高亮时根据布尔数组的状态切换添加高亮标签，避免标签嵌套混乱

终端运行版本（ANSI颜色高亮） #

import re
from docx import Document

# 配置参数
PEPTIDE_LIST = ['MKWVTFISLLLLFSSAYSRGV', 'SSAYSRGVFRRDTHKSEIAH', 'KPKATEEQLKTVMENFVAFVDKCCA']
FILE_PATH = 'BSA.docx'  # 支持docx或txt格式
HIGHLIGHT_START = '\033[93;43m'  # 黄色背景高亮，可自行修改
HIGHLIGHT_END = '\033[0m'

def read_long_sequence(file_path):
    """读取长序列，跳过首行标识"""
    if file_path.endswith('.docx'):
        doc = Document(file_path)
        # 从第二段开始拼接序列，跳过首行标识
        full_text = '\n'.join([p.text for p in doc.paragraphs[1:]])
    elif file_path.endswith('.txt'):
        with open(file_path, 'r', encoding='utf-8') as f:
            lines = f.readlines()
            full_text = ''.join(lines[1:])
    # 去掉序列中的空白字符
    return re.sub(r'\s+', '', full_text)

def get_match_mask(long_seq, peptides):
    """获取匹配位置的掩码数组"""
    seq_len = len(long_seq)
    mask = [False] * seq_len
    for pep in peptides:
        pep_len = len(pep)
        if pep_len == 0 or pep_len > seq_len:
            continue
        # 遍历所有可能的起始位置，找所有匹配（支持重叠）
        for i in range(seq_len - pep_len + 1):
            if long_seq[i:i+pep_len] == pep:
                # 标记匹配到的所有位置
                for j in range(i, i+pep_len):
                    mask[j] = True
    return mask

def highlight_sequence(long_seq, mask):
    """根据掩码生成高亮后的序列"""
    highlighted = []
    in_highlight = False
    for char, is_match in zip(long_seq, mask):
        if is_match and not in_highlight:
            highlighted.append(HIGHLIGHT_START)
            in_highlight = True
        elif not is_match and in_highlight:
            highlighted.append(HIGHLIGHT_END)
            in_highlight = False
        highlighted.append(char)
    # 收尾：如果最后一个字符是高亮状态，加结束标签
    if in_highlight:
        highlighted.append(HIGHLIGHT_END)
    return ''.join(highlighted)

if __name__ == '__main__':
    # 读取长序列
    long_seq = read_long_sequence(FILE_PATH)
    total_len = len(long_seq)
    if total_len == 0:
        print("未读取到有效序列")
        exit()
    # 获取匹配掩码
    match_mask = get_match_mask(long_seq, PEPTIDE_LIST)
    # 计算覆盖度
    highlight_count = sum(match_mask)
    coverage = (highlight_count / total_len) * 100
    # 输出结果
    print("高亮后的序列：")
    print(highlight_sequence(long_seq, match_mask))
    print(f"\n序列总长度：{total_len}")
    print(f"高亮字符数：{highlight_count}")
    print(f"序列覆盖度：{coverage:.2f}%")

Jupyter Notebook运行版本（HTML高亮） #

如果在Jupyter中运行，希望更美观的高亮效果，可以替换高亮函数部分：

from IPython.display import display, HTML

HIGHLIGHT_START = '<span style="background-color: yellow; font-weight: bold;">'
HIGHLIGHT_END = '</span>'

# 输出时替换为HTML显示
display(HTML(highlight_sequence(long_seq, match_mask)))

可选：保存高亮后的docx文件 #

如果需要把高亮结果保存为新的docx文件，可以用python-docx的高亮属性实现：

from docx.shared import Pt
from docx.enum.text import WD_COLOR_INDEX

def save_highlighted_docx(long_seq, mask, output_path='高亮结果.docx'):
    doc = Document()
    p = doc.add_paragraph()
    in_highlight = False
    current_text = ''
    for char, is_match in zip(long_seq, mask):
        if is_match == in_highlight:
            current_text += char
        else:
            run = p.add_run(current_text)
            if in_highlight:
                run.font.highlight_color = WD_COLOR_INDEX.YELLOW
            current_text = char
            in_highlight = is_match
    # 处理最后一段文本
    if current_text:
        run = p.add_run(current_text)
        if in_highlight:
            run.font.highlight_color = WD_COLOR_INDEX.YELLOW
    doc.save(output_path)

效果说明 #

• 支持任意数量的重叠肽段匹配，不会漏标
• 覆盖度统计自动去重，重叠区域的字符只会被计算一次
• 高亮样式可自行修改配置参数调整

内容的提问来源于stack exchange，提问作者RuMa