I'm interested in differential exon usage between genotypes in my RNAseq experiment. I've already aligned reads using HISAT2 on the Galaxy web platform, and now I'm trying to generate exon count tables for analysis using DEXSeq. As I'm a bioinformatics novice but have some very limited R experience, I've been using the work flows described in Appendix B of the Bioconductor DEXSeq documentation (http://bioconductor.org/packages/release/bioc/vignettes/DEXSeq/inst/doc/DEXSeq.pdf) and in the parathyroidGenesSE vignette (http://bioconductor.org/packages/release/data/experiment/vignettes/parathyroidSE/inst/doc/parathyroidSE.pdf) to preprocess my data - i.e. generate the count tables. This has been successful through creation of the summarized experiment object using summarizeOverlaps. However, when I then try to generate a DEXSeq dataset from this summarized experiment using the DEXSeqDataSetFromSE command, I receive an error saying that "coercing an AtomicList object to an atomic vector is supported only for objects with top-level elements of length <= 1". I'm not sure what this means, or how to fix it!
Here's what I've run, including the error traceback output:
> hse = makeTxDbFromBiomart(biomart="ensembl",dataset="hsapiens_gene_ensembl") > exonicParts=disjointExons(hse, aggregateGenes=TRUE) > fls=list.files("G:/x/galaxy hisat2 output from trimmed files", pattern=".bam",full=TRUE) > bamlst=BamFileList(fls, yieldSize=100000,obeyQname=TRUE) > SE=summarizeOverlaps(exonicParts,bamlst,mode="Union",singleEnd=FALSE,ignore.strand=TRUE,inter.feature=FALSE,fragments=TRUE) > colData(SE)$condition=c("A","A","A","B","B","B","nc","WT","WT","WT") > DEXSeqDataSetFromSE( SE, design= ~ sample + exon + condition:exon) Error in as.vector(x, mode) : coercing an AtomicList object to an atomic vector is supported only for objects with top-level elements of length <= 1 > traceback() 8: stop("coercing an AtomicList object to an atomic vector ", "is supported only for\n", " objects with top-level elements of length <= 1") 7: as.vector(x, mode) 6: as.vector(x, mode) 5: as.vector(x, mode = "character") 4: .local(x, ...) 3: as.character(mcols(SE)$gene_id) 2: as.character(mcols(SE)$gene_id) 1: DEXSeqDataSetFromSE(SE, design = ~sample + exon + condition:exon)
Here's the session info output:
R version 3.4.2 (2017-09-28) Platform: x86_64-w64-mingw32/x64 (64-bit) Running under: Windows 7 x64 (build 7601) Service Pack 1 Matrix products: default locale:  LC_COLLATE=English_United Kingdom.1252 LC_CTYPE=English_United Kingdom.1252 LC_MONETARY=English_United Kingdom.1252 LC_NUMERIC=C  LC_TIME=English_United Kingdom.1252 attached base packages:  parallel stats4 stats graphics grDevices utils datasets methods base other attached packages:  DEXSeq_1.24.1 RColorBrewer_1.1-2 DESeq2_1.18.1 BiocParallel_1.12.0 GenomicAlignments_1.14.1 Rsamtools_1.30.0  Biostrings_2.46.0 XVector_0.18.0 SummarizedExperiment_1.8.0 DelayedArray_0.4.1 matrixStats_0.52.2 GenomicFeatures_1.30.0  AnnotationDbi_1.40.0 Biobase_2.38.0 GenomicRanges_1.30.0 GenomeInfoDb_1.14.0 IRanges_2.12.0 S4Vectors_0.16.0  BiocGenerics_0.24.0 biomaRt_2.34.0 loaded via a namespace (and not attached):  RMySQL_0.10.13 tidyr_0.7.2 bit64_0.9-7 splines_3.4.2 Formula_1.2-2 assertthat_0.2.0  statmod_1.4.30 latticeExtra_0.6-28 blob_1.1.0 GenomeInfoDbData_0.99.1 progress_1.1.2 RSQLite_2.0  backports_1.1.1 lattice_0.20-35 glue_1.2.0 digest_0.6.12 checkmate_1.8.5 colorspace_1.3-2  htmltools_0.3.6 Matrix_1.2-12 plyr_1.8.4 XML_3.98-1.9 pkgconfig_2.0.1 genefilter_1.60.0  zlibbioc_1.24.0 xtable_1.8-2 purrr_0.2.4 snow_0.4-2 scales_0.5.0 annotate_1.56.1  htmlTable_1.11.0 tibble_1.3.4 ggplot2_2.2.1 nnet_7.3-12 lazyeval_0.2.1 survival_2.41-3  magrittr_1.5 memoise_1.1.0 hwriter_1.3.2 foreign_0.8-69 data.table_1.10.4-3 tools_3.4.2  prettyunits_1.0.2 stringr_1.2.0 locfit_1.5-9.1 munsell_0.4.3 cluster_2.0.6 bindrcpp_0.2  compiler_3.4.2 rlang_0.1.4 grid_3.4.2 RCurl_1.95-4.8 rstudioapi_0.7 htmlwidgets_0.9  bitops_1.0-6 base64enc_0.1-3 gtable_0.2.0 DBI_0.7 R6_2.2.2 gridExtra_2.3  knitr_1.17 dplyr_0.7.4 rtracklayer_1.38.1 bit_1.1-12 bindr_0.1 Hmisc_4.0-3  stringi_1.1.6 Rcpp_0.12.14 geneplotter_1.56.0 rpart_4.1-11 acepack_1.4.1
Do you have any ideas where I might be going wrong? Thanks in advance.