这是RNA-Seq 上游分析的大致流程，比对+定量。当然实验目的若只需要定量已知基因，也可以选择free-alignment 的流程工具如kallisto/Salmon/Sailfish，其优点是可用于RNA-seq的基因表达的快速定量，但是对于小RNA和表达量低的基因分析效果并不好（2018年刚发表的一篇文章对free-alignment 的工具进行了质量评估，doi: https://doi.org/10.1101/246967 ）。基于比对的流程，比对工具也有很多选择，如Hisat，STAR,Tophat(hista可以替代tophat),bowtie等, 还有据说速度超快的Subread。根据比对工具的选择，定量软件也可以配套选择，如
STAR/bowtie—RSEM；
Hisat —featureCounts/HTseq；
Subread—featureCounts。
由于后续要鉴定新的转录本，因此需要对转录本组装，组装可以选择cufflinks,或者stringtie。文章是用的Tophat+cufflinks。我的后续分析打算用Hisat2比对，stringtie组装，featureCounts定量，同时打算尝试下Subread+featureCounts的流程。
Hisat2+stringtie+featureCounts；
Subread+featureCounts 。
sra格式的数据需要先用fastq-dump转换, --split-3 表示双端测序，--gzip将生成的fastq文件压缩。
质控检测采用fastqc和multiqc联合使用， multiqc有利于多个样本同时展示质量检测结果。FastQC的检测报告左侧会给出测序质量的一个summary,通常并不是所有的检测项都会是绿色通过，会有一些警告和fail的项目。主要可以从Per base sequence quality 看一下测序碱基质量，Per sequence GC content 看一下GC含量，如果实际的GC含量（红线）出现双峰，且导致后期的序列比对很低时，需要引起注意，有可能是存在样品污染；再者就是看一下Adapter是否去除及去除的是否干净。这里的Adapter虽然显示通过，但是去除的并不干净，所以后续还会进一步去除Adapter。
Trimmomatic 支持多线程，处理数据速度快，主要用来去除 Illumina 平台的 fastq 序列中的接头，并根据碱基质量值对 fastq 进行修剪。软件有两种过滤模式，分别对应 SE 和 PE 测序数据，同时支持 gzip 和 bzip2 压缩文件。另外也支持 phred-33 和 phred-64 格式互相转化。
Trimmomatic 过滤数据的步骤与命令行中过滤参数的顺序有关，通常的过滤步骤如下：
ILLUMINACLIP : 过滤 reads 中的 Illumina 测序接头和引物序列，并决定是否去除反向互补的 R1/R2 中的 R2, 该引物序列可以在Trimmomatic软件的安装目录下找到，双端通常选择TruSeq3-PE-2。
SLIDINGWINDOW: 从 reads 的 5’ 端开始，进行滑窗质量过滤，切掉碱基质量平均值低于阈值的滑窗。
MAXINFO : 一个自动调整的过滤选项，在保证 reads 长度的情况下尽量降低测序错误率，最大化 reads 的使用价值。
LEADING : 从 reads 的开头切除质量值低于阈值的碱基。
TRAILING : 从 reads 的末尾开始切除质量值低于阈值的碱基。
CROP : 从 reads 的末尾切掉部分碱基使得 reads 达到指定长度。
HEADCROP : 从 reads 的开头切掉指定数量的碱基。
MINLEN : 如果经过剪切后 reads 的长度低于阈值则丢弃这条 reads。
AVGQUAL : 如果 reads 的平均碱基质量值低于阈值则丢弃这条 reads。
TOPHRED33 : 将 reads 的碱基质量值体系转为 phred-33。
TOPHRED64 : 将 reads 的碱基质量值体系转为 phred-64。
参考 ： NGS 数据过滤之 Trimmomatic 详细说明
Hisat2 比对，首先用hisat2-build 构建index，然后比对，输出sam格式，再用samtools将sam转为bam格式，并排序构建索引。
featureCounts 目前已经整合在subread，subread是一个综合的软件，还包括比对的工具和对应的R包， featrueCounts的定量效果和HTSeq差不多，但是featureCounts的优点非常快！
二进制版本的软件解压即可使用
（PS：存储不够了，后续选择两组数据做流程分析。）
SRR4042230和SRR4042231的比对率分别为88.02%和88.81%，总比对时间将近5小时。
hisat_counts.txt.summary 包含一些基本的统计信息。
可以看出subjunc比对不到2小时，比hisat2快3个多小时。
根据第1列是Geneid，第7,8列是counts数，用awk提取出geneID和counts。
一个植物转录组项目的实战
史上最快的转录组流程-subread
转录组定量-FEATURECOUNTS