怎么使用rmats进行可变剪切的分析

引言

可变剪切（Alternative Splicing, AS）是基因表达调控的重要机制之一，它使得一个基因能够产生多个不同的mRNA转录本，从而增加蛋白质的多样性。可变剪切在生物体的发育、细胞分化、疾病发生等过程中起着关键作用。因此，研究可变剪切事件对于理解基因表达的复杂性和疾病机制具有重要意义。

rMATS（Replicate Multivariate Analysis of Transcript Splicing）是一款用于分析RNA-seq数据中可变剪切事件的工具。它能够检测五种常见的可变剪切事件类型：外显子跳跃（Exon Skipping, SE）、内含子保留（Intron Retention, IR）、5’端可变剪切（Alternative 5’ Splice Site, A5SS）、3’端可变剪切（Alternative 3’ Splice Site, A3SS）和互斥外显子（Mutually Exclusive Exons, MXE）。rMATS通过比较不同条件下的RNA-seq数据，识别出显著差异的可变剪切事件。

本文将详细介绍如何使用rMATS进行可变剪切的分析，包括软件安装、数据准备、运行分析以及结果解读。

1. 软件安装

1.1 系统要求

rMATS支持在Linux和macOS系统上运行。为了确保rMATS的正常运行，系统需要满足以下要求：

• Python 2.7或3.x
• NumPy
• SciPy
• PySam
• R（用于生成统计图表）

1.2 安装步骤

1. 下载rMATS：首先，从rMATS的GitHub仓库（https://github.com/Xinglab/rmats-turbo）下载最新版本的源代码。

   git clone https://github.com/Xinglab/rmats-turbo.git
   cd rmats-turbo

1. 安装依赖：确保系统中安装了所需的Python库和R。

   pip install numpy scipy pysam

如果系统中没有安装R，可以通过以下命令安装：

   sudo apt-get install r-base

1. 编译rMATS：在rMATS的根目录下运行以下命令进行编译。

   make

编译完成后，rMATS的可执行文件将生成在rmats-turbo目录下。

2. 数据准备

2.1 RNA-seq数据

rMATS需要输入的RNA-seq数据为BAM格式。BAM文件是经过比对后的序列数据，通常由比对工具（如STAR、HISAT2等）生成。确保BAM文件已经按照样本条件分组，并且每个样本的BAM文件都包含对应的索引文件（.bai）。

2.2 参考基因组和注释文件

rMATS需要参考基因组序列文件和基因注释文件（GTF格式）。这些文件可以从公共数据库（如Ensembl、UCSC等）下载。

• 参考基因组文件：通常为FASTA格式（如hg38.fa）。
• 基因注释文件：通常为GTF格式（如hg38.gtf）。

3. 运行rMATS

3.1 基本命令

rMATS的基本命令格式如下：

python rmats.py --b1 b1.txt --b2 b2.txt --gtf hg38.gtf --od output_dir -t paired --readLength 150 --nthread 10

• --b1 和 --b2：分别指定两个条件下的BAM文件列表文件。b1.txt和b2.txt是文本文件，每行包含一个BAM文件的路径。
• --gtf：指定基因注释文件（GTF格式）。
• --od：指定输出目录。
• -t：指定测序类型，paired表示双端测序，single表示单端测序。
• --readLength：指定测序读长。
• --nthread：指定使用的线程数。

3.2 参数说明

• --b1 和 --b2：这两个参数分别指定两个条件下的BAM文件列表。每个列表文件包含多个BAM文件的路径，每行一个路径。例如：

  # b1.txt
  /path/to/sample1_rep1.bam
  /path/to/sample1_rep2.bam

  # b2.txt
  /path/to/sample2_rep1.bam
  /path/to/sample2_rep2.bam

• --gtf：指定基因注释文件（GTF格式），用于识别外显子和内含子区域。

• --od：指定输出目录，rMATS将在此目录下生成分析结果。

• -t：指定测序类型，paired表示双端测序，single表示单端测序。

• --readLength：指定测序读长，通常为150或100。

• --nthread：指定使用的线程数，加快分析速度。

3.3 运行示例

假设我们有两个条件下的RNA-seq数据，分别为condition1和condition2，每个条件有两个重复样本。BAM文件列表如下：

# condition1.txt
/path/to/condition1_rep1.bam
/path/to/condition1_rep2.bam

# condition2.txt
/path/to/condition2_rep1.bam
/path/to/condition2_rep2.bam

运行rMATS的命令如下：

python rmats.py --b1 condition1.txt --b2 condition2.txt --gtf hg38.gtf --od output_dir -t paired --readLength 150 --nthread 10

4. 结果解读

rMATS运行完成后，将在指定的输出目录下生成多个文件，主要包括：

• AS_events.MATS.JC.txt：包含所有检测到的可变剪切事件的详细信息，包括事件类型、基因名称、外显子/内含子位置、剪切变化量（Inclusion Level Difference, ILD）等。

• AS_events.MATS.JCEC.txt：与上述文件类似，但使用了更严格的过滤条件。

• SE.MATS.JC.txt：仅包含外显子跳跃（SE）事件的结果。

• RI.MATS.JC.txt：仅包含内含子保留（RI）事件的结果。

• A5SS.MATS.JC.txt：仅包含5’端可变剪切（A5SS）事件的结果。

• A3SS.MATS.JC.txt：仅包含3’端可变剪切（A3SS）事件的结果。

• MXE.MATS.JC.txt：仅包含互斥外显子（MXE）事件的结果。

4.1 结果文件格式

以SE.MATS.JC.txt为例，文件内容如下：

ID	GeneID	geneSymbol	chr	strand	exonStart_0base	exonEnd	upstreamES	upstreamEE	downstreamES	downstreamEE	ID	IJC_SAMPLE_1	SJC_SAMPLE_1	IJC_SAMPLE_2	SJC_SAMPLE_2	IncFormLen	SkipFormLen	PValue	FDR	IncLevel1	IncLevel2	IncLevelDifference
SE_1	ENSG001	GeneA	chr1	+	1000	2000	800	900	2100	2200	SE_1	10,20	5,10	30,40	15,20	100	50	0.01	0.05	0.6	0.7	-0.1

• ID：可变剪切事件的唯一标识符。
• GeneID 和 geneSymbol：基因的Ensembl ID和符号。
• chr 和 strand：染色体和链信息。
• exonStart_0base 和 exonEnd：外显子的起始和结束位置（0-based）。
• upstreamES 和 upstreamEE：上游外显子的起始和结束位置。
• downstreamES 和 downstreamEE：下游外显子的起始和结束位置。
• IJC_SAMPLE_1 和 SJC_SAMPLE_1：条件1中包含和跳过该外显子的junction reads数。
• IJC_SAMPLE_2 和 SJC_SAMPLE_2：条件2中包含和跳过该外显子的junction reads数。
• IncFormLen 和 SkipFormLen：包含和跳过该外显子的转录本长度。
• PValue 和 FDR：差异剪切事件的P值和FDR（False Discovery Rate）。
• IncLevel1 和 IncLevel2：条件1和条件2中的外显子包含水平（Inclusion Level）。
• IncLevelDifference：条件1和条件2之间的外显子包含水平差异。

4.2 结果可视化

rMATS还提供了R脚本用于生成可视化图表。在输出目录下，可以找到rmats_plot.R脚本。运行该脚本可以生成差异剪切事件的可视化图表。

Rscript rmats_plot.R --input-dir output_dir --output-dir plots

生成的图表将保存在plots目录下，包括外显子跳跃、内含子保留等事件的可视化结果。

5. 总结

rMATS是一款功能强大的工具，能够从RNA-seq数据中检测和分析可变剪切事件。通过本文的介绍，读者可以掌握rMATS的安装、数据准备、运行分析以及结果解读的基本流程。希望本文能够帮助研究人员更好地利用rMATS进行可变剪切的分析，从而深入理解基因表达的调控机制。