RNA-seq数据分析完全指北-03：去除奇怪的RNA

1、GC双峰和重复序列来自哪里？

不同于我们常见的polyA富集方法，号称全转录组测序的rRNA depletion建库对于实验的要求更高，并且在建库过程中引入的我们并不想分析的序列也更多。

一个真正意义上的全转录组，包括哪些内容呢？

1. 编码RNA
2. 非编码RNA
   1. 非编码大RNA：lncRNA和rRNA
   2. 非编码小RNA：
      • tRNA
      • 核酶
      • 细胞器的RNA组分
      • 小分子RNA：
        • miRNA
        • piRNA、snRNA、snoRNA等等

但是在实际上，我们关心哪些RNA呢？主要是mRNA、lncRNA、miRNA以及circRNA，或者一些人会去关注snRNA等等。其余RNA对于普通研究者来说都不会去接触，但实际上这些我们不想去关注的RNA才是占总RNA比例最大的部分，也是引入重复序列和GC偏好最严重的序列。所以，在进行数据分析时，要对这些RNA进行去除。

当然，如果这就是你想研究的内容，那么就根据自己的课题，进行更加个性化的分析吧！

2、下载rRNA序列

jimmy曾在一篇推文中提到，去除rRNA可以去除GC双峰的右峰

2.1、进入NCBI的Nucleotide，输入txid9606[Organism:exp]

2.2、勾选上图中的rRNA，并按下图方式下载FASTA序列

3、Hisat2构建索引并输出未比对的fastq序列

3.1、构建索引

hisat2-build -p 4 rRNA.fasta rRNA

3.2、输出没有比对到rRNA的序列

for i in {48..53}
do
a0="hisat2 -x ~/reference/linux/hisat2/otherRNA/rRNA "
a1="-1 SRR111783${i}_1.fastq.gz "
a2="-2 SRR111783${i}_2.fastq.gz "
a3="--un-conc-gz ../2.rrRNA/SRR111783${i}_rmr_%.fq.gz -p 16 -S ../2.rrRNA/SRR111783${i}.sam"
echo $a0$a1$a2$a3 
done > rmRNA.sh

nohup bash rmRNA.sh &

3.3、fastq文件比较

可以看到，去除rRNA序列之后，fastq文件大小都减少了大约20%，也可以通过查看nohup.out查看细节。

3.4、再次质控并与初始质控文件比较

4、左峰是什么？

既然右峰是rRNA，那么左峰有没有可能是tRNA呢？

具体操作和前面类似，不再赘述，只看最后的结果

4.1、查看nohup.out文件

我惊了，不是tRNA，那是啥？

5、其他序列

其实还有一些结构性RNA需要去除，包括scRNA、SRP RNA还有Ribonuclease P RNA Component H1等，获得这些序列的方法类似，但是过程要更加繁琐一些，这里就不具体介绍了。我把整理好的otherRNA.fa（包括rRNA、tRNA和otherRNA）上传到了百度网盘，需要的读者可以自取。

链接：https://pan.baidu.com/s/11odAtO-tqWex4nIDxmxffg 提取码：lvgv

5.1、fastq文件对比

再次减少了10%~20%

5.2、再次质控并与初始质控文件比较

可以看到前5对fastq文件现在质量已经可以勉强使用了，但是最后一个文件仍然有很大的问题。这么奇怪的GC含量，会不会是有其他物种污染呢？