1. 课程概述

各位同学,咱们今天正式开课。先聊聊功能注释这件事。

你拿到一个陌生的基因或蛋白质序列,第一反应是什么?我猜大概是——跑个BLAST,看看它像谁。这很自然,BLAST就像生物信息学里的“瑞士军刀”,简单粗暴,谁都用过。但干了十几年这行,我得说一句:光靠BLAST,容易翻车。尤其是当你面对的是非模式生物、宏基因组数据,或者那些“看起来很靠谱但实际是假阳性”的比对结果时。

那怎么办?我的答案是:把Pfam请进来,让它们俩打配合。

1.1 为什么BLAST一个人扛不住?

BLAST的核心思路,说白了就是“找长得像的”。它比较你给我的query序列和数据库里每条序列的局部相似性。速度快,上手简单,这没得说。

但问题也出在这——它只看“长得像”,不看“是不是一家子”。

  • 同源 vs. 相似:BLAST告诉你A和B有80%的相似性,但它没法保证A和B有共同的祖先(同源)。你想想看,趋同进化、低复杂度区域“假匹配”,都会让BLAST给出一个漂亮的E-value,然后把你带沟里去。
  • 结构域识别短板:很多蛋白质是多结构域的。BLAST只能给你一个粗粒度的全长比对,很难说清楚“这个蛋白的哪个片段在干哪个活”。
  • 敏感度天花板:对于进化距离很远的同源蛋白,序列可能已经面目全非。BLAST靠的是成对比对(pairwise alignment),在这种情况下捡出信号的能力有限。
我曾经接手过一个项目,同事用BLAST把一段宏基因组的contig注释成了“新型抗生素合成酶”,兴奋得不得了。我一看,嗯,确实和数据库里某条序列有70%的相似性。但当我用Pfam一查,真正匹配上的Pfam结构域是个转运蛋白家族——功能和抗生素合成八竿子打不着。那条“高相似度”的序列,恰好是一个保守的跨膜区造成的假象。

这个坑,我踩过。所以现在,BLAST只排第一关,绝不轻信。

1.2 Pfam带来了什么?

Pfam不是看“序列像不像”,而是看“你属于哪个家族”。它维护的是蛋白质结构域和家族的多序列比对(MSA)以及对应的隐马尔可夫模型(HMM)。

用HMM来搜索,相当于拿着整个家族的“共识画像”去比对,而不是拿着某个成员的“大头照”去匹配。这带来了两个巨大的好处:

  • 灵敏度高:即使序列相似度掉到20%,如果它保留了家族的保守模式,HMM仍能把它捞出来。这一点BLAST很难做到。
  • 结构域级别的分辨率:Pfam告诉你“第120到300位氨基酸是一个蛋白激酶结构域,第400到500位是一个SH3结构域”。这种粒度的注释,BLAST给不了。

我个人习惯是,做功能注释之前先想清楚:我需要的是“这个蛋白像谁”(BLAST),还是“这个蛋白属于哪个家族,具备什么结构域架构”(Pfam)。

1.3 联合注释:1+1 > 2

既然BLAST和Pfam各有千秋,那把它们串起来,效果就上去了。我推荐的策略是这样的:

  1. 先用BLAST快速过一遍:对query序列做BLASTP或BLASTX,拿到最初的top hit和物种分布,快速评估数据质量。这一步快,性价比高。
  2. 再用Pfam精确扫描:用hmmscan(或者HMMER的phmmer)对同一条序列搜索Pfam数据库。这一步会给你结构域边界、E-value、家族ID。
  3. 交叉验证:把BLAST的top hit手动比对到Pfam的结果上。如果BLAST说是“激酶”,Pfam也出“激酶结构域”,那基本稳了。如果BLAST说是“转录因子”,Pfam只匹配到一个“锌指”,那你就要警惕——可能BLAST只抓到了那个保守的DNA结合域,而全长的同源性是假的。

说白了,BLAST负责“广撒网”,Pfam负责“精确定位”。

记住一个原则:BLAST结果引导方向,Pfam结果确认功能。两者缺一不可。

1.4 实战实例:一条序列的注释之旅

为了让你更直观地理解,我们拿一个简单的场景来看。假设我有一条query序列,长度大概450个氨基酸。这是实际项目中遇到的一个案例。

第一步:BLASTP搜索。

$ blastp -query query.fa -db nr -out blast_out.txt -evalue 1e-5 -outfmt 6
输出摘要:
seq_id_001   95% identity  | E-value: 0.0  | 假设的环化酶 [某细菌]
seq_id_002   60% identity  | E-value: 1e-120 | 未表征蛋白 [另一细菌]
seq_id_003   30% identity  | E-value: 1e-15  | 环化酶类似物 [真核生物]

一看,top hit是“环化酶”,95%的相似性,E-value完美。很多人到这一步就停了,直接写论文“这是一个环化酶”。

我建议你别急。咱们上Pfam看看更细的:

第二步:hmmscan搜索。

$ hmmscan --domtblout pfam_out.txt Pfam-A.hmm query.fa
输出摘要:
#target name    accession  tlen  query name   accession   qlen  E-value  score  bias
Cyclase_fam     PF01234.20  350   query_seq    -           450   2.3e-50  180.2  0.0
ADH_zinc_fam    PF00107.25  140   query_seq    -           450   4.5e-12   45.1  0.1

看到了吗?这里Pfam找到了两个结构域:一个是环化酶家族(覆盖第10-350位),另一个是含锌的醇脱氢酶家族(覆盖第360-450位)。

第三步:联合解读。

BLAST只给了我一个“环化酶”的标签,因为它和那条参考序列的全长相似度足够高。但Pfam告诉我:这个蛋白是个双结构域融合蛋白,可能执行着环化反应+氧化还原反应偶联的功能。如果我只信了BLAST,就会错过它C端那个关键的氧化还原结构域。

这种“结构域重组”现象在自然界非常普遍,尤其在微生物次级代谢通路里。不联合注释,你真的会错过很多重要信息。

1.5 课程知识体系一览

好,概念铺垫到这。我们把这一章以及整个课程的核心逻辑,用一张图串起来。

联合注释策略知识体系 未知蛋白质序列 BLAST(快速初筛) 成对比对·速度优先·粗粒度 Pfam(精确扫描) HMM模型·高灵敏度·结构域级 交叉验证与联合解读 BLAST引导方向 · Pfam确认功能 高质量功能注释输出 课程覆盖:BLAST配置与优化 · Pfam数据库解析 · HMMER使用 · 联合策略开发 · 实战案例 · 流程自动化

这张图清晰展示了我们整个课程的脉络:从一条未知序列出发,经BLAST和Pfam双通道比对,最终汇聚到交叉验证环节,得到高质量的功能注释。后续的29章,都会围绕这个框架展开细节。

1.6 你将从这门课带走什么

这门课我们不谈虚的。每一章都会给你可执行的方法:

  • BLAST怎么跑才能兼顾速度和精度? —— 不是无脑-nr,也不是无脑-evalue 1e-5。
  • Pfam数据库的HMM模型怎么选? —— Pfam-A和Pfam-B有什么区别?什么时候用full,什么时候用seed?
  • 联合策略具体怎么写脚本? —— 我会给出现成的shell/Python流水线,你可以直接改改就用。
  • 遇到冲突怎么判? —— BLAST说是A,Pfam说是B,你听谁的?我有实战案例教你逐级排查。

嗯,大致就是这样。第一章既是开胃菜,也是定调。记住:BLAST和Pfam是互补的,不是替换的。只有把两者拧成一股绳,你的功能注释才算真的做扎实了。

小贴士:如果你刚接触这些工具,别急着布置大数据集。先拿一条自己熟悉的蛋白序列跑一遍BLAST和hmmscan,对比着看输出。这个“练兵”的过程,比读十篇教程都有用。