[BioC] edgeR: Likelihood ratio test error
Gordon K Smyth
smyth at wehi.EDU.AU
Wed Aug 8 09:13:17 CEST 2012
Dear KJ Lim,
The problem is that you are following an example from the User's Guide
that goes with the developmental version of edgeR, a newer version of egeR
which has not yet been released. The syntax for the glmLRT() function has
been changed compared with the official release version of edgeR that you
are using.
For your version of edgeR you need
glmLRT(y, hl.fit, contrast=hl.contrasts[,"HvsL.00"])
instead of
glmLRT(hl.fit, contrast=hl.contrasts[,"HvsL.00"])
In other words, you can follow the examples from the developmental User's
Guide, but you will need to add an extra argument to each glmLRT call. The
first argument of glmLRT() needs to be the data object used to create the
fit object for glmLRT() in your version of edgeR.
Don't forget that you need to look at the help pages for the functions you
are using, not just the User's Guide, because the help pages explain the
correct syntax for each function.
Best wishes
Gordon
----------------- original message ------------------
[BioC] edgeR: Likelihood ratio test error
KJ Lim jinkeanlim at gmail.com
Tue Aug 7 15:58:10 CEST 2012
Dear edgeR users community,
Good day.
I'm analysing my RNA-Seq data which has 2 different genotypes (treeHS,
treeLS) and time of treatment (0H, 3H,
24H, 96H) with edgeR. I would like to learn the differential
expression of treeHS vs treeLS at specific time points i.e. 24H.
Thanks to Prof Gordon and his colleagues, the latest edgeR user's
guide (Chapter 3) is very useful for my case.
I encountered an error message when I carried out the likelihood ratio
test with either makeContrasts or coefficients approach. The error
message:
> hl.contrasts <- makeContrasts(HvsL.00 = H.00H-L.00H, HvsL.03 =
(H.03H-H.00H)-(L.03H-L.00H),LvsH.00 = L.00H-H.00H,LvsH.03 =
(L.03H-L.00H)-(H.03H-H.00H), levels=hl.design)
> lrt.HvsL00 <- glmLRT(hl.fit, contrast=hl.contrasts[,"HvsL.00"])
Error in array(x, c(length(x), 1L), if (!is.null(names(x)))
list(names(x), :
dims [product 0] do not match the length of object [11]
---------------------------------------------------------------------------------------
> colnames(hl.fit)
[1] "H.00H" "H.03H" "H.24H" "H.96H" "L.00H" "L.03H" "L.24H" "L.96H"
> lrt.HvsL00 <- glmLRT(hl.fit,contrast=c(1,0,0,0,-1,0,0,0))
Error in array(x, c(length(x), 1L), if (!is.null(names(x)))
list(names(x), :
dims [product 0] do not match the length of object [11]
Could the community kindly please light me on this matter? What
correction I should make in this case?
Thank you very much for your help.
Best regards,
KJ Lim
> sessionInfo()
R version 2.15.1 (2012-06-22)
Platform: i486-pc-linux-gnu (32-bit)
locale:
[1] LC_CTYPE=en_US.utf8 LC_NUMERIC=C
[3] LC_TIME=en_US.utf8 LC_COLLATE=en_US.utf8
[5] LC_MONETARY=en_US.utf8 LC_MESSAGES=en_US.utf8
[7] LC_PAPER=C LC_NAME=C
[9] LC_ADDRESS=C LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.utf8 LC_IDENTIFICATION=C
attached base packages:
[1] splines stats graphics grDevices utils datasets methods
[8] base
other attached packages:
[1] edgeR_2.6.10 limma_3.12.1
loaded via a namespace (and not attached):
[1] tools_2.15.1
On 6 August 2012 10:11, Gordon K Smyth <smyth at wehi.edu.au> wrote:
> Dear KJ Lim,
>
> Well, your "tree" column has 16 different entries, a different entry for
> every library. Naturally this causes a problem.
>
> So either you have no replication of any experimental condition, or you
> have mistakingly used different labels for the same condition.
>
> In your earlier email below, you claimed to have just two genotypes
> (tree=HS and tree=LS), but in your latest targets file you have 16
> different genotypes. I'm guessing that the genotypes are actually H and
> L, so the entries in the targets file should be just H and L.
>
> That's as much advice as I have time to give you, and I hope that others
> will help you further. I feel that I've given you good general advice
> that is applicable to your experiment.
>
>
> Best wishes Gordon
>
> ---------------------------------------------
> Professor Gordon K Smyth,
> Bioinformatics Division,
> Walter and Eliza Hall Institute of Medical Research,
> 1G Royal Parade, Parkville, Vic 3052, Australia.
> Tel: (03) 9345 2326, Fax (03) 9347 0852,
> http://www.statsci.org/smyth
>
> On Mon, 6 Aug 2012, KJ Lim wrote:
>
>> Dear Prof Gordon,
>>
>> Good day. Thanks for you prompt replied.
>>
>> Below is the output of cbind(targets,Group=Group)
>>
>>> cbind(targets,Group=Group)
>>
>> files treat tree X Group
>> 1 ./rawData/HS0H01.txt 00H H1 NA 00H.H1
>> 2 ./rawData/HS0H02.txt 00H H2 NA 00H.H2
>> 3 ./rawData/HS3H01.txt 03H H3 NA 03H.H3
>> 4 ./rawData/HS3H02.txt 03H H4 NA 03H.H4
>> 5 ./rawData/HS1D01.txt 24H H5 NA 24H.H5
>> 6 ./rawData/HS1D02.txt 24H H6 NA 24H.H6
>> 7 ./rawData/HS4D01.txt 96H H7 NA 96H.H7
>> 8 ./rawData/HS4D02.txt 96H H8 NA 96H.H8
>> 9 ./rawData/LS0H01.txt 00H L1 NA 00H.L1
>> 10 ./rawData/LS0H02.txt 00H L2 NA 00H.L2
>> 11 ./rawData/LS3H01.txt 03H L3 NA 03H.L3
>> 12 ./rawData/LS3H02.txt 03H L4 NA 03H.L4
>> 13 ./rawData/LS1D01.txt 24H L5 NA 24H.L5
>> 14 ./rawData/LS1D02.txt 24H L6 NA 24H.L6
>> 15 ./rawData/LS4D01.txt 96H L7 NA 96H.L7
>> 16 ./rawData/LS4D02.txt 96H L8 NA 96H.L8
>>
>> Thank you for your time.
>>
>> Best regards,
>> KJ Lim
>>
>>
>> On 6 August 2012 07:33, Gordon K Smyth <smyth at wehi.edu.au> wrote:
>>>
>>> The error message is pretty self-explanatory. Apparently your design
>>> matrix has as many columns as there are libraries, so there are no
>>> degrees of freedom left from which to estimate variability.
>>> According to the code and output you have given, this message should
>>> be impossible, so you must have changed the data in some way that I
>>> cannot see.
>>>
>>> If you had shown the output from cbind(targets,Group=Group), then I
>>> might have been able to say something useful.
>>>
>>> Best wishes
>>> Gordon
>>>
>>>
>>> On Mon, 6 Aug 2012, KJ Lim wrote:
>>>
>>>> Dear Prof Gordon,
>>>>
>>>> Good day.
>>>>
>>>> Thanks for your time to update the edgeR User's Guide. It is useful.
>>>>
>>>> I combined all my experiment factors into one combined factor and
>>>> described the experiment matrix like the example in the User's guide:
>>>>
>>>> > Group <- factor(paste(targets$treat, targets$tree,sep="."))
>>>> > cbind(targets,Group=Group)
>>>> > hl.design <- model.matrix(~0+Group)
>>>>
>>>> But, I encountered an error when I performed the estimate dispersion
>>>> for my data
>>>>
>>>> > hl <- estimateGLMCommonDisp(hl, hl.design)
>>>> Warning message:
>>>> In estimateGLMCommonDisp.default(y = y$counts, design = design, :
>>>> No residual df: setting dispersion to NA
>>>>
>>>> I tried to figure out what I have done wrong, unfortunately, I have
>>>> no luck on that. Could you or the community kindly please light me
>>>> for this matter?
>>>>
>>>> Thank you very much for your time.
>>>>
>>>>> sessionInfo()
>>>>
>>>>
>>>> R version 2.15.1 (2012-06-22)
>>>> Platform: i486-pc-linux-gnu (32-bit)
>>>>
>>>> locale:
>>>> [1] LC_CTYPE=en_US.utf8 LC_NUMERIC=C
>>>> [3] LC_TIME=en_US.utf8 LC_COLLATE=en_US.utf8
>>>> [5] LC_MONETARY=en_US.utf8 LC_MESSAGES=en_US.utf8
>>>> [7] LC_PAPER=C LC_NAME=C
>>>> [9] LC_ADDRESS=C LC_TELEPHONE=C
>>>> [11] LC_MEASUREMENT=en_US.utf8 LC_IDENTIFICATION=C
>>>>
>>>> attached base packages:
>>>> [1] stats graphics grDevices utils datasets methods base
>>>>
>>>> other attached packages:
>>>> [1] edgeR_2.6.10 limma_3.12.1
>>>>
>>>> loaded via a namespace (and not attached):
>>>> [1] tcltk_2.15.1 tools_2.15.1
>>>>
>>>> Best regards,
>>>> KJ Lim
______________________________________________________________________
The information in this email is confidential and intend...{{dropped:4}}
More information about the Bioconductor
mailing list