高通量测序数据分析：RNA-seq

本文围绕RNA-seq学习路线进行生信入门，主要内容有：

☆ RNA-seq方法原理
☆ RNA-seq的生物信息分析
1.数据获取
测序数据下载与处理（SRA Toolkit）
测序数据质控与过滤（fastp）
2.序列比对（SAMtools、HISAT2）
3.序列组装（StringTie、TACO）
4.表达定量和差异表达分析（Salmon、DESeq2）
5.GO和KEGG富集分析（clusterProfiler）

---

☆ RNA-seq方法原理

目的是要给mRNA测序，得到样本的基因表达信息。

• llumina的Truseq RNA建库方法：

带Poly（T）探针的磁珠与总RNA进行杂交，吸附其中的带Poly（A）尾巴的mRNA
Mg”离子溶液处理RNA，把RNA打成短片段 被打断的mRNA片段，用随机引物逆转出第一链的cDNA，再合成双链cDNA
在双链CDNA的两端加“A"碱基，并连上"Y“型的接头
经过PCR扩增，成为可以上机的文库

起始总RNA质量控制：用电泳方法。rRNA占有总RNA的大部分，形成的峰越高/尖，RIN（RNA完整度评分值）越高，8以上质量比较好。
测到的RNA片段 mapping到基因组上，进行样品的reads在参考基因上的分布均匀性(Gene coverage)统计。两端平衡的时候表示mRNA降解少（3’高降解多）。

---

☆ RNA-seq的生物信息分析

一、深度测序数据获取

和EBI、DDBJ组成INSDC，数据内容相同所以找NCBI就行。

（一）NCBI常用数据库

GenBank：遗传序列数据库，收集了所有公开的DNA序列及其注释 GEO (Gene Expression Omnibus)
：收集整理各种表达芯片数据，后来加入了甲基化、lncRNA、miRNA、CNV等其他芯片，还有高通量测序数据。文献中常见GSM和GSE开头的编号，分别是GEO
Sample和GEO Series的数据 PubMed / PMC (PubMed
Central)：前者把测序数据和文章联系起来，后者可以进行全文检索，无法访问校园网时可以替代Web of Knowledge
RefSeq：为所有常见生物提供非冗余、人工挑选过的参考序列，通常包含：参考基因组、参考转录组、参考蛋白序列、参考SNP信息、参考CNV信息等等

（二）测序数据的下载和处理：SRA Toolkit

1. 测序数据序列格式
   （1）FASTA：表示生物序列的文本格式，基因组和EST序列常常采用
   
   （2）FASTQ格式：表示生物序列及其质量的文本格式
   
   （3）ncbi SRA (Sequence Read Archive) ：存储高通量测序原始数据和比对信息，把FASTQ格式文件压缩为SRA格式
   
   绝大多数分析工具不支持SRA，需要使用配套工具包SRA Toolkit先行处理

1. SRA toolkit软件下载

在官网选择适合自己的版本下载。

#我选的ubuntu版本，其他一样，把下载链接修改一下就好了
wget http://ftp-trace.ncbi.nlm.nih.gov/sra/sdk/2.10.5/sratoolkit.2.10.5-ubuntu64.tar.gz  

用conda install sra-tools失败，只好用wget方法或者手动下载到linux盘符下。把安装包下载后用tar xzvf 解压，再配置完PATH就安装好了。
检查配置：

prefetch -V

2.用SRAtoolkit下载并处理NCBI数据

将 .sra文件转换为 .fstaq.gz文件的工具。用NCBI的SRR数据测试一下。
（1）下载
理论上下载东西都可以用wget，但是太慢了。单个数据下载还好，批量下载

prefetch SRRxxxxxxx -O .  #-O . 指定到当前路径，否则默认路径难找

一个数据下了好久，大概1个多小时。不知道怎么优化。

（2）解压

fastq-dump SRRxxxxxxx.sra #解压后从sra文件变为fastq文件

双端测序数据要加–split-files，否则解压后两端的数据不会分开，难以被其他软件读取 如果所用分析软件支持读取gzip，建议加上–gzip，将解压后的数据用gzip压缩，避免占用过多空间

fastq-dump --split-files --gzip xxx.sra

（三）测序数据质控与过滤： fastp

输出HTML和JSON报告，前者方便阅读，后者方便软件读取
单端：fastp -i raw.fq -o clean.fq
双端：fastp -i raw_1.fq -I raw_2.fq -o clean_1.fq -O clean_2.fq
有必要附加的参数：-l 36 -j xxx.json -h xxx.html

默认报告文件名 fastp.json 和 fastp.html，处理多个样本时极易互相覆盖，建议改为样本名称

fastp参数设置

 # I/O options 输入输出序列文件
  -i <单端-输入文件名>
  -o <单端-输出文件名>
  -I <双端-输入文件名>
  -O <双端-输出文件名>
  
#过滤后的最短序列长度
  -l 36  #默认15，建议设为36或40

# reporting options 报告参数
  -j   <the json format report file name >
  -h   <the html format report file name >
  -R   "report_title"

二、序列比对：HISAT2

• 注释格式介绍
  （1）GFF/GTF格式：一般用于基因组和基因注释
  （2）