SyRI --- 从全基因组组装中识别基因组重排和局部序列差异

发表于 更新于 分类于 ， 阅读次数： 本文字数： 4.3k 阅读时长 ≈ 4 分钟

Introduction

首先，了解什么是结构变异？

名词解释：

  • Structural Variations (SVs)：结构变异
  • Translocation：易位
  • Inversions：倒位
  • Duplication：重复
  • Deletion：缺失
  • Insertion：插入
  • InDel：插入缺失
  • Copy Number Variants（CNVs）：拷贝数变异
  • Copy Number Polymorphism（CNP）：拷贝数多态性
  • Presence absence Variation（PAV）：获得与缺失变异
  • Genomic Imbalances：基因组失衡
  • Single nucleotide polymorphism（SNP）：单核苷酸多态性
  • whole-genome alignments（WGA）：全基因组比对

SVs 结构变异，包括 长度在 50bp 以上的长片段序列的插入或缺失（Indel）、染色体倒位、序列串联倍增、染色体内部或染色体之间的序列易位、拷贝数变异（CNV）以及一些形式更为复杂的变异。

Indel，插入缺失，指的是在基因组的某个位置上所发生的小片段序列的插入或者删除，其长度通常在 50bp 以下。

CNV 拷贝数变异，一般是指长度为由 kb 到 Mb 级别组成大片段序列的拷贝数增加或减少。

SNP，单核苷酸多态性，由单个核苷酸 A、T、G、C 的改变而引起的 DNA 序列的改变，造成个体之间基因组的多样性。SNP 位点的分布是不均匀的，在非编码区比在编码区更为常见；一般来说，自然选择倾向于保留最有利于遗传适应性的 SNP 位点。在人类的遗传变异中，约 90% 为 SNP 变异，也就是说在人类基因组中，每隔 100 至 300 个碱基就会存在一个 SNP 位点。

局部序列差异包括：SNP、Indel、SVs。

预测基因组差异的许多先进的方法都是利用针对参考序列的短读或长读比对，这样的比对可以高精度的找到局部序列差异，但准确预测结构差异（structural differences）仍然具有挑战性。全基因组组装能够识别复杂的重排，因为与原始序列 reads 相比，组装的重叠群通常更长且质量更高。

SyRI 全称为 “Synteny and Rearrangement Identifier”；一种用于染色体水平组装的成对全基因组比较工具，可以识别两个全基因组组装之间的结构和序列差异；首先寻找重排区域，然后搜索序列中的差异，这些差异因位于同线（syntenic）或重排（rearranged）区域而不同。

基因组可以在结构和序列上有所不同，如果高度相似的区域在不同基因组之间具有不同的拷贝数、位置或方向，则会出现结构差异。在这里，这些区域被称为重排区域，而所有保守区域都称为同线区域。相反，序列差异 是导致 SNP、Indel 等的核苷酸序列变化。

Installation

很幸运的是，官方将 SyRI 在 conda 提供支持，所以我们可以通过 conda 很轻松的安装环境：

1
2
3
4
5
6
7
8
9
10
11
conda create -n syri_env python=3.9
conda activate syri_env
conda install -c bioconda syri -y
conda install -c bioconda plotsr -y

$ python --version
Python 3.9.12
$ syri --version
1.5.4
$ plotsr --version
0.5.3

需要注意的是，syri 包并不携带 plotsr，需要通过 conda 而外 install。

更多详情：SyRI-github

Getting Started

根据官方提供的例子，将 SyRI 的整个分析流程走一遍，更多内容：SyRI-Pipeline。

设置变量

1
2
3
cwd="."     # Change to working directory
PATH_TO_SYRI="../syri/bin/syri" #Change the path to point to syri executable
PATH_TO_PLOTSR="../syri/bin/plotsr" #Change the path to point to plotsr executable

获取示例数据

1
2
3
4
5
6
cd $cwd
## Get Yeast Reference genome
wget ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCA/000/146/045/GCA_000146045.2_R64/GCA_000146045.2_R64_genomic.fna.gz
## Get Query genome
wget ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCA/000/977/955/GCA_000977955.2_Sc_YJM1447_v1/GCA_000977955.2_Sc_YJM1447_v1_genomic.fna.gz

解压并重链接

1
2
3
4
5
6
7
gzip -df GCA_000146045.2_R64_genomic.fna.gz
gzip -df GCA_000977955.2_Sc_YJM1447_v1_genomic.fna.gz
## Remove mitochondrial DNA
head -151797 GCA_000977955.2_Sc_YJM1447_v1_genomic.fna > GCA_000977955.2_Sc_YJM1447_v1_genomic.fna.filtered
ln -sf GCA_000146045.2_R64_genomic.fna refgenome
ln -sf GCA_000977955.2_Sc_YJM1447_v1_genomic.fna.filtered qrygenome

需要注意的是，SyRI 处理的是染色体水平的两个基因组组装，所以这里也期望输入两个染色体级别的组装，并要求同源染色体的编号（染色体 id）完全相同。如果情况并非如此，syri 会尝试使用全基因组比对来寻找同源基因组，但这种方法是启发式的，可能会导致次优结果。此外，建议两个基因组（fasta 文件）应具有相同数量的染色体。

运行全基因组比对（WGA）

1
2
# Using minimap2 for generating alignment. Any other whole genome alignment tool can also be used.
minimap2 -ax asm5 --eqx refgenome qrygenome > out.sam

建议用户测试不同的对齐设置，以找到适合其生物学问题的对齐分辨率。一些对齐工具会找到更长的对齐（有很多间隙），而其他对齐工具会找到更小更碎片化的对齐。较小的比对通常具有较高的比对同一性分数，并且更有助于识别较小的基因组结构重排。但它们也可能导致冗余对齐显着增加，从而导致对齐工具和 SyRI 的运行时间增加。

或者采用 MUMmer 获得 WGA。

Run SyRI

1
2
3
4
python3 $PATH_TO_SYRI -c out.sam -r refgenome -q qrygenome -k -F S
## OR
samtools view -b out.sam > out.bam
python3 $PATH_TO_SYRI -c out.bam -r refgenome -q qrygenome -k -F B

最终得到两个含有基因组结构差异的文件：syri.out、syri.vcf。

绘制 SyRI 预测的基因组结构 --- Plotsr

1
python3 $PATH_TO_PLOTSR syri.out refgenome qrygenome -H 8 -W 5

不幸的是，通过 conda 安装的 plotsr 运行上面命令进行可视化将会出校 \(\textcolor{red}{ERROR}\)。

更多详细内容：plotsr-github

首先，你需要准备一个 genomes.txt 文件，在有 WGA 的情况下运行 plotsr --sr syri.out -H 8 -W 5 --genomes genomes.txt -o output_plotsr.png。

1
2
3
4
5
$ cat genomes.txt
#file	name
<path>/refgenome	refgenome
<path>/qrygenome	qrygenome

根据官方提供的信息，可以在第三列添加 tags 标签，可以定义绘图的一些参数，其中第一列为基因组 fasta 文件的路径，第二列为对应的名称（呈现在图中）。这是一个 genomes.txt 示例：genomes.txt。

至此，整个分析流程算是简单完成了，但是绘图可以根据官方提供的参数设置，生成更加全面的可视化。

对于全基因组比对，可以选择 MUMmer

1
2
3
4
5
nucmer --maxmatch -c 100 -b 500 -l 50 refgenome qrygenome       # Whole genome alignment. Any other alignment can also be used.
delta-filter -m -i 90 -l 100 out.delta > out.filtered.delta     # Remove small and lower quality alignments
show-coords -THrd out.filtered.delta > out.filtered.coords      # Convert alignment information to a .TSV format as required by SyRI
python3 $PATH_TO_SYRI -c out.filtered.coords -d out.filtered.delta -r refgenome -q qrygenome
python3 $PATH_TO_PLOTSR syri.out refgenome qrygenome -H 8 -W 5

References

[1] https://zhuanlan.zhihu.com/p/161472051

[2] https://genomebiology.biomedcentral.com/articles/10.1186/s13059-019-1911-0

[3] https://github.com/schneebergerlab/syri

[4] https://github.com/schneebergerlab/plotsr

[5] https://schneebergerlab.github.io/syri/pipeline.html