10分钟用ANARCI给抗体编号：从一条序列到千条FASTA的实战笔记

10分钟用ANARCI给抗体编号：从一条序列到千条FASTA的实战笔记

【免费下载链接】ANARCIAntibody Numbering and Antigen Receptor ClassIfication项目地址: https://gitcode.com/gh_mirrors/an/ANARCI

开篇：审稿人一句话，逼你学会编号

我第一回接触抗体编号，是被审稿意见逼的。当时论文里有个小鼠单抗的序列，审稿人只回了一句："Please provide the IMGT numbering of the antibody variable domains and annotate the CDRs." 大意是：把 CDR 位置按 IMGT 标准标出来。

我盯着那条一百多个氨基酸的字符串看了十分钟，然后翻开 IMGT 编号表开始手抄：FRI 的残基对齐 1 号，CDR1 从 27 号开始……抄到第 50 个残基就乱了，因为 CDR1 长一点的序列会插进 A、B、C 这样的插入码。最后我抄出来的编号，自己都不敢信是对的。

如果你也有过类似经历——手里攥着一条抗体序列，却不知道它该叫"1、2、3"还是"26、27、28"，不知道怎么把 CDR1、CDR2、CDR3 从一串字母里精确切出来——那这篇笔记就是给你写的。我要带你把牛津蛋白信息学组（OPIG）开发的ANARCI用起来。这个名字是 Antibody Numbering and Antigen Receptor ClassIfication 的缩写，它干的事一句话说清：给你一条抗体或 TCR 序列，它用隐马尔可夫模型自动识别物种和链型，然后按 IMGT、Kabat、Chothia 等 6 种编号方案给出标准编号，顺带把 CDR 边界、种系基因（germline）归属全部算好。

接下来我们不做"安装→功能介绍→FAQ"式的流水账。我们直接演一遍：从你拿到一条序列、到把上千条序列跑完导出 CSV 的完整过程。每一步都是真实场景，每一段代码你都能直接跑。

魔法时刻：一条命令，编号和CDR全出来

先给你看 ANARCI 最爽的地方。假设你手里有一条小鼠重链可变区序列（就是我从项目自带示例里摘的，来自抗体 12E8 的重链，真实存在于项目Example_scripts_and_sequences/12e8.fasta中）：

ANARCI -i EVQLQQSGAEVVRSGASVKLSCTASGFNIKDYYIHWVKQRPEKGLEWIGWIDPEIGDTEYVPKFQGKATMTADTSSNTAYLQLSSLTSEDTAVYYCNAGHDYDRGRFPYWGQGTLVTVSA

一秒之内，屏幕上会吐出类似这样的结果：

# 12E8:H # ANARCI numbered # Domain 1 of 1 # Most significant HMM hit #|species|chain_type|e-value|score| #|mouse|H|8.6e-58|184.9| # Scheme = imgt H 1 Q H 2 V H 3 Q H 4 L H 5 Q ... H 26 G H 27 F <- CDR1 起点 ... H 105 A <- CDR3 区域 ... //

注意那行#|mouse|H|8.6e-58|184.9|——它告诉我们：ANARCI 自动识别出这是小鼠（mouse）的重链（H），e-value 8.6e-58 表示这个判断的可信度极高。你不需要告诉它物种，不需要告诉它链型，它自己用 HMMER 把序列和人类、小鼠、大鼠、兔等物种的 V/J 种系 HMM 数据库比对一遍，谁的分最高、最显著，就用谁的。每条氨基酸都有标准位置编号，CDR 在哪里自然一目了然。

这就是魔法时刻：从一串"字母面条"到结构化、可比对的编号结果，你只敲了一行命令。现在让我们把魔法用在一个真实的项目流程里。

案例主线：把 2000 条序列变成一张编号表

假设你现在的任务真实一点：从实验室的杂交瘤筛选结果里拿到了一批候选抗体序列，存在一个 FASTA 文件里，有 2000 条。老板说："把每条序列的链型分出来，按 IMGT 编号，标出 CDR，导出成表格。"手动做？不可能。下面是完整的四步。

第一步：把环境搭好，10 分钟内能跑通

安装没有想象中复杂。ANARCI 是 Python 包，但核心的序列比对靠的是 HMMER 这个外部程序，所以要先把它请进来。用 conda 是最省事的：

conda install -c conda-forge biopython -y conda install -c bioconda hmmer=3.3.2 -y git clone https://gitcode.com/gh_mirrors/an/ANARCI cd ANARCI python setup.py install

装完之后验证一下：

ANARCI --help

能看到用法说明就说明成功了。这里有个小细节要提前告诉你：setup.py install的过程中会从 IMGT/GENE-DB 下载种系基因数据并构建 HMM 模型（项目里的INSTALL文件明确写了这一步），所以第一次安装需要联网，也请耐心等它跑完。

验证完，先把项目自带的一对示例序列跑一遍，确认环境 OK：

ANARCI -i Example_scripts_and_sequences/12e8.fasta

12e8.fasta里其实藏了四条链（H、L、M、P），M 是轻链 L 的重复、P 是重链 H 的重复。你会看到它们被正确分成了重链和轻链，这就是第一步的成果：你能得到什么——一条条带物种、链型、编号注释的清晰输出。

第二步：处理真正的 FASTA 文件

项目里给你准备了一个很好的练习素材：Example_scripts_and_sequences/antibody_sequences.fasta，里面有近 2000 条真实 PDB 抗体链序列。直接喂给 ANARCI：

ANARCI -i Example_scripts_and_sequences/antibody_sequences.fasta -o my_numbering.txt

-o指定输出文件，不加的话结果会刷满整个终端。每条序列的结果以//分隔，没被识别成抗体的序列（比如某个溶菌酶混进去了）会原样输出名字，并告诉你没有显著比对结果。

但垂直格式的编号，人读起来舒服，做下游分析却不好处理。所以更常用的是 CSV 模式——这是我最推荐你用的姿势：

ANARCI -i Example_scripts_and_sequences/antibody_sequences.fasta -s i --csv -o my_imgt --ncpu 4 --assign_germline

逐块解释这段命令在干什么：

• -s i：选择编号方案，i是 IMGT 的缩写。项目支持的 6 种方案缩写分别是i（IMGT）、k（Kabat）、c（Chothia）、m（Martin，即增强版 Chothia）、a（AHo）、w（Wolfguy）。等下我们会专门对比它们。
• --csv：把编号结果按链型分开，写成横向排列的 CSV 文件。所有序列按编号方案对齐，同一位点排在同一列——这正是你想要的"表格"。
• -o my_imgt：CSV 模式必须给一个输出文件名前缀，ANARCI 会生成my_imgt_H.csv、my_imgt_L.csv这样的文件，按链型自动分类。
• --ncpu 4：让 HMMER 用 4 个核并行。项目里的性能脚本Example_scripts_and_sequences/run_numbering_benchmark.sh用 4 进程处理约一万条序列，注释写着"不到 5 分钟"，你可以拿 2000 条序列感受一下这个速度。
• --assign_germline：额外输出"与该序列最相似的种系基因"注释，比如IGHV1-12*01和序列同一性分数，这个信息对分析抗体的成熟程度很有用。

跑完后打开my_imgt_H.csv，你会看到每一行是一条序列、每一列是 IMGT 的一个位点，CDR 区域一目了然。你能得到什么：一份干净的、可直接丢进 Excel 或 pandas 的抗体编号对齐表。

第三步：把 ANARCI 嵌进你的 Python 流程

命令行适合一次性任务，但如果你要写脚本批量分析、或者把编号结果接进自己的管线，就该用 Python API 了。项目里有两个公开函数，先看最省事的number：

from anarci import number seq = "EVQLQQSGAEVVRSGASVKLSCTASGFNIKDYYIHWVKQRPEKGLEWIGWIDPEIGDTEYVPKFQGKATMTADTSSNTAYLQLSSLTSEDTAVYYCNAGHDYDRGRFPYWGQGTLVTVSA" numbering, chain_type = number(seq, scheme="kabat") print("链型：", chain_type) # 输出 'H' print("第一个结构域的 Kabat 编号：", numbering)

number很贴心：传一条序列进去，返回编号结果和链型，适合快速验证"这条序列是不是抗体、是重链还是轻链"。

如果你想要完整信息——物种、e-value、比对位点、HMMER 的全部命中表——用anarci主函数。示例脚本Example_scripts_and_sequences/anarci_API_example.py写得非常清楚，我把骨架拆给你看：

from anarci import anarci # 输入是 (名字, 序列) 的列表，可以一次塞很多条 sequences = [ ("12e8:H", "EVQLQQSGAEVVRSGASVKLSCTASGFNIKDYYIHWVKQRPEKGLEWIGWIDPEIGDTEYVPKFQGKATMTADTSSNTAYLQLSSLTSEDTAVYYCNAGHDYDRGRFPYWGQGTLVTVSA"), ("lysozyme:A", "KVFGRCELAAAMKRHGLDNYRGYSLGNWVCAAKFESNFNTQATNRNTDGSTDYGILQINSRWWCNDGRTPGSRNLCNIPCSALLSSDITASVNCAKKIVSDGNGMNAWVAWRNRCKGTDVQAWIRGCRL"), ] numbering, alignment_details, hit_tables = anarci(sequences, scheme="imgt", output=False) for i in range(len(sequences)): if numbering[i] is None: print(sequences[i][0], "没有识别出抗体/TCR 结构域") else: print(sequences[i][0], "识别出", len(numbering[i]), "个结构域") print(alignment_details[i][0]["chain_type"], alignment_details[i][0]["species"])

注意我把一条溶菌酶序列混了进去——这不是抗体。ANARCI 会老老实实告诉你它没有找到结构域（numbering[i] is None）。这个特性恰好可以用来做一件事：把你的序列库先过一遍 ANARCI，快速筛出哪些是抗体，非抗体的原样跳过。

函数返回的三个列表，和输入的序列一一对应：

你能得到什么：一个可编程的抗体编号引擎，几十行代码就能接入自己的分析管线。

第四步：上万条序列的批处理

当序列量再上一个数量级——比如你要给整个 PDB 数据库里的免疫球蛋白重编号——项目也给了现成的答案。打开Example_scripts_and_sequences/run_numbering_benchmark.sh看看：

ANARCI -i pdb_sequences.fa.txt.gz -s i --csv -o /tmp/pdb_imgt --ncpu 4 --assign_germline ANARCI -i pdb_sequences.fa.txt.gz -s k --csv -o /tmp/pdb_kabat --ncpu 4 --assign_germline

看到了吗，这个脚本一口气用 6 种编号方案各跑了一遍同一批数据（那个.gz结尾的输入文件 ANARCI 直接支持，不用手动解压）。你可以照葫芦画瓢：同一份 FASTA，换不同的-s参数，就能得到不同编号方案下的对齐结果，用于对比方案间的差异。

而在 Python 侧，如果你要处理的序列实在太多，就用run_anarci而不是裸的anarci——后者走的是 HMMER 内部的并行，前者是 ANARCI 在 Python 层面的原生多进程，任务量大的时候效率更稳。

你能得到什么：一条可扩展的批处理链路，从 2000 条到几万条，只是时间问题，不是能力问题。

能力地图：你现在站在哪一级

把刚才走过的东西整理成一张能力地图，方便你对号入座，也方便"跳级"：

新手把前两级跑通，就能解决 90% 的日常需求；第三级是遇到瓶颈时再回来翻的锦囊。

过来人的避坑笔记

下面这些坑，都是我在实际使用中踩过或从项目代码注释里读出来的，比任何 FAQ 都实在。

坑一：别拿 ANARCI 当物种注释工具。项目 README 第一页就明确写了：ANARCI 用物种的 V/J 种系比对来决定编号时该用哪套 HMM，但官方明确不推荐把它当作主要的物种鉴定工具。e-value 低只能说明"它最像小鼠重链"，不等于"它就是小鼠的"。跨物种嵌合抗体、人源化抗体、驼类 VHH 都可能给出偏颇的物种判断。要报物种，请用专门的工具，ANARCI 的结果只当参考。

坑二：非抗体序列不会被硬编号。我们前面混入的溶菌酶就是例子——ANARCI 不会硬凑一个编号出来，没有显著比对结果就返回None。这其实是好事：天然的筛选器。但反过来，如果你拿到的序列太短（比如只有 50 个残基），或者经过了大量工程改造、偏离天然抗体太远，它也可能"漏判"。这时候可以调低bit_score_threshold（默认 80）让识别更宽松，但小心：调太低会把其他免疫球蛋白样分子误认成抗体。

坑三：scFv 会识别出两个结构域。示例脚本里有个scfv:A的序列，VH 和 VL 通过 linker 连在同一条链上。ANARCI 会诚实地报告识别出两个结构域，numbering[i]是个长度为 2 的列表。写代码时如果默认"每条序列只有一个结构域"，遇到 scFv 就会下标越界或取错域——记得先判断len(numbering[i])。

坑四：不同编号方案，同一段序列，CDR 边界不一样。这是新手最容易懵的地方。同一个 CDR-H1，在 Kabat 里可能是 31-35，在 Chothia 里却算 26-32，在 IMGT 里又变成 27-38。所以：

• 论文里写编号，务必同时声明用了哪个方案，不然审稿人（没错，就是你遇到的那位）会再打回来。
• Kabat、Chothia、Martin 只适用于抗体免疫球蛋白（IG），TCR 得用 IMGT 或 AHo。
• Wolfguy 方案的 CDRL1 长度分类依赖硬编码的序列 motif，处理少见长度时要留意它的注释逻辑（这部分逻辑在lib/python/anarci/schemes.py的_get_wolfguy_L1里）。

坑五：插入码不是 Bug，是特性。长 CDR 会引入字母插入码，比如 IMGT 的 CDR3 插入是围绕 111/112 对称编号的，Kabat 则用 100A、100B 这样的后缀。这不是编号错了，而是方案本身的设计。解析结果时，把"位置编号"和"插入码"当两个字段处理，别拼成一个整数。

坑六：想给 PDB 结构重编号？有现成脚本。如果你的序列不是 FASTA 而是 PDB 结构文件，项目里的Example_scripts_and_sequences/ImmunoPDB.py能用 ANARCI 给结构中的每条链按指定方案重编号，还能标注 CDR 区域、按界面半胱氨酸距离做 H/L 配对，甚至把 scFv 拆成两条链。用法很简单：

python ImmunoPDB.py -i your.pdb -o renumbered.pdb -s imgt

处理 TCR 结构就把--receptor tr加上。注意它依赖 Biopython 和 muscle，缺了会报错。

收尾：现在就把手弄脏

回到开头那个场景。审稿人要 IMGT 编号和 CDR 标注，如果当时的我会用 ANARCI，整个过程只需要两步：把序列丢进去，把输出贴进论文。10 分钟，完事。

现在轮到你了。动手路径我给到最具体的程度：

1. 装好环境，先跑Example_scripts_and_sequences/12e8.fasta，把输出里的物种、链型、编号逐行读一遍。
2. 用--csv处理antibody_sequences.fasta，打开生成的 CSV，肉眼扫一遍对齐效果。
3. 打开Example_scripts_and_sequences/anarci_API_example.py，把它跑通，然后改成处理你自己的序列。
4. 找一条自己的序列，分别用-s i和-s k跑一次，对比 CDR 边界的差异——这个对比做完，你对"编号方案"的理解就超过大多数只会背概念的人了。

ANARCI 的核心价值，是把"抗体编号"这件本该是体力活的事，变成一行命令的事。用顺了之后你会发现，真正值钱的时间省下来，应该花在理解你的抗体为什么这么编号、以及这些位置背后的结构含义上——那才是研究该去的地方。

【免费下载链接】ANARCIAntibody Numbering and Antigen Receptor ClassIfication项目地址: https://gitcode.com/gh_mirrors/an/ANARCI