2
$\begingroup$

This question was also asked on Biostars and the Bioconductor Forum

I want to find differential accessibility between conditions (WS, WF, MS, MF). I can compare any pair conditions without error except this pair MS vs MF. Data is the bam and broadpeak files from ATAC seq pipeline and the code is like a pipeline of the DiffBind R package. The sample is human and uses Illumina sequencing technology.

Would you please suggest what wrong in this case? I pretty much run the same code as I used for another data set and didn't have this error. Thank you so much!

result <- dba(sampleSheet=samples)
result <- dba.blacklist(result) 
result <- dba.count(result, bParallel=FALSE)
result <- dba.count(result)
result <- dba.normalize(result)
result <- dba.contrast(result, minMembers=2)
result <- dba.analyze(result)
result <- dba.plotProfile(result)

Generating report-based DBA object... Error: No valid contrasts/methods specified.

There is post about this error but I still got error when trying those solutions. Such as the order of the function, the path to the bam and broad peak files in sample sheet.

dba.show(result, bContrasts=TRUE)
     Factor    Group Samples  Group2 Samples2 DB.DESeq2
1 Condition diseased       2 control        2         0

I don't know why the value of DB.DESeq2 in this case was 0 which maybe the reason of the error.

This is the result after running dba.contrast:

4 Samples, 89585 sites in matrix:
   ID Tissue   Factor Condition Replicate    Reads FRiP
1  MS    non ATAC-seq   control         1 23731398 0.17
2 MSR    non ATAC-seq   control         2 10499000 0.16
3  MF    non ATAC-seq  diseased         1 31535510 0.18
4 MFR    non ATAC-seq  diseased         2 18488559 0.18

Design: [~Condition] | 1 Contrast:
     Factor    Group Samples  Group2 Samples2
1 Condition diseased       2 control        2

If I try to compare MS vs MF, I will get that error but if I remove MS and replace by other sample, I don't have that error anymore.

I don't have that error if I add more samples and the value of DB.SESeq2 at the first row from dba.show() is not 0.

Full samplesheet:

SampleID,Tissue,Factor,Condition,Treatment,Replicate,bamReads,ControlID,bamControl,Peaks,PeakCallers
ws1,non,ATAC-seq,ws,,1,ATAC_seq_output/bwa/merged_library/ws_REP1.mLb.clN.sorted.bam,,,ATAC_seq_output/bwa/merged_library/macs2/broad_peak/ws_REP1.mLb.clN_peaks.broadPeak,bed
ws2,non,ATAC-seq,ws,,2,ATAC_seq_output/bwa/merged_library/ws_REP2.mLb.clN.sorted.bam,,,ATAC_seq_output/bwa/merged_library/macs2/broad_peak/ws_REP2.mLb.clN_peaks.broadPeak,bed
wf1,non,ATAC-seq,wf,,1,ATAC_seq_output/bwa/merged_library/wf_REP1.mLb.clN.sorted.bam,,,ATAC_seq_output/bwa/merged_library/macs2/broad_peak/wf_REP1.mLb.clN_peaks.broadPeak,bed
wf2,non,ATAC-seq,wf,,2,ATAC_seq_output/bwa/merged_library/wf_REP2.mLb.clN.sorted.bam,,,ATAC_seq_output/bwa/merged_library/macs2/broad_peak/wf_REP2.mLb.clN_peaks.broadPeak,bed
ms1,non,ATAC-seq,ms,,1,ATAC_seq_output/bwa/merged_library/ms_REP1.mLb.clN.sorted.bam,,,ATAC_seq_output/bwa/merged_library/macs2/broad_peak/ms_REP1.mLb.clN_peaks.broadPeak,bed
ms2,non,ATAC-seq,ms,,2,ATAC_seq_output/bwa/merged_library/ms_REP2.mLb.clN.sorted.bam,,,ATAC_seq_output/bwa/merged_library/macs2/broad_peak/ms_REP2.mLb.clN_peaks.broadPeak,bed
mf1,non,ATAC-seq,mf,,1,ATAC_seq_output/bwa/merged_library/ws_REP1.mLb.clN.sorted.bam,,,ATAC_seq_output/bwa/merged_library/macs2/broad_peak/ws_REP1.mLb.clN_peaks.broadPeak,bed
mf2,non,ATAC-seq,mf,,2,ATAC_seq_output/bwa/merged_library/ws_REP2.mLb.clN.sorted.bam,,,ATAC_seq_output/bwa/merged_library/macs2/broad_peak/ws_REP2.mLb.clN_peaks.broadPeak,bed

enter image description here

$\endgroup$
5
  • $\begingroup$ Hi @gringer, I added the sample sheet in the question. Any other pair comparisons have non zero DB.SESeq2 value. $\endgroup$
    – Chris
    Jul 19, 2023 at 20:43
  • 1
    $\begingroup$ The conditions in that sample sheet (ws / wf / ms / mf) don't match the dba.contrast result (control / diseased). Can you show the result of dba.contrast after doing another pair comparison? $\endgroup$
    – gringer
    Jul 20, 2023 at 0:17
  • $\begingroup$ It is the full sample sheet. Because I don't want to compare ws and wf so I delete them. Maybe there is a way to manipulate the output but I still haven't known so I delete ws and wf. $\endgroup$
    – Chris
    Jul 20, 2023 at 6:20
  • 1
    $\begingroup$ How did "diseased" and "control" get into the first result after dba.contrast? They don't appear on the sample sheet at all. Why does the ID column look different from the sample sheet? I'm being pedantic about this because it seems important for the specific error that you're having. $\endgroup$
    – gringer
    Jul 20, 2023 at 10:16
  • $\begingroup$ I replaced ms in condition column by control and mf by diseased. Thank you. Hope it is clearer now. $\endgroup$
    – Chris
    Jul 20, 2023 at 22:14

1 Answer 1

2
$\begingroup$

I have very little idea what's going on here, but in the absence of any other answers or help, I'll try to work through this from a no-knowledge start point to see if I can hopefully help shift this problem closer towards a solution.

The specific error you have encountered is No valid contrasts/methods specified. Looking at the DiffBind Vignette, it provides this information about defining contrasts:

Before running the differential analysis, we need to tell DiffBind how to model the data, including which comparison(s) we are interested in. This is done using the dba.contrast function, as follows:

> tamoxifen <- dba.contrast(tamoxifen, reorderMeta=list(Condition="Responsive"))  
> tamoxifen
11 Samples, 2845 sites in matrix:
       ID Tissue Factor  Condition  Treatment Replicate   Reads FRiP
 1 BT4741  BT474     ER  Resistant Full-Media         1  652697 0.16
 2 BT4742  BT474     ER  Resistant Full-Media         2  663370 0.15
 3  MCF71   MCF7     ER Responsive Full-Media         1  346429 0.31
 4  MCF72   MCF7     ER Responsive Full-Media         2  368052 0.19
 5  MCF73   MCF7     ER Responsive Full-Media         3  466273 0.25
 6  T47D1   T47D     ER Responsive Full-Media         1  399879 0.11
 7  T47D2   T47D     ER Responsive Full-Media         2 1475415 0.06
 8 MCF7r1   MCF7     ER  Resistant Full-Media         1  616630 0.22
 9 MCF7r2   MCF7     ER  Resistant Full-Media         2  593224 0.14
10  ZR751   ZR75     ER Responsive Full-Media         1  706836 0.33
11  ZR752   ZR75     ER Responsive Full-Media         2 2575408 0.22

After showing your contrast result, you explained a partial removal of the error:

If I try to compare MS vs MF, I will get that error but if I remove MS and replace by other sample, I don't have that error anymore.

Specifically relating to MS vs MF, I can see why there would be an issue with the condition contrast, especially if it's using DESeq2 in the background. The short answer is that there is not enough variation within the groups.

DESeq2 cares about making sure that conditions added to the model are separable from samples. This is further explained in the DESeq2 manual here. The issue is that there's no way to distinguish between the control condition and the MS sample (or, similarly, between the diseased condition and the MF sample). In order to get past this error, you need both control and diseased in your MS group or in your MF group (or both, if you want to be really rigorous), and ideally a few more replicates of each condition.

$\endgroup$
0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.