[BioC] Lumi methylation analysis errors in lmFit
Sean Davis
sdavis2 at mail.nih.gov
Thu Jun 2 23:41:00 CEST 2011
On Thu, Jun 2, 2011 at 5:13 PM, Jasreet <jovialj86 at gmail.com> wrote:
> Hi Sean,
> So I tried using only 2 samples by subsetting as well as using all 12
> samples to find differentially expressed genes.
>
> It doesn't throw out an error now in the lmFit step, rather in the next one:
>> dim(selDataMatrix)
> [1] 485574 12
>> dim(design)
> [1] 12 12
> >fit <- lmFit(selDataMatrix, design)
>> fit <- eBayes(fit)
> *Error in ebayes(fit = fit, proportion = proportion, stdev.coef.lim =
> stdev.coef.lim, :
> No residual degrees of freedom in linear model fits*
Hi, Jasreet.
This is expected. As I mentioned, it is not possible to perform a
statistical test of differences between two samples and the error
message just says that but uses statistical terminology to do so.
This limitation is not specific to your data or to limma.
You can calculate a fold-change and then rank CpG sites based on that
if you like.
Sean
>> traceback()
> 3: stop("No residual degrees of freedom in linear model fits")
> 2: ebayes(fit = fit, proportion = proportion, stdev.coef.lim =
> stdev.coef.lim,
> trend = trend)
> 1: eBayes(fit)
>
> On Thu, Jun 2, 2011 at 2:06 PM, Sean Davis <sdavis2 at mail.nih.gov> wrote:
>
>> On Thu, Jun 2, 2011 at 2:51 PM, Jasreet <jovialj86 at gmail.com> wrote:
>> > Hi Sean,
>> >
>> > So cant I select just two samples to see the differential expression? Do
>> I
>> > have to do it for all 12 together?
>>
>> Well, you will not be able to perform differential methylation between
>> two samples using limma, no. The simplest way to do that is to make a
>> ratio between the two samples and then use arbitrary thresholds to
>> determine differential methylation. This is not a limma or even
>> programming issue but a statistical one.
>>
>> As for subsetting samples, the design matrix is not useful for doing
>> that. You would need to subset your data first. You can use normal
>> subsetting behavior of R to get that effect. For example, to subset
>> the first two samples, you could do:
>>
>> mldat[,1:2]
>>
>> You'll see that object then includes only two samples.
>>
>> Sean
>>
>>
>> > On Thu, Jun 2, 2011 at 1:23 PM, Sean Davis <sdavis2 at mail.nih.gov> wrote:
>> >
>> >> Hi, Jasreet. See below.
>> >>
>> >> On Thu, Jun 2, 2011 at 2:07 PM, Jasreet <jovialj86 at gmail.com> wrote:
>> >> > Hi Sean,
>> >> >
>> >> > Here is a brief overview of each object and some steps .. I basically
>> >> want
>> >> > differentially methylated genes between two samples (which I am
>> selecting
>> >> > through sampleType) out of total 12 samples.
>> >> >
>> >> >> mldat
>> >> > MethyLumiM (storageMode: lockedEnvironment)
>> >> > assayData: 485577 features, 12 samples
>> >> > element names: detection, exprs, methylated, unmethylated
>> >> > protocolData: none
>> >> > phenoData
>> >> > sampleNames: 3y0eY1 3y0eY2 ... 3y0eXY (12 total)
>> >> > varLabels: sampleID label
>> >> > varMetadata: labelDescription
>> >> > featureData
>> >> > featureNames: cg00000029 cg00000108 ... rs9839873 (485577 total)
>> >> > fvarLabels: TargetID ProbeID_A ProbeID_B COLOR_CHANNEL
>> >> > fvarMetadata: labelDescription
>> >> > experimentData: use 'experimentData(object)'
>> >> > Annotation:
>> >> >
>> >> >
>> >> >
>> >> >>presentCount <- detectionCall(mldat)
>> >> >
>> >> >> presentCount[1:3]
>> >> > cg00000029 cg00000108 cg00000109
>> >> > 12 12 12
>> >> >
>> >> >> dim(dataMatrix)
>> >> > [1] 485577 12
>> >> >
>> >> >>selDataMatrix <- dataMatrix[presentCount > 0,]
>> >> >
>> >> >> dim(presentCount)
>> >> > NULL
>> >> >
>> >> >> dim(selDataMatrix)
>> >> > [1] 485574 12
>> >> >
>> >> >>probeList <- rownames(selDataMatrix)
>> >> >
>> >> >> sampleType
>> >> > [1] "3y0eY2" "3y0eY0"
>> >> >
>> >> >>design <- model.matrix(~ factor(sampleType))
>> >>
>> >> If I am not mistaken, sampleType is length 2 while you have 12
>> >> samples? Not sure how that will work.
>> >>
>> >> Sean
>> >>
>> >> > On Thu, Jun 2, 2011 at 12:13 PM, Sean Davis <sdavis2 at mail.nih.gov>
>> >> wrote:
>> >> >
>> >> >> Hi, Jasreet. It looks like the design matrix and the selDataMatrix
>> >> >> object are not the same dimensions. You may need to post some code
>> if
>> >> >> you can't find the issue yourself. In particular, how did you
>> >> >> construct selDataMatrix and "design".
>> >> >>
>> >> >> Sean
>> >> >>
>> >> >>
>> >> >> On Thu, Jun 2, 2011 at 11:17 AM, Jasreet <jovialj86 at gmail.com>
>> wrote:
>> >> >> > Hi all
>> >> >> >
>> >> >> > I am new to R and have been using it for Illumina Infinium
>> Methylation
>> >> >> 450 k
>> >> >> > Array analysis using the lumi package.
>> >> >> >
>> >> >> > After normalization, I tried to identify differentially expressed
>> >> genes
>> >> >> > using the example illusrated in the lumi PDF but look like during
>> >> lmFit,
>> >> >> I
>> >> >> > keep getting the same error and hence, I cant move forward.
>> >> >> >
>> >> >> >
>> >> >> > *Error Message : *
>> >> >> >> fit <- lmFit(selDataMatrix, design)
>> >> >> > Error in lm.fit(design, t(M)) : incompatible dimension
>> >> >> > *
>> >> >> > *
>> >> >> > I have also copied sessioninfo and traceback.Any help would be
>> highly
>> >> >> > appreciated.
>> >> >> >
>> >> >> > Thanks
>> >> >> > -jess
>> >> >> > *
>> >> >> >
>> >> >> > Session Info :*
>> >> >> >> sessionInfo()
>> >> >> > R version 2.13.0 (2011-04-13)
>> >> >> > Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit)
>> >> >> >
>> >> >> > locale:
>> >> >> > [1] en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8
>> >> >> >
>> >> >> > attached base packages:
>> >> >> > [1] stats graphics grDevices utils datasets methods
>> base
>> >> >> >
>> >> >> > other attached packages:
>> >> >> > [1] limma_3.8.2
>> >> >> > MASS_7.3-13
>> >> >> IlluminaHumanMethylation450k.db_1.0.7
>> >> >> > org.Hs.eg.db_2.5.0 RSQLite_0.9-4
>> >> >> > [6] DBI_0.2-5
>> >> >> > AnnotationDbi_1.14.1 lumi_2.4.0
>> >> >> > nleqslv_1.8.5 Biobase_2.12.1
>> >> >> >
>> >> >> > loaded via a namespace (and not attached):
>> >> >> > [1] affy_1.30.0 affyio_1.20.0 annotate_1.30.0
>> >> >> > grid_2.13.0 hdrcde_2.15 KernSmooth_2.23-5
>> >> >> > lattice_0.19-26 Matrix_0.999375-50
>> >> >> > [9] methylumi_1.8.0 mgcv_1.7-6 nlme_3.1-101
>> >> >> > preprocessCore_1.14.0 tools_2.13.0 xtable_1.5-6
>> >> >> >
>> >> >> > [[alternative HTML version deleted]]
>> >> >> >
>> >> >> > _______________________________________________
>> >> >> > 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
>> >> >> >
>> >> >>
>> >> >
>> >> > [[alternative HTML version deleted]]
>> >> >
>> >> > _______________________________________________
>> >> > 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
>> >> >
>> >>
>> >
>> > [[alternative HTML version deleted]]
>> >
>> > _______________________________________________
>> > 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
>> >
>>
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> 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
>
More information about the Bioconductor
mailing list