Entering edit mode
My script was working a month ago, but now it's not. I'm struggling with getting my DESeq2 analysis working. Here's what I'm trying to do:
factorNames <- c("Experiment", "Line");
dds <-
DESeqDataSetFromMatrix(count.mat, meta.df,
as.formula(paste0("~ ", paste(factorNames,
collapse=" + "))));
dds <- DESeq(dds);
results(dds); # error happens after this command
Problem
What I would like to appear is a results table. I get the same error when I only specify a single variable in the formula, i.e.:
dds <- DESeqDataSetFromMatrix(count.mat, meta.df, ~ Line);
dds <- DESeq(dds);
results(dds); # error happens after this command
And this is the output:
Error in .wrap_in_length_one_list_like_object(value, name, x) :
failed to coerce 'list(value)' to a DESeqResults object of length 1
What else I've Tried (which still produce the error)
Reordering my library declarations at the start of my script:
library(DESeq2); library(tidyverse); library(magrittr);
## vs
library(tidyverse); library(magrittr); library(DESeq2);
Updating DESeq2 via BiocManager::install("DESeq2")
Possibly Useful Things
Here's the count matrix layout:
> str(count.mat)
num [1:9934, 1:18] 1 0 7 0 0 0 0 45 379 392 ...
- attr(*, "dimnames")=List of 2
..$ : chr [1:9934] "ENSMUST00000192271_rev" "ENSMUST00000127957_rev" "ENSMUST00000073783_fwd" "ENSMUST00000122876_rev" ...
..$ : chr [1:18] "CGJan19_006.BC08" "AJ_01-03_Jun19.BC01" "AJ_01-03_Jun19.BC02" "AJ_01-03_Jun19.BC03" ...
And my metadata data frame:
> meta.df
# A tibble: 18 x 5
SampleID Label Experiment Line Replicate
<chr> <chr> <chr> <fct> <fct>
1 CGJan19_006.BC08 CG006.GL261ρ0_D106#A CG006 p0D106 A
2 AJ_01-03_Jun19.BC01 AJJun19.GL261ρ0_D106#B AJJun19 p0D106 B
3 AJ_01-03_Jun19.BC02 AJJun19.GL261ρ0_D106#C AJJun19 p0D106 C
4 AJ_01-03_Jun19.BC03 AJJun19.GL261ρ0_D106#D AJJun19 p0D106 D
5 CGNov19_011.BC06 CGNov19.GL261ρ0_V1#A CG011 p0V1 A
6 AJ_01-03_Jun19.BC09 AJJun19.GL261ρ0_V1#B AJJun19 p0V1 B
7 AJ_01-03_Jun19.BC10 AJJun19.GL261ρ0_V1#C AJJun19 p0V1 C
8 CGJan19_006.BC07 CG006.GL261ρ0_V3#A CG006 p0V3 A
9 AJ_01-03_Jun19.BC04 AJJun19.GL261ρ0_V3#B AJJun19 p0V3 B
10 AJ_01-03_Jun19.BC05 AJJun19.GL261ρ0_V3#C AJJun19 p0V3 C
11 AJ_01-03_Jun19.BC06 AJJun19.GL261ρ0_V3#D AJJun19 p0V3 D
12 CGNov19_011.BC07 CGNov19.GL261.TMZR#A CG011 TMZR A
13 CGNov19_011.BC08 CGNov19.GL261.TMZR#B CG011 TMZR B
14 CGNov19_011.BC09 CGNov19.GL261.TMZR#C CG011 TMZR C
15 CGJan19_006.BC06 CG006.GL261.WT#A CG006 WT A
16 AJ_01-03_Jun19.BC07 AJJun19.GL261.WT#C AJJun19 WT C
17 AJ_01-03_Jun19.BC08 AJJun19.GL261.WT#D AJJun19 WT D
18 CGNov19_011.BC05 CGNov19.GL261.WT#D CG011 WT D
And finally, sessionInfo():
> sessionInfo()
R version 3.6.1 (2019-07-05)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Debian GNU/Linux bullseye/sid
Matrix products: default
BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.9.0
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.9.0
locale:
[1] LC_CTYPE=en_NZ.UTF-8 LC_NUMERIC=C
[3] LC_TIME=en_NZ.UTF-8 LC_COLLATE=en_NZ.UTF-8
[5] LC_MONETARY=en_NZ.UTF-8 LC_MESSAGES=en_NZ.UTF-8
[7] LC_PAPER=en_NZ.UTF-8 LC_NAME=C
[9] LC_ADDRESS=C LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_NZ.UTF-8 LC_IDENTIFICATION=C
attached base packages:
[1] parallel stats4 stats graphics grDevices utils datasets
[8] methods base
other attached packages:
[1] DESeq2_1.24.0 SummarizedExperiment_1.16.0
[3] DelayedArray_0.12.0 BiocParallel_1.20.0
[5] matrixStats_0.55.0 Biobase_2.46.0
[7] GenomicRanges_1.38.0 GenomeInfoDb_1.22.0
[9] IRanges_2.20.1 S4Vectors_0.24.1
[11] BiocGenerics_0.32.0 magrittr_1.5
[13] forcats_0.4.0 stringr_1.4.0
[15] dplyr_0.8.3 purrr_0.3.3
[17] readr_1.3.1 tidyr_1.0.0
[19] tibble_2.1.3 ggplot2_3.2.1
[21] tidyverse_1.3.0
loaded via a namespace (and not attached):
[1] nlme_3.1-142 bitops_1.0-6 fs_1.3.1
[4] bit64_0.9-7 lubridate_1.7.4 RColorBrewer_1.1-2
[7] httr_1.4.1 tools_3.6.1 backports_1.1.5
[10] utf8_1.1.4 R6_2.4.1 rpart_4.1-15
[13] Hmisc_4.3-0 DBI_1.0.0 lazyeval_0.2.2
[16] colorspace_1.4-1 nnet_7.3-12 withr_2.1.2
[19] tidyselect_0.2.5 gridExtra_2.3 bit_1.1-14
[22] compiler_3.6.1 cli_1.1.0 rvest_0.3.5
[25] htmlTable_1.13.3 xml2_1.2.2 scales_1.1.0
[28] checkmate_1.9.4 genefilter_1.68.0 digest_0.6.23
[31] foreign_0.8-72 XVector_0.26.0 base64enc_0.1-3
[34] pkgconfig_2.0.3 htmltools_0.4.0 dbplyr_1.4.2
[37] htmlwidgets_1.5.1 rlang_0.4.2 readxl_1.3.1
[40] RSQLite_2.1.2 rstudioapi_0.10 generics_0.0.2
[43] jsonlite_1.6 acepack_1.4.1 RCurl_1.95-4.12
[46] GenomeInfoDbData_1.2.2 Formula_1.2-3 Matrix_1.2-18
[49] fansi_0.4.0 Rcpp_1.0.3 munsell_0.5.0
[52] lifecycle_0.1.0 stringi_1.4.3 zlibbioc_1.32.0
[55] blob_1.2.0 grid_3.6.1 crayon_1.3.4
[58] lattice_0.20-38 haven_2.2.0 splines_3.6.1
[61] annotate_1.64.0 hms_0.5.2 locfit_1.5-9.1
[64] zeallot_0.1.0 knitr_1.26 pillar_1.4.2
[67] geneplotter_1.62.0 XML_3.98-1.20 reprex_0.3.0
[70] glue_1.3.1 latticeExtra_0.6-28 data.table_1.12.6
[73] BiocManager_1.30.10 modelr_0.1.5 vctrs_0.2.1
[76] cellranger_1.1.0 gtable_0.3.0 assertthat_0.2.1
[79] xfun_0.11 xtable_1.8-4 broom_0.5.3
[82] survival_3.1-7 memoise_1.1.0 AnnotationDbi_1.48.0
[85] cluster_2.1.0
Update [2019-Dec-20]
After comments from Mike Love, I've now uninstalled some packages from the base Debian system, and reinstalled packages into my local folders, and I'm getting a different error at a different step:
> dds <- DESeqDataSetFromMatrix(count.mat, meta.df, ~ Line);
converting counts to integer mode
Error in validObject(.Object) :
invalid class “DESeqDataSet” object: superclass "vector_OR_Vector" not defined in the environment of the object's class
So this seems to be an error due to a broken installation, and it's probably going to be difficult to get help for a problem relating to the specific combination of packages that my system decided to install. I'll have a go at wiping out the local R installation (i.e. rm -rf ~/R/x86_64-pc-linux-gnu-library), and reinstall DESeq2, and see where that gets me....
... and now it works. Thank you for bearing with me while I do the software equivalent of, "Have you tried turning it off and on again?"
> dds <-
+ DESeqDataSetFromMatrix(count.mat, meta.df,
+ as.formula(paste0("~ ", paste(factorNames,
+ collapse=" + "))));
converting counts to integer mode
Warning message:
In DESeqDataSet(se, design = design, ignoreRank) :
some variables in design formula are characters, converting to factors
> dds <- DESeq(dds);
estimating size factors
estimating dispersions
gene-wise dispersion estimates
mean-dispersion relationship
final dispersion estimates
fitting model and testing
> results(dds);
log2 fold change (MLE): Line WT vs p0D106
Wald test p-value: Line WT vs p0D106
DataFrame with 9934 rows and 6 columns
baseMean log2FoldChange lfcSE
<numeric> <numeric> <numeric>
ENSMUST00000192271_rev 26563.2466352323 -1.09277749904499 0.377854584351656
ENSMUST00000127957_rev 9026.17510382961 -0.168498805001082 0.429131257538863
ENSMUST00000073783_fwd 7867.65681244628 -0.969776170520552 0.415097461404838
ENSMUST00000122876_rev 4950.21480608967 -0.906867087205864 0.316912053075736
ENSMUST00000154732_rev 3787.28318276015 -0.808128426231181 0.354406394204153
... ... ... ...
ENSMUST00000201831_fwd 0.755790685509435 -0.923018792648103 2.86138231601317
ENSMUST00000168645_fwd 0.670469974454824 -1.58950422743093 2.98741494666397
ENSMUST00000044373_fwd 1.49097034760506 -0.584856475816726 2.38605077355759
ENSMUST00000034452_fwd 0.588833538409114 -0.934693043540383 3.00417556365216
ENSMUST00000102758_rev 2.42279945004559 -0.0280350877040124 2.24570756265536
stat pvalue
<numeric> <numeric>
ENSMUST00000192271_rev -2.89205833222863 0.00382726880343014
ENSMUST00000127957_rev -0.392650971097864 0.694577287407033
ENSMUST00000073783_fwd -2.33626138603325 0.0194776227007531
ENSMUST00000122876_rev -2.86157335577621 0.00421543940158844
ENSMUST00000154732_rev -2.2802309423505 0.0225939950836063
... ... ...
ENSMUST00000201831_fwd -0.322577932869232 0.747014908051351
ENSMUST00000168645_fwd -0.532066772045149 0.594679750039667
ENSMUST00000044373_fwd -0.245114849314253 0.806367505494896
ENSMUST00000034452_fwd -0.311131298333338 0.755700808493421
ENSMUST00000102758_rev -0.0124838550531768 0.990039583510677
padj
<numeric>
ENSMUST00000192271_rev 0.0826304581229773
ENSMUST00000127957_rev 0.9449149448463
ENSMUST00000073783_fwd 0.196809244865185
ENSMUST00000122876_rev 0.0884652702389295
ENSMUST00000154732_rev 0.21487097859105
... ...
ENSMUST00000201831_fwd NA
ENSMUST00000168645_fwd NA
ENSMUST00000044373_fwd NA
ENSMUST00000034452_fwd NA
ENSMUST00000102758_rev NA
The very last step, producing the results table
Can you check if you have a valid installation of Bioc? There is a valid() function in BiocManager.
No...
19 packages out-of-date; 18 packages too new. I'll see if running the recommended installations works for me:I have now uninstalled as many packages from the Debian system as it will let me get away with (i.e. without demanding that I also remove
r-base), so that they could be updated in my user folder, andBiocManager::valid()now returns true:But
BiocManager::install()still complains about non-updateable packages:Unfortunately things were broken in doing that uninstall. I also ran
install.packages("pkgconfig", "png", "htmlwidgets", "xfun", "tidyverse", "xtable", "bit64", "blob", "memoise")to get thelibrary(DESeq2)command happy. Presumably these were packages that were installed in the Debian base system, andBiocManager::valid()glossed over them.* *sigh **
Debian....
I also can't guarantee that it will resolve this bug you saw, but hopefully it may.
Having a combination of packages across releases can lead to mysterious errors.