Question: Problems running JunctionSeq
0
gravatar for jbono
3.0 years ago by
jbono0
United States
jbono0 wrote:

Hi,

I am trying to run the JunctionSeq package, but am running into a problem with the runJunctionsSeqAnalyses function that I can't seem to solve.  I entered the following:

> jscs <- runJunctionSeqAnalyses(
+     sample.files = countFiles,
+     sample.names = decoder$sample.ID,
+     condition = decoder$group.ID,
+     flat.gff.file = "outputData/countTables/withNovel.forJunctionSeq.gff.gz",
+     nCores = 1,
+     verbose=TRUE,
+     debug.mode = TRUE

It appears there might be an issue with DESeq (see below). Furthermore,the function does not seem to produce the jscs object because when I try the writeCompleteResults function it says it cannot find the jscs object. Any help on what might be causing the problem would be greatly appreciated.

This is what it returns from runJunctionSeqAnalyses (DEseq part is at the end):

> STARTING runJunctionSeqAnalyses (v1.3.5) (Thu Nov 17 14:47:55 2016)
> rJSA: sample.files:  outputData/countTables/C1/QC.spliceJunctionAndExonCounts.withNovel.forJunctionSeq.txt.gz, outputData/countTables/C2/QC.spliceJunctionAndExonCounts.withNovel.forJunctionSeq.txt.gz, outputData/countTables/C3/QC.spliceJunctionAndExonCounts.withNovel.forJunctionSeq.txt.gz, outputData/countTables/O1/QC.spliceJunctionAndExonCounts.withNovel.forJunctionSeq.txt.gz, outputData/countTables/O2/QC.spliceJunctionAndExonCounts.withNovel.forJunctionSeq.txt.gz, outputData/countTables/O3/QC.spliceJunctionAndExonCounts.withNovel.forJunctionSeq.txt.gz
> rJSA: sample.names:  C1, C2, C3, O1, O2, O3
> rJSA: condition:  caper, caper, caper, control, control, control
> rJSA: analysis.type:  junctionsAndExons
> rJSA: use.junctions:  TRUE
> rJSA: use.novel.junctions:  TRUE
> rJSA: use.exons:  TRUE
> rJSA: nCores:  1
> rJSA: use.covars:  
> rJSA: test.formula0:  ~ sample + countbin
> rJSA: test.formula1:  ~ sample + countbin + condition:countbin
> rJSA: use.multigene.aggregates:  FALSE
> rJSA: Reading Count files... Thu Nov 17 14:47:55 2016.
-> STARTING readJunctionSeqCounts (Thu Nov 17 14:47:55 2016)
---> RJSC; (v1.3.5)
---> RJSC: samplenames: C1,C2,C3,O1,O2,O3
---> RJSC: flat.gff.file: outputData/countTables/withNovel.forJunctionSeq.gff.gz
---> RJSC: use.exons:TRUE
---> RJSC: use.junctions:TRUE
---> RJSC: use.novel.junctions:TRUE
---> File read complete.
---> Extracted counts. Found 160953 features so far.
---> Extracted gene-level counts. Found: 16101 genes and aggregate-genes.
---> Removed gene features. Found: 144852 features to be included so far.
---> Note: 108353 counting bins from overlapping genes
--->          There are 7248 multigene aggregates.
--->          There are 25080 genes that are part of an aggregate.
---> Removed multigene-aggregate features. Found: 36499 features to be included so far.
---> Final feature count: 36499 features to be included in the analysis.
---> Extracted feature counts.
---> counts complete.
-----> reading annotation...
-----> formatting annotation...
-----> initial generation...
-----> creating jscs...
-----> generating count vectors... (Thu Nov 17 14:48:10 2016)
> Using single-core execution.
    getAllJunctionSeqCountVectors: dim(counts) = 36499,6 (Thu Nov 17 14:48:10 2016)
    getAllJunctionSeqCountVectors: dim(gct) = 16101,6
    getAllJunctionSeqCountVectors: out generated. dim = 36499,12 (Thu Nov 17 14:48:10 2016)
-----> count vectors generated (Thu Nov 17 14:48:10 2016)
-----> generating DESeqDataSet... (Thu Nov 17 14:48:10 2016)
Error in DESeqDataSet(se, design = design, ignoreRank) : 
  all samples have 0 counts for all genes. check the counting script.

This is returned from the writeCompleteResults function:

> STARTING writeCompleteResults (v1.3.5) (Thu Nov 17 14:57:36 2016)
> wcr: outfile.prefix:  outputData/junctionseq
> wcr: FDR.threshold:  0.01
> wcr: save.allGenes:  TRUE
> wcr: save.sigGenes:  TRUE
> wcr: save.fit:  FALSE
> wcr: save.VST:  FALSE
> wcr: bedtrack.format:  BED
Error in save(jscs, file = paste0(outfile.prefix, "jscs.RData")) : 
  object 'jscs' not found

Thank you, I appreciate it!

Jeremy

The output of the sessioninfo command is below:

R version 3.2.3 (2015-12-10)
Platform: x86_64-apple-darwin13.4.0 (64-bit)
Running under: OS X 10.11.3 (El Capitan)

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] stats4    parallel  stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] geneplotter_1.48.0         annotate_1.48.0            XML_3.98-1.5               AnnotationDbi_1.32.3      
 [5] genefilter_1.52.1          BiocParallel_1.4.3         Hmisc_4.0-0                ggplot2_2.2.0             
 [9] Formula_1.2-1              survival_2.40-1            lattice_0.20-34            locfit_1.5-9.1            
[13] stringr_1.1.0              plotrix_3.6-3              statmod_1.4.26             DESeq2_1.10.1             
[17] RcppArmadillo_0.7.500.0.0  Rcpp_0.12.7                JunctionSeq_1.3.5          SummarizedExperiment_1.0.2
[21] Biobase_2.30.0             GenomicRanges_1.22.4       GenomeInfoDb_1.6.3         IRanges_2.4.8             
[25] S4Vectors_0.8.11           BiocGenerics_0.16.1        BiocInstaller_1.20.3      

loaded via a namespace (and not attached):
 [1] splines_3.2.3        colorspace_1.3-0     htmltools_0.3.5      chron_2.3-47         DBI_0.5-1           
 [6] foreign_0.8-67       RColorBrewer_1.1-2   lambda.r_1.1.9       plyr_1.8.4           zlibbioc_1.16.0     
[11] munsell_0.4.3        gtable_0.2.0         futile.logger_1.4.3  latticeExtra_0.6-28  knitr_1.15          
[16] htmlTable_1.7        acepack_1.4.1        xtable_1.8-2         scales_0.4.1         XVector_0.10.0      
[21] gridExtra_2.2.1      digest_0.6.10        stringi_1.1.2        grid_3.2.3           tools_3.2.3         
[26] magrittr_1.5         RSQLite_1.0.0        lazyeval_0.2.0       tibble_1.2           cluster_2.0.5       
[31] futile.options_1.0.0 Matrix_1.2-7.1       data.table_1.9.6     assertthat_0.1       rpart_4.1-10        
[36] nnet_7.3-12 

junctionseq • 391 views
ADD COMMENTlink written 3.0 years ago by jbono0
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 16.09
Traffic: 353 users visited in the last hour