Dear list, My name is Eva Benito, I'm trying to implement DEXSeq for my RNASeq data. I aligned with tophat2 and counted with the corresponding scripts for dexseq in HTSeq. No problems there. Then I read the data into R and because I have quite some samples (18), I decided to give the TRT version of the DEXSeq functions a try, since they were specifically designed to optimize running time for a large number of samples (as I found them in this post here: https://stat.ethz.ch/pipermail/bioconductor/2013-January/050252.html). I keep getting a warning after I fit the dispersion function (see below) which does not appear when I run the regular version of estimateDispersions() and I'm wondering whether there's something I'm missing and whether there's something wrong with this. I found these two posts (https://stat.ethz.ch/pipermail/bioconductor/2012-February/043644.html and http://grokbase.com/t/r/bioconductor/11bx58k9hk/bioc-error-in- dexseq-fitdispersionfunction) and it seems that back then the solution was to update DEXSeq to the newest version, but I am running the latest version as far as I know and I even installed the development version and still got the warning. Any input in this respect would be greatly appreciated! Thanks a lot in advance and best regards, Eva Here's the warning I get: > subset = fitDispersionFunction(subset) Warning messages: 1: In glmgam.fit(mm, disps[good], coef.start = coefs) : Too much damping - convergence tolerance not achievable 2: In glmgam.fit(mm, disps[good], coef.start = coefs) : Too much damping - convergence tolerance not achievable Here's the code I'm using: library(DEXSeq) library(parallel) files = list.files(pattern = 'CA_[0-9]') metadata = read.delim('metadata_CA.txt', stringsAsFactors = TRUE) cds.exons = read.HTSeqCounts(files, metadata, flattenedfile = '/usr/users/ebenito/mm9_ensembl_dexseq.gff') subset = cds.exons[1:1000,] subset = estimateSizeFactors(subset) subset = estimateDispersionsTRT(subset) subset = fitDispersionFunction(subset) subset = testForDEUTRT(subset) res = DEUresultTable(subset) And here's the sessionInfo(): R version 3.0.1 (2013-05-16) Platform: x86_64-unknown-linux-gnu (64-bit) locale: [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8 [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8 [7] LC_PAPER=C LC_NAME=C [9] LC_ADDRESS=C LC_TELEPHONE=C [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C attached base packages: [1] parallel stats graphics grDevices utils datasets methods [8] base other attached packages: [1] DEXSeq_1.7.3 Biobase_2.20.0 BiocGenerics_0.6.0 loaded via a namespace (and not attached): [1] biomaRt_2.16.0 Biostrings_2.28.0 bitops_1.0-5 [4] GenomicRanges_1.12.4 hwriter_1.3 IRanges_1.18.1 [7] RCurl_1.95-4.1 Rsamtools_1.12.3 statmod_1.4.17 [10] stats4_3.0.1 stringr_0.6.2 XML_3.98-1.1 [13] zlibbioc_1.2.0 -- output of sessionInfo(): R version 3.0.1 (2013-05-16) Platform: x86_64-unknown-linux-gnu (64-bit) locale: [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8 [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8 [7] LC_PAPER=C LC_NAME=C [9] LC_ADDRESS=C LC_TELEPHONE=C [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C attached base packages: [1] parallel stats graphics grDevices utils datasets methods [8] base other attached packages: [1] DEXSeq_1.7.3 Biobase_2.20.0 BiocGenerics_0.6.0 loaded via a namespace (and not attached): [1] biomaRt_2.16.0 Biostrings_2.28.0 bitops_1.0-5 [4] GenomicRanges_1.12.4 hwriter_1.3 IRanges_1.18.1 [7] RCurl_1.95-4.1 Rsamtools_1.12.3 statmod_1.4.17 [10] stats4_3.0.1 stringr_0.6.2 XML_3.98-1.1 [13] zlibbioc_1.2.0 -- Sent via the guest posting facility at bioconductor.org.

Hi Eva > Here's the warning I get: > >> subset = fitDispersionFunction(subset) > Warning messages: > 1: In glmgam.fit(mm, disps[good], coef.start = coefs) : > Too much damping - convergence tolerance not achievable > 2: In glmgam.fit(mm, disps[good], coef.start = coefs) : > Too much damping - convergence tolerance not achievable The warning means that the dispersion fit may not describe well the meausred versions, and hence, you should check it. To do so, make a plot such as Figure 1 in the vignette. (In the vignette, we use a function called "plotDispEst" for this purpose which is defined in the Sweave file. Copy and paste it from the DEXSeq.Rnw file -- and yes, we will move this to a more convenient place soon.) Check whether the red line fits the estimated dispersions well, and maybe compare the plot obtained form using the old function with that from the new TRT functions. Maybe post the plots if you are unsure. Also, you seem to work on a small subset of the exons: > subset = cds.exons[1:1000,] Subsetting is useful for testing purposes, but going down to only 1000 exons might thin out information too much. Also, etter chose them at random, because often, the first items in a list of things are not typical. Try subset = cds.exons[ sample.int( nrow(cds.exons), 10000 ), ] Simon