1. BLAST 简介：发展历史、基本原理与常用变体
说到 BLAST，这玩意儿在生物信息学圈子里，基本就是「瑞士军刀」级别的存在。我个人做序列分析这么多年，跟它打交道的频率比跟家人吃饭还高。今天咱们就从根上捋一捋，BLAST 到底是什么，它怎么来的，以及我们什么时候该用哪种版本。
1.1 它从哪来？—— 一段简短的发展史
上世纪 80 年代，做序列比对是一件相当痛苦的事。那时候主流算法是 Needleman-Wunsch 和 Smith-Waterman，虽然结果准，但算起来慢得要命。你想比对一条 1000 bp 的序列跟整个数据库去比？嗯，我当年用工作站跑过，没几个小时下不来。
直到 1990 年，NCBI 的一群狠人（Altschul 等人）提出了 BLAST 算法。这篇经典文章发表在哪？Journal of Molecular Biology，我至今还存着 PDF。它最大的贡献就是——把比对速度提升了几个数量级。
后来 BLAST 一直在演进：Gapped BLAST 引入了空位，PSI-BLAST 能迭代搜索，还有专门针对短序列的 MegaBLAST。我个人最常用的就是 MegaBLAST，尤其是在处理 NGS 数据的时候，真的快。
关键里程碑：
- 1990 — BLAST 原版算法发布
- 1997 — Gapped BLAST + PSI-BLAST 登场
- 2002 — MegaBLAST 上线，专治短 reads
- 2009 — BLAST+ 版本发布（彻底重写了底层代码）
1.2 核心原理：BLAST 到底是怎么找相似序列的？
说白了，BLAST 不是一上来就笨比着全局比对。它的策略很聪明——先找「种子」，再扩展。
核心三步走：
- 建立索引（k-mer 词典）
将查询序列切分成短片段（k-mer，默认长度 11aa 或 28nt）。数据库也提前建好索引。 - 种子检测（Seed-and-extend）
在数据库里找到与查询 k-mer 高度匹配的「热点」位置。这些位置就是潜在的比对起点。 - 扩展与打分
从种子位置向两侧延伸，用动态规划（有限范围）拼出精确比对。采用空位罚分 + 替换矩阵（BLOSUM62 或 PAM 族）计算得分。
避坑指南：我曾经遇到一个新手，直接用默认参数比对病毒短片段，结果漏了好多相似序列。为什么？因为 k 值太大，短序列切不出有效种子。这种情况我建议调小 -word_size。你想想看，默认 28 的 word size 对 50 bp 的 reads 来说，基本是灾难。
整个逻辑可以总结成下面这张图：
1.3 BLAST 家族：到底该用 blastn 还是 blastp？
这个我几乎每周都会被问到。其实选哪个，核心就一条规则：你的序列类型决定了变体选择。下面这张表是我自己整理出来的，记住了它，你就能少绕很多弯路。
| 变体 | 查 询 | 数据库 | 典型用途 |
|---|---|---|---|
| blastn | DNA/RNA | 核苷酸 | 物种鉴定、同源基因搜寻 |
| blastp | 蛋白质 | 蛋白质 | 功能注释、保守域分析 |
| blastx | DNA（翻译） | 蛋白质 | EST/RNA-seq 的蛋白功能标注 |
| tblastn | 蛋白质 | DNA（6框翻译） | 未知基因的跨物种搜索 |
| tblastx | DNA（翻译） | DNA（翻译） | 远源序列的敏感性比对 |
⚠️ 注意区分 blastx 和 tblastx。 blastx 是把你的 DNA 翻译成蛋白去比对蛋白库；而 tblastx 是两边都翻译，两边都比对。显然后者慢很多，但有时能发现跨物种的微弱相似性。我当年做鱼类免疫基因鉴定时，就靠 tblastx 找到了跟人类 IL-1β 同源的一个碎片——用单纯 blastn 完全没戏。
1.4 命令行实战：手把手跑一次 BLAST
说了这么多，来点实操。BLAST+ 的命令行其实很直白，我个人习惯用这种格式：
# blastn 示例：用一条未知菌16S序列比对ncbi_nt库
blastn -db nt -query unknown_16S.fa -out result.txt -outfmt 6 -evalue 1e-5我来拆解读一下：
- -db nt — 数据库，这里是 NCBI 的核苷酸总库
- -query — 输入 fasta 文件
- -outfmt 6 — 输出表格格式（后续章节我们主要解这个格式）
- -evalue 1e-5 — 期望值阈值，大于此值的 hits 会被丢弃
如果比对蛋白质序列，改成 blastp 就成：
blastp -db nr -query my_protein.fa -out protein_blast.out -outfmt 6个人习惯：我通常在比较任务开始时先把输出格式固定为 -outfmt 6，后面用 Python 或 R 来解析。这比看默认的 pairwise 文本舒服多了。当然，做教学演示时我也会开个 -outfmt 0 来看具体比对详情。
1.5 如何设定参数？几个关键选项
参数这东西，不同需求差距很大。但有几个核心参数你必须理解：
| 参数 | 含义 | 我的建议 |
|---|---|---|
| -word_size | 种子长度 | 短序列调小（7-11），长序列用默认 |
| -evalue | 期望阈值 | 严格搜索用 1e-10，宽松用 1e-3 |
| -num_alignments / -max_target_seqs | 最多返回条数 | 日常我设 100，探索性设 500 |
| -gapopen / -gapextend | 空位罚分（开/扩） | 软件默认即可，除非你做特殊优化 |
| -matrix | 替换矩阵 | 蛋白常用 BLOSUM62，近源用 BLOSUM80 |
为什么我会特意提到 BLOSUM62？因为它是 BLAST 默认矩阵，覆盖了从 62% 到 90% 相似度的蛋白质。如果你手头是高度保守的序列，可以考虑用 80，如果有远源同源就试试 45。
1.6 常见误区：新手容易踩的坑
- 选了错误变体 — 用 blastn 查蛋白序列，这种比对基本废了
- 忽略了 E-value 设置 — 设成 1e-300 这种极小数，会滤掉一切真实结果
- 没建本地索引 — 每次单独建库其实很慢，建议跑
makeblastdb建一次 - 忽略序列质量 — N 多的序列比对结果会非常差，我建议先清理一下
避坑指南：我曾经给一个合作者处理数据，他的 blastx 结果显示 E-value 全是 0 或极大值。后来发现他把蛋白当成核酸输进去了——BLAST 不会报错，但数据完全不能看。所以，跑之前一定要确认文件格式是核酸还是蛋白。
好了，BLAST 的基础就讲到这里。这章的核心就是把工具的原理和家族成员摸透，后面的章节我们会用它真正干活——解析比对结果图谱，从那些表格和点状图中看出生物学意义来。