SyRI：一款从组装的基因组中检测结构变异的实用软件
简介
SyRI，全称Synteny and Rearrangement Identifier，是一个使用 whole-genome assemblies（WGA）预测基因组之间差异的工具。全基因组比对结果用作SyRI的输入。SyRI可识别共线性区域，结构重排（倒置，易位和重复），局部变异（SNP，indels，CNV等）以及未比对上的区域。
SyRI使用一种空前的方法，它从识别最长的共线性区域开始。由于所有非共线性区域都对应于在两个基因组之间重排的区域，因此鉴定共线区域也会同时鉴定了所有结构重排。在该鉴定步骤之后，将所有非共线性区域分类为倒位，易位或重复。这种方法将SR识别的挑战性问题转换为SR分类相对容易的问题。SR : Synteny Rearrangement
此外，SyRI还可以识别所有共线性区域和结构重排区域内的局部变异。Local variations 包括小变异（如SNP和小indel），structural variations（如大indel，CNV（拷贝数变异）和HDR）。从比对结果中解析出短的变异，通过比较共线性区域或重排区域的连续排列之间的 overlaps 和 gaps 来预测结构变异。
一、安装指南
1. 安装SyRI的前提条件
- C / C ++编译器：g ++
- Python3.5和以下软件包：Cython，numpy，scipy，pandas（版本：0.23.4），python-igraph，biopython，psutil，pysam和matplotlib
所有这些软件包都可以通过anaconda获得，并且可以使用以下方式安装：
conda install cython numpy scipy pandas=0.23.4 biopython psutil matplotlib=3.0.0
conda install -c conda-forge python-igraph
conda install -c bioconda pysam
2. Whole genome aligner
SyRI使用整个基因组比对作为输入。用户可以使用任何比对软件。SyRI接受SAM / BAM格式或带有CIGAR字符串格式的tab-separated value比对输入。如果用户要使用MUMmer，则可以使用.delta文件代替CIGAR字符串。有关更多信息，请参见文件格式。
3. 安装指南
- 下载（或克隆）SyRI的存储库
git clone https://github.com/schneebergerlab/syri.git
- 打开目录
setup.py（我们将其称为cwd），然后在终端运行：
python3 setup.py install # Install syri
chmod +x syri/bin/syri syri/bin/chroder syri/bin/plotsr # Make files executable
所有可执行文件都在中cwd/syri/bin/。
二、工作实例
设定变量
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
转到工作目录并下载酵母 reference genome and a query assembly
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
解压缩基因组的预处理
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。因此，建议用户对fasta文件进行预处理，以确保同源染色体在对应于两个基因组的两个fasta文件中具有完全相同的id。如果不是这种情况，syri会尝试使用全基因组比对来找到同源基因组，但是这种方法可能导致得不到最佳结果。另外，建议两个基因组（fasta文件）应具有相同数量的染色体。
全基因组比对
# Using minimap2 for generating alignment. Any other whole genome alignment tool can also be used.
minimap2 -ax asm5 --eqx refgenome qrygenome > out.sam
注意：建议用户测试不同的比对设置，以找到适合其生物学问题的比对分辨率。一些比对工具会找到更长的比对（有很多gaps），而其他工具则会找到更小的零碎比对。较小的比对通常具有较高的比对 identity scores，并且在识别较小的基因组结构重排中更有用。但是它们也可能导致冗余比对的显著增加，从而导致 aligner 和SyRI的运行时间增加。
Ryn SyRI
python3 $PATH_TO_SYRI -c out.sam -r refgenome -q qrygenome -k -F S
或者
samtools view -b out.sam > out.bam
python3 $PATH_TO_SYRI -c out.bam -r refgenome -q qrygenome -k -F B
SyRI将输出基因组结构差异文件syri.out和syri.vcf。
绘制由SyRI预测的基因组结构
python3 $PATH_TO_PLOTSR syri.out refgenome qrygenome -H 8 -W 5
使用SyRI从MUMmer产生的全基因组比对中识别基因组重排
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
三、文件格式
输入文件格式
SyRI采用TSV格式的完整基因组比对坐标，并包含以下几列：
| 列号 | 值 | 类型 |
|---|---|---|
| 1 | 基因组A的起始位置（从1开始，包括起始位置） | 整型 |
| 2 | 基因组A的终止位置（从1开始，包括末端位置） | 整型 |
| 3 | 基因组B的起始位置（从1开始。包括起始位置。） | 整型 |
| 4 | 基因组B的终止位置（从1开始。包括结束位置。） | 整型 |
| 5 | 基因组A上的比对长度 | 整型 |
| 6 | 基因组B的比对长度 | 整型 |
| 7 | alingment identity （百分比，0-100） | 浮动 |
| 8 | 基因组A中的比对方向 aDir（总是1） | 整型 |
| 9 | B基因组中的比对方向 bDir（1表示定向比对，-1表示反向比对） | 整型 |
| 10 | 基因组A的染色体ID | 字符串 |
| 11 | 基因组B的染色体ID | 字符串 |
| 12 | 比对对应的CIGAR字符串（可选；“ =”表示match，“ X”表示mismatch，“ D”表示deletion，“ I”表示insertion） | 字符串 |
需要以multi-fasta格式提供基因组。或者，nucmer生成的.delta也代替CIGAR字符串来call SNP。
输出文件格式
SyRI以TSV格式和VCF文件格式输出结果。
TSV格式规格
| 列号 | 值 | 类型 |
|---|---|---|
| 1个 | 基因组A的染色体ID | 字符串 |
| 2 | 基因组A的起始位置（从1开始，包括起始位置） | 整型 |
| 3 | 基因组A的终止位置（从1开始，包括末端位置） | 整型 |
| 4 | 基因组A的序列（仅适用于SNP和indel） | 字符串 |
| 5 | 基因组B的序列（仅适用于SNP和indel） | 字符串 |
| 6 | 基因组B的染色体ID | 字符串 |
| 7 | 基因组B的起始位置（从1开始，包括起始位置） | 整型 |
| 8 | 基因组B的终止位置（从1开始，包括末端位置） | 整型 |
| 9 | Unique ID （注释类型+数字） | 字符串 |
| 10 | Parent ID （注释类型+数字） | 字符串 |
| 11 | 注释类型 | 字符串 |
| 12 | 拷贝状态（重复） | 字符串 |
在这里，注释类型含义如下：
| 注解 | 含义 | 注解 | 含义 |
|---|---|---|---|
| SYN | Syntenic region | SYNAL | Alignment in syntenic region |
| INV | Inverted region | INVAL | Alignment in inverted region |
| TRANS | Translocated region | TRANSAL | Alignment in translocated region |
| INVTR | Inverted translocated region | INVTRAL | Alignment in inverted translocated region |
| DUP | Duplicated region | DUPAL | Alignment in duplicated region |
| INVDP | Inverted duplicated region | INVDPAL | Alignment in inverted duplicated region |
| NOTAL | Un-aligned region | SNP | Single nucleotide polymorphism |
| CPG | Copy gain in genome B | CPL | Copy loss in genome B |
| HDR | Highly diverged regions | TDM | Tandem repeat |
| INS | Insertion in genome B | DEL | Deletion in genome B |
1 ##fileformat=VCFv4.3
2 ##fileDate=20200930
3 ##source=syri
4 ##ALT=<ID=SYN,Description="Syntenic region">
5 ##ALT=<ID=INV,Description="Inversion">
6 ##ALT=<ID=TRANS,Description="Translocation">
7 ##ALT=<ID=INVTR,Description="Inverted Translocation">
8 ##ALT=<ID=DUP,Description="Duplication">
9 ##ALT=<ID=INVDP,Description="Inverted Duplication">
10 ##ALT=<ID=SYNAL,Description="Syntenic alignment">
11 ##ALT=<ID=INVAL,Description="Inversion alignment">
12 ##ALT=<ID=TRANSAL,Description="Translocation alignment">
13 ##ALT=<ID=INVTRAL,Description="Inverted Translocation alignment">
14 ##ALT=<ID=DUPAL,Description="Duplication alignment">
15 ##ALT=<ID=INVDPAL,Description="Inverted Duplication alignment">
16 ##ALT=<ID=HDR,Description="Highly diverged regions">
17 ##ALT=<ID=INS,Description="Insertion in non-reference genome">
18 ##ALT=<ID=DEL,Description="Deletion in non-reference genome">
19 ##ALT=<ID=CPG,Description="Copy gain in non-reference genome">
20 ##ALT=<ID=CPL,Description="Copy loss in non-reference genome">
21 ##ALT=<ID=SNP,Description="Single nucleotide polymorphism">
22 ##ALT=<ID=TDM,Description="Tandem repeat">
23 ##ALT=<ID=NOTAL,Description="Not Aligned region">
24 ##INFO=<ID=END,Number=1,Type=Integer,Description="End position on reference genome">
25 ##INFO=<ID=ChrB,Number=1,Type=String,Description="Chromoosme ID on the non-reference genome">
26 ##INFO=<ID=StartB,Number=1,Type=Integer,Description="Start position on non-reference genome">
27 ##INFO=<ID=EndB,Number=1,Type=Integer,Description="End position on non-reference genome">
28 ##INFO=<ID=Parent,Number=1,Type=String,Description="ID of the parent SR">
29 ##INFO=<ID=VarType,Number=1,Type=String,Description="Start position on non-reference genome">
30 ##INFO=<ID=DupType,Number=1,Type=String,Description="Copy gain or loss in the non-reference genome">
Copy状态描述了基因组B（copygain，即基因组B具有额外的拷贝）或者基因组A的（copyloss，即基因组A具有额外的拷贝）的重复区域。
Parent ID对应于其中存在alignment或 local variation 的注释块中（共线性区域或结构重排）的unique ID。因此，如果在基因组A的Chr1：10和基因组B的Chr2：542有一个易位区域（unique ID TRANS1）存在A-> T SNP（unique ID SNP1），则相应的条目将为：
Chr1 10 10 A T Chr2 542 542 SNP1 TRANS1 SNP -
VCF格式
以上信息被转换为VCF（v4.3）文件格式，其中基因组A被视为参考基因组。但是，由于VCF基于参考基因组位置，因此我们不会在VCF文件中输出基因组B中未比对上的区域。
四、Identifying genomic differences using SyRI
SyRI要求染色体水平的基因组才能准确鉴定SR。如果没有染色体级别的，则可以使用chroder实用程序创建伪染色体级程序集。
全基因组比对
SyRI使用全基因组比对作为输入。这些可以使用用户选择的全基因组比对工具生成（请参阅安装指南。在这里，我们将使用MUMmer3软件包。
首先，使用NUCmer对基因组（多fasta格式）进行比对。
nucmer --maxmatch -c 500 -b 500 -l 100 refgenome qrygenome;
这里-c，-b和-l参数是控制比对分辨率的，需要根据基因组的大小和复杂程度进行调整。更多详细信息在这里。
NUCmer将生成一个out.delta文件作为输出。使用delta-filter对比对结果进行过滤，然后使用show-coords转换成一个制表符分隔的格式。
delta-filter -m -i 90 -l 100 out.delta > out_m_i90_l100.delta;
show-coords -THrd out_m_i90_l100.delta > out_m_i90_l100.coords;
为了适合特定的基因组和科学问题的值，可以更改-i和-l的值。可在此处获得更多信息。
为了识别结构重排（包括重复），不过滤overlapping alignments。在上面的示例中，--maxmatch（对于nucmer）可识别所有比对。-m参数（delta-filter）删除冗余比对，虽然它不是必要的，但它有助于减少比对数这和由SyRI所需时间和存储资源。最后，-THrd（对于show-coords）将.delta比对格式转换为SyRI所需的比对坐标组成的.tsv格式。
对于使用MUMmer3生成的比对，不需要CIGAR字符串。对于其他比对软件（ aligners），以SAM / BAM或CIGAR字符串存储的.tsv格式的全基因组比对结果都可以解析，用来鉴定SNP和小的indel（但是，无需CIGAR字符串即可鉴定结构重排和结构变异）。
使用 syri识别SR
SyRI将基因组比对坐标作为输入。另外，如果还需要鉴定结构变异，则还需要两个基因组的fasta文件。此外，对于 short variation identification，当没有CIGAR字符串时，则需要.delta文件（由NUCmer生成）。
用法和参数是：
usage: syri [-h] -c INFILE [-r REF] [-q QRY] [-d DELTA] [-F {T,S,B}] [-o FOUT]
[-k] [--log {DEBUG,INFO,WARN}] [--lf LOG_FIN] [--dir DIR]
[--prefix PREFIX] [--seed SEED] [--nc NCORES] [--novcf] [--nosr]
[-b BRUTERUNTIME] [--unic TRANSUNICOUNT] [--unip TRANSUNIPERCENT]
[--inc INCREASEBY] [--no-chrmatch] [--nosv] [--nosnp] [--all]
[--allow-offset OFFSET] [--cigar] [-s SSPATH]
Input Files:
-c INFILE File containing alignment coordinates in a tsv format
(default: None)
-r REF Genome A (which is considered as reference for the
alignments). Required for local variation (large
indels, CNVs) identification. (default: None)
-q QRY Genome B (which is considered as query for the
alignments). Required for local variation (large
indels, CNVs) identification. (default: None)
-d DELTA .delta file from mummer. Required for short variation
(SNPs/indels) identification when CIGAR string is not
available (default: None)
optional arguments:
-h, --help show this help message and exit
-F {T,S,B} Input file type. T: Table, S: SAM, B: BAM (default: T)
-o FOUT Output file name (default: syri)
-k Keep internediate output files (default: False)
--log {DEBUG,INFO,WARN}
log level (default: INFO)
--lf LOG_FIN Name of log file (default: syri.log)
--dir DIR path to working directory (if not current directory).
All files must be in this directory. (default: None)
--prefix PREFIX Prefix to add before the output file Names (default: )
--seed SEED seed for generating random numbers (default: 1)
--nc NCORES number of cores to use in parallel (max is number of
chromosomes) (default: 1)
--novcf Do not combine all files into one output file
(default: False)
SR identification:
--nosr Set to skip structural rearrangement identification
(default: False)
-b BRUTERUNTIME Cutoff to restrict brute force methods to take too
much time (in seconds). Smaller values would make
algorithm faster, but could have marginal effects on
accuracy. In general case, would not be required.
(default: 60)
--unic TRANSUNICOUNT Number of uniques bps for selecting translocation.
Smaller values would select smaller TLs better, but
may increase time and decrease accuracy. (default:
1000)
--unip TRANSUNIPERCENT
Percent of unique region requried to select
translocation. Value should be in range (0,1]. Smaller
values would selection of translocation which are more
overlapped with other regions. (default: 0.5)
--inc INCREASEBY Minimum score increase required to add another
alignment to translocation cluster solution (default:
1000)
--no-chrmatch Do not allow SyRI to automatically match chromosome
ids between the two genomes if they are not equal
(default: False)
ShV identification:
--nosv Set to skip structural variation identification
(default: False)
--nosnp Set to skip SNP/Indel (within alignment)
identification (default: False)
--all Use duplications too for variant identification
(default: False)
--allow-offset OFFSET
BPs allowed to overlap (default: 5)
--cigar Find SNPs/indels using CIGAR string. Necessary for
alignment generated using aligners other than nucmers
(default: False)
-s SSPATH path to show-snps from mummer (default: show-snps)
SR识别
如果两个assemblies的染色体ID不相同，SyRI会尝试寻找同源染色体，然后将其ID映射为相同。可以使用--no-chrmatch参数关闭此行为。
本节中的其他参数规定了如何识别易位和重复（TD）。对于重叠候选TD的小型网络，SyRI使用蛮力方法来找到TD的最佳集合。可以使用-b参数限制此方法运行的时间。如果对于网络，暴力破解方法花费的时间超过了分配的时间，则它将自动切换为随机贪婪方法。--unic和--unip参数说明how unique a candidate TD need to be。与共线性和inversions高度重叠且未通过这些阈值的candidate将被滤除。从候选TD网络，可以选择不同的候选集。使用--inc 阈值来确定一组新的候选者是否比当前候选者更好，因此可以选择它作为解决方案。
识别局部变异的参数
--allow-offset参数用于定义阈值，以决定在一个注释块内的两个连续alignment是否重叠。Alignments中，重叠碱基的数量多于--allow-offset就被鉴定为CNV（copyloss / copygain）。通过使用CIGAR字符串进行比对或使用MUMmer程序包中的show-snps（需要.delta文件），可以识别出小的变异（SNP，小的indel ）。用户可以在使用CIGAR时设置--cigar。如果show-snps不在环境路径中，则-s可以用来提供路径。默认情况下，SyRI不会报告重复区域内的小变异，因为它们缺乏区域之间的一对一映射，从而使小变异变得模棱两可。但是，用户可以设置--all 它将返回所有 annotated alignments中的所有的小变异。
五、Pseudo-genome assembly generation using chroder
从生物学上讲， syntenic regions are meaningful only for homologous regions.。因此，要找到结构重排（取决于对共线性区域的识别），两个基因组中的DNA分子必须彼此同源，并且可能在两个基因组分子之间进行一对一的映射。通常，这将要求基因组处于染色体级别。如果一个或两个assemblies不在染色体水平，则可以利用两个assemblies之间的同源性来生成假染色体。在这种情况下，可以使用chroder工具来生成假染色体。它可以使用不完整基因组与参考基因（染色体水平组装的）的同源性，以及两个不完全 assemblies之间的同源性来生成pseudo-genomes。因此，也可以分析非染色体级组装的物种。
chroder分析了两个不完整基因组组装之间的全基因组比对。当其中一个处于染色体水平时，则将来自另一个的 scaffolds 相对于染色体锚定。 scaffolds 的顺序是根据其在染色体中排列的位置确定的。类似地，假染色体中的scaffold方向取决于与染色体的对齐方向（scaffold 是否与正义序列或反义序列对齐）。具有反向排列的scaffold是反向互补的，因此所有scaffold相对于染色体都处于正向。
当两个assemblies 都处于scaffold水平时，chroder将与另一个assemblies中的 scaffold比对上的assemblies中的 scaffold group到一起，（相互之间都是这样来处理，可以结合下图理解）。这将可能来自染色体的相邻区域的序列连成了scaffolds链。然后对scaffolds进行排序和定向，然后用Ns连接以创建伪染色体样分子。
将基因组A的scaffolds （绿色）与基因组B的scaffolds （蓝色）对齐。与同一个scaffolds 对齐的scaffolds 组合在一起。例如，Scaf_B_1和Scaf_B_2都与Scaf_A_1对齐，因此认为它们来自于染色体的相邻区域。类似地，可以将Scaf_A_1和Scaf_A_2视为来自相邻区域。重复该过程以找到所有可能来自同一染色体的scaffolds，然后将它们连接起来生成假染色体
usage: chroder [-h] [-n NCOUNT] [-o OUT] [-noref] coords ref qry
positional arguments:
coords Alignment coordinates in a tsv format
ref Assembly of genome A in multi-fasta format
qry Assembly of genome B in multi-fasta format
optional arguments:
-h, --help show this help message and exit
-n NCOUNT number of N's to be inserted
-o OUT output file prefix
-noref Use this parameter when no assembly is at chromosome level
此方法具有高度的启发性，不应视为scaffolding 的替代方法。此外，生成基于同源性的假染色体将导致更高的假阴性，因为在假染色体生成期间将去除 assemblies 之间的一些结构重排。
六、绘制基因组结构
使用以下方法绘制基因组结构 plotsr
plotsr使用SyRI的输出，运用以下命令生成基因组结构图：
python /path/to/plotsr syri.out /path/to/refgenome /path/to/qrygenome
用法和可用参数：
usage: Arguments for plotting SRs predicted by SyRI [-h] [-s S] [-f F] [-H H]
[-W W] [-o {pdf,png,svg}]
[-d D]
reg r q
positional arguments:
reg syri.out file generated by SyRI
r path to reference genome
q path to query genome
optional arguments:
-h, --help show this help message and exit
-s S minimum size of a SR to be plotted
-f F font size
-H H height of the plot
-W W width of the plot
-o {pdf,png,svg} output file format (pdf, png, svg)
-d D DPI for the final image
输出图示例：
使用plotsr生成的样例图
当前限制：
-
两个基因组中的同源染色体需要来自同一条链。如果染色体来自不同的链，则染色体之间的比对将被颠倒。由于SyRI找共线性区域然后检查同向的比对，因此它将无法找到共线性区域，并且可能会崩溃。
该问题的当前解决方案是手动检查alignments。如果同源染色体之间的大多数比对是反向的，则查询基因组中的染色体需要反向互补。然后，需要将校正后的查询基因组与参考基因组进行比对。我们正在研究一种可以生成点图以自动识别和反向互补此类反向染色体的方法。 -
大的异位和复制(consisting of multiple alignments) 可能导致较高的内存使用和CPU运行时间。
问题及解决方案
如直接使用python的pip命令下载pandas包
pip3.5 install pandas==0.23.4 -i http://pypi.douban.com/simple --trusted-host pypi.douban.com
cython scipy biopython
python-igraph
Ref
https://github.com/schneebergerlab/syri
https://schneebergerlab.github.io/syri/
https://pubmed.ncbi.nlm.nih.gov/31842948/
Goel, M., Sun, H., Jiao, W. et al. SyRI: finding genomic rearrangements and local sequence differences from whole-genome assemblies. Genome Biol 20, 277 (2019) doi:10.1186/s13059-019-1911-0
