I used diffbind to analysis the data with 5 time point. I wanted to find out the unique peaks in different time point. so I run the code:
dbObj.OL <- dba.overlap(dbObj,dbObj$masks$time)
but it return:
Warning message:
In pv.overlap(DBA, mask = mask) : Too many peaksets in mask.
when I used four of the five time point, it worked very well.
what should I do?
As far a returning all the combinations of overlapping peaksets goes, the default mode of dba.overlap() is limited to overlapping four sets (it underlies dba.plotVenn(), which is likewise limited).
You may want to look at the results from using dba.overlap() with mode=DBA_OLAP_ALL. This will show you the results of all pair-wise comparisons. There's an example on the help page.
The $called matrix does indeed have a row for every consensus peak, with a 1 if the sample had a called peak that overlaps that consensus peak.
I ran into a similar problem with getting a consensus matrix by having 5-7 replicates per group. I used the called matrix to pull out the union and unique peaks for 2 groups of interest. Is this the right approach and the correct use of the called matrix?
CON.LPS.called <- dat$called[,c(3,4,8,10,29,32,22)] # subset group
MIA.LPS.called <- dat$called[,c(6,7,9,26,30)] # subset group
CON.LPS.sum <- (rowSums(CON.LPS.called))/7 # percent overlap
MIA.LPS.sum <- (rowSums(MIA.LPS.called))/5 # percent overlap
CON.LPS.consensus <- CON.LPS.sum > 0.5 # keep 50% overlap
MIA.LPS.consensus <- MIA.LPS.sum > 0.5 # keep 50% overlap
union <- CON.LPS.consensus==TRUE & MIA.LPS.consensus==TRUE
MIA.unique.peaks <- MIA.LPS.consensus==TRUE & CON.LPS.consensus==FALSE
CON.unique.peaks <- CON.LPS.consensus==TRUE & MIA.LPS.consensus==FALSE
neither <- CON.LPS.consensus==FALSE & MIA.LPS.consensus==FALSE
I think what you are doing look like it should work
As you are only working with to groups here, not five, you should be able to do what you want using the documented interface instead of using the internal $called matrix.
If, for example, you had the Condition factor set to CON and MIA (where CON was set for samples 3,4,8,10,29,32, and 22, and MIA was set for samples 6,7,9,26, and 30):
dat.consensus <- dba.peakset(dat, consensus=DBA_CONDITION, minOverlap=.5)
dat.consensus <- dba(dat.consensus, mask=dat.consensus$Consensus, minOverlap=1)
dat.overlaps <- dba.overlap(dat.consensus, mask=1:2)
CON.unique.peaks <- dat.overlaps$onlyA
MIA.unique.peaks <- dat.overlaps$onlyB
union <- dat.overlaps$inAll
getting neither as well is a little tricker, but you could get the complete set of peaks and then identify the ones that don't overlap
all.peaks <- dba.peakset(dba(dat.min, Overlap=1), bRetrieve=TRUE)
dat.consensus <- dba.peakset(dat.consensus, all.peaks)
neither <- dba.overlap(dat.consensus, mask=1:3)$onlyC
You can call dba,plotVenn() instead of dba.overlap() to visualize the results as you compute them.
A problem is that the full dat object has 6 groups, and for this analysis I'm looking at the consensus of only 2 of those groups. I changed DBA_TISSUE to CON, MIA or NA for the groups I want to ignore. I already have a consensus of all the samples for the previous DAR analysis across multiple contrasts. I deleted the consensus mask but that didn't work.
It does not like NA in the grouping variable. I remade the dba object with just the 2 groups samples and it all worked fine. Thanks!
Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
version 2.3.6
I found some information in
dbObj2$called. Are they the peak information which used to merge for consensus peak set?