RNA-seq分析：Step4（比对）

目录

前记

一、建立基因组索引

1、gtf文件转换

2、提取可变剪切和外显子信息

3、建立索引

二、reads比对

三、SAM文件处理

四、比对结果可视化

1、对bam文件建立索引

2、使用IGV可视化bam文件

​编辑 五、比对结果质量评估

六、后记

---

前记

预处理之后的fastq文件，可以使用HISAT2软件将reads比对至参考基因组上。

HISAT2软件将reads比对至参考基因组上的步骤如下：

1. 对参考基因组进行建索引。对参考基因组进行预处理，生成索引文件，包括基因组序列的FM索引和SA索引。

2. 将原始的reads进行预处理。包括去除接头序列、去除低质量的碱基，以及进行过滤，得到质量较高的、可靠的reads。

3. 利用建立的FM索引和SA索引，将处理后的reads与参考基因组进行比对。首先，通过FM索引，快速找到与reads匹配的可能位置；然后，通过SA索引，在找到的可能位置中，比对精确定位reads的位置。

4. 在reads比对完成后，进行后处理和优化。对比对的结果进行排序和去重，以及对那些多次比对的reads进行处理，得到最终的比对结果。

5. 对比对结果进行注释和分析。将reads比对的结果与基因组注释信息进行结合，得到基因的表达信息、SNP位点信息等。

这些步骤综合运用，可以高效地将reads比对至参考基因组上，并得到比对结果。

一、建立基因组索引

1、gtf文件转换

#将tair10.gff3转换位tair10.gtf文件
gffread tair10.gff3 -T -o tair10.gtf &

2、提取可变剪切和外显子信息

#提取可变剪切信息
hisat2_extract_splice_sites.py tair10.gtf > tair10.ss &

#提取外显子信息
hisat2_extract_exons.py tair10.gtf > tair10.exon &

3、建立索引

#hisat2-bulid建立索引，添加可变剪切和外显子信息，速度慢吃内存
hisat2-build --ss tair10.ss --exon tair10.exon  tair10.fa ~/rnaseq/tair10_genome/index/tair10 &

#直接建立索引，速度较快（本例采用）
hisat2-build tair10.fa ~/rnaseq/tair10_genome/index/tair10 &

将索引文件输出到index文件夹下，产生8个ht2索引文件。

二、reads比对

基因组索引完成之后，就可以进行比对了。

#使用hisat2将过滤后的fastq文件比对至参考基因组上
hisat2 -p -2 --dta -x ~/rnaseq/tair10_genome/index/tair10 -U ~/rnaseq/fastx_results/SRR3418005_filtered.fastq -S SRR3418005.sam &

输出得到比对后的sam文件。

三、SAM文件处理

使用samtools对sam文件进行排序并转换为bam文件。

samtools sort -@2 -m 200M -o SRR3418005.bam SRR3418005.sam &
#额外线程数为2，每个线程最大内存为200M

最终得到的sam和bam文件，为减少内存占用，可只保留bam文件用于后续分析。

四、比对结果可视化

1、对bam文件建立索引

#使用samtools建立索引
samtools index SRR3418005.bam SRR3418005.bai &

2、使用IGV可视化bam文件

#打开igv浏览器
igv.sh

可视化结果：

 五、比对结果质量评估

#检测SRR3418005的bam文件质量
qualimap rnaseq -bam SRR3418005.bam -outdir /home/cxgg/rnaseq/hisat2_results -outformat PDF:HTML --java-mem-size=1G 

基因组质量评估：（五）mapping法：5. 用软件QualiMap统计BAM文件 - 知乎 (zhihu.com)

六、后记

忘了一个重要事情，本文使用的转录组原始数据除以上SRR3418005之外，还有SRR3418006、SRR3418019和SRR3418020三个数据。其中，SRR3418005和SRR3418019是两个生物学重复，为实验组；SRR3418006和SRR3418020是两个生物学重复，为对照组。其他数据操作与此相同，在此不做展示。

先写到这里了，后续讲述转录本的拼接和整合的内容。

2023.8.22

----CXGG