Hi,

I am trying to compute the coverage histogram of my bedpe file (converted from paired-end bam file) by following the HelloRanges tutorial, but it gives me error while using eval:

>code <- bedtools_coverage("-a 1.bedpe -hist -b org.genome")
​>code

{
    genome <- Seqinfo(genome = NA_character_)
    gr_a <- import("1.bedpe", genome = genome)
    gr_b <- import("org.genome", genome = genome)
    cov <- unname(coverage(gr_b)[gr_a])
    tab <- t(table(cov))
    tab <- cbind(tab, all = rowSums(tab))
    covhist <- DataFrame(as.table(tab))
    colnames(covhist) <- c("coverage", "a", "count")
    len <- c(lengths(cov, use.names = FALSE), sum(lengths(cov)))
    covhist$len <- rep(len, each = nrow(tab))
    covhist <- subset(covhist, count > 0L)
    covhist$fraction <- with(covhist, count/len)
    ans <- gr_a
    covhistList <- split(covhist, ~a)[, -2L]
    mcols(ans)$coverage <- head(covhistList, -1L)
    metadata(ans)$coverage <- covhistList$all
    ans
}
>ans <- eval(code)
Error: logical subscript contains NAs

Kindly guide!

Thanks

You can just run the code trunk between { }  like your own code.

I would run it line by line and see what happens.

I am guessing a variable in the code trunk is having wrong value or wrong format.