Question

News:sorted the results of function dba.report in DiffBind R package

0

Entering edit mode

bioqianyang • 0

@bioqianyang-18743

Last seen 5.4 years ago

Hello,
 
I am using your R package these days, and I got the results of .dba.report () function. I want to how to sort the results ? For example, I ran the code:
 
bptH3K27AcPeak.DB<-dba.report(bptH3K27AcPeak_a, th=0.05, bUsePval=TRUE, fold=2)
 
And I got the bptH3K27AcPeak.DB
 
Then I want to use sorted bptH3K27AcPeak.DB files, for example, sorted by P-value. And do the dba.plotHeatmap function using top 1000 peaks.
 
Thank you very much !

diffbind News • 1.4k views

ADD COMMENT • link updated 5.4 years ago by Rory Stark ★ 5.2k • written 5.4 years ago by bioqianyang • 0

score 2 · Answer 1 · 2018-12-08

2

Entering edit mode

Rory Stark ★ 5.2k

@rory-stark-5741

Last seen 4 days ago

Cambridge, UK

The returned report should already be sorted by p-value (lowest to highest).

You can change how the report is sorted based on any of the metadata columns. For example, to sort by absolute value of the fold change (greatest absolute fold change to lowest):

> bptH3K27AcPeak.DB[order(abs(bptH3K27AcPeak.DB$Fold),decreasing=FALSE),]

To control what sites are used in dba.plotHeatmap(), you can specify report=bptH3K27AcPeak.DB. You can set the value for maxSites to the number of sites in the report (and make sure the report has the sites you want to plot). For example, to plot the 100 sites in the report with the highest absolute fold changes:

> highestAbsFold <- bptH3K27AcPeak.DB[order(abs(bptH3K27AcPeak.DB$Fold), decreasing=TRUE),][1:100,]
> dba.plotHeatmap(bptH3K27AcPeak_a, correlations=FALSE,
               report=highestAbsFold, maxSites=100)

ADD COMMENT • link 5.4 years ago Rory Stark ★ 5.2k

0

Entering edit mode

Thank you very much !

ADD REPLY • link 5.4 years ago bioqianyang • 0

0

Entering edit mode

Hello, Rory Stark

I run the dba.plotHeatmap function described below, and I got the same heatmap with different condition. I am not sure if the parameter report does not work when the correlation=TRUE ? Or other reasons ?

Thank you !

> highestAbsFold1 <- bptH3K27AcPeak.DB[order(abs(bptH3K27AcPeak.DB$Fold), decreasing=TRUE),][1:10]
> dba.plotHeatmap(bptH3K27AcPeak_a, correlations=TRUE,report=highestAbsFold1,maxSites=100)
> dba.plotHeatmap(bptH3K27AcPeak_a, correlations=TRUE,report=highestAbsFold1,maxSites=10)
> highestAbsFold1 <- bptH3K27AcPeak.DB[order(abs(bptH3K27AcPeak.DB$Fold), decreasing=TRUE),][1:5]
> highestAbsFold1
GRanges object with 5 ranges and 6 metadata columns:
seqnames ranges strand | Conc Conc_normal Conc_tumor Fold p-value FDR
<Rle> <IRanges> <Rle> | <numeric> <numeric> <numeric> <numeric> <numeric> <numeric>
49335 chrX 115004581-115005081 * | 7.54 -0.41 7.95 -8.36 4.16e-05 0.0595
19931 chr17 46805721-46806221 * | 5.15 -0.41 5.56 -5.97 3.57e-05 0.0595
45457 chr8 95220198-95220698 * | 4.96 -0.41 5.37 -5.78 0.000649 0.146
13038 chr13 91519165-91519665 * | 4.84 -0.41 5.24 -5.65 1.06e-05 0.0393
42705 chr7 48183233-48183733 * | 5.3 0.06 5.7 -5.64 2.02e-05 0.0545
-------
seqinfo: 35 sequences from an unspecified genome; no seqlengths
> dba.plotHeatmap(bptH3K27AcPeak_a, correlations=TRUE,report=highestAbsFold1,maxSites=5)
> dba.plotHeatmap(bptH3K27AcPeak_a, correlations=TRUE,report=highestAbsFold1)
> dba.plotHeatmap(bptH3K27AcPeak_a, correlations=TRUE)

ADD REPLY • link 5.4 years ago bioqianyang • 0

1

Entering edit mode

You still need to specify a contrast number (the contrast the report was generated from). In your original call to dba.report() you used the default contrast=1, so you should add contrast=1 to the call to dba.plotHeatmap().

When no contrast is specified, dba.plotHeatmap() always uses the global binding matrix.