[BioC] DEXseq: fitDispersionFunction error
Alejandro Reyes
alejandro.reyes at embl.de
Fri Feb 17 16:15:49 CET 2012
Dear Elena,
You can always use the svn to get the latests versions:
http://www.bioconductor.org/developers/source-control/
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
>> >
>>
>
More information about the Bioconductor
mailing list