[BioC] DEXseq: fitDispersionFunction error

Sun Feb 19 23:46:52 CET 2012

Dear Elena,

Feb/17/12 4:15 PM, Alejandro Reyes scripsit::
>
> You can always use the svn to get the latests versions:
> http://www.bioconductor.org/developers/source-control/
>

Or here:
http://www.bioconductor.org/packages/devel/bioc/html/DEXSeq.html

	Best wishes
	Wolfgang

> Bests,
> Alejandro
>
>
>
>> Dear Alejandro,
>>
>> I can't seem to find the development version that you were talking
>> about. I am running the version available for BioC 2.9:
>> http://watson.nci.nih.gov/bioc_mirror/packages/2.9/bioc/html/DEXSeq.html
>>
>> ...Could you tell me where to download v 1.1.6 of your software?
>>
>> Thanks,
>> Elena
>>
>> On 2/17/2012 6:11 AM, Alejandro Reyes wrote:
>>> Dear Elena,
>>>
>>> I am glad you like the package!
>>> We recently fixed this bug, it should be fine if you use the most
>>> recent development
>>> version of DEXSeq (1.1.6)
>>>
>>> This should be solved by then, if not, could you send me your ecs
>>> object to give it a
>>> closer look?
>>>
>>> Thanks!
>>>
>>> Alejandro
>>>
>>>
>>> Greetings all,
>>>
>>> I ran into another problem with the DEXseq program, and I was hoping
>>> that somebody could help me once more. I successfully ran through a
>>> full analysis of one time point of my dataset, and was hoping ideally to
>>> extend the analysis to my full time course, using the GLM test. However,
>>> I ran into this error when running fitDispersionFunction:
>>>
>>> > ecs<- fitDispersionFunction(ecs)
>>> Warning messages:
>>> 1: In glmgam.fit(mm, disps[good], start = coefs) :
>>> Too much damping - convergence tolerance not achievable
>>> 2: In glmgam.fit(mm, disps[good], start = coefs) :
>>> Too much damping - convergence tolerance not achievable
>>>
>>> I couldn't find the answer in the list archive. Do the DEXseq authors or
>>> anyone else know what this error means?
>>>
>>> I'd be truly grateful for any advice about this - I really like the
>>> DEXseq package, and I'm hoping to use/reference it in my manuscript if I
>>> can get it to work! My full script is below.
>>>
>>> Sincerely,
>>> Elena
>>> -------------------------------------
>>> > library(DEXSeq)
>>> > annotationfile = file.path("../DEXSeq_annotations.gff")
>>> > annotationfile
>>> [1] "../DEXSeq_annotations.gff"
>>> > samples<- read.table("PLdesign.txt",header=T,row.names=1)
>>> > samples
>>> treatment time
>>> PL_D1hr_1 vehicle one
>>> PL_D1hr_2 vehicle one
>>> PL_U1hr_1 drug one
>>> PL_U1hr_2 drug one
>>> PL_D2hr_1 vehicle two
>>> PL_D2hr_2 vehicle two
>>> PL_U2hr_1 drug two
>>> PL_U2hr_2 drug two
>>> PL_D4hr_1 vehicle four
>>> PL_D4hr_2 vehicle four
>>> PL_U4hr_1 drug four
>>> PL_U4hr_2 drug four
>>> PL_D6hr_1 vehicle six
>>> PL_D6hr_2 vehicle six
>>> PL_U6hr_1 drug six
>>> PL_U6hr_2 drug six
>>> > fullFilenames<- list.files(full.names=TRUE,pattern="exon.txt")
>>> > fullFilenames
>>> [1] "./PL_1hrD1_exon.txt" "./PL_1hrD2_exon.txt" "./PL_1hrU1_exon.txt"
>>> [4] "./PL_1hrU2_exon.txt" "./PL_2hrD1_exon.txt" "./PL_2hrD2_exon.txt"
>>> [7] "./PL_2hrU1_exon.txt" "./PL_2hrU2_exon.txt" "./PL_4hrD1_exon.txt"
>>> [10] "./PL_4hrD2_exon.txt" "./PL_4hrU1_exon.txt" "./PL_4hrU2_exon.txt"
>>> [13] "./PL_6hrD1_exon.txt" "./PL_6hrD2_exon.txt" "./PL_6hrU1_exon.txt"
>>> [16] "./PL_6hrU2_exon.txt"
>>> > ecs<- read.HTSeqCounts(countfiles = fullFilenames,design =
>>> samples,flattenedfile = annotationfile)
>>> # Check out ecs object
>>> > head(counts(ecs))
>>> ./PL_1hrD1_exon.txt ./PL_1hrD2_exon.txt ./PL_1hrU1_exon.txt
>>> 2L52.1:001 12 16 7
>>> 2L52.1:002 32 22 19
>>> 2L52.1:003 18 16 14
>>> 2L52.1:004 19 16 20
>>> 2L52.1:005 28 22 9
>>> 2L52.1:006 19 15 22
>>> ./PL_1hrU2_exon.txt ./PL_2hrD1_exon.txt ./PL_2hrD2_exon.txt
>>> 2L52.1:001 4 7 4
>>> 2L52.1:002 19 23 16
>>> 2L52.1:003 20 15 17
>>> 2L52.1:004 11 8 17
>>> 2L52.1:005 17 23 27
>>> 2L52.1:006 20 19 23
>>> ./PL_2hrU1_exon.txt ./PL_2hrU2_exon.txt ./PL_4hrD1_exon.txt
>>> 2L52.1:001 2 4 0
>>> 2L52.1:002 18 13 13
>>> 2L52.1:003 15 5 15
>>> 2L52.1:004 12 4 8
>>> 2L52.1:005 18 14 13
>>> 2L52.1:006 15 10 13
>>> ./PL_4hrD2_exon.txt ./PL_4hrU1_exon.txt ./PL_4hrU2_exon.txt
>>> 2L52.1:001 2 1 1
>>> 2L52.1:002 9 6 8
>>> 2L52.1:003 5 3 17
>>> 2L52.1:004 2 3 7
>>> 2L52.1:005 17 4 5
>>> 2L52.1:006 16 4 6
>>> ./PL_6hrD1_exon.txt ./PL_6hrD2_exon.txt ./PL_6hrU1_exon.txt
>>> 2L52.1:001 4 2 0
>>> 2L52.1:002 20 13 4
>>> 2L52.1:003 8 7 7
>>> 2L52.1:004 7 4 6
>>> 2L52.1:005 14 16 9
>>> 2L52.1:006 12 13 6
>>> ./PL_6hrU2_exon.txt
>>> 2L52.1:001 2
>>> 2L52.1:002 12
>>> 2L52.1:003 8
>>> 2L52.1:004 4
>>> 2L52.1:005 6
>>> 2L52.1:006 4
>>> > design(ecs)
>>> treatment time
>>> ./PL_1hrD1_exon.txt vehicle one
>>> ./PL_1hrD2_exon.txt vehicle one
>>> ./PL_1hrU1_exon.txt drug one
>>> ./PL_1hrU2_exon.txt drug one
>>> ./PL_2hrD1_exon.txt vehicle two
>>> ./PL_2hrD2_exon.txt vehicle two
>>> ./PL_2hrU1_exon.txt drug two
>>> ./PL_2hrU2_exon.txt drug two
>>> ./PL_4hrD1_exon.txt vehicle four
>>> ./PL_4hrD2_exon.txt vehicle four
>>> ./PL_4hrU1_exon.txt drug four
>>> ./PL_4hrU2_exon.txt drug four
>>> ./PL_6hrD1_exon.txt vehicle six
>>> ./PL_6hrD2_exon.txt vehicle six
>>> ./PL_6hrU1_exon.txt drug six
>>> ./PL_6hrU2_exon.txt drug six
>>> > head(fData(ecs))
>>> geneID exonID testable dispBeforeSharing dispFitted dispersion
>>> 2L52.1:001 2L52.1 E001 NA NA NA NA
>>> 2L52.1:002 2L52.1 E002 NA NA NA NA
>>> 2L52.1:003 2L52.1 E003 NA NA NA NA
>>> 2L52.1:004 2L52.1 E004 NA NA NA NA
>>> 2L52.1:005 2L52.1 E005 NA NA NA NA
>>> 2L52.1:006 2L52.1 E006 NA NA NA NA
>>> pvalue padjust chr start end strand transcripts
>>> 2L52.1:001 NA NA chrII 1867 1911 + 2L52.1
>>> 2L52.1:002 NA NA chrII 2506 2694 + 2L52.1
>>> 2L52.1:003 NA NA chrII 2738 2888 + 2L52.1
>>> 2L52.1:004 NA NA chrII 2931 3036 + 2L52.1
>>> 2L52.1:005 NA NA chrII 3406 3552 + 2L52.1
>>> 2L52.1:006 NA NA chrII 3802 3984 + 2L52.1
>>> # All in all, Exon Count Set looks good
>>> > ecs<- estimateSizeFactors(ecs)
>>> > sizeFactors(ecs)
>>> ./PL_1hrD1_exon.txt ./PL_1hrD2_exon.txt ./PL_1hrU1_exon.txt
>>> ./PL_1hrU2_exon.txt
>>> 1.7543472 1.6354950 1.6969669
>>> 1.5496025
>>> ./PL_2hrD1_exon.txt ./PL_2hrD2_exon.txt ./PL_2hrU1_exon.txt
>>> ./PL_2hrU2_exon.txt
>>> 1.4932349 1.3392481 1.4051111
>>> 1.1780430
>>> ./PL_4hrD1_exon.txt ./PL_4hrD2_exon.txt ./PL_4hrU1_exon.txt
>>> ./PL_4hrU2_exon.txt
>>> 0.9026051 0.5879963 0.5206059
>>> 0.5726255
>>> ./PL_6hrD1_exon.txt ./PL_6hrD2_exon.txt ./PL_6hrU1_exon.txt
>>> ./PL_6hrU2_exon.txt
>>> 0.7345296 0.7536120 0.7763199
>>> 0.7529934
>>> > formuladispersion<- count ~ sample + ( exon + treatment ) * time
>>> > ecs<- estimateDispersions( ecs, formula = formuladispersion,
>>> nCores=11)
>>> Estimating Cox-Reid exon dispersion estimates using 11 cores. (Progress
>>> report: one dot per 100 genes)
>>> ..........................................................................................................................................................>
>>>
>>>
>>> > /usr/local/lib64/R/library
>>> > ecs<- fitDispersionFunction(ecs)
>>> Warning messages:
>>> 1: In glmgam.fit(mm, disps[good], start = coefs) :
>>> Too much damping - convergence tolerance not achievable
>>> 2: In glmgam.fit(mm, disps[good], start = coefs) :
>>> Too much damping - convergence tolerance not achievable
>>>
>>> > sessionInfo()
>>> R version 2.14.0 (2011-10-31)
>>> 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] stats graphics grDevices utils datasets methods base
>>>
>>> other attached packages:
>>> [1] multicore_0.1-7 DEXSeq_1.0.2 Biobase_2.14.0
>>>
>>> loaded via a namespace (and not attached):
>>> [1] hwriter_1.3 plyr_1.7.1 statmod_1.4.14 stringr_0.6
>>> tools_2.14.0
>>> >
>>>
>>
>
> _______________________________________________
> Bioconductor mailing list
> Bioconductor at r-project.org
> https://stat.ethz.ch/mailman/listinfo/bioconductor
> Search the archives:
> http://news.gmane.org/gmane.science.biology.informatics.conductor

-- 
Best wishes
	Wolfgang

Wolfgang Huber
EMBL
http://www.embl.de/research/units/genome_biology/huber