Hi All,
I am interested in using the ChIPseqR package for the MNASE data. I am new to using R, and I quite don't understand the error i am receiving.
My input is paired end BAM file, and I executed the following as per the reference manual.
aln1<- readAligned(path,pattern,type="BAM")
sp <- strandPileup(aln1,chrLen=chr,extend=1,plot=FALSE)
Error: UserArgumentMismatch if 'shift' is not 0L, then it must be a vector of integer values see ?"AlignedRead-class"
I get the same error while using simpleNucCall straight. I can only find this shift argument in the "coverage" function, which strandPileup is dependent I guess.
I appreciate any suggestion or help on this.
PS: My BAM file contains reads from all chromosomes, so the chr here is named integer vector.
Thanks,
RT
traceback()
16: stop(cond)
15: .throw(SRError("UserArgumentMismatch", "if '%s' is not 0L, then it must be a vector of integer values\n see %s",
"shift", "?\"AlignedRead-class\""))
14: .throw(SRError("UserArgumentMismatch", "if '%s' is not 0L, then it must be a vector of integer values\n see %s",
"shift", "?\"AlignedRead-class\""))
13: .local(x, shift, width, weight, ...)
12: coverage(fwd, extend = ext1, coords = coords, ...)
11: coverage(fwd, extend = ext1, coords = coords, ...)
10: (function (filter, len, extend, aligned, ...)
{
sfilter1 <- filter & strand(aligned) == "+"
sfilter2 <- filter & strand(aligned) == "-"
fwd <- aligned[sfilter1]
rev <- aligned[sfilter2]
ext1 <- -width(fwd) + as.integer(extend)
ext2 <- -width(rev) + as.integer(extend)
counts <- list(coverage(fwd, extend = ext1, coords = coords,
...)[[1]], coverage(coverage(rev, extend = ext2, coords = coords,
...))[[1]])
counts <- .fixCounts(counts, len)
if (!compress)
counts <- decompress(counts)
counts
})(dots[[1L]][[1L]], dots[[2L]][[1L]], dots[[3L]][[1L]], aligned = <S4 object of class "AlignedRead">,
list())
9: .Method(..., FUN = FUN, MoreArgs = MoreArgs, SIMPLIFY = SIMPLIFY,
USE.NAMES = USE.NAMES)
8: eval(expr, envir, enclos)
7: eval(.dotsCall, env)
6: eval(.dotsCall, env)
5: standardGeneric("mapply")
4: mapply(function(filter, len, extend, aligned, ...) {
sfilter1 <- filter & strand(aligned) == "+"
sfilter2 <- filter & strand(aligned) == "-"
fwd <- aligned[sfilter1]
rev <- aligned[sfilter2]
ext1 <- -width(fwd) + as.integer(extend)
ext2 <- -width(rev) + as.integer(extend)
counts <- list(coverage(fwd, extend = ext1, coords = coords,
...)[[1]], coverage(coverage(rev, extend = ext2, coords = coords,
...))[[1]])
counts <- .fixCounts(counts, len)
if (!compress)
counts <- decompress(counts)
counts
}, chrFilter, chrLen, extend, MoreArgs = list(aligned = aligned,
list(...)), SIMPLIFY = FALSE)
3: .local(aligned, chrLen, ...)
2: strandPileup(aln1, chrLen = chr, extend = 1, plot = FALSE)
1: strandPileup(aln1, chrLen = chr, extend = 1, plot = FALSE)
sessioninfo()
R version 3.1.1 (2014-07-10) Platform: x86_64-w64-mingw32/x64 (64-bit)
locale: [1] LC_COLLATE=English_Singapore.1252 LC_CTYPE=English_Singapore.1252 [3] LC_MONETARY=English_Singapore.1252 LC_NUMERIC=C [5] LC_TIME=English_Singapore.1252
attached base packages: [1] grid parallel stats graphics grDevices utils datasets methods [9] base
other attached packages: [1] ChIPseqR_1.18.0 nucleR_1.12.0 Biobase_2.24.0 [4] Gviz_1.8.4 rtracklayer_1.24.2 PICS_2.8.0 [7] PING_2.8.0 chipseq_1.14.0 ShortRead_1.22.0 [10] GenomicAlignments_1.0.6 Rsamtools_1.16.1 BiocParallel_0.6.1 [13] BSgenome_1.32.0 Biostrings_2.32.1 XVector_0.4.0 [16] GenomicRanges_1.16.4 GenomeInfoDb_1.0.2 IRanges_1.22.10 [19] BiocGenerics_0.10.0 BiocInstaller_1.14.2
loaded via a namespace (and not attached): [1] acepack_1.3-3.3 AnnotationDbi_1.26.1 base64enc_0.1-2 [4] BatchJobs_1.4 BBmisc_1.7 biomaRt_2.20.0 [7] biovizBase_1.12.3 bitops_1.0-6 brew_1.0-6 [10] checkmate_1.4 cluster_1.15.2 codetools_0.2-8 [13] colorspace_1.2-4 DBI_0.3.1 dichromat_2.0-0 [16] digest_0.6.4 fail_1.2 fBasics_3010.86 [19] foreach_1.4.2 foreign_0.8-61 Formula_1.1-2 [22] GenomicFeatures_1.16.3 HilbertVis_1.22.0 Hmisc_3.14-5 [25] hwriter_1.3.2 iterators_1.0.7 lattice_0.20-29 [28] latticeExtra_0.6-26 MASS_7.3-33 matrixStats_0.10.0 [31] munsell_0.4.2 nnet_7.3-8 plyr_1.8.1 [34] R.methodsS3_1.6.1 RColorBrewer_1.0-5 Rcpp_0.11.3 [37] RCurl_1.95-4.3 rpart_4.1-8 RSQLite_0.11.4 [40] scales_0.2.4 sendmailR_1.2-1 splines_3.1.1 [43] stabledist_0.6-6 stats4_3.1.1 stringr_0.6.2 [46] survival_2.37-7 timeDate_3010.98 timeSeries_3010.97 [49] timsac_1.3.3 tools_3.1.1 VariantAnnotation_1.10.5 [52] XML_3.98-1.1 zlibbioc_1.10.0
I'm not a user of ChIPseqr, but it might be that
extend
is supposed to be a integer (rather than floating-point numeric). Force it to be an integer by appending 'L' or usingas.integer(1)
, e.g.,strandPileup(aln1,chrLen=chr,extend=1L,plot=FALSE).
Is there particular functionality that ChIPseqr provides? I ask becausereadAligned
is not really a modern function for reading BAM files -- useGenomicAlignments::readGAlignments
or related functions instead.Thanks for the suggestion Martin. I made it integer, but I still receive the same error! I am using readAligned , because the manual suggests that the input data to strandPileup or the SimpleNucCall be in objects of class AlignedRead or data.frame. readGAlignments gives me a GAlignmentPairs object, which needs some conversion and I will give it a try.
Apologies for not responding earlier, I only just saw this. I'm investigating the problem now (the good news is that I can reproduce the issue). I should be able to commit a fix in the next few days.