DiffBind analysis - minCount parameter in dba.count()
1
0
Entering edit mode
8 weeks ago
vartika • 0

Why does changing minCount from 0 to 1 in dba.count() (with summits set to false and control reads supplied) change the called status of a peak? With minCount set to 1, a peak is being called in all samples, but with minCount set to 0, the same peak is being called in a few of them?

So, if minCount is set to 1, it will show the peak is called anyway? But that is not the case, cause there are peaks that are not called in half of the samples but called in the other half with minCount set to 1?

Am I missing something basic? What might be happening in my case?

ChIP-Seq DiffBind • 166 views
ADD COMMENT
1
Entering edit mode

If you can supply a bit more info, I can see if I can reproduce this.

How are you determining how many samples a peak is called in? Can you show a bit of code demonstrating:

  • your first call to dba.count()
  • show a peak and which samples it is called in
  • your second call to dba.count() changing minCount
  • show how the peak is appearing to have been called in more/fewer samples

    Also, can you confirm which version of DiffBind you are using?

ADD REPLY
0
Entering edit mode

Hi Rory! Thank you for your response. I was referring to the columns "Called1" and "Called2" in the report, they change when I change the minCount.

Diffbind version 3.0.13

Line from the Code:

minCount = 1

exp <- dba.count(exp, bSubControl = TRUE, minCount = 1, summits = FALSE)

minCount = 0

exp <- dba.count(exp, bSubControl = TRUE, minCount = 0, summits = FALSE)

I made no other changes for the two, the following were set:

exp$config$design <- FALSE

exp$config$doBlacklist <- FALSE

exp$config$doGreylist <- FALSE

After count:

exp <- dba.normalize(exp,normalize = DBA_NORM_LIB, library =DBA_LIBSIZE_FULL )

exp <- dba.contrast(exp,categories = DBA_CONDITION)

exp <-dba.analyze(exp)

exp.DB <- dba.report(exp, contrast =1, bCalled = TRUE, bCalledDetail = TRUE, bCounts = TRUE, bNormalized = TRUE)

Please find the link for a google sheet below which has three examples from the report.

Three examples from the report generated

ADD REPLY
0
Entering edit mode

I've tried this on a couple of my own data sets but can not reproduce it (the Called values are always the same when I compare sum(rep0$Called1 != rep1$Called1) and sum(rep0$Called2 != rep1$Called2)).

I'd like to get to the bottom of this. If you could send me a link to your two exp objects immediately after calling dba.count() I can try to see what is going on.

ADD REPLY
0
Entering edit mode

Okay, I'm not sure where I have gone wrong. Here you go, exp objects after dba.count

ADD REPLY
2
Entering edit mode
8 weeks ago
Rory Stark ★ 1.0k

I am able to reproduce this and have identified, and fixed, a bug in how the Called status is maintained after filtering (when minCounts is lowered to 0, some counts are lessened, and additional sites are filtered based on the requirement that at least one sample have an RPKM value of 1 or more).

I have checked in a fix that will appear as DiffBind_3.0.14, available via BiocManager::install() in the next couple of days.

ADD COMMENT
0
Entering edit mode

oh! Thank you for answering to my question and fixing the bug so quickly! Also, being a newbie to bioinformatics, it is very inspiring to see the time and efforts put into developing such packages for all of us to use. Great work! Thanks, again.

ADD REPLY
0
Entering edit mode

This was actually a non-trivial bug, so I appreciate you giving me the information I needed to reproduce and resolve it!

ADD REPLY

Login before adding your answer.

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