Question

How to get all possible comparisons in DiffBind

1

Entering edit mode

bony.dekumar ▴ 10

@bonydekumar-19283

Last seen 22 months ago

United States

I was trying to analyze differential binding for four different set of experiments done in duplicate

SampleI Tissue  Factor  Condition   Replicate   bamReads    ControlID    bamControl Peaks    PeakCaller
NS_1    NS  Un  NS  1   15-011-001.sorted.bam           ns_s1_peaks.narrowPeak  narrowPeak
NS_2    NS  UN  NS  2   15-011-005.sorted.bam           ns_s2_peaks.narrowPeak  narrowPeak
HSF_1   HSF treated HSF 1   15-011-002.sorted.bam           hsf_s1_peaks.narrowPeak narrowPeak
HSF_2   HSF treated HSF 2   15-011-006.sorted.bam           hsf_s2_peaks.narrowPeak narrowPeak
LPS_1   LPS treated LPS 1   15-011-003.sorted.bam           lps_s1_peaks.narrowPeak narrowPeak
LPS_2   LPS treated LPS 2   15-011-007.sorted.bam           lps_s2_peaks.narrowPeak narrowPeak
LPS_HSF_1   HSF+LPS treated HSF+LPS 1   15-011-004.sorted.bam           hsf_lpf_s1_peaks.narrowPeak narrowPeak
LPS_HSF_1   HSF+LPS treated HSF+LPS 2   15-011-008.sorted.bam           hsf_lpf_s2_peaks.narrowPeak narrowPeak

hsf=dba(sampleSheet="samples_1.csv")
hsf = dba.count(hsf)
hsf<- dba.contrast(hsf, minMembers=2,categories=DBA_CONDITION)
hsf<- dba.analyze(hsf)
hsf.DB <- dba.report(hsf)
hsf.DB
GRanges object with 1297 ranges and 6 metadata columns:
        seqnames              ranges strand |      Conc   Conc_NS  Conc_HSF
           <Rle>           <IRanges>  <Rle> | <numeric> <numeric> <numeric>
   4525    chr10   19009922-19011726      * |      9.83      9.41     10.15
  34430     chr2 118112535-118114224      * |      9.29      9.64      8.83
  12272    chr12     3873249-3875190      * |      9.25      8.77      9.61
   4017     chr1 191266448-191269505      * |     10.02      9.64     10.32
  14287    chr12 108637001-108639115      * |     10.45     10.09     10.74
    ...      ...                 ...    ... .       ...       ...       ...
  18950    chr14   63125596-63125926      * |      8.57      8.33      8.77
  54141     chr7 139247666-139248397      * |      8.19      7.92      8.41
  32468     chr2   30441240-30441852      * |      7.87      7.57      8.12
  35762     chr2 156178342-156180177      * |      9.96      9.77     10.12
   9741    chr11   75463189-75463574      * |      7.46      7.13      7.72
             Fold   p-value       FDR
        <numeric> <numeric> <numeric>
   4525     -0.74  3.65e-11   3.3e-07
  34430      0.81  4.07e-11   3.3e-07
  12272     -0.84  5.27e-11   3.3e-07
   4017     -0.68  1.74e-10  8.15e-07
  14287     -0.66  4.78e-10   1.8e-06
    ...       ...       ...       ...
  18950     -0.44   0.00343    0.0499
  54141     -0.49   0.00344      0.05
  32468     -0.55   0.00345      0.05
  35762     -0.35   0.00345      0.05
   9741     -0.59   0.00345      0.05

It seems that i am getting only i comparison ( NS Vs HSF). How to get other comparisons also .

Thanks

diffbind • 4.0k views

ADD COMMENT • link updated 2.8 years ago by Rory Stark ★ 5.2k • written 5.9 years ago by bony.dekumar ▴ 10

score 3 · Answer 1 · 2019-01-09

3

Entering edit mode

Rory Stark ★ 5.2k

@rory-stark-5741

Last seen 5 weeks ago

Cambridge, UK

It looks like you have 4 conditions: NS, HSF, LPS, and HSF+LPS, is that right? And you'd like to see the 6 possible two-way contrasts?

If this is the case, you have actually computed them all already! To see them, you need to specify which contrasts you want in the call to dba.report() by setting the contrast parameter. In the code snippet you provided, you didn't specify this, so you got the default contrast (contrast=1).

You can see all the contrasts and their corresponding contrast numbers by printing out the DBA object:

> hsf

You can also retrieve the contrasts programmatically:

> contrasts <- dba.show(hsf, bContrasts=TRUE)

then cycle through the contrasts.

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

0

Entering edit mode

Thanks . It helped to create other contrast.

ADD REPLY • link 5.9 years ago bony.dekumar ▴ 10

0

Entering edit mode

Hi Rory, I am having similar situation. If we take the example above. Then, to get all possible comparisons, I have to write code chunks, as following?

#Lets assume "result=our DBA object having samples for multiple pairwise comparison"

dba.count(result) # for counting the reads for all conditions.

##for first pair
result1<- dba.contrast(result, for first pair)
result1<- dba.analyze(result1)
result1.db<- dba.report(result1)

##for second pair
result2<- dba.contrast(result, code for second pair)
result2<- dba.analyze(result2)
result2.db<- dba.report(result2)

#For 3rd, 4th..etc....

And so on...

Thanks again for this amazing package and all your help!

ADD REPLY • link 5.9 years ago bioinfouser2 ▴ 10

0

Entering edit mode

You can include all the contrasts in a single DBA object. You can call dba.contrast() multiple times, adding all the contrasts, then use a single call to dba.analyze() to run all the analyses (in parallel by default). Finally, to get the results, call dba.report() with contrast=n, changing the contrast number to get each specific contrast. Likewise for the plotting functions, set contrast=n to plot the results for each contrast.

ADD REPLY • link 5.9 years ago Rory Stark ★ 5.2k

0

Entering edit mode

Hi Rory,

Is DiffBind able to somehow run all these comparisons simultaneously? If I had, for instance, three conditions and wanted to obtain the set of differential peaks, would I be able to obtain these peaks without having to run all 3-way pairwise comparisons?

Thanks, Pedro

ADD REPLY • link 4.6 years ago Pedro Baldoni ▴ 40

score 0 · Answer 2 · 2020-06-09

0

Entering edit mode

Rory Stark ★ 5.2k

@rory-stark-5741

Last seen 5 weeks ago

Cambridge, UK

Currently, DiffBindonly performs pair-wise comparisons.

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

0

Entering edit mode

Hi, Dr Stark. If I want to some interactions work, Can I directly output diffbind object to DESeq2?

ADD REPLY • link 4.5 years ago Guandong Shang ▴ 40

0

Entering edit mode

The upcoming release of DiffBind will support the ability to include arbitrary designs and complex contrasts, made up of any of the metadata factors Tissue, Factor, Condition, Treatment, Replicate, and Caller. It also has a simple interface to retrieve the DESeq2 DESeqDataSet object.

If you're interested in trying if before the release in a few months, it should be available within 24 hours to download and install here: https://bioconductor.org/packages/3.12/bioc/html/DiffBind.htm

Any version after 2.99.1 has the ability to handle arbitrary design formulae. I am available to help as I am interested in having it tested by some users before release.

ADD REPLY • link 4.4 years ago Rory Stark ★ 5.2k

score 0 · Answer 3 · 2022-01-19

0

Entering edit mode

Rory Stark ★ 5.2k

@rory-stark-5741

Last seen 5 weeks ago

Cambridge, UK

Note that contrasts are handled differently since DiffBind_3.0. You can still calculate results for multiple contrasts simultaneously, but now you can include an arbitrary design formula, and more flexible contrast specifications, to allow any valid GLM design.

See the documentation for dba.contrast().

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

0

Entering edit mode

Good day,

Thanks for the helpful comments here. I have been able to get all my contrasts and see them. However, I can not out put anything other than the first pair. To this point all is good:

dba_report_filtered_df <- dbaObj_report_df %>% dplyr::filter (Fold >= 2 | Fold <= -2)

dba_report_final <- dba_report_filtered_df

dba_report_final <- dba_report_filtered_df

dba.show(dbObj, bContrasts=T) NULL

dba.show(dbObj_Contrasts, bContrasts=T) Factor Group Samples Group2 Samples2

1 Factor Day0 2 Day3_lysis 2

2 Factor Day0 2 Day3_nolysis 2

3 Factor Day0 2 Day10 2

4 Factor Day3_lysis 2 Day3_nolysis 2

5 Factor Day3_lysis 2 Day10 2

6 Factor Day10 2 Day3_nolysis 2

I am using the command: write_tsv (dba_report_filtered_df, "Results/Day0vsDay3lysis.tsv")

I have not been able to find information on how to change the output away from the default of contrast=1

Thanks in advance Chris

ADD REPLY • link 2.8 years ago c.eskiw • 0

score 0 · Answer 4 · 2022-02-02

0

Entering edit mode

Rory Stark ★ 5.2k

@rory-stark-5741

Last seen 5 weeks ago

Cambridge, UK

The dba.report() function, as well as the plotting functions, take a contrast= parameter (which defaults to contrast=1). So you just need to specify which contrast you want in each of those functions.

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