Diffbind: Does Heatmap use count or score? and for which intervals ?
1
0
Entering edit mode
ZheFrench ▴ 50
@zhefrench-11689
Last seen 16 months ago
France

Hey,

I saved a first chipQC object using consensus=TRUE,bCounts=TRUE,summits=250 options.

Then I reused it latter inside DiffBind.

First I was wondering if it uses score or count to make heatmap ? Documentation it misleading, it is said it uses scores but if you have add also your Input, necessarely it uses count from all region defined in input control  because we don't have done peaks calling for input ... But when you don't retrieve a chipQC object which used consensus do you use score or count ? Not clear for me. 

Second can you confirm , intervals that are used to make heatmap in the following workflow ?

experiment =

ChIPQC(samples,annotation="hg38",consensus=TRUE,bCounts=TRUE,summits=250,chromosomes=c("chr1","chr2","chr3","chr4","chr5","chr6","chr7","chr8","chr9","chr10","chr11","chr12","chr13","chr14","chr15","chr16","chr17","chr18","chr19","chr20","chr21","chr22","chrX","chrY"),blacklist="/home/jean-philippe.villemin/mount_archive2/commun.luco/ref/genes/GRCh38_PRIM_GENCODE_R25/hg38.blacklist.bed.gz") 

dbasavefile <- dba.save(experiment, file=name_file_to_save, dir='.', pre='dba_', ext='RData', bMinimize=FALSE) 

#' Create Report.
ChIPQCreport(experiment,facet=T,lineBy=c("Replicate"),facetBy=c("Condition","Factor"),reportName=paste0(opt$name,"Marks",collapse = "_"), reportFolder=opt$name,colourBy=c("Replicate"))

Then when I load my dba.RDATA object and try to plot heatmap, I get also control from Input. So necessarely it uses count from all region defined in control  because we don't have done peaks calling for input ...

dba.show(dba_chip)

          ID Tissue  Factor Condition Treatment Replicate Caller Intervals

1      L13S13  MCF10   K27AC        T1 Tamoxifen         1   macs     94888
2     L142S14  MCF10   K27AC        T1 Tamoxifen         3   macs    100680
3      L17S17  MCF10   K27AC        T7 Tamoxifen         1   macs    106629
4     L182S18  MCF10   K27AC        T7 Tamoxifen         3   macs     87323
5      L15S15  MCF10   K27AC      unT7 unTreated         1   macs    106930
6     L162S16  MCF10   K27AC      unT7 unTreated         3   macs    109774
7        L7S7  MCF10  K27ME3        T1 Tamoxifen         1   macs    159118
8        L8S8  MCF10  K27ME3        T1 Tamoxifen         3   macs    140369
9      L11S11  MCF10  K27ME3        T7 Tamoxifen         1   macs    150973
10     L12S12  MCF10  K27ME3        T7 Tamoxifen         3   macs    131621
11       L9S9  MCF10  K27ME3      unT7 unTreated         1   macs    143701
12     L10S10  MCF10  K27ME3      unT7 unTreated         3   macs    157131
13       L1S1  MCF10   K4ME1        T1 Tamoxifen         1   macs    173416
14       L2S2  MCF10   K4ME1        T1 Tamoxifen         2   macs    169116
15       L5S5  MCF10   K4ME1        T7 Tamoxifen         1   macs    164891
16       L6S6  MCF10   K4ME1        T7 Tamoxifen         2   macs    179542
17       L3S3  MCF10   K4ME1      unT7 unTreated         1   macs    189781
18       L4S4  MCF10   K4ME1      unT7 unTreated         2   macs    170599
19   T1_INPUT  MCF10 Control        T1 Tamoxifen        c1           269818
20   T7_INPUT  MCF10 Control        T7 Tamoxifen        c2           269818
21 unT7_INPUT  MCF10 Control      unT7 unTreated        c3           269818


dba.plotHeatmap(dba_chip_init)

So here it plots every sample for the 269818 intervals based on counts  ?

dba_chip <- dba(dba_chip_init, mask=dba_chip_init@DBA$masks[[opt$name]])
dba_chip$config$fragmentSize <- dba_chip$config$fragmentSize[dba_chip_init@DBA$masks[[opt$name]]]
dba_chip$config$fragmentSize

print("Mask")
dba.show(dba_chip)

[1] "Mask"
    ID Tissue Factor Condition Treatment Replicate Caller Intervals
1 L1S1  MCF10  K4ME1        T1 Tamoxifen         1   macs    173416
2 L2S2  MCF10  K4ME1        T1 Tamoxifen         2   macs    169116
3 L5S5  MCF10  K4ME1        T7 Tamoxifen         1   macs    164891
4 L6S6  MCF10  K4ME1        T7 Tamoxifen         2   macs    179542
5 L3S3  MCF10  K4ME1      unT7 unTreated         1   macs    189781
6 L4S4  MCF10  K4ME1      unT7 unTreated         2   macs    170599

 

So it plots heatmap for 269818 intervals based on counts but with a mask on K4ME1 ?

dba.plotHeatmap(dba_chip)

dba_chip_consensus <- dba.peakset(dba_chip, consensus=c(type),minOverlap=0.99)
dba.show(dba_chip_consensus)

    ID Tissue Factor Condition Treatment Replicate Caller Intervals
1 L1S1  MCF10  K4ME1        T1 Tamoxifen         1   macs    173416
2 L2S2  MCF10  K4ME1        T1 Tamoxifen         2   macs    169116
3 L5S5  MCF10  K4ME1        T7 Tamoxifen         1   macs    164891
4 L6S6  MCF10  K4ME1        T7 Tamoxifen         2   macs    179542
5 L3S3  MCF10  K4ME1      unT7 unTreated         1   macs    189781
6 L4S4  MCF10  K4ME1      unT7 unTreated         2   macs    170599
7   T1  MCF10  K4ME1        T1 Tamoxifen       1-2   macs    139069
8   T7  MCF10  K4ME1        T7 Tamoxifen       1-2   macs    138384
9 unT7  MCF10  K4ME1      unT7 unTreated       1-2   macs    118478

consensus_peaks <- dba.peakset(dba_chip_consensus, bRetrieve=TRUE)

dba_chip <- dba.count(dba_chip, peaks=consensus_peaks)

dba.show(dba_chip)

    ID Tissue Factor Condition Treatment Replicate Caller Intervals FRiP
1 L1S1  MCF10  K4ME1        T1 Tamoxifen         1 counts    182104 0.57
2 L2S2  MCF10  K4ME1        T1 Tamoxifen         2 counts    182104 0.57
3 L5S5  MCF10  K4ME1        T7 Tamoxifen         1 counts    182104 0.58
4 L6S6  MCF10  K4ME1        T7 Tamoxifen         2 counts    182104 0.57
5 L3S3  MCF10  K4ME1      unT7 unTreated         1 counts    182104 0.57
6 L4S4  MCF10  K4ME1      unT7 unTreated         2 counts    182104 0.55

dba.plotHeatmap(dba_chip)

Does it plot heatmap for 182104 intervals ?

 

diffbind • 1.0k views
ADD COMMENT
0
Entering edit mode
Rory Stark ★ 5.2k
@rory-stark-5741
Last seen 14 days ago
Cambridge, UK

All heatmaps are generated using a consensus peak set (same set of merged peaks for all samples in the plot). If you have done the counting step, as in these examples, it plots the correlation values using the scores for each sample. 

If, instead of calling dba.show(), you just type the name of the object, it will tell you how many sites are in the consensus binding matrix.

-Rory

ADD COMMENT

Login before adding your answer.

Traffic: 670 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6