Entering edit mode
Iddo Ben-dov
▴
20
@iddo-ben-dov-6603
Last seen 10.3 years ago
hello,
i am testing the change of miRNA expression between 2 time periods in a set of matched patient samples (100 ‘before', 100 ‘after’).
after filtering very low coverage miRNA i am left with ~500 miRNA
i then define:
countData <- as.matrix(round(raw_reads)) # raw_reads is the dataframe of miRNA counts (rows) by sample (columns) colData <- data.frame (patient, period) # period is a ‘before' and ‘after’ factor; patient is 1:100 ddsE <- DESeqDataSetFromMatrix(countData = countData, colData = colData, design = ~ period+patient) ddsE<-DESeq(ddsE)
the following messages appear:
> ddsE<-DESeq(ddsE) estimating size factors estimating dispersions gene-wise dispersion estimates Error in simpleLoess(y, x, w, span, degree, parametric, drop.square, normalize, : NA/NaN/Inf in foreign function call (arg 1)
this error does not occur if i further diminish my gene list by filtering additional low coverage genes.
can anyone suggest a solution that will allow me to test all the miRNA of interrest using DESeq2?
thanks,
iddo
> traceback() 8: simpleLoess(y, x, w, span, degree, parametric, drop.square, normalize, control$statistics, control$surface, control$cell, iterations, control$trace.hat) 7: loess(lpo ~ s, span = 0.1) 6: fitDispInR(y = counts(objectNZ)[refitDisp, , drop = FALSE], x = modelMatrix, mu = mu[refitDisp, , drop = FALSE], logAlphaPriorMean = rep(0, sum(refitDisp)), logAlphaPriorSigmaSq = 1, usePrior = FALSE) 5: estimateDispersionsGeneEst(object, maxit = maxit, quiet = quiet) 4: .local(object, ...) 3: estimateDispersions(object, fitType = fitType, quiet = quiet) 2: estimateDispersions(object, fitType = fitType, quiet = quiet) 1: DESeq(ddsE) > sessionInfo() R version 3.1.1 (2014-07-10) Platform: x86_64-apple-darwin13.1.0 (64-bit) locale: [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8 attached base packages: [1] parallel stats graphics grDevices utils datasets methods base other attached packages: [1] ggplot2_1.0.0 RColorBrewer_1.0-5 gplots_2.14.2 ReportingTools_2.4.0 [5] AnnotationDbi_1.26.1 Biobase_2.24.0 RSQLite_0.11.4 DBI_0.3.1 [9] knitr_1.6 DESeq2_1.4.5 RcppArmadillo_0.4.450.1.0 Rcpp_0.11.3 [13] GenomicRanges_1.16.4 GenomeInfoDb_1.0.2 IRanges_1.22.10 BiocGenerics_0.10.0 [17] plyr_1.8.1 genefilter_1.46.1 matrixStats_0.10.0 pheatmap_0.7.7 [21] scatterplot3d_0.3-35 edgeR_3.6.8 limma_3.20.9 loaded via a namespace (and not attached): [1] acepack_1.3-3.3 annotate_1.42.1 AnnotationForge_1.6.1 base64enc_0.1-2 [5] BatchJobs_1.4 BBmisc_1.7 BiocParallel_0.6.1 biomaRt_2.20.0 [9] Biostrings_2.32.1 biovizBase_1.12.3 bitops_1.0-6 brew_1.0-6 [13] BSgenome_1.32.0 Category_2.30.0 caTools_1.17.1 checkmate_1.4 [17] cluster_1.15.3 codetools_0.2-9 colorspace_1.2-4 dichromat_2.0-0 [21] digest_0.6.4 evaluate_0.5.5 fail_1.2 foreach_1.4.2 [25] foreign_0.8-61 formatR_1.0 Formula_1.1-2 gdata_2.13.3 [29] geneplotter_1.42.0 GenomicAlignments_1.0.6 GenomicFeatures_1.16.3 ggbio_1.12.10 [33] GO.db_2.14.0 GOstats_2.30.0 graph_1.42.0 grid_3.1.1 [37] gridExtra_0.9.1 GSEABase_1.26.0 gtable_0.1.2 gtools_3.4.1 [41] Hmisc_3.14-5 hwriter_1.3.2 iterators_1.0.7 KernSmooth_2.23-13 [45] labeling_0.3 lattice_0.20-29 latticeExtra_0.6-26 locfit_1.5-9.1 [49] MASS_7.3-35 Matrix_1.1-4 munsell_0.4.2 nnet_7.3-8 [53] PFAM.db_2.14.0 proto_0.3-10 R.methodsS3_1.6.1 R.oo_1.18.0 [57] R.utils_1.34.0 RBGL_1.40.1 RCurl_1.95-4.3 reshape2_1.4 [61] rpart_4.1-8 Rsamtools_1.16.1 rtracklayer_1.24.2 scales_0.2.4 [65] sendmailR_1.2-1 splines_3.1.1 stats4_3.1.1 stringr_0.6.2 [69] survival_2.37-7 tools_3.1.1 VariantAnnotation_1.10.5 XML_3.98-1.1 [73] xtable_1.7-4 XVector_0.4.0 zlibbioc_1.10.0
Why are you rounding your read counts? ie, your
countData <- as.matrix(round(raw_reads))
. These should already be integers, or are you sending data intoDESeq2
that's not *actual* counts? (ie. RSEM,for instance(?))these are actual counts, but occasionally a read maps with distance 1 to multiple locations and these reads are split