deepTools 使用指南

deepTools

deepTools 是一套基于python开发的工具,适用于有效处理分析高通量测序数据,可用于ChIP-seq, RNA-seq 或 MNase-seq。


#1. deepTools 系列工具

deepTools workflow

##1.1 deepTools 系列工具信息汇总

tool type input files main output file(s) application
multiBamSummary data integration 2 or more BAM interval-based table of values perform cross-sample analyses of read counts –> plotCorrelation, plotPCA
multiBigwigSummary data integration 2 or more bigWig interval-based table of values perform cross-sample analyses of genome-wide scores –> plotCorrelation, plotPCA
plotCorrelation visualization bam/multiBigwigSummary output clustered heatmap visualize the Pearson/Spearman correlation
plotPCA visualization bam/multiBigwigSummary output 2 PCA plots visualize the principal component analysis
plotFingerprint QC 2 BAM 1 diagnostic plot assess enrichment strength of a ChIP sample
computeGCBias QC 1 BAM 2 diagnostic plots calculate the exp. and obs. GC distribution of reads
correctGCBias QC 1 BAM, output from computeGCbias 1 GC-corrected BAM obtain a BAM file with reads distributed according to the genome’s GC content
bamCoverage normalization BAM bedGraph or bigWig obtain the normalized read coverage of a single BAM file
bamCompare normalization 2 BAM bedGraph or bigWig normalize 2 files to each other (e.g. log2ratio, difference)
computeMatrix data integration 1 or more bigWig, 1 or more BED zipped file for plotHeatmap or plotProfile compute the values needed for heatmaps and summary plots
estimateReadFiltering information 1 or more BAM files table of values estimate the number of reads filtered from a BAM file or files
alignmentSieve QC 1 BAM file 1 filtered BAM or BEDPE file filters a BAM file based on one or more criteria
plotHeatmap visualization computeMatrix output heatmap of read coverages visualize the read coverages for genomic regions
plotProfile visualization computeMatrix output summary plot (“meta-profile”) visualize the average read coverages over a group of genomic regions
plotCoverage visualization 1 or more BAM 2 diagnostic plots visualize the average read coverages over sampled genomic positions
bamPEFragmentSize information 1 BAM text with paired-end fragment length obtain the average fragment length from paired ends
plotEnrichment visualization 1 or more BAM and 1 or more BED/GTF A diagnostic plot plots the fraction of alignments overlapping the given features
computeMatrixOperations miscellaneous 1 or more BAM and 1 or more BED/GTF A diagnostic plot plots the fraction of alignments overlapping the given features

##1.2 BAM 和bigWig文件处理工具

​ 利用两个或多个bam文件计算基因组区段reads覆盖度;BED-file 指定基因组区域,bins模式可用于全基因组范围分析;产生的结果(.npz)可用于plotCorrelation进行相关性分析和用于plotPCA进行主成分分析。

$ deepTools2.0/bin/multiBamSummary bins \
 --bamfiles testFiles/*bam \ # using all BAM files in the folder
 --minMappingQuality 30 \
 --region 19 \ # limiting the binning of the genome to chromosome 19
 --labels H3K27me3 H3K4me1 H3K4me3 HeK9me3 input \
 -out readCounts.npz --outRawCounts readCounts.tab

 $ head readCounts.tab
 #'chr'     'start' 'end'   'H3K27me3'      'H3K4me1'       'H3K4me3'       'HeK9me3'       'input'
 19 10000   20000   0.0     0.0     0.0     0.0     0.0
 19 20000   30000   0.0     0.0     0.0     0.0     0.0
 19 30000   40000   0.0     0.0     0.0     0.0     0.0
 19 40000   50000   0.0     0.0     0.0     0.0     0.0
 19 50000   60000   0.0     0.0     0.0     0.0     0.0
 19 60000   70000   1.0     1.0     0.0     0.0     1.0
 19 70000   80000   0.0     1.0     7.0     0.0     1.0
 19 80000   90000   15.0    0.0     0.0     6.0     4.0
 19 90000   100000  73.0    7.0     4.0     16.0    5.0
  • multiBigwigSummary
    ​ 与multiBamSummary相比,输入文件格式是bigWig 。

  • correctGCBias
    ​ 矫正GC-bias;

  • bamCoverage

    norm_IGVsnapshot_indFiles

    ​ bamCoverage 利用测序数据比对结果转换为基因组区域reads覆盖度结果。可以自行设定覆盖度计算的窗口大小(bin);bamCoverage 内置了各种标准化方法:scaling factor, Reads Per Kilobase per Million mapped reads (RPKM), counts per million (CPM), bins per million mapped reads (BPM) and 1x depth (reads per genome coverage, RPGC).

Example : bamCoverage 用于ChIPseq分析

bamCoverage --bam a.bam -o a.SeqDepthNorm.bw \
    --binSize 10
    --normalizeUsing RPGC
    --effectiveGenomeSize 2150570000
    --ignoreForNormalization chrX
    --extendReads
    --outFileFormat bedgraph
  • bamCompare
    ​ 两个BAM 文件相比较,计算二者之间窗口中的reads丰度比率。
usage:  bamCompare -b1 treatment.bam -b2 control.bam -o log2ratio.bw
  • bigwigCompare
  • computeMatrix
    ​ 给基因组区段打分,产生的文件可用于plotHeatmapplotProfiles作图;基因组区段可以是基因或其他区域,使用BED格式文件定义即可。

computeMatrix 有两种不同的模式

computeMatrix two modes

  • reference-point(relative to a point): 计算某个点的信号丰度
  • scale-regions(over a set of regions): 把所有基因组区段缩放至同样大小,然后计算其信号丰度
    如下命令查看帮助:
$ computeMatrix scale-regions –help
$ computeMatrix scale-regions -S <biwig file(s)> -R <bed file> -b 1000
$ computeMatrix reference-point –help
$ computeMatrix reference-point -S <biwig file(s)> -R <bed file> -a 3000 -b 3000

Example 1:单个输入文件 (reference-point mode)

$ computeMatrix reference-point \ # choose the mode
       --referencePoint TSS \ # alternatives: TSS, TES, center
       -b 3000 -a 10000 \ # define the region you are interested in
       -R testFiles/genes.bed \
       -S testFiles/log2ratio_H3K4Me3_chr19.bw  \
       --skipZeros \
       -o matrix1_H3K4me3_l2r_TSS.gz \ # to be used with plotHeatmap and plotProfile
       --outFileSortedRegions regions1_H3K4me3_l2r_genes.bed

​ 注:point-BED文件指定基因组区段的起始位置

Example 2:多个输入文件 (scale-regions mode)

$ deepTools2.0/bin/computeMatrix scale-regions \
  -R genes_chr19_firstHalf.bed genes_chr19_secondHalf.bed \ # separate multiple files with spaces
  -S testFiles/log2ratio_*.bw  \ or use the wild card approach
  -b 3000 -a 3000 \
  --regionBodyLength 5000 \
  --skipZeros -o matrix2_multipleBW_l2r_twoGroups_scaled.gz \
  --outFileNameMatrix matrix2_multipleBW_l2r_twoGroups_scaled.tab \
  --outFileSortedRegions regions2_multipleBW_l2r_twoGroups_genes.bed
Note that the reported regions will have the same coordinates as the ones in 

##1.3 质控工具

  • plotCorrelation
    ​ 基于multiBamSummary 或multiBigwigSummary结果计算样品间的相关性。并且还可以通过ScatterplotHeatmap进行展示。

Example 1:Scatterplot

$ deepTools2.0/bin/plotCorrelation \
-in scores_per_transcript.npz \
--corMethod pearson --skipZeros \
--plotTitle "Pearson Correlation of Average Scores Per Transcript" \
--whatToPlot scatterplot \
-o scatterplot_PearsonCorr_bigwigScores.png   \
--outFileCorMatrix PearsonCorr_bigwigScores.tab
scatterplot_PearsonCorr_bigwigScores

Example 2:Heatmap

$ deepTools2.0/bin/plotCorrelation \
    -in readCounts.npz \
    --corMethod spearman --skipZeros \
    --plotTitle "Spearman Correlation of Read Counts" \
    --whatToPlot heatmap --colorMap RdYlBu --plotNumbers \
    -o heatmap_SpearmanCorr_readCounts.png   \
    --outFileCorMatrix SpearmanCorr_readCounts.tab
heatmap_SpearmanCorr_readCounts1
  • plotPCA
    ​ 基于multiBamSummary 或multiBigwigSummary结果进行主成分分析,并作出基于两个主成分的图和前五个特征代表性的图。

Example

$ deepTools2.0/bin/plotPCA -in readCounts.npz \
-o PCA_readCounts.png \
-T "PCA of read counts"
PCA_readCounts
  • plotFingerprint
    ​ 对样本比对结果reads累积情况进行展示。一定长度窗口(bin)上reads数进行计数,然后排序,再依次累加画图。input (能测到90%DNA片段)在基因组理论上是均匀分布,随着测序深度增加趋近于直线,实验组在排序越高的窗口处reads累积速度越快,说明这些区域富集的越特异。
    QC_fingerprint

Example

$ deepTools2.0/bin/plotFingerprint \
 -b testFiles/*bam \
--labels H3K27me3 H3K4me1 H3K4me3 H3K9me3 input \
--minMappingQuality 30 --skipZeros \
--region 19 --numberOfSamples 50000 \
-T "Fingerprints of different samples"  \
--plotFile fingerprints.png \
--outRawCounts fingerprints.tab
fingerprints1
  • bam PEFragmentSize
    ​ 计算bam文件中双端reads的fragment size长度。
  • compute GCBias
    ​ 计算GC-bias
  • plot Coverage
    ​ 计算样品测序深度。随机抽取1 million bp ,计算reads数,统计碱基覆盖率和覆盖次数。

##1.4 热图和总结图

Example 1: 根据computeMatrix结果画热图

# run compute matrix to collect the data needed for plotting
$ computeMatrix scale-regions -S H3K27Me3-input.bigWig \
                                 H3K4Me1-Input.bigWig  \
                                 H3K4Me3-Input.bigWig \
                              -R genes19.bed genesX.bed \
                              --beforeRegionStartLength 3000 \
                              --regionBodyLength 5000 \
                              --afterRegionStartLength 3000
                              --skipZeros -o matrix.mat.gz
$ plotHeatmap -m matrix.mat.gz \
      -out ExampleHeatmap1.png \
plot Heatmap

Example 2: plotHeatmap还可以进行聚类分析

$ plotHeatmap -m matrix_two_groups.gz \
     -out ExampleHeatmap2.png \
     --colorMap RdBu \
     --whatToShow 'heatmap and colorbar' \
     --zMin -3 --zMax 3 \
     --kmeans 4  #聚类参数

[图片上传失败...(image-d28547-1556115525034)]
其他参数

颜色自定义:--colorList 'black, yellow' 'white,blue' '#ffffff,orange,#000000'
去掉热图边框:--boxAroundHeatmaps no

Example 1: 根据样本画图

# run compute matrix to collect the data needed for plotting
$ computeMatrix scale-regions -S H3K27Me3-input.bigWig \
                                 H3K4Me1-Input.bigWig  \
                                 H3K4Me3-Input.bigWig \
                              -R genes19.bed genesX.bed \
                              --beforeRegionStartLength 3000 \
                              --regionBodyLength 5000 \
                              --afterRegionStartLength 3000
                              --skipZeros -o matrix.mat.gz

$ plotProfile -m matrix.mat.gz \
              -out ExampleProfile1.png \
              --numPlotsPerRow 2 \
              --plotTitle "Test data profile"
Example Profile1

Example 2: 根据基因画图

$ plotProfile -m matrix.mat.gz \
     -out ExampleProfile2.png \
     --plotType=fill \ # add color between the x axis and the lines
     --perGroup \ # make one image per BED file instead of per bigWig file
     --colors red yellow blue \
     --plotTitle "Test data profile"
Example Profile2

Example 3: 聚类画图

$ plotProfile -m matrix.mat.gz \
      --perGroup \
      --kmeans 2 \
      -out ExampleProfile3.png
Example Profile3

Example 4: 画热图

$ plotProfile -m matrix.mat.gz \
      --perGroup \
      --kmeans 2 \
      -plotType heatmap \
      -out ExampleProfile3.png

Example Profile4

plotEnrichment
​ 统计样本BED文件中peak或者GTF文件中feature 在chipseq结果中富集情况

Example

$ plotEnrichment -b Input.bam H3K4Me1.bam H3K4Me3.bam \
--BED up.bed down.bed \
--regionLabels "up regulated" "down regulated" \
-o enrichment.png
plot Enrichment

参考:

The tools of deeptools

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 207,113评论 6 481
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 88,644评论 2 381
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 153,340评论 0 344
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 55,449评论 1 279
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 64,445评论 5 374
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 49,166评论 1 284
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 38,442评论 3 401
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 37,105评论 0 261
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 43,601评论 1 300
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,066评论 2 325
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,161评论 1 334
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,792评论 4 323
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 39,351评论 3 307
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,352评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,584评论 1 261
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 45,618评论 2 355
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,916评论 2 344