[BioC] ChIPpeakAnno: makeVennDiagram and sampling peaks

Zhu, Lihua (Julie) Julie.Zhu at umassmed.edu
Thu Mar 15 15:41:31 CET 2012


If you assume 10% total histone are available for modification, then you
would set totalTest = 3 * 10^9 bp / 146 * 0.1  which is about 2 million.

For your new question, here is some code snippet that might address your

    t1 =findOverlappingPeaks(peaks1, peaks2, maxgap=0,
NameOfPeaks1="TF", NameOfPeaks2="Histone", select="First")

AllPeaks = c(peaks1[!rownames(peaks1) %in% rownames(t1$Peaks1withOverlap),],
peaks2[!rownames(peaks2) %in% rownames(t1$Peaks2withOverlap),],

Totals = rownames(AllPeaks)
Sample.n = dim(t1$MergedPeaks)[1]

##### put the following code snippets in a loop allow you to sample from the
total peak population multiple times
s1 = AllPeaks[rownames(AllPeaks) %in% sample(Totals, Sample.n),]
go.s1 = getEnrichedGO(annotatePeakInBatch(s1,....), ....)

Please let me know if you come up with more elegant ways to do this. Thanks!

Best regards,


On 3/14/12 5:35 PM, "Ron Hart" <rhart at rci.rutgers.edu> wrote:

> Julie,
> In response to my last question and your phone call, I tried several values of
> totalTest based on the recommendations but I could only get either a 0 or a 1
> value.  For histone marks, I used as the largest estimate the total number of
> possible histone overlaps (3 x 10^9 bp / 146 bp per nucleosome).  Then I tried
> the sum of the two sets of marks,  but nothing made sense for me.  So I gave
> up trying to get a p-value.  It¹s really not important for my study.
> New question.  I¹m using the overlap function to extract intersecting peaks in
> common between two marks.  Everything is working great.  But I¹d like to
> compare the result to a random sampling of the same number of peaks from the
> union set of both marks.  I think this sort of a bootstrapping approach would
> be convincing that my enriched GO list was unique to the actual intersection
> of the two sets of peaks.
> Ideally, I¹d like to merge two annotated peak objects and then sample them for
> the number I observed in the intersection set (which I know). Since I¹m not
> that familiar with working with IRanges-based objects, can you suggest a code
> snippet that would work for me?
> Does that make sense?
> Ron

More information about the Bioconductor mailing list