📚 课程目录第4章基因结构错误

4. 基因结构错误：外显子边界错误、内含子错误

基因注释里，出问题最多的就是基因结构。说白了，就是外显子和内含子的位置没标对。

你想想看，外显子边界差一个碱基，CDS翻译就可能移码，蛋白质序列全错。内含子边界不守GT-AG规则，剪接机器根本认不出来。这事我遇到过不止一次——看似小误差，到下游分析全炸了。

4.1 外显子边界错误

外显子边界错误，本质上是注释的外显子起止位置与真实剪接位点对不上。我在处理一批水稻RNA-seq数据时，发现某个基因的转录本注释普遍偏移了3个碱基。一查原因，是建库时用了旧版参考，边界标记系统性地偏了。

常见的外显子边界错误有3类：

边界偏移：外显子起点或终点多了/少了几bp。偏移量如果是3的倍数，CDS还能保住完整读码框；如果不是，直接移码。
边界不一致：同一基因的不同转录本，对同一个外显子的边界标注打架。这种情况在可变剪接基因里特别常见。
CDS边界与转录边界冲突：CDS区域的起止位置和外显子的边界不匹配，比如CDS起始不在外显子第一个碱基上。

看个具体例子——下面这条GFF记录，外显子边界标错了3个碱基：

# 错误：外显子起点偏移了3 bp（应该从1003开始）
chr1    .   exon    1000    1200    .   +   .   transcript_id "rna1"; gene_id "gene1"
chr1    .   CDS     1000    1200    .   +   0   transcript_id "rna1"; gene_id "gene1"

# 正确：边界对齐真实剪接位点
chr1    .   exon    1003    1200    .   +   .   transcript_id "rna1"; gene_id "gene1"
chr1    .   CDS     1003    1200    .   +   0   transcript_id "rna1"; gene_id "gene1"

你看，1000和1003就差3个碱基。但如果你做变异注释，这3 bp可能让一个错义突变变成同义突变——结论完全反了。

我的习惯：拿到新注释，第一件事就是随机抽50个外显子，和RNA-seq比对齐位置。偏移超过2 bp就说明建库或注释流程有问题，得整体排查。

4.2 内含子错误

内含子错误比外显子边界错误更隐蔽，因为它不直接影响CDS序列，但会影响剪接识别和转录本结构。

三大常见类型：

GT-AG规则违反：绝大多数内含子以GT开头、AG结尾（少数是GC-AG）。如果边界序列不符合这个规则，剪接体根本不会识别它。
内含子长度异常：典型真核生物内含子长度在50 bp到10 kb之间。少于20 bp或大于100 kb都值得怀疑。我在检查一个真菌基因组时发现过一条25 kb的内含子——后来确认是组装错误导致的假内含子。
内含子相位错误：内含子相位（0/1/2）指的是内含子位于密码子的哪个位置。相位错误会导致翻译产物在剪接后发生移码或提前终止。

来看一个GT-AG违规的例子：

# 错误：内含子边界序列是GT-AT，不是GT-AG
chr2    .   intron  5000    5200    .   +   .   transcript_id "rna2"; gene_id "gene2"

# 实际边界序列：
# 5'端：GTAT...（应该为GT）
# 3'端：...AT（应该为AG）

为什么会这样？很多情况是注释工具把内含子区域直接标成"外显子之间的区间"，根本没检查边界序列。我建议各位千万别轻信自动注释，一定要拿边界序列做二次验证。

注意：不是所有非GT-AG的内含子都是错误。少数可变剪接事件使用GC-AG边界。但这种情况占比不到1%，遇到时要先确认是否有实验证据支持。

4.3 基因结构错误的检测方法

检测基因结构错误，我有一套固定的检查流程。下面这张图把错误类型和检测方法串起来了，你可以对着它搭自己的质控管线：

对应的检测手段我整理了一个速查表：

错误类型	检测方法	推荐工具 / 命令
外显子边界偏移	与RNA-seq比对，检查边界覆盖断点	gff3validator, bedtools
外显子边界不一致	同基因跨转录本比较外显子起止	AGAT, custom script
CDS/转录边界冲突	检查CDS相位的帧与exon边界是否匹配	gff3validator
GT-AG规则违反	提取内含子5'/3'端序列，检查dinucleotide	seqkit, grep, awk
内含子长度异常	统计所有内含子长度分布，标记离群值	AGAT, awk, R
内含子相位错误	检查相位与密码子位置的一致性	gff3validator, AGAT

具体操作上，我常用下面这组命令快速筛查内含子边界问题：

# 从GFF提取所有内含子坐标
grep -P "\tintron\t" annotation.gff | \
  awk '{print $1, $4, $5, $7}' > introns.bed

# 提取边界序列并检查GT-AG
cat introns.bed | while read chr start end strand; do
  seqkit subseq -r $start:$end -s $chr genome.fa | \
    grep -P "^GT.*AG$"
done

效率建议：内含子边界检查用并行加速。我在检查一个4万条内含子的注释时，用xargs -P 8把时间从20分钟压到了3分钟。

4.4 实战案例与修复策略

说一个我实际修过的案例。某植物基因组注释里，一个转录因子基因的5号外显子边界标注有问题。原始GFF长这样：

# 问题GFF记录
chr3    .   transcript  50000   55000   .   +   .   transcript_id "tf1"
chr3    .   exon        50000   50700   .   +   .   transcript_id "tf1"
chr3    .   exon        50800   51200   .   +   .   transcript_id "tf1"
chr3    .   exon        51300   55000   .   +   .   transcript_id "tf1"
chr3    .   CDS         50000   50700   .   +   0   transcript_id "tf1"
chr3    .   CDS         50800   51200   .   +   1   transcript_id "tf1"
chr3    .   CDS         51300   55000   .   +   2   transcript_id "tf1"

肉眼一看没啥大问题。但我用RNA-seq数据做sashimi图时发现，第二外显子（50800–51200）的5'端实际剪接点在50803处，差了3 bp。

检查了边界序列后确认：50800处不是剪接供体位点，真正的GT在50803。这说明注释工具把外显子起点标早了3 bp。

修复后的记录：

# 修复后
chr3    .   exon        50000   50700   .   +   .   transcript_id "tf1"
chr3    .   exon        50803   51200   .   +   .   transcript_id "tf1"
chr3    .   exon        51300   55000   .   +   .   transcript_id "tf1"
chr3    .   CDS         50000   50700   .   +   0   transcript_id "tf1"
chr3    .   CDS         50803   51200   .   +   1   transcript_id "tf1"
chr3    .   CDS         51300   55000   .   +   2   transcript_id "tf1"

修复后翻译的蛋白质序列多了1个亮氨酸，功能域完整了。如果不修，这个转录因子在下游分析里会一直被当成假基因或截短蛋白。

核心经验：外显子边界错误靠"人眼+规则校验"很难全发现。我的做法是——批量提取所有外显子–内含子连接处的序列，用k-mer分布做异常检测。连接处序列保守性异常的区域，十有八九是边界标错了。

避坑提醒：我曾经对一个组装质量一般的基因组做注释质控，结果发现大量"内含子长度异常"——不是注释错了，是组装把内含子区域拼出了gap。这种情况要先用组装完整性评估（如BUSCO、QUAST）做个判断，别上来就改注释，改完也是白改。

总结一句话：基因结构错误是质控里最常见的硬伤。外显子边界查偏移、内含子边界查GT-AG、长度查分布。这三刀砍下去，80%的结构问题都能揪出来。