共线性分析MCScanX¶

一句话概述：MCScanX是检测基因组间和基因组内共线性区块（synteny blocks）的标准工具，用于揭示基因组结构变异、全基因组复制事件和基因组进化历史。2024年Nature Protocols发布了最新使用指南。

核心知识点速查表¶

概念	说明
共线性(Synteny)	不同基因组区域中基因的排列顺序保持一致（白话：基因的排列顺序像排队一样对得上）
共线性区块	一组基因在两个基因组区域中顺序一致的片段
全基因组复制(WGD)	整个基因组加倍事件，植物中很常见
Anchor gene	共线性区块中的锚定基因对
Ks	同义替换率，用于估算复制时间
基因复制类型	WGD/串联/近端/分散/单体复制
Dot plot	点阵图，直观展示两基因组间的共线性
Circos图	环形图，展示染色体间的共线性关系

一、安装与输入准备¶

# === 安装MCScanX ===
conda install -c bioconda mcscanx    # conda安装
# 或从GitHub编译
git clone https://github.com/wyp1125/MCScanX.git   # 克隆源码
cd MCScanX && make                                    # 编译

# === 准备输入文件（2个文件） ===

# 文件1: xyz.gff - 基因位置文件
# 格式: chr gene_id start end
# 注意: 前缀必须与BLAST数据库中的基因ID匹配
# 示例:
# sp1_chr1  sp1_gene001  1000  2000
# sp1_chr1  sp1_gene002  3000  5000
# sp2_chr1  sp2_gene001  1500  2500

# === 从GFF3文件提取基因位置 ===
# 提取基因坐标（MCScanX格式）
awk '$3=="gene" {
    split($9, a, "[=;]");           # 按=和;分割第9列
    gsub(/.*ID=/, "", a[2]);        # 提取基因ID
    print $1, a[2], $4, $5          # 输出: 染色体 基因ID 起始 终止
}' species1.gff3 > sp1_genes.txt

awk '$3=="gene" {
    split($9, a, "[=;]");
    gsub(/.*ID=/, "", a[2]);
    print $1, a[2], $4, $5
}' species2.gff3 > sp2_genes.txt

# 合并为MCScanX输入格式
cat sp1_genes.txt sp2_genes.txt > xyz.gff    # 合并两个物种

# === 文件2: xyz.blast - 全对全BLAST结果 ===
# 蛋白质全对全比对（BLASTP或DIAMOND）

# 合并蛋白质序列
cat species1.pep.fa species2.pep.fa > all_proteins.fa    # 合并

# 使用DIAMOND进行全对全比对（比BLAST快500倍）
diamond makedb \
    --in all_proteins.fa \        # 输入蛋白质
    --db all_proteins              # 数据库名

diamond blastp \
    --db all_proteins \            # 数据库
    --query all_proteins.fa \      # 查询序列（自身比对）
    --outfmt 6 \                   # BLAST表格输出格式
    --evalue 1e-10 \               # E-value阈值
    --max-target-seqs 5 \          # 每个query保留前5个hit
    --threads 16 \                 # 线程数
    --out xyz.blast                 # 输出文件

二、运行MCScanX¶

# === 运行MCScanX主程序 ===
# 输入文件必须放在同一目录，且前缀相同（如xyz.gff和xyz.blast）
MCScanX xyz                        # 前缀名（会自动找xyz.gff和xyz.blast）

# === 主要参数 ===
MCScanX \
    -s 5 \              # 最少anchor基因对数（默认5，越小越宽松）
    -m 25 \             # anchor对之间最大间隔基因数（默认25）
    -e 1e-10 \          # E-value阈值
    -b 2 \              # 对比的对象数量（2=成对比较）
    xyz                  # 前缀名

# === 输出文件 ===
# xyz.collinearity     - ★共线性区块详细信息
# xyz.html             - 网页可视化结果
# xyz.tandem           - 串联复制基因对

# 查看共线性区块
head -30 xyz.collinearity
# 输出示例:
# ## Alignment 0: score=1234 e_value=0 N=25 sp1_chr1&sp2_chr3 plus
#   0-  0:  sp1_gene001  sp2_gene201  1e-100
#   0-  1:  sp1_gene002  sp2_gene202  1e-95
#   ...

三、基因复制类型分类¶

# === 使用duplicate_gene_classifier ===
# 将基因分为5种复制类型
duplicate_gene_classifier xyz    # 输入前缀

# 输出: xyz.pairs.class
# 复制类型分类:
# 1. WGD/segmental  - 全基因组复制/片段复制（共线性区块中的）
# 2. tandem         - 串联复制（相邻位置的重复基因）
# 3. proximal       - 近端复制（间隔<10个基因）
# 4. transposed     - 转座复制
# 5. dispersed      - 分散复制（不属于以上任何类型）

四、可视化¶

# === 使用MCScanX内置可视化 ===

# Dot plot（点阵图）
java -jar MCScanX/downstream_analyses/dot_plotter.jar \
    -g xyz.gff \           # 基因位置文件
    -s xyz.collinearity \  # 共线性结果
    -c dot.ctl \           # 控制文件
    -o dot_plot.png         # 输出图片

# Circle plot（环形图）
java -jar MCScanX/downstream_analyses/circle_plotter.jar \
    -g xyz.gff \
    -s xyz.collinearity \
    -c circle.ctl \
    -o circle_plot.png

# === 使用Python JCVI/MCscan可视化（更美观） ===
# pip install jcvi

# 第1步：准备bed和cds文件
python -m jcvi.formats.gff bed \
    --type=mRNA \
    --key=Name \
    species1.gff3 -o species1.bed    # 转换GFF3为bed格式

python -m jcvi.formats.gff bed \
    --type=mRNA \
    --key=Name \
    species2.gff3 -o species2.bed

# 第2步：提取CDS序列
python -m jcvi.formats.gff cds \
    species1.gff3 genome1.fasta -o species1.cds   # 提取CDS

# 第3步：成对共线性分析
python -m jcvi.compara.catalog ortholog \
    species1 species2 \
    --no_strip_names                                # 运行共线性分析

# 第4步：绘制共线性图
python -m jcvi.graphics.dotplot \
    species1.species2.anchors                       # 绘制dot plot

python -m jcvi.graphics.karyotype \
    seqids layout \
    species1.species2.anchors                       # 绘制染色体共线性图

五、Ks分析与WGD检测¶

# === 计算共线性基因对的Ks值 ===
# 用于检测全基因组复制事件

# 使用wgd工具
pip install wgd                     # 安装wgd

wgd dmd \
    species.pep.fa \                # 蛋白质序列
    -o wgd_output                    # 输出目录

wgd ksd \
    wgd_output/species.pep.fa.tsv \ # 同源基因对
    species.cds.fa \                 # CDS序列
    -o ks_output                     # 输出目录

# Ks分布图解读：
# 一个峰 = 一次WGD事件
# 峰的位置越小 = WGD事件越近期
# 多个峰 = 经历了多次WGD

六、面试高频考点¶

Q1: 共线性分析能告诉我们什么？¶

全基因组复制(WGD)：通过物种内共线性区块检测
基因组进化：不同物种基因组结构的保守与变异
染色体重排：倒位、易位、融合等结构变异
基因复制机制：区分WGD复制、串联复制等

Q2: MCScanX和JCVI的区别？¶

MCScanX：经典工具，功能全面，有基因复制分类功能
JCVI(MCscan)：Python实现，可视化更美观，适合发表图片
实际工作中经常配合使用：MCScanX做分析，JCVI做可视化

Q3: 如何区分WGD和串联复制？¶

WGD复制：基因对分布在不同染色体的共线性区块中，Ks值集中在某个峰
串联复制：基因紧密相邻排列在同一染色体上，通常功能相同或相似

常见报错与解决¶

报错	原因	解决方案
`No collinear blocks`	BLAST结果太少或参数太严	放宽-s参数到3，降低E-value
`Gene IDs don't match`	gff和blast中的基因ID不一致	检查ID前缀是否匹配
`Java heap space`	可视化时内存不足	`java -Xmx4g -jar ...`
`Too few anchors`	物种太远或基因组太碎	检查基因组质量/换更近缘物种

速查表¶

# === MCScanX 完整流程 ===
# 1. 准备gff和blast文件
# 2. 运行MCScanX
MCScanX -s 5 -m 25 xyz
# 3. 基因复制分类
duplicate_gene_classifier xyz
# 4. 可视化（JCVI更推荐）
python -m jcvi.compara.catalog ortholog sp1 sp2
python -m jcvi.graphics.dotplot sp1.sp2.anchors