day26 ChIP-seq deeptools可视化查看rea
在chip_seq数据分析中,通常会对peak区域在基因组上的分布进行探究,查看其分布是否存在规律,比如是否在转录起始位点,或者转录终止位点附近存在富集,此时我们可以通过deeptools这个工具来实现。教程参考。
http://www.weinformatics.cn/4820741406-2/
https://zhuanlan.zhihu.com/p/166282791?utm_source=wechat_session
https://www.plob.org/article/11392.html
http://www.360doc.com/content/17/1218/22/19913717_714335774.shtml
1. computeMatrix 此工具依赖于之前步骤生成的bw文件,再给出bed或者gtf文件的区域(例如TSS),计算每个基因区域的结合得分,生成中间文件用以给plotHeatmap和plotProfiles作图。
reference-point mode:针对一个点,计算其周围信号分布。比如转录起始点。如果针对的是TSS转录起始点,那么对应的bed里面第二列的位置。bed里第一列是染色体号,第二列是转录起点,第三列是转录终点。
scale- regions mode: 针对一个区间,包含起点,终点。计算范围内的信号分布。
2. 单一bam文件从bam到bw,再用computeMatrix,最后用plotheatmap展示热图,用plotprofile展示折线图的脚本。
project=/data/zds209/ssresult
ls $project | if [ ! -d visual ]; then
mkdir -p $project/visual;
fi
ls $project/clean/*.rm.bam | while read id
do
file=$(basename $id)
sample=${file%%.*}
bamCoverage -b $id -o $project/visual/$sample.bw
computeMatrix reference-point -p 5 \
--referencePoint TSS \ #也可以是TES, centre
-b 3000 -a 3000\ # b上游,a下游,上游和下游多少 bp 之间的区间
-S $sample.bw \ #提供bigwig文件
-R /data/zds209/database/ucsc.gh38.bed\ #提供基因的注释信息(必须有)
--binSize 10 \
-o $project/visual/$sample.gz\
--outFileSortedRegions $project/visual/$sample.bed
plotHeatmap -m $project/visual/$sample.gz # -m是前面computeMatrix产生的.gz文件。
-out $project/visual/$sample.png \
--colorMap RdBu\
--whatToShow 'plot, heatmap and colorbar'\
done
名词:转录起始位点(TSS)和转录终止位点(TES)
3.获取基因组的bed、gtf注释文件
由于之前做单细胞cellranger,需要注释文件夹,查到里面有gene的gtf格式文件。但是这次需要bed格式的,要了解各种格式的意义用途。参考教程:
https://zhuanlan.zhihu.com/p/209019545
基因组注释文件gff和gtf
gff全称General featureformat,主要是用来注释基因组。gtf全称Gene transfer format,主要是用来对基因进行注释。
从ensemble下载的gtf文件前5行一般是以#开头的注释信息,而bed没有前面的注释,上来第一行就是数据。。
可以从UCSC上下载之,得到一个物种的所有基因的TST区域的bed文件。
简言之,http://genome.ucsc.edu/ 网站,tools里面选择table browser
human hg38的有28M,小鼠mm10有20M。
4.多个文件从bam,到bw,再用computeMatrix,最后用plotheatmap展示热图,用plotprofile展示折线图的脚本。
project=/data/zds209/ChIP-seqtest
ls $project | if [ ! -d visual ]; then
mkdir -p $project/visual;
fi
ls $project/bam/*.sorted.bam | while read id
do
file=$(basename $id)
sample=${file%%.*}
#bamCoverage -b $id -o $project/visual/$sample.bw
computeMatrix reference-point -p 5 --referencePoint TSS -b 3000 -a 3000 -S $project/visual/*.bw -R /data/zds209/database/mm10.refseq.bed -o $project/visual/ChIP-test.gz --outFileSortedRegions $project/visual/ChIP-test.bed --skipZeros
plotHeatmap -m $project/visual/ChIP-test.gz -out $project/visual/ChIP-test.png --colorMap RdBu --whatToShow 'plot, heatmap and colorbar' --kmeans 3
plotProfile -m ChIP-test.gz -out $ChIP-testProfile.png --plotType=fill --perGroup --numPlotsPerRow 2 --colors red yellow blue green --plotTitle "ChIP-test220612"
done
5.在Shell脚本中想要将Linux语句执行的结果赋值给某个变量,一定要将语句用反引号(Esc键下面的键)包裹起来,用单引号、双引号或者不用引号,该变量都得不到该语句执行的结果。
image.png
这是一种把某个文件夹下的文件列出来,并且放到一个变量中,可以用于后面的脚本编程。
再比如,可以把一个txt文本里面的内容放到一个变量中。
image.png