1. 从“大海捞针”到“精准制导”：新抗原预测为何成为肿瘤免疫治疗的关键

如果你在肿瘤免疫治疗或者生物信息学领域待过一段时间，一定会对“新抗原”这个词不陌生。它听起来有点学术，但背后的逻辑其实很直接：我们的免疫系统就像一支精锐部队，而癌细胞是伪装潜入的叛军。这支精锐部队识别叛军的主要方式，就是看它们身上有没有特殊的“身份牌”——也就是抗原。大多数癌细胞和正常细胞长得太像，免疫系统容易“脸盲”。但有一类抗原，是癌细胞独有的，由基因突变产生，正常细胞里压根没有，这就是“新抗原”。找到它，就等于拿到了叛军的精准画像，可以指导免疫系统进行“定点清除”。

这就是新抗原预测的核心价值所在。它不是一个单一的软件或算法，而是一整套从原始数据到最终候选名单的计算工作流程。想象一下，你手头有一份肿瘤样本和一份正常样本的基因测序数据，里面是数十亿个碱基对的序列。新抗原预测工作流要做的，就是从这片数据的汪洋大海里，捞出那些由非同义突变产生的、能被MHC分子呈递的、并且能激活T细胞的短肽序列。这个过程，无异于一场精密的“分子考古”和“免疫学推演”。

我接触这个领域快十年了，从最早的手动拼接各种脚本，到如今相对成熟的流程化工具，踩过的坑不计其数。很多人以为有了现成的流程，点点鼠标就能出结果，但实际情况远非如此。一个可靠的新抗原预测，其计算工作流程至少涉及突变识别、HLA分型、肽段-MHC结合亲和力预测、免疫原性评估等多个关键环节，每个环节的工具选择、参数设置、结果解读都充满了学问。今天，我就结合自己的实战经验，把这套工作流程里里外外拆解一遍，不仅告诉你“怎么做”，更重点分享“为什么这么做”以及“哪些地方最容易翻车”。

2. 工作流程全景图：四大核心模块的协同作战

在深入每个模块之前，我们必须先建立起一个全局视野。一个新抗原预测的计算工作流程，本质上是一个多步骤的过滤和优先级排序管道。原始数据（通常是肿瘤和配对的正常组织的全外显子组或全基因组测序数据）经过这个管道，层层筛选，最终输出一个按可能性排序的新抗原候选列表。

整个流程可以清晰地划分为四个核心模块，它们环环相扣：

1. 体细胞突变检测模块：这是整个流程的基石。目标是从肿瘤测序数据中，准确地找出那些在正常组织中不存在的基因突变。这里的输出是一份突变列表，包括单核苷酸变异（SNV）、小的插入缺失（InDel）等。
2. HLA分型模块：免疫系统通过主要组织相容性复合体（MHC，在人类中称为HLA）分子来呈递抗原肽。不同个体的HLA基因型别（分型）千差万别，这直接决定了哪些肽段能被呈递。因此，我们必须从测序数据中推断出样本的HLA分型。
3. 肽段-MHC结合预测模块：对于突变位点，我们需要将其翻译成可能的多肽序列（通常是8-11个氨基酸长）。然后，利用计算模型预测这些肽段与患者特定HLA型别结合的强弱。结合亲和力高的肽段，才有机会被呈递到细胞表面。
4. 免疫原性预测与优先级排序模块：能结合不等于能激活免疫。最后一步，我们需要评估这些候选肽段触发T细胞免疫反应（即免疫原性）的潜力。同时，结合突变的功能重要性、表达水平等多维度信息，对候选新抗原进行综合打分和排序。

下面这个表格概括了每个模块的核心任务、常用工具及其输出，方便你快速建立认知框架：

注意：这个流程主要针对由基因编码区突变产生的新抗原。对于由基因融合、非编码区突变或转录后修饰产生的新抗原，需要额外或不同的分析模块，本文聚焦于最主流、最经典的流程。

3. 模块一：体细胞突变检测——一切分析的起点与最大误差源

突变检测的准确性直接决定了后续所有分析的可靠性。这里最大的误区是认为“用默认参数跑一遍流程就行了”。实际上，这是最容易引入系统性偏差的环节。

3.1 工具选型：没有银弹，只有组合拳

目前主流的三款工具是GATK Mutect2、Strelka2和VarScan2。它们基于不同的统计模型，各有优劣：

• Mutect2：基于贝叶斯模型，在肿瘤纯度较高、测序深度足够时表现稳健，对Indel检测较好，是许多流程的首选。
• Strelka2：采用基于等位基因频率的联合分型模型，在低频率突变检测上可能更敏感。
• VarScan2：基于启发式过滤和统计检验，速度较快，配置相对简单。

我的经验是，不要只依赖一个工具。在实际项目中，我通常会采用至少两种工具进行交叉验证。例如，用Mutect2进行主呼叫，然后用Strelka2的结果作为补充，只取两者交集或经过严格过滤的并集。这样可以大幅降低假阳性率。一个常见的组合是：Mutect2 -> 与Strelka2取交集 -> 使用bcftools filter或自定义脚本进行深度、频率、质量等过滤。

3.2 输入数据准备：BAM文件处理的魔鬼细节

突变检测工具的输入是比对后的BAM文件。这里有几个关键点常被忽略：

1. 标记重复序列：必须使用如Picard或GATK的MarkDuplicates工具处理。重复读取来自PCR扩增或光学重复，不标记会导致等位基因频率估算错误。
2. 碱基质量重校准：测序仪本身的系统误差会导致碱基质量分数不准。GATK的BaseRecalibrator和ApplyBQSR步骤可以校正这一点，这对依赖质量分数的突变检测器（如Mutect2）尤为重要。
3. 肿瘤与正常样本的协调性：确保肿瘤和正常样本的BAM文件使用了完全相同的参考基因组版本、比对工具和参数。一个细微的不一致都可能导致大量假阳性“突变”。

3.3 过滤策略：平衡敏感性与特异性的艺术

工具输出的原始VCF文件包含大量假阳性突变。过滤是门艺术，过严会丢失真实信号（假阴性），过松则后续分析负担沉重且不可靠。我通常会实施一个多层次过滤流水线：

1. 工具内置过滤：首先通过工具自带的过滤器（如Mutect2的FilterMutectCalls）。
2. 通用硬过滤：使用bcftools filter或SnpSift，设置阈值，例如：
   • QUAL > 20（变异质量）
   • DP > 10（总深度，肿瘤和正常均需考虑）
   • AF > 0.05（等位基因频率，在肿瘤中）
   • 在正常样本中AF < 0.02（确保是体细胞突变）
3. 数据库过滤：利用gnomAD等人群频率数据库，过滤掉在健康人群中频率较高的常见多态性位点（通常设定频率>0.1%或1%）。
4. 功能注释与筛选：使用ANNOVAR、SnpEff或VEP对突变进行注释。我们只关心编码区的非同义突变，所以会过滤掉同义突变、内含子区、UTR区等。这一步会淘汰掉约90%的突变，但它们是无效信号。

实操心得：过滤参数没有金标准。我建议先用一个已知突变的数据集（如细胞系混合样本）测试你的流程，根据召回率和精确度微调阈值。另外，务必保存过滤前后的突变数量统计，这是衡量数据质量和流程严格度的重要指标。

4. 模块二：HLA分型——决定呈递格局的“锁芯”

如果把新抗原比作“钥匙”，那么HLA分子就是“锁芯”。钥匙能不能用，首先得看锁芯的型号。HLA分型的准确性至关重要，错误的分型会导致后续所有结合预测完全偏离方向。

4.1 分型原理与工具选择：从序列到型别

HLA分型工具（如OptiType、Polysolver）的工作原理，是从测序数据中提取覆盖HLA基因区域的读取，与已知的HLA等位基因序列数据库进行比对，通过统计模型或枚举法，推断出最有可能的HLA型别。它们通常支持从RNA-seq或WES/WGS数据中分型。

• OptiType：采用整数线性规划模型，将分型问题转化为优化问题，结果非常精确，尤其对HLA-I类基因（A, B, C）支持很好，是目前很多流程的默认选择。
• Polysolver：较早的工具，基于外显子捕获和贝叶斯模型，也较为常用。
• HLA-HD：对二代测序数据支持较好，能给出高分辨率的4位数字分型（如A*02:01）。

我的建议是：首选OptiType进行HLA-I类分型。如果条件允许，可以用另一个工具（如HLA-HD）进行验证。对于HLA-II类基因（DP, DQ, DR），分型更复杂，准确度相对低一些，工具支持也不如I类完善，需要特别注意。

4.2 关键注意事项与验证

1. 输入数据质量：分型工具对覆盖HLA基因区域的测序深度很敏感。WES数据由于捕获探针设计，在HLA区域覆盖可能不均匀。WGS或RNA-seq数据通常更好。务必检查工具输出的覆盖度报告。
2. 分辨率问题：临床或研究需求可能需要高分辨率分型（4位或8位数字）。OptiType通常输出2位或4位。要明确你的下游预测工具（如NetMHC）需要什么分辨率。大多数工具需要4位分型。
3. 验证：如果样本来源允许，强烈建议用湿实验方法（如PCR-SSO、测序）对关键样本的HLA分型进行验证，至少是随机抽查。计算分型在极端罕见等位基因或数据质量差时可能出错。
4. 杂合子与纯合子：正确识别HLA位点是杂合（两个不同的等位基因）还是纯合（两个相同的等位基因）很重要。纯合子意味着该位点只有一种“锁芯”，预测时需要考虑到这一点。

5. 模块三：pMHC结合预测——计算免疫学的核心战场

这是新抗原预测中最具计算生物学特色的部分。我们需要为每个过滤后的非同义突变，生成可能的多肽序列，并预测它们与患者HLA分子的结合强度。

5.1 肽段生成：从突变到候选肽

对于一个点突变，我们需要考虑该突变在蛋白质序列上下文中的所有可能肽段。通常，我们会以突变氨基酸为中心，分别向上游和下游延伸，生成一系列固定长度的肽段（如8、9、10、11聚体）。例如，对于一个SNV，我们会生成包含该突变位置的所有可能长度的肽段。

这里的一个关键决策是：是否包含野生型（正常）序列的肽段？答案是必须包含。因为我们需要对比突变肽段和野生型肽段与HLA的结合差异。一个强结合的突变肽段，如果其对应的野生型肽段结合也很强，那么它被免疫系统有效识别的可能性就会降低（因为胸腺阴性选择可能已清除了针对该野生型肽段的T细胞）。

5.2 预测算法与工具实战

当前主流的预测工具（NetMHC、NetMHCpan、MHCflurry）都基于机器学习模型，它们通过大量已知的结合数据训练，预测新肽段的结合亲和力。

• NetMHC/NetMHCpan：经典工具，基于人工神经网络。NetMHC针对已知等位基因，NetMHCpan可预测任何等位基因。输出指标包括IC50（半最大抑制浓度，数值越小结合越强）和%Rank（与随机肽段库相比的排名，<0.5%或<2%常作为强结合阈值）。
• MHCflurry：基于深度神经网络（长短期记忆网络LSTM），在某些等位基因上表现可能优于NetMHCpan，且支持本地化部署和训练，速度较快。

如何选择？对于常规分析，NetMHCpan 4.0或更新版本是一个稳健的起点。MHCflurry 2.0因其易用性和性能也越来越受欢迎。我个人的流程中会并行运行两者，观察结果的一致性。

运行预测时，你需要准备：

1. 一个FASTA文件，包含所有需要预测的肽段序列（突变型和野生型）。
2. 患者的HLA等位基因列表。
3. 选择肽段长度（通常为8-11）。

一个典型的MHCflurry命令如下：

# 假设 peptides.fasta 包含肽段， alleles.txt 包含HLA型别 mhcflurry-predict peptides.fasta --alleles-file alleles.txt --out predictions.csv

输出文件会包含每个肽段-等位基因对的预测结合亲和力（IC50）和%Rank。

5.3 阈值设定：强结合不等于新抗原

这是最大的误解之一。预测工具给出的“强结合剂”（strong binder）只是一个必要条件，远非充分条件。通常的初筛阈值是：

• IC50 < 50 nM或%Rank < 0.5（对于NetMHCpan，强结合剂）
• IC50 < 500 nM或%Rank < 2（弱结合剂）

我们首先会筛选出突变肽段为强/弱结合剂，而对应的野生型肽段结合很弱（如IC50 > 5000 nM）或不是结合剂的肽段。这一步筛选能排除大量“假性新抗原”。

踩坑实录：不要盲目迷信低IC50值。我曾遇到一个突变肽段，预测IC50值极低（<10nM），看起来是绝佳候选。但进一步检查发现，该突变位于一个表达量极低的基因上（FPKM < 1）。没有表达，蛋白都无法合成，何谈呈递？因此，必须结合表达量数据。

6. 模块四：免疫原性预测与综合排序——从“能结合”到“能激活”

经过前三步，我们得到了一份“可能被呈递”的肽段列表。但最终目标是找到那些“能激活T细胞免疫反应”的新抗原。这一步不确定性最大，也是当前研究的焦点。

6.1 免疫原性预测模型

免疫原性（Immunogenicity）指肽段引发T细胞反应的能力。它不仅仅取决于pMHC结合强度，还涉及：

• T细胞受体（TCR）识别：肽段-MHC复合物的结构是否容易被TCR识别。
• 抗原加工：肽段能否被蛋白酶体切割、被TAP转运。
• 同源野生型肽段的耐受性：如前所述，如果野生型肽段也存在，可能导致中枢耐受。

一些工具尝试整合这些因素：

• NetCTLpan：在结合亲和力基础上，加入了蛋白酶体切割和TAP转运效率的预测。
• PRIME：一个更复杂的模型，整合了结合亲和力、序列特征等来预测免疫原性。
• pVACtools等流程：会计算突变肽段与野生型肽段的差异，差异越大，逃避免疫耐受的可能性越高。

然而，必须清醒认识到，目前所有的计算免疫原性预测模型，其准确率都远未达到临床可靠的水平。它们更多是提供一种相对的优先级排序参考。

6.2 多维度信息整合排序

因此，一个更务实的策略是构建一个综合评分体系，整合多个维度的证据，对候选新抗原进行排序。我常用的评分维度包括：

1. 结合亲和力（权重最高）：突变肽段的%Rank或IC50值，标准化为0-1的分数。
2. 突变肽与野生型肽的结合差异：计算两者结合分数的比值或差值，差异越大，分数越高。
3. 基因表达水平：从RNA-seq数据获取突变基因的TPM或FPKM值。没有表达，一切归零。可以设置一个表达阈值（如TPM > 1），并进行对数转换后评分。
4. 突变等位基因频率：在肿瘤样本中的突变频率（VAF）。VAF越高，意味着表达该突变蛋白的癌细胞比例越高，作为靶点的价值可能越大。
5. 突变的功能重要性：通过注释（如SIFT, PolyPhen）判断突变是否可能影响蛋白功能。有害突变可能更关键，但这不是绝对标准。
6. 肽段-MHC复合物稳定性预测：除了结合亲和力，复合物在细胞表面的半衰期（稳定性）也影响免疫原性。工具有NetMHCstabpan。
7. 克隆性：来自体细胞进化树分析，是否为克隆性突变（存在于所有癌细胞中）。亚克隆突变可能因免疫编辑而丢失。

你可以为每个维度分配权重，计算一个加权总分。例如：综合得分 = 0.4*结合分数 + 0.2*表达分数 + 0.2*VAF分数 + 0.1*差异分数 + 0.1*其他...

这个权重需要根据你的研究目标和经验调整。最终输出是一个按综合得分降序排列的表格，前20-50个候选者通常会被优先考虑进行下游实验验证（如肽段合成、T细胞功能实验）。

6.3 流程自动化与可重复性

手动执行以上所有步骤是不现实的。因此，需要使用工作流管理系统将其自动化。我的选择是：

• Snakemake：基于Python，规则定义清晰，非常适合生物信息学流程。它能优雅地处理复杂的依赖关系和集群投递。
• Nextflow：基于Groovy/DSL，声明式语法，内置强大的容错和重试机制，对Docker/容器支持极好。
• pVACtools：如果你想要一个开箱即用、集成了多个预测工具和排序方法的专用套件，这是一个很好的选择。它封装了从VCF到最终排名的许多步骤。

无论选择哪种，核心是保证流程的可重复性。必须使用容器（Docker/Singularity）封装所有工具和依赖，并严格管理软件版本。每次运行都应记录完整的参数和版本信息。

7. 验证、局限性与未来展望

计算预测的终点是湿实验验证。最常见的验证方法是合成排名靠前的候选肽段，在体外用患者的免疫细胞（如外周血单个核细胞PBMC）进行刺激，通过ELISpot、细胞内因子染色等技术检测抗原特异性T细胞反应。预测排名与实验验证结果的相关性，是评估你整个工作流程效能的黄金标准。

必须承认，当前的计算新抗原预测流程存在明显局限：

1. 假阴性/假阳性：算法不完美，会漏掉真实新抗原，也会推荐无效候选。
2. HLA-II类预测不成熟：对CD4+ T细胞重要的HLA-II类新抗原预测准确度远低于I类。
3. 忽略转录后修饰：如磷酸化、糖基化产生的新抗原无法被当前基于序列的流程捕获。
4. 肿瘤异质性：活检样本可能无法代表所有肿瘤克隆，预测的新抗原可能只针对部分癌细胞。

未来的方向包括整合多组学数据（如质谱检测到的实际呈递肽段数据）、应用更复杂的深度学习模型、以及结合单细胞测序来理解肿瘤微环境中的免疫状态。但无论如何进化，一个严谨、透明、可重复的计算工作流程，始终是连接基因组数据与免疫治疗应用的桥梁。构建这个流程的过程，本身就是对肿瘤免疫学深刻理解的一次实践。