知识中心 - 北京概普生物科技有限公司(GapTech)

QuickRNASeq——RNA-Seq数据处理及可视化工具

生信干货 sxr ·2018年8月28日 21:49

一、背景

RNA-Seq测序已被广泛应用于细胞转录组和表达水平的研究，同时RNA-Seq数据也可应用于鉴定某些癌症相关的融合基因等。但RNA-Seq数据量庞大，因此急需一款既能处理大量数据又能生成报告的工具，于是QuickRNASeq应运而生。

QuickRNASeq是一款用于大规模RNA-Seq数据分析和复杂数据集的实时交互式可视化应用。该应用程序可通过多种开源工具，以高效的为用户提供友好、易于共享和生成报告的使用体验。在此，我们介绍了在用户环境中运行QuickRNASeq所需的步骤，并演示了高级交互式可视化模块的一些基本但功能强大的实用程序。

二、软件

QuickRNASeq软件的运行基于5款开源软件。STAR软件将reads比对到参考基因组或组装转录组；Subread包里的FeatureCounts软件用于对比对到基因、外显子、启动子、基因组特征区域的reads计数；VarScan用于检测变异；RSeQC用于控制RNA-seq测序质量；Samtools用于处理SAM格式的比对结果。

软件下载地址

STAR：https://github.com/alexdobin/STAR/releases

Subread：http://subread.sourceforge.net/

VarScan：http://varscan.sourceforge.net/

RSeQC：http://rseqc.sourceforge.net/

Samtools：http://sourceforge.net/projects/samtools/files/

QuickRNASeq v1.2版：https://sourceforge.net/projects/quickrnaseq

三、方法

QuickRNASeq流程分为3个主要步骤：

1)准备RNA-seq输入数据，包括样品配置文件

2)分别处理样品数据，包括比对、计数和QC分析

3)整合所有样品的分析结果，最后生成报告。

1、以Gencode release23版本的GRCh38基因组为例，为RNA-seq分析建立索引。

1.1、下载基因组文件

ftp://ftp.sanger.ac.uk/pub/gencode/Gencode_human/release_23/GRCh38.primary_assembly.genome.fa.gz

1.2、下载GTF格式的基因注释文件

ftp://ftp.sanger.ac.uk/pub/gencode/Gencode_human/release_23/gencode.v23.annotation.gtf.gz

1.3、genome和gtf文件解压并重命名

解压并重命名为 GRCh38.primary.genome.fa和GRCh38.gencode.v23.gtf；此例中，基因组文件放在/opt/fasta目录下，GTF文件放在 /opt/gencode 目录下

1.4、利用QuickRNASeq生成注释文件和BED文件，进入/opt/gencode目录，执行命令

/opt/ngsapp/QuickRNASeq/gtf2bed.plGRCh38.gencode.v23.gtf > GRCh38.gencode.v23.bed/opt/ngsapp/QuickRNASeq/gtf2annot.pl GRCh38.gencode.v23.gtf > GRCh38.gencode.v23.annot/opt/ngsapp/bin/samtools faidx GRCh38.primary.genome.fa

1.5、建立基因组索引

在/opt/STAR/ 目录下创建文件夹GRCh38_gencode23_100 ，进入该目录，建立基因组索引。其中sjdbOverhang 参数通常设置为read长度-1，而此处设置为99，意为只给长度在100bp以上的reads建立索引，命令：

/opt/ngsapp/STAR_2.4.0k/bin/Linux_x86_64/STAR--runThreadN 32 --runMode genomeGenerate -- genomeDir/opt/STAR/GRCh38_gencode23_100 –genomeFastaFiles/opt/fasta/GRCh38.primary.genome.fa –sjdbGTFfile /opt/gencode/GRCh38.gencode.v23.gtf-- sjdbOverhang 99

1.6、找到MHC基因所在的染色体位置。MHC基因具有高度多态性，所以这段区域非常适合用于检验SNP一致性。STAR生成的索引目录中有一个坐标文件“chrNameLength.txt”，可以从中查找到6号染色体的坐标。此例中，坐标位置为1-170805979。以上步骤完成后，就可以在配置文件中设置参考基因的相关参数。

2、安装包中含有一个示例目录“TestIyRun”可供参考，共涉及到48个样品。

2.1、准备样品的注释文件和ID文件

2.1.1、注释文件

注释文件应为tab分割的text格式文件。前两列信息建议为”sample_id”和”subject_id”。示例的sample.annotation.txt文件包含“Run”, “subject_id”, “histological_type”和 “sex” 等信息。

2.1.2、ID文件

每个样品ID占一行。示例中allIDs.txt 文件列出了48个样品ID，如第一个ID为“SRR607214”。

2.2、准备每个样品的fastq文件

双端测序有两个fastq文件。命名格式为sample_id_1.fastq.gz 和sample_id_2.fastq.gz。例如，SRR607214样品的两个fastq文件分别为：SRR607214_1.fastq.gz和SRR607214_2.fastq.gz；

单端测序则只有一个fastq文件。

随后将此步文件存放于fastq文件夹中。

2.3、设置配置文件

run.config配置文件包含测序、基因组和软件相关信息。基因组和软件参数只有在软件更新或基因组、索引、注释文件路径发生改变时修改即可。测序信息相关的参数需要每次根据实际情况进行修改。

2.3.1、创建FASTQ_DIR目录存放fastq数据

2.3.2、设置fastq后缀

QuickRNASeq会自动为allID.txt中每个样品对应的fastq文件添加后缀：“_1.FASTQ_SUFFIX” 和“_2. FASTQ_SUFFIX”，因此allID.txt中的样品ID要与fastq文件名对应起来。示例中设置参数FASTQ_SUFFIX=fastq.gz，程序就会在FASTQ_DIR目录下寻找以fastq.gz为后缀的文件。

2.3.3、配置序列方向

STRAND=0代表无方向的RNA-seq；STRAND=1代表正向；STRAND=2代表反向。

2.3.4、设置测序深度

有两个测序深度的选项。当测序reads数在40-80 Mb之间时设置为“regular”，当高于100 Mb时设置为“deep”。例如，SEQUENCE_DEPTH=regular。

2.3.5、设置测序类型

设置是双端还是单端测序，如SEQUENCE_TYPE=pair

2.3.6、配置物种基因组索引和GTF文件路径

GENOME_FASTA=/opt/fasta/GRCh38.primary.genome.faGENOME_INDEX=/opt/STAR/GRCh38_gencode23_100GENOME_ANNOTATION=/opt/gencode/GRCh38.gencode.v23.annotGTF_FILE=/opt/gencode/GRCh38.gencode.v23.gtfBEDFILE=/opt/gencode/GRCh38.gencode.v23.bed CHR_REGION=chr6:1-170805979

2.4、配置环境变量

STAR_RNA=/opt/ngsapp/STAR_2.4.0k/bin/Linux_x86_64

FEATURECOUNTS=/opt/ngsapp/subread-1.4.6/binRSeQC=/opt/ngsapp/anaconda/bin VARSCAN_JAR=/opt/ngsapp/bin/VarScan.v2.4.0.jarSAMTOOLS=/opt/ngsapp/bin

2.5、运行QuickRNASeq

在运行软件之前，需要确认已安装R 3.1或更高版本，以及以下R包：ggplot2, edgeR, scales, and reshape2。并确保$PATH环境变量中包含了QuickRNASeq_1.2的路径。

执行star-fc-qc.sh运行比对、计数、QC和SNP检测分析。此步计算量大，建议在HPC集群上利用LSF调度进程。最终会为每个样品单独生成一个结果目录。如果没有HPC集群，可以在Linux平台中执行star-fc-qc.ws.sh脚本。

运行示例：

#配置环境变量
export QuickRNASeq={QuickRNASeq_installation_Directory}

#如exportQuickRNASeq=/opt/ngsapp/QuickRNASeq_1.2 export PATH=$QuickRNASeq:$PATH
star-fc-qc.sh allIDs.txt run.config

#若是单计算节点则执行以下命令

star-fc-qc.ws.sh allIDs.txt run.config

2.6、整合分析结果并生成报告

整合sample.annotation.txt文件中所有样品的分析结果。执行Rscript $QuickRNASeq /QuickRNASeq_html/getEnsemblAnno.R可生成GENE_ANNOTATION参数所需的基因注释文件。执行以下命令生成报告：

exportGENOME_ANNOTATION=/opt/gencode/hg19.gencode.v19.annot
export GENE_ANNOTATION=/opt/gencode/Ensembl_v75_hg19_Gencode.v19_human.txt.gznohup star-fc-qc.summary.sh sample.annotation.txt &> Results.log
#若是单计算节点则执行以下命令
star-fc.summary.sh sample.annotation.txt

2.7、输出文件

输出结果存放在Results目录下，包括7个html文件（gex.html, index.html*, longitudinal.html,qc_fc-counting-summary.html, qc_expr_count_RPKM.html, qc_overview.html, qc_star-mapping-summary.html）和3个文件夹（package,QC, summary）。打开Results目录下的index.html文件即可查看报告。Index.html报告中包含所有数据和图表的下载地址。

3、交互式绘图，示例图可见fig1

3.1、绘制单个基因的箱线图详见fig2图解。

3.2、绘制热图详见fig3图解。

3.3、将报告开放到github.com

3.3.1、按照fig4在网站https://github.com创建GitHub账号。

3.3.2、按以下命令将报告上传至GitHub知识库，其中红色字体需要根据实际情况更改

git clone https://github.com/username/RNASeq_1.git

cd RNASeq_1
git checkout --orphan gh-pages
cp -R path2result/Results/*

git add
git commit -a -m "Adding RNASeq_1 resultsfrom QuickRNASeq" git push origin gh-pages

访问网址查看报告：http://username.github.io/RNASeq_1

示例数据集在http://baohongz.github.io/QuickRNASeq

有生信分析请留言

TCGA | 小工具 | 数据库 |组装| 注释 | 基因家族 | Pvalue

基因预测 |bestorf | sci | NAR | 在线工具 | 生存分析 | 热图

舞台|基因组 | 黄金测序 | 套路 | 杂谈组装 | 进化 | 测序简史