Search
Question: Error for plotDEXSeq and warning for estimateDispersions
0
gravatar for Emmanouela Repapi
3.2 years ago by
United Kingdom
Emmanouela Repapi20 wrote:

Hello,

I have been trying to use DEXSeq to do a DEU analysis and plot some of the genes that are coming up in the results. The analysis seems to be working fine (even though I am getting a warning when I estimate the dispersion) and I have about 2000 hits, but when I try to visualise the differences using the plot I get this error:

 plotDEXSeq(dxr1, "ENSG00000163660", displayTranscripts=T,legend=TRUE,  cex.axis=1.2, cex=1.3, lwd=1)
Error in start(unlist(genomicData)) : 
  error in evaluating the argument 'x' in selecting a method for function 'start': Error in unlist(genomicData) : argument not a list

Could this error have something to do with the warning that Im getting: 

dxd <- estimateDispersions( dxd )
using supplied model matrix
using supplied model matrix
Warning message:
In estimateDispersionsFit(object, fitType = fitType, quiet = quiet) :
  the parametric fit of dispersion estimates over the mean of counts
failed, which occurs when the trend is not well captured by the
function y = a/x + b. A local regression fit is automatically performed,
and the analysis can continue. You can specify fitType='local' or 'mean'
to avoid this message if re-running the same data.
When using local regression fit, the user should examine plotDispEsts(dds)
to make sure the fitted line is not sharply curving up or down based on
the position of individual points.

 

I tried plotting the dispersion (using plotDispEsts) and it doesn't look too alarming so the error sounds unrelated to me. (let me know if it would be helpful to upload it) 

Any help would be very much appreciated! Many thanks in advance for your reply!

Best wishes,

Emmanouela

Here is my code: 

dxd <- DEXSeqDataSetFromHTSeq(my_samples, sampleData=sample_ids, flattenedfile=my_gft, design= ~ sample + exon + Tissue:exon ) # I am looking at DEU in different tissues (and my sample_ids have a column named Tissue instead of condition)
dxd <- estimateSizeFactors( dxd )
dxd <- estimateDispersions( dxd )
plotDispEsts( dxd )

dxd <- testForDEU( dxd )
dxd <- estimateExonFoldChanges( dxd, fitExpToVar="Tissue")
dxr1 <- DEXSeqResults( dxd )
plotDEXSeq(dxr1, "ENSG00000163660", displayTranscripts=T,legend=TRUE,  cex.axis=1.2, cex=1.3, lwd=1)

And here is my sessionInfo():

sessionInfo()
R version 3.0.1 (2013-05-16)
Platform: x86_64-unknown-linux-gnu (64-bit)

locale:
 [1] LC_CTYPE=en_GB.ISO-8859-1       LC_NUMERIC=C                    LC_TIME=en_GB.ISO-8859-1        LC_COLLATE=en_GB.ISO-8859-1    
 [5] LC_MONETARY=en_GB.ISO-8859-1    LC_MESSAGES=en_GB.ISO-8859-1    LC_PAPER=C                      LC_NAME=C                      
 [9] LC_ADDRESS=C                    LC_TELEPHONE=C                  LC_MEASUREMENT=en_GB.ISO-8859-1 LC_IDENTIFICATION=C            

attached base packages:
[1] parallel  stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] DEXSeq_1.10.6             BiocParallel_0.4.1        DESeq2_1.4.5              RcppArmadillo_0.4.450.1.0 Rcpp_0.11.3              
 [6] GenomicRanges_1.14.4      XVector_0.2.0             IRanges_1.20.7            Biobase_2.22.0            BiocGenerics_0.8.0       

loaded via a namespace (and not attached):
 [1] annotate_1.40.1      AnnotationDbi_1.24.0 BatchJobs_1.4        BBmisc_1.7           biomaRt_2.18.0       Biostrings_2.30.1   
 [7] bitops_1.0-6         brew_1.0-6           checkmate_1.4        codetools_0.2-9      DBI_0.3.1            digest_0.6.4        
[13] fail_1.2             foreach_1.4.2        genefilter_1.44.0    geneplotter_1.40.0   grid_3.0.1           hwriter_1.3.2       
[19] iterators_1.0.7      lattice_0.20-29      locfit_1.5-9.1       RColorBrewer_1.0-5   RCurl_1.95-4.3       Rsamtools_1.14.3    
[25] RSQLite_0.11.4       sendmailR_1.2-1      splines_3.0.1        statmod_1.4.20       stats4_3.0.1         stringr_0.6.2       
[31] survival_2.37-7      tools_3.0.1          XML_3.98-1.1         xtable_1.7-4         zlibbioc_1.8.0      

 

ADD COMMENTlink modified 3.2 years ago by Alejandro Reyes1.5k • written 3.2 years ago by Emmanouela Repapi20

I think the error is unrelated to the warning (which is coming from DESeq2). Also I've changed the warning to a message in the next release, because there's nothing the user needs to do. It's just a notification.

ADD REPLYlink written 3.2 years ago by Michael Love15k

That's what I thought too. Thanks for clarifying. 

Emmanouela

ADD REPLYlink written 3.2 years ago by Emmanouela Repapi20
0
gravatar for Alejandro Reyes
3.2 years ago by
Alejandro Reyes1.5k
Dana-Farber Cancer Institute, Boston, USA
Alejandro Reyes1.5k wrote:

Hi Emmanouela,

Could you include the code you used to generate your DEXSeqDataSet object?

Alejandro

ADD COMMENTlink written 3.2 years ago by Alejandro Reyes1.5k

Hi Alejandro,

Here is the code that I used to get the counts objects:

module load HTSeq
python /Filers/package/R/3.0.1/lib64/R/library/DEXSeq/python_scripts/dexseq_prepare_annotation.py ensembl_hg19_genes.gtf ensembl_hg19_genes_flattened.gtf
python /Filers/package/R/3.0.1/lib64/R/library/DEXSeq/python_scripts/dexseq_count.py -f bam -p yes -r pos -s no ensembl_hg19_genes_flattened.gtf A1_final_filtered.bam A1_final_filtered_counts

And then in R:

my_samples <- list.files(pattern="filtered_counts")
dxd <- DEXSeqDataSetFromHTSeq(my_samples, sampleData=sample_ids, flattenedfile="ensembl_hg19_genes_flattened.gtf", design= ~ sample + exon + Tissue:exon )

and the rest as above. 

Having looked at the error and the code in more detail the problem seems to be that plotDEXSeq expects dxr1$genomicData to be a list, when its a GRanges object, and so not a list. 

> class(dxr1$genomicData)
[1] "GRanges"
attr(,"package")
[1] "GenomicRanges"
>  head(unlist(dxr1$genomicData))
Error in head(unlist(dxr1$genomicData)) : 
  error in evaluating the argument 'x' in selecting a method for function 'head': Error in unlist(dxr1$genomicData) : argument not a list

Could that be a bug of some sort? 

Many thanks,

Emmanouela

ADD REPLYlink written 3.2 years ago by Emmanouela Repapi20
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 2.2.0
Traffic: 349 users visited in the last hour