News:sorted the results of function dba.report in DiffBind R package
1
0
Entering edit mode
@bioqianyang-18743
Last seen 2.9 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 • 623 views
ADD COMMENT
2
Entering edit mode
Rory Stark ★ 4.1k
@rory-stark-5741
Last seen 21 days ago
CRUK, 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
0
Entering edit mode

Thank you very much ! 

ADD REPLY
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
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.

ADD REPLY
0
Entering edit mode

Oh, I see. Thank you so much !

ADD REPLY

Login before adding your answer.

Traffic: 460 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