可变剪接（Alternative splicing；又称“选择性剪接”）是一种在真核生物中非常普遍的基因表达方式，具体表现为一个基因的外显子以不同的组合方式剪接形成不同的成熟RNA，从而在不同的时空环境和状态下形成不同的蛋白质，执行不同的生物学功能。常见的可变剪接软件包括rMATS，Asprofile以及miso等。本文主要介绍rMATS软件的使用，并对结果利用rmats2sashimiplot可视化。
rMATS是一个从RNA-Seq数据中检测差异选择性剪接事件的计算工具，根据RNA-Seq数据，rMATS可以自动检测和分析与所有主要类型的可变剪接模式相对应的可变剪接事件。rMATS可识别的可变剪切事件有5种（Figure1，http://rnaseq-mats.sourceforge.net/）:
1）Skippedexon (SE) 外显子跳跃
2）Alternative5' splice site (A5SS) 5’端可变剪切
3）Alternative3' splice site (A3SS) 3’端可变剪切
4）Mutuallyexclusive exons (MXE) 互斥可变外显子
5）Retainedintron (RI) 内含子保留
Figure 1：rMATS可识别的可变剪切类型
1.分析所需软件
1）Linux操作系统
2）rMATS
3）rmats2sashimiplot
安装使用rMATS需预先安装如下软件：
4）InstallPython 2.7.x and corresponding versions of NumPy and SciPy
5）Downloadand install pysam (rMATS was tested with v0.9.1.4)
6）Downloadand install samtools (version 1.2 or later)
7）Downloadand install STAR (version 2.5 or later)
2.软件下载及安装
1）Python 安装
代码语言：javascript复制下载：https://www.python.org/
下载后得到文件：Python-2.7.12.tgz，安装过程如下：
tar-zxvf Python-2.7.12.tgz
cd /software/python/Python-2.7.12
./configure--prefix=/software/python/python2712install
make
makeinstall
vim ~/.bashrc
exportPATH="$PATH:/software/python/python2712install/bin";
source ~/.bashrcPython安装完成后，需要安装如下python模块，python模块一般均可从https://pypi.org/网站下载得到。
NumPy，SciPy，pysam
下面以python setup.py install安装的方式介绍Python模块的安装过程（需要预先安装setuptools模块）。NumPy和pysam模块可以选择该安装方式。
首先下载相关模块文件，以NumPy为例：
代码语言：javascript复制下载：https://pypi.org/project/scipy/#files.下载后得到文件numpy-1.16.4.zip，安装过程如下：unzip numpy-1.16.4.zip
cd numpy-1.16.4/
python2.7 setup.py install在安装SciPy时，由于其依赖模块较多，建议下载.whl格式文件，并通过pip install的方式安装（需预先安装pip模块）。
代码语言：javascript复制下载：scipy-0.19.0-cp27-cp27m-manylinux1_x86_64.whl
pip2.7 install scipy-0.19.0-cp27-cp27m-manylinux1_x86_64.whl2）samtools
代码语言：javascript复制下载：https://sourceforge.net/projects/samtools/files/安装：
代码语言：javascript复制tar -xjf samtools-1.3.1.tar.bz2
cd samtools-1.3.1/
./configure--prefix=/software/samtools/samtools-1.3.1-install
make
makeinstall
vim ~/.bashrc
exportPATH="$PATH:/software/samtools/samtools-1.3.1-install/bin";
source ~/.bashrc3）STAR
代码语言：javascript复制下载：https://github.com/alexdobin/STAR安装：
代码语言：javascript复制unzip STAR-master.zip
vim ~/.bashrc
exportPATH="$PATH:/software/STAR/STAR-master/bin/Linux_x86_64";
source ~/.bashrc4）rMATS
代码语言：javascript复制下载：http://rnaseq-mats.sourceforge.net/下载后得到rMATS.4.0.2.tgz文件
安装：
代码语言：javascript复制cd /software/rMATS
tar -zxcf /software/rMATS
cd rMATS.4.0.2可以发现该文件夹下出现了两个文件（rMATS-turbo-Linux-UCS2 和rMATS-turbo-Linux-UCS4），如何确定使用哪个文件夹下的文件呢？首先进入python控制台，然后输入以下内容：
代码语言：javascript复制>>>import sys
>>>print sys.maxunicode
1114111
Thisoutput indicates that your python is built with --enable-unicode=ucs4, and youshould use rMATS-turbo-xxx-UCS4.
>>>import sys
>>>print sys.maxunicode
65535This output indicates that your python is built with --enable-unicode=ucs2, and youshould use rMATS-turbo-xxx-UCS2.
笔者在操作过程中，操作系统出现的结果是1114111，因此使用rMATS-turbo-Linux-UCS4下的文件。测试一下是否可以运行：
代码语言：javascript复制python/software/rMATS/rMATS.4.0.2/rMATS-turbo-Linux-UCS4/rmats.py --help5）rmats2sashimiplot
代码语言：javascript复制下载：https://github.com/Xinglab/rmats2sashimiplot/安装：
代码语言：javascript复制unzip rmats2sashimiplot-master.zip
cd rmats2sashimiplot-master/
python2.7 setup.py install在运行过程中发现还需要安装matplotlib模块，该模块的安装可参照上述.whl格式文件安装过程，在此不在进一步描述。
3.运行rMATS
一般情况下，新一代测序数据文件都非常大，本次以官网测试数据为例。
http://rnaseq-mats.sourceforge.net/rmats4.0.2/下载测试数据testData.tgz（Figure 2）
Figure 2：testData.tgz数据文件
Fastq和bam文件均为二代测序过程中常见的文件，此处不做过多描述，主要看一下txt格式文件中的信息。
代码语言：javascript复制cat s1.txt
231ESRP.25K.rep-1.R1.fastq:231ESRP.25K.rep-1.R2.fastq,231ESRP.25K.rep-2.R1.fastq:231ESRP.25K.rep-2.R2.fastqcat s2.txt
231EV.25K.rep-1.R1.fastq:231EV.25K.rep-1.R2.fastq,231EV.25K.rep-2.R1.fastq:231EV.25K.rep-2.R2.fastq1）通过fastq文件开始运行：
代码语言：javascript复制python2.7 /software/rMATS/rMATS.4.0.2/rMATS-turbo-Linux-UCS4/rmats.py --s1 s1.txt --s2 s2.txt --gtf Homo_sapiens.GRCh37.75.gtf --bi STARindexFolder -od outDir -tpaired --readLength 50 --libType fr-unstranded --nthread 102）通过bam文件开始运行：
代码语言：javascript复制python2.7 /software/rMATS/rMATS.4.0.2/rMATS-turbo-Linux-UCS4/rmats.py --b1 b1.txt --b2 b2.txt --gtf Homo_sapiens.GRCh37.75.gtf --od outDir -t paired --readLength 50--libType fr-unstranded --nthread 10参数详解：
代码语言：javascript复制--s1 s1.txt A text file contains FASTQfile(s) for the sample_1.文件内以逗号分隔重复样本的FASTQ文件名(Onlyif using fastq)
--s2 s2.txt A text file contains FASTQfile(s) for the sample_2.文件内以逗号分隔重复样本的FASTQ文件名(Onlyif using fastq)
--b1 b1.txt A text file records mappingresults for the sample_1 in bam format. 文件内以逗号分隔重复样本的bam文件名(Only if using bam)
--b2 b2.txt A text file records mappingresults for the sample_2 in bam format. 文件内以逗号分隔重复样本的bam文件名(Only if using bam)
-treadType Type of read used in theanalysis. readType is either 'paired' or 'single'. 'paired' is for paired-enddata and 'single' is for single-end data
--readLength 测序reads的长度
--gtf gtfFile 需要输入的gtf格式文件
--bi STARIndexFolderThe folder name of the STARbinary indexes (i.e., the name of the folder that contains SA file). Forexample, use ~/STARindex/hg19 for hg19. (Only if using fastq)
--odoutDir 所有输出文件的路径
--tophatAnchorThe "anchor length" or "overhang length" used in thealigner. At least “anchor length” NT must be mapped to each end of a givenjunction. The default is 6. (This parameter applies only if using fastq)
--nthreadThe number of thread.
--cstat The cutoff splicing difference. The cutoffused in the null hypothesis test for differential splicing. The default is0.0001 for 0.01% difference. Valid: 0 ≤ cutoff < 1
--tstat The number of thread for statistical model.
--statoff Turn statistics part off.
-libTypelibraryTypeLibrary type. Default isunstranded (fr-unstranded). Use fr-firststrand or fr-secondstrand forstrand-specific data.软件运行后主要产生如下结果文件（Figure 3）：rMATS的结果文件是以各个可变剪切事件的分布的，主要由以下几类构成，详细可参考http://rnaseq-mats.sourceforge.net/rmats4.0.2/user_guide.htm.
AS_Event.MATS.JC.txtevaluates splicing with only reads that span splicing junctions.
AS_Event.MATS.JCEC.txtevaluates splicing with reads that span splicing junctions and reads on target(striped regions on Figure 1).
fromGTF.AS_Event.txtall possible alternative splicing (AS) events derived from GTF and RNA.
JC.raw.input.AS_Event.txtevaluates splicing with only reads that span splicing junctions.
JCEC.raw.input.AS_Event.txtevaluates splicing with reads that span splicing junctions and reads on target(striped regions on Figure 1).
其中JC和JCEC的区别在于JC考虑跨越剪切位点的reads，而JCEC不仅考虑前者的reads还考虑到比对到Figure 1中条纹的区域（也就是说没有跨越剪切位点的reads），一般情况使用JC的结果就够了（如果只是单纯的比较两组样品间可变剪切的差异的话）
Figure 3：rMATS结果文件
这几类文件中我们需要重点关注的是AS_Event.MATS.JC.txt，因为其他文件的信息大多数都整合在该文件中，以SE.MATS.JC.txt为例，该文件主要包含如下各列信息：
代码语言：javascript复制1）ID
2）GeneID
3）geneSymbol
4）chr
5）strand
6）exonStart_0base
7）exonEnd
8）upstreamES
9）upstreamEE
10）downstreamES
11）downstreamEE
12）ID
13）IJC_SAMPLE_1：样品在inclusion junction下的count数，重复样本的结果以逗号分隔
14）SJC_SAMPLE_1：样品在skipping junction counts（SJC）下的count数，重复样本的结果以逗号分隔
15）IJC_SAMPLE_2：样品在inclusion junction下的count数，重复样本的结果以逗号分隔
16）SJC_SAMPLE_2：样品在skipping junction counts（SJC）下的count数，重复样本的结果以逗号分隔
17）IncFormLen：lengthof inclusion form, used for normalization
18）SkipFormLen：lengthof skipping form, used for normalization
19）PValue：Significanceof splicing difference between two sample groups（两组样品可变剪切的统计学显著差异指标）
20）FDR：FalseDiscovery Rate calculated from p-value（校正后的PValue）
21）IncLevel1：inclusionlevel for SAMPLE_1 replicates (comma separated) calculated from normalizedcounts
22）IncLevel2：inclusionlevel for SAMPLE_2 replicates (comma separated) calculated from normalizedcounts
23）IncLevelDifference：average(IncLevel1)– average(IncLevel2)4.利用rmats2sashimiplot可视化
代码语言：javascript复制rmats2sashimiplot --b1 231ESRP.25K.rep-1.bam,231ESRP.25K.rep-2.bam --b2 231EV.25K.rep-1.bam,231EV.25K.rep-2.bam -t SE -e /software/rMATS/test/outDir/SE.MATS.JC.txt --l1 231ESRP.25K --l2 231EV.25K -o SE_plot注意：各个文件中染色体个格式均需保持一致，要么使用1,2,3......21,22,X,Y，要么使用chr1,chr2,chr3......chr21,chr22,chrX,chrY
参数详解：
代码语言：javascript复制--b1 B1 sample_1 in bam format(s1_rep1.bam[,s1_rep2.bam])
--b2 B2 sample_2 in bam format(s2_rep1.bam[,s2_rep2.bam])
-t rMATS结果中产生的可变剪切类型{SE,A5SS,A3SS,MXE,RI}
-e EVENTS_FILE The rMATS output event file (Onlyif using rMATSformat result as event file).
--l1 L1 The label for first sample.
--l2 L2 The label for second sample.
-o OUT_DIR The output directory.程序产生的图片如下所示（Figure 4，Figure5，Figure6）：
Figure 4：HRAS基因可变剪切
Figure 5：RPLP2基因可变剪切
Figure 6：RNH1基因可变剪切
好了，可变剪切分析和可视化就这样完成了，貌似有些复杂，不过只要你不抛弃不放弃，还是可以实现的。当然谁都不可能一次就成功，生信分析的坑还是不少的，这次不行，就再试一次，一步步的跨过这些坑才会成为生信大神！
