在上一期中我们得到hisat2的比对结果，现在一起来看下比对情况。
一，比对结果解读
比对结果文件如下，每一个样本一个排序后的bam文件，bam的索引bai以及存有比对信息的log日志。
image
选取一个样本的log查看比对详细情况：
cat SRR866991.log
58186286 reads; of these:
58186286 (100.00%) were unpaired; of these:
4153900 (7.14%) aligned 0 times
35063517 (60.26%) aligned exactly 1 time
18968869 (32.60%) aligned >1 times
92.86% overall alignment rate
我们可以看到这个样本的总比对率为92.86%，这是非常重要的一个指标，一般来说，样本的总比对率需要达到80%以上，70%多也能接受。如果你的总比对率太低如20%以下，有可能有以下几种情况：
- 1.参考基因组序列与fastq文件中的reads不是来源于同一个物种，要么是参考基因组用错了，要么是fastq文件搞错了。
- 2.序列中含有其他序列如细菌等污染。
这个时候可以对你的fastq序列抽取其中一部分reads做NCBI的nr/nt库的blast比对来确定fastq中的序列来源。
比对log中其他信息：
- 第一行：58186286 reads; of these，统计了这个样本的总reads数
- 第二行：58186286 (100.00%) were unpaired; of these:，read是单端
- 第三行：4153900 (7.14%) aligned 0 times，没有比对上参考基因组序列的read数和比例
- 第四行：35063517 (60.26%) aligned exactly 1 time，比对上参考基因组并且只比对上一个位置的序列。
- 第五行：18968869 (32.60%) aligned >1 times，比对上参考基因组且比对到多个位置。
- 第六行：92.86% overall alignment rate，总比对率
清楚了每一行的信息之后，我们来算一下总比对率是如何算的：
>((58186286-4153900)/58186286)*100 = 92.86%
也就是说总比对率其实包含了多比对read的部分。
将我们的log日志进行整理得到总比对的表格或者一个html文档
># multiqc可以将log日志中的数据进行可视化方便一起查看多个样本的比对结果。
multiqc *log
image
所有样本的总比对，有一个样本的比对率有些低，可以回去查查这个样的的质控等情况。
image
详细比对结果：
image
综合以上的情况，SRR866993号样本在这里的结果不是特别好，后面分析的有效数据会比较少。
二，bam转bw文件并绘热图
现在我们来将bam文件转成bw文件，可以使用IGV来查看比对结果，以及peak可视化（IGV教程参考隔壁公众号：保姆级 IGV 基因组浏览器使用指南（图文详解）https://mp.weixin.qq.com/s/Dtin-r85QRBexR2Rik9p2g）。
bam转bw有代码资料的那篇文献中代码使用的是最简单的参数：
image
这里我们加上--centerReads：可获得更陡峭的峰
># bam2bw
ls ../mapping/*bam | while read id
do
sam=${id##*/};sam1=${sam%%.*};
nohup bamCoverage -p 6 -b ${id} --centerReads --ignoreDuplicates -o ${sam1}.bw &
done
# 制作gene的bed文件
cat Mus_musculus.GRCm38.104.gtf |awk 'BEGIN{ FS=" ";OFS="\t" }{if($3 ~/gene/) print $1,$4-1,$5-1,$10,$6,$7}'|sed -e s'/"//'g -e s'/;//'g >Mus_musculus.GRCm38.104.bed
# computeMatrix
computeMatrix scale-regions -p 8 \
-b 3000 -a 3000 \
-R Mus_musculus.GRCm38.101.bed \
--regionBodyLength 6000 \
-S SRR*.bw \
-o All_Sample.TSS-body-TES_distribute.gz
# plot heatmap
plotHeatmap --dpi 200 --heatmapHeight 15 \
--matrixFile All_Sample.TSS-body-TES_distribute.gz \
--outFileName All_Sample.TSS-body-TES_Heatmap.png
image
三，IGV可视化查看文献中得到的peak情况
根据样本的表型信息，我们知道这个数据集总共有六个样本，每个样本有一个IP和一个Input。KO组三个生物学重复和WT组三个生物学重复，详细信息如下：
image
先去GEO上下载作者的两个peak结果，为bed格式：KO的peak结果
image
另外一个：WT的peak结果
image
在KO样本中随便找一个peak：
image
IGV可视化后，在这个大致的区间内确实是有peak的。
image
网友评论