基因组重复序列分析¶

一句话概述：基因组中约50%的序列是重复序列（转座子、卫星DNA等），RepeatMasker通过同源搜索屏蔽已知重复，RepeatModeler2从头发现新的转座元件家族，2025年Terrier深度学习分类器进一步提升了TE分类的准确性和覆盖率。

核心知识点速览¶

概念	白话解释
重复序列	基因组中出现多次的DNA序列（人类基因组~50%是重复）
转座元件（TE）	能在基因组中"跳来跳去"的DNA序列
散在重复	分散在基因组各处的TE（SINE/LINE/LTR/DNA转座子）
串联重复	连续排列的短序列重复（微卫星/卫星DNA）
LINE	长散在核元件（如L1，~6kb，自主转座）
SINE	短散在核元件（如Alu，~300bp，非自主转座）
LTR逆转座子	两端有长末端重复序列的逆转座子（类似内源性逆转录病毒）
RepeatMasker	基于同源搜索的重复序列屏蔽工具
RepeatModeler2	从头发现新TE家族的工具
Dfam	转座元件家族的开放数据库

一、重复序列基础¶

1.1 重复序列分类¶

基因组重复序列的层级分类：

重复序列
├── 散在重复（Interspersed Repeats）= 转座元件
│   ├── I类（逆转座子，copy-and-paste机制）
│   │   ├── LTR逆转座子（有长末端重复）
│   │   │   ├── Gypsy/Ty3
│   │   │   ├── Copia/Ty1
│   │   │   └── ERV（内源性逆转录病毒）
│   │   └── 非LTR逆转座子（无长末端重复）
│   │       ├── LINE（长散在核元件）
│   │       │   ├── L1（人类最活跃，~6kb）
│   │       │   └── L2
│   │       └── SINE（短散在核元件）
│   │           ├── Alu（灵长类特有，~300bp）
│   │           └── MIR
│   └── II类（DNA转座子，cut-and-paste机制）
│       ├── TIR转座子（有末端反向重复）
│       │   ├── hAT
│       │   ├── Mariner/Tc1
│       │   └── MULE
│       └── Helitron（滚环复制）
│
├── 串联重复（Tandem Repeats）
│   ├── 卫星DNA（大块着丝粒/端粒区）
│   ├── 小卫星（10-60bp重复单元）
│   └── 微卫星/SSR（1-6bp重复单元，如ACACAC）
│
└── 低复杂度序列（如poly-A、AT-rich）

人类基因组重复序列组成：
  LINE: ~20%（其中L1最多）
  SINE: ~13%（主要是Alu）
  LTR: ~8%（主要是ERV）
  DNA转座子: ~3%
  其他重复: ~6%
  总计: ~50%

1.2 为什么要分析重复序列？¶

重复序列分析的重要性：

1. 基因组注释必需步骤
   基因预测前必须先屏蔽重复序列
   否则基因预测工具会把TE误判为基因

2. 进化研究
   TE插入/扩增可以追踪物种进化历史
   不同物种的TE组成差异很大

3. 功能研究
   TE可以作为调控元件（增强子/启动子）
   Alu元件与RNA编辑密切相关
   L1仍在人类基因组中活跃转座

4. 疾病相关
   新的TE插入可以导致基因突变
   TE激活与癌症和衰老相关
   内源性逆转录病毒与自身免疫疾病相关

5. 基因组组装
   重复序列是基因组组装的最大障碍
   了解重复结构有助于评估组装质量

二、RepeatMasker分析¶

2.1 安装¶

# 方法1：Conda安装（推荐）
conda install -c bioconda repeatmasker  # 安装RepeatMasker

# 方法2：Docker/Singularity容器
# Dfam-TETools容器包含所有工具
singularity pull docker://dfam/tetools:latest  # 下载容器

# 方法3：从源码安装
wget https://www.repeatmasker.org/RepeatMasker/RepeatMasker-4.1.9.tar.gz  # 下载
tar -xzf RepeatMasker-4.1.9.tar.gz  # 解压
cd RepeatMasker
perl ./configure  # 配置（需要搜索引擎如RMBlast）

# 安装搜索引擎RMBlast
conda install -c bioconda rmblast  # 安装RMBlast

# 检查安装
RepeatMasker -v  # 查看版本

2.2 基本使用¶

# RepeatMasker基本用法：用已知TE库搜索基因组中的重复序列

# ===== 标准运行 =====
RepeatMasker \
  -species human \              # 物种（自动选择相应TE库）
  -pa 8 \                       # 并行线程数
  -gff \                        # 输出GFF格式注释
  -dir output_dir/ \            # 输出目录
  -xsmall \                     # 用小写字母标记重复（而非N或X）
  genome.fa                     # 输入基因组FASTA

# 常用参数说明：
# -species: 物种名（human/mouse/rice/drosophila等）
# -lib: 使用自定义重复库（优先于-species）
# -pa: 并行线程数
# -gff: 输出GFF格式
# -xsmall: 小写屏蔽（便于后续分析仍能看到序列）
# -nolow: 不屏蔽低复杂度序列
# -no_is: 不屏蔽散在重复
# -s: 慢速但灵敏模式
# -q: 快速但不太灵敏模式

# ===== 输出文件说明 =====
# genome.fa.out: 主输出文件（每个重复的详细信息）
# genome.fa.out.gff: GFF格式注释
# genome.fa.masked: 屏蔽后的基因组（重复序列被替换为N）
# genome.fa.tbl: 统计摘要表

# ===== 使用自定义重复库 =====
# 如果分析非模式生物，先用RepeatModeler建库
RepeatMasker \
  -lib custom_repeats.fa \      # 使用自定义TE库
  -pa 8 \
  -gff \
  -xsmall \
  genome.fa

2.3 结果解读¶

# 查看统计摘要
cat output_dir/genome.fa.tbl

# 典型输出（人类基因组示例）：
# ==================================================
# file name: genome.fa
# sequences:             24
# total length:  3088269832 bp
# GC level:         40.89 %
# bases masked:  1515524231 bp ( 49.07 %)
# ==================================================
#                 number of      length    percentage
#                 elements*      occupied   of sequence
# --------------------------------------------------
# SINEs:            1728424    228855370 bp    7.41 %
#   ALUs            1181287    166523652 bp    5.39 %
#   MIRs             547137     62331718 bp    2.02 %
# LINEs:             919526    573282985 bp   18.56 %
#   LINE1            509699    467498994 bp   15.14 %
#   LINE2            349832     82919367 bp    2.69 %
# LTR elements:      477759    252783477 bp    8.19 %
#   ERVL             113795     42575456 bp    1.38 %
#   ERVL-MaLRs       194765     75003786 bp    2.43 %
#   ERV_classI       119069    105413329 bp    3.41 %
# DNA elements:      466879    102523277 bp    3.32 %
# Unclassified:            0            0 bp    0.00 %
# Total interspersed repeats:  1157445109 bp   37.49 %
# Small RNA:          10949      1701340 bp    0.06 %
# Satellites:         30665     73785752 bp    2.39 %
# Simple repeats:    816225     37802665 bp    1.22 %
# Low complexity:    126282      6451610 bp    0.21 %
# ==================================================

2.4 R语言分析RepeatMasker输出¶

# 解析RepeatMasker输出
library(tidyverse)

# 读取.out文件（跳过头部3行）
rm_output <- read.table("genome.fa.out",
                         skip = 3,            # 跳过头部
                         fill = TRUE,          # 填充缺失列
                         stringsAsFactors = FALSE)

# 设置列名
colnames(rm_output) <- c("SW_score", "perc_div", "perc_del",
                          "perc_ins", "query_seq", "query_begin",
                          "query_end", "query_left", "strand",
                          "repeat_name", "repeat_class", "repeat_begin",
                          "repeat_end", "repeat_left", "ID")

# 清理数据
rm_clean <- rm_output %>%
  mutate(
    length = as.numeric(query_end) - as.numeric(query_begin) + 1,  # 长度
    divergence = as.numeric(perc_div)  # 序列散度
  ) %>%
  filter(!is.na(length))  # 去除无效行

# 1. TE类别统计
class_summary <- rm_clean %>%
  group_by(repeat_class) %>%
  summarise(
    count = n(),                           # 数量
    total_bp = sum(length),                # 总碱基数
    mean_div = mean(divergence, na.rm = TRUE)  # 平均散度
  ) %>%
  arrange(desc(total_bp))

print(class_summary)

# 2. TE组成饼图
ggplot(class_summary %>% head(10),
       aes(x = "", y = total_bp, fill = repeat_class)) +
  geom_bar(stat = "identity", width = 1) +
  coord_polar("y") +
  labs(title = "Repeat Element Composition") +
  theme_void()

# 3. TE散度分布（反映TE年龄）
ggplot(rm_clean, aes(x = divergence, fill = repeat_class)) +
  geom_histogram(binwidth = 1, alpha = 0.6) +  # 直方图
  facet_wrap(~ repeat_class, scales = "free_y") +  # 按类分面
  labs(x = "Sequence Divergence (%)", y = "Count",
       title = "TE Age Distribution (Kimura Distance)") +
  theme_minimal()
# 散度低 = 最近插入的（年轻的TE）
# 散度高 = 古老的TE

# 4. 染色体上TE密度分布
te_density <- rm_clean %>%
  mutate(bin = as.numeric(query_begin) %/% 1e6) %>%  # 1Mb窗口
  group_by(query_seq, bin) %>%
  summarise(te_count = n(),
            te_bp = sum(length))

ggplot(te_density %>% filter(query_seq == "chr1"),
       aes(x = bin, y = te_bp / 1e6)) +
  geom_line(color = "steelblue") +
  labs(x = "Position (Mb)", y = "TE Coverage (Mb)",
       title = "TE Density along Chr1") +
  theme_minimal()

三、RepeatModeler2从头发现¶

3.1 运行RepeatModeler2¶

# RepeatModeler2用于从新基因组中发现未知的TE家族
# 特别适合非模式生物

# 安装
conda install -c bioconda repeatmodeler  # Conda安装

# 步骤1：建立基因组数据库
BuildDatabase \
  -name genome_db \             # 数据库名
  -engine rmblast \             # 搜索引擎
  genome.fa                     # 输入基因组

# 步骤2：运行RepeatModeler2
RepeatModeler \
  -database genome_db \         # 步骤1创建的数据库
  -pa 16 \                      # 并行线程
  -engine rmblast \             # 搜索引擎
  -LTRStruct                    # 启用LTR结构搜索（RepeatModeler2新功能）

# RepeatModeler2包含3个从头发现算法：
# 1. RECON — 基于序列自比对和聚类
# 2. RepeatScout — 基于种子扩展
# 3. LtrHarvest/Ltr_retriever — 结构化LTR逆转座子搜索（新增）

# 输出文件：
# genome_db-families.fa — TE家族的共识序列
# genome_db-families.stk — 种子比对文件
# 可以直接作为RepeatMasker的-lib输入

# 步骤3：用发现的库运行RepeatMasker
RepeatMasker \
  -lib genome_db-families.fa \  # 使用RepeatModeler的输出
  -pa 8 \
  -gff \
  -xsmall \
  genome.fa

3.2 TE分类¶

# RepeatModeler2自带分类器（RepeatClassifier）
# 但可能有些家族标记为"Unknown"

# 2025年新工具Terrier——深度学习TE分类器
# 比DeepTE、TERL、TEclass2分类准确率更高
# 覆盖97.1%的Repbase序列

# 使用Terrier对Unknown TE进行分类
python terrier.py \
  --input unknown_families.fa \    # Unknown TE序列
  --output classified_families.fa \ # 分类后的TE序列
  --model terrier_model.h5         # 预训练模型

# TEtrimmer（2025）——自动化TE边界修正
# 比EDTA和RepeatModeler2更快、更准确地注释TE
TEtrimmer \
  --input genome_db-families.fa \  # RepeatModeler输出
  --genome genome.fa \             # 参考基因组
  --output trimmed_families.fa     # 修正后的TE库

四、完整注释流水线¶

# 基因组重复序列注释的标准流程

# ===== 步骤1：从头发现新TE家族 =====
BuildDatabase -name genome_db -engine rmblast genome.fa
RepeatModeler -database genome_db -pa 16 -LTRStruct

# ===== 步骤2：合并自定义库和已知库 =====
# 把RepeatModeler发现的TE和Dfam/RepBase的已知TE合并
cat genome_db-families.fa Dfam_consensus.fa > combined_lib.fa

# ===== 步骤3：RepeatMasker屏蔽 =====
RepeatMasker \
  -lib combined_lib.fa \
  -pa 16 \
  -gff \
  -xsmall \
  -dir repeat_annotation/ \
  genome.fa

# ===== 步骤4：分析和可视化 =====
# 使用R脚本（见上面的分析代码）

# ===== 步骤5：用于基因注释 =====
# 把屏蔽后的基因组用于基因预测
# Augustus/BRAKER/MAKER等工具需要屏蔽后的基因组

常见报错与解决¶

报错信息	原因	解决方案
`RepeatMasker: search engine not found`	没安装RMBlast	conda install rmblast
`No repeat library found`	物种名不识别	用-lib指定自定义库
`RepeatModeler: out of memory`	基因组太大	增加内存或分染色体处理
`Unknown TE families`	RepeatClassifier无法分类	用Terrier深度学习重分类
`BuildDatabase: fasta error`	序列名有特殊字符	清理序列名（去除空格和特殊符号）
`Low masking percentage`	物种TE含量少或库不全	先用RepeatModeler建自定义库

速查表¶

# 重复序列注释流程
RepeatModeler2(从头发现) → 合并已知TE库(Dfam)
→ RepeatMasker(同源屏蔽) → 统计分析 → 基因注释

# RepeatMasker关键参数
-species: 物种名（使用Dfam/RepBase库）
-lib: 自定义TE库
-pa: 并行线程数
-gff: 输出GFF注释
-xsmall: 小写屏蔽
-s/-q: 灵敏/快速模式

# RepeatModeler2关键参数
-database: BuildDatabase创建的数据库名
-pa: 并行线程数
-LTRStruct: 启用LTR结构搜索（推荐）
-engine: rmblast或abblast

# TE分类方法
RepeatClassifier: RepeatModeler自带（基本）
Terrier(2025): 深度学习分类器（覆盖率97.1%）
TEtrimmer(2025): 自动化边界修正

# 人类基因组TE组成
LINE(~20%) > SINE(~13%) > LTR(~8%) > DNA(~3%)
总散在重复: ~45%
串联重复+低复杂度: ~5%
总重复: ~50%

# 数据库
Dfam: 开放TE家族数据库（推荐）
RepBase: 商业TE数据库（需要订阅）
NCBI RepeatDB: NCBI的重复序列资源

面试高频问题¶

Q1：转座元件I类和II类的区别是什么？¶

答：I类（逆转座子）通过"复制-粘贴"（copy-and-paste）机制转座——先转录为RNA，逆转录为cDNA，再插入基因组新位置，原位置的拷贝保留。包括LTR逆转座子（如ERV）和非LTR逆转座子（如LINE/SINE）。II类（DNA转座子）通过"剪切-粘贴"（cut-and-paste）机制——转座酶直接将DNA从原位置切出并插入新位置。I类每次转座会增加一个拷贝，所以I类TE在基因组中数量多；II类转座不增加拷贝数。

Q2：RepeatMasker和RepeatModeler的关系是什么？¶

答：两者互补。RepeatModeler2是从头发现（de novo）工具——从新基因组序列中发现未知的TE家族，生成TE共识序列库。RepeatMasker是同源搜索工具——用已知的TE库（Dfam数据库或RepeatModeler生成的库）搜索基因组中的重复序列并屏蔽。标准流程是先用RepeatModeler2发现物种特异的TE，再合并已知库后用RepeatMasker屏蔽。RepeatMasker只能找到库里有的TE，RepeatModeler2可以发现全新的TE家族。

Q3：为什么基因预测前要屏蔽重复序列？¶

答：①TE中含有的ORF和剪接信号会被基因预测工具误判为真实基因——例如LINE-1编码的ORF1p和ORF2p会被预测为蛋白编码基因；②串联重复和低复杂度序列会干扰序列比对——导致同源基因搜索产生虚假比对；③重复序列中的错配会影响基因模型的准确性。屏蔽方法有两种：硬屏蔽（用N替换）适合ab initio预测，软屏蔽（用小写字母）适合同源搜索（BLAST仍然可以匹配但会降低权重）。

Q4：TE散度分析能告诉我们什么？¶

答：TE散度（Kimura距离）反映TE插入基因组的时间。新插入的TE与共识序列（参考序列）几乎完全一致（散度低），随着时间推移积累突变导致散度增加。因此散度分布图相当于TE的"年龄谱"：①低散度峰代表近期活跃的TE爆发事件；②高散度的TE代表古老的插入。比较不同物种的TE散度谱可以揭示物种特异的TE扩增事件——如人类基因组中Alu元件有一个低散度峰（仍在活跃转座），而L2在高散度区（古老的已不活跃的TE）。

Q5：2025年TE分析有什么新进展？¶

答：①Terrier（2025年Briefings in Bioinformatics）——深度学习TE分类器，覆盖97.1%的Repbase序列，在非模式生物中表现优于DeepTE和TEclass2，解决了RepeatModeler输出中大量"Unknown"分类的问题；②TEtrimmer（2025年）——自动化TE边界修正工具，比EDTA和RepeatModeler2更快速和准确地注释TE的精确边界；③长读长测序（PacBio HiFi/ONT）大幅提升了重复区域的组装质量，使得全长TE的结构分析成为可能；④EDTA（Extensive de novo TE Annotator）整合了多种从头发现和同源方法，提供更全面的TE注释流水线。