Hi:
I have set of genomic interval in GRanges objects. However I intend to filter them conditionally, while export them as bed format files. I have implemented function that possibly accomplish this task, but I have an error. I have used dplyr
packages for facilitating output more efficiently. I simulated same data set as data.frame and it worked perfectly. How can I fix this issue? Any quick solution ?
mini example:
grs <- GRangesList( foo = GRanges( seqnames=Rle("chr1", 4),ranges=IRanges(c(3,33,54,91), c(26,42,71,107)), rangeName=c("a1", "a2", "a3", "a4"), score=c(22, 6,13, 7)), bar = GRanges(seqnames=Rle("chr1", 6),ranges=IRanges(c(5,12,16,21,37,78), c(9,14,19,29,45,84)), rangeName=c("b1", "b2", "b3", "b4", "b5", "b6"), score=c(3, 5, 3, 9, 4, 3)), cat = GRanges(seqnames=Rle("chr1", 7),ranges=IRanges(c(1,8,18,35,42,59,81), c(6,13,27,40,46,63,114)), rangeName=c("c1", "c2", "c3", "c4","c5","c6","c7"), score= c(2.1, 3, 5.1, 3.5, 7, 2, 10)) )
helper function :
.cast.Pvl<- function(x, pvalueBase = 1L, ...) { if(is.null(x$pvalue)){ x$p.value <- 10^(score(x)/(- pvalueBase)) colnames(mcols(x))[3] <- "p.value" } else { x } return(x) }
I intend to implement this sketch function :
mylist2 <- lapply(seq_along(grs), function(i) { if(is.null(i$p.value)) { i <- .cast.Pvl(i, 1L) } else { splitter <- function(i, tau.w) { require(dplyr) DF <- grs[[i]] DF <- as(DF, "data.frame") DF %>% filter(p.value >= tau.w) %>% write.csv(., sprintf("dropped.%s.csv", i), row.names = FALSE) total.ERs <- filter(DF, p.value <= tau.w) return(total.ERs) } } res <- splitter(i, 1.0E-04) })
I have an error like this:
Error in i$p.value : $ operator is invalid for atomic vectors
can't figure out what is possible reason. Can anyone tell me what's going on? How can I fix this bug ? How can I facilitate the output of myList as I expected ? Thanks a lot
Thanks for your respond. The main purpose of I implement above sketch function, given list of genomic interval, first I need to filter them by its p.value, but I let my function only return genomic regions whose p.value less then given threshold, while I need also export dropped features as Bed files (just statistical showed which one dropped). Michael suggest do not use writing my own bed exported. What's the feasible approach ? Any suggestion to make more code more clear and clean ? Thank you
The rtracklayer package contains the function
export()
which can take aGRanges
object and output a number of different formats including bed files.