I am attempting to use DESeq2 to model differential abundance in microbiome count data. I have been performing my other analyses in phyloseq, and I used the phyloseq_to_deseq2() conversion prior to attempting differential abundance analysis:
diffAbund = subset_samples(realData, CaptureNumber == 10)
diffAbund = phyloseq_to_deseq2(diffAbund, ~ Sex)
No error messages so far. But when I run DESeq, I get the following message:
diffAbund = DESeq(diffAbund, test = "Wald", fitType = "parametric")
estimating size factors
Error in estimateSizeFactorsForMatrix(counts(object), locfunc = locfunc, : every gene contains at least one zero, cannot compute log geometric means
A google search revealed that using estimateSizeFactors with type = "iterate" could solve the problem. However, when I attempt that strategy I get the following error message:
estimateSizeFactors(diffAbund, type = "iterate") Error in estimateSizeFactorsIterate(object) : iterative size factor normalization did not converge
I have checked the data for NAs to make sure that isn't the source of the problem, and there are none in the dataset. Any suggestions for how to get around this issue would be much appreciated.
Here is my session info:
R version 3.4.1 (2017-06-30) Platform: x86_64-w64-mingw32/x64 (64-bit) Running under: Windows >= 8 x64 (build 9200) Matrix products: default locale: [1] LC_COLLATE=English_United States.1252 [2] LC_CTYPE=English_United States.1252 [3] LC_MONETARY=English_United States.1252 [4] LC_NUMERIC=C [5] LC_TIME=English_United States.1252 attached base packages: [1] stats4 parallel stats graphics grDevices utils [7] datasets methods base other attached packages: [1] DESeq2_1.16.1 SummarizedExperiment_1.6.3 [3] DelayedArray_0.2.7 matrixStats_0.52.2 [5] Biobase_2.36.2 GenomicRanges_1.28.3 [7] GenomeInfoDb_1.12.2 IRanges_2.10.2 [9] S4Vectors_0.14.3 BiocGenerics_0.22.0 [11] lmerTest_2.0-33 lme4_1.1-13 [13] Matrix_1.2-10 reshape2_1.4.2 [15] vegan_2.4-3 lattice_0.20-35 [17] permute_0.9-4 dplyr_0.7.1 [19] ggplot2_2.2.1 phyloseq_1.20.0 [21] ancom.R_1.1-3 doParallel_1.0.10 [23] iterators_1.0.8 shiny_1.0.5 [25] Rcpp_0.12.12 foreach_1.4.3 loaded via a namespace (and not attached): [1] nlme_3.1-131 bitops_1.0-6 [3] bit64_0.9-7 RColorBrewer_1.1-2 [5] tools_3.4.1 backports_1.1.0 [7] R6_2.2.2 DT_0.2 [9] rpart_4.1-11 DBI_0.7 [11] Hmisc_4.0-3 lazyeval_0.2.0 [13] mgcv_1.8-17 colorspace_1.3-2 [15] ade4_1.7-6 nnet_7.3-12 [17] gridExtra_2.2.1 bit_1.1-12 [19] compiler_3.4.1 htmlTable_1.9 [21] exactRankTests_0.8-29 sandwich_2.4-0 [23] labeling_0.3 scales_0.4.1 [25] checkmate_1.8.3 mvtnorm_1.0-6 [27] genefilter_1.58.1 stringr_1.2.0 [29] digest_0.6.12 foreign_0.8-69 [31] minqa_1.2.4 XVector_0.16.0 [33] base64enc_0.1-3 pkgconfig_2.0.1 [35] htmltools_0.3.6 htmlwidgets_0.9 [37] rlang_0.1.1 RSQLite_2.0 [39] bindr_0.1 zoo_1.8-0 [41] jsonlite_1.5 BiocParallel_1.10.1 [43] acepack_1.4.1 RCurl_1.95-4.8 [45] magrittr_1.5 modeltools_0.2-21 [47] GenomeInfoDbData_0.99.0 Formula_1.2-2 [49] biomformat_1.4.0 munsell_0.4.3 [51] ape_4.1 stringi_1.1.5 [53] multcomp_1.4-6 MASS_7.3-47 [55] zlibbioc_1.22.0 rhdf5_2.20.0 [57] plyr_1.8.4 blob_1.1.0 [59] grid_3.4.1 Biostrings_2.44.1 [61] splines_3.4.1 annotate_1.54.0 [63] multtest_2.32.0 locfit_1.5-9.1 [65] knitr_1.16 igraph_1.0.1 [67] geneplotter_1.54.0 codetools_0.2-15 [69] XML_3.98-1.9 glue_1.1.1 [71] latticeExtra_0.6-28 data.table_1.10.4 [73] nloptr_1.0.4 httpuv_1.3.5 [75] gtable_0.2.0 assertthat_0.2.0 [77] mime_0.5 coin_1.2-1 [79] xtable_1.8-2 survival_2.41-3 [81] tibble_1.3.3 memoise_1.1.0 [83] AnnotationDbi_1.38.1 bindrcpp_0.2 [85] cluster_2.0.6 TH.data_1.0-8
Thank you Michael, that worked perfectly!