This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Problem with finding overlapping chips peaks with ChIPpeakAnno

I have 2 datasets 1. control and other is 2. RA and I am trying to find overlapping chip seq peaks in theses 2 datasets with the script below:

#To read peak tables#
macspeaks<- read.table("MYCN_ControlPeaks.txt", header=TRUE, sep="\t", stringsAsFactors=FALSE)
head (macspeaks)
macspeaks.RA<- read.table("MYCN_RAPeaks.txt", header=TRUE, sep="\t", stringsAsFactors=FALSE)
head (macspeaks.RA)
--------------------
# function to convert table into a IRanges#`

Control.peaks = RangedData(IRanges(start=macspeaks$Start, end=macspeaks$End, names=paste(macspeaks$Chr,macspeaks$Start,sep=":")),strand="*")
head (Control.peaks)

RA.peaks = RangedData(IRanges(start=macspeaks.RA$Start, end=macspeaks.RA$End, names=paste(macspeaks.RA$Chr,macspeaks.RA$Start,sep=":")),strand="*")
head (RA.peaks)
-----------------------
#to find overlapping peaks in control and RA#

t1 =findOverlapsOfPeaks (Control.peaks, RA.peaks, maxgap=0L, minoverlap=1L, NameOfPeaks="Control","RA", select="all", annotate=1)
r = t1$OverlappingPeaks
pie(table(r$overlapFeature))
as.data.frame(t1$MergedPeaks)

-------------------------

While running it I am getting the following error in the find overlapping peaks

+ t1 =findOverlapsOfPeaks (Control.peaks, RA.peaks, maxgap=0L, minoverlap=1L, NameOfPeaks="Control","RA", select="all", annotate=1)
Error in findOverlapsOfPeaks(Control.peaks, RA.peaks, maxgap = 0L, minoverlap = 1L,  : 
  The length of input peaks list should no more than 5

> r = t1$OverlappingPeaks
Error in t1$OverlappingPeaks :
  object of type 'closure' is not subsettable
> pie(table(r$overlapFeature))
Error in eval(expr, envir, enclos) : object 'r' not found
> as.data.frame(t1$MergedPeaks)
Error in as.data.frame(t1$MergedPeaks) :
  error in evaluating the argument 'x' in selecting a method for function 'as.data.frame': Error in t1$MergedPeaks : object of type 'closure' is not subsettable

I am still in the learning stage of R and have very limited knowledge.

Kindly help to solve the problem.

Thank you

chip-seq r

1 answer

One problem is that you have an error in your function call: NameOfPeaks="Control","RA"

This will cause an error because the comma separates the arguments to the function, so the "RA" string is all alone as an argument. If you want to hand NameOfPeaks more than one thing you have to use c(), as in NameOfPeaks=c("Control","RA"). Either way, I don't know what purpose this argument serves for this function. To debug your problem, two approaches you can take are (1) it can be helpful to create some toy data, to see where things are failing, and (2) simplify your function call to the bare minimum.

# load the library
library(ChIPpeakAnno)

# create 5 random segments on an interval
# that are at least 50 bases wide
peakStarts <- sample(1:5000,5)
peakEnds <- peakStarts + 50 + sample(1:300,1)

# populate Control with fake data
Control.peaks = RangedData(IRanges(start=peakStarts,
                  end=peakEnds,
                  names=paste(rep("Chr1",ntimes=length(peakStarts)),peakStarts,sep=":")),
                  strand="*")

# get 5 new random locations
peakStarts <- sample(1:5000,5)
peakEnds <- peakStarts + 50 + sample(1:300,1)

# populate the RA with fake data
RA.peaks = RangedData(IRanges(start=peakStarts,
             end=peakEnds,
             names=paste(rep("Chr1",ntimes=length(peakStarts)),peakStarts,sep=":")),
             strand="*")

# check to see that we have a legit data structure
head (RA.peaks)

# RangedData with 5 rows and 1 value column across 1 space
#            space       ranges |      strand
#          <factor>    <IRanges> | <character>
# Chr1:1675        1 [1675, 1802] |           *
# Chr1:4661        1 [4661, 4788] |           *
# Chr1:3516        1 [3516, 3643] |           *
# Chr1:1959        1 [1959, 2086] |           *
# Chr1:3031        1 [3031, 3158] |           *

# find overlaps
t1 <- findOverlapsOfPeaks(Control.peaks, RA.peaks, maxgap=0, minoverlap=1)

You might check the help for ?findOverlapsOfPeaks

Thanks will try and get back to you.

Thanks,
Debbie

Log in to answer this question.