Error in last step of Limma analysis (eBayes)
1
0
Entering edit mode
peirinl • 0
@peirinl-14078
Last seen 15 months ago

A running limma analysis of a dataset, and I got an error in the final step, eBayes(fit2)

This is my code

i <- pData(eset)$Sample_title
Part <- ifelse(pData(eset)$Part == "Tn", "nT", "T")
Grade <- ifelse(pData(eset)$Grade == "Low", "nH", "H")
sample_byGrade <- factor(paste(Grade, Part, sep="."))

design_Grade <- model.matrix(~ 0  + sample_byGrade)
colnames(design_Grade) <- levels(sample_byGrade)

corfit <- duplicateCorrelation(eset,design_Grade,block = i)

contrast_matrix_Grade<- makeContrasts(HighVsLowinT = H.T-nH.T, HighVsLowinnT = H.nT-H.nT, 
                                      TVsnTinHigh = H.T-H.nT, TVsnTinLow = nH.T-nH.nT,
                                      levels = design_Grade)

fit <- lmFit(eset,design = design_Grade,block = i, correlation = corfit$consensus)
fit2<- contrasts.fit(fit,contrast_matrix_Grade)
fit2<- eBayes(fit2)

The error message:

   Error in eigen(cor.matrix, symmetric = TRUE) : 
  infinite or missing values in 'x'
Besides: Warning messages:
1: In ebayes(fit = fit, proportion = proportion, stdev.coef.lim = stdev.coef.lim,  :
  Estimation of var.prior failed - set to default value
2: In cov2cor(object$cov.coefficients) :
  diag(.) had 0 or NA entries; non-finite result is doubtful

I see previous reply and some said that it is possible due to NA in the data. I do anyNA(exprs(eset)) and it is FALSE, and I did do filtering before running LIMMMA so I believe my data got no NA values.

But some of my column in pData(eset) did have NA value, but these columns are not used during this analysis. Where is my problem?

Thanks everyone.

ebayes limma • 1.0k views
ADD COMMENT
0
Entering edit mode

Please show us pData(eset).

Is it possible that Sample_title takes a different value for each sample? In that case, setting block=i wouldn't make any sense and would lead to an error such as you have.

Have you checked that corfit$consensus is between -1 and 1 (not equal to -1 or 1)?

ADD REPLY
0
Entering edit mode

The pData is here. Not sure if this work out. 

https://www.dropbox.com/s/35hqpm11umynx19/eset.txt?dl=0

Sample_title, i, is the subject name

i
 [1] "E005" "E005" "E012" "E012" "E021" "E021" "E028" "E028" "E033" "E033" "E041" "E041" "E042" "E042"
[15] "E045" "E045" "E048" "E048" "E052" "E052" "E058" "E058" "E064" "E064" "E073" "E073" "E082" "E082"
[29] "E086" "E086" "E088" "E088" "E094" "E094" "E109" "E109" "E119" "E119" "E120" "E120" "E123" "E123"
[43] "E130" "E130" "E133" "E133" "E136" "E136" "E138" "E138" "E139" "E139" "E152" "E152" "E162" "E162"
[57] "E165" "E165" "E178" "E178" "E179" "E179" "E181" "E181" "E182" "E182" "E186" "E186" "E203" "E203"
[71] "E208" "E208" "E209" "E209" "E218" "E218" "E220" "E220" "E229" "E229" "E234" "E234" "E245" "E245"
[85] "E246" "E246" "E249" "E249" "E250" "E250"

> corfit$consensus
[1] 0.05972559

 

ADD REPLY
0
Entering edit mode

The pData file on drop box doesn't agree with the output you give for Sample_title, i.

In the pData file, the Sample_title values are "E005-Tc8", "E005-Tn", "E012-Tn" etc. In the pData file, every sample has a different title.

ADD REPLY
0
Entering edit mode

Ok.

I modified i <- pData(eset)$Sample_title into this: 

i <- sub("-.*", "", pData(eset)$Sample_title), then the data become as what i show up there.

ADD REPLY
0
Entering edit mode

Sorry for wasting your time, but I make a stupid mistake. 

HighVsLowinnT = H.nT-H.nT     should be     HighVsLowinnT = H.nT-nH.nT.

Thank you for helping me. 

ADD REPLY
0
Entering edit mode
peirinl • 0
@peirinl-14078
Last seen 15 months ago

Sorry, I think I make a huge STUPID mistake. 

contrast_matrix_Grade<- makeContrasts(HighVsLowinT = H.T-nH.T, HighVsLowinnT = H.nT-H.nT, 
                                      TVsnTinHigh = H.T-H.nT, TVsnTinLow = nH.T-nH.nT,
                                      levels = design_Grade)

 

should be 

HighVsLowinnT = H.nT-nH.nT

 

I apologize for wasting everyone's time.

 

ADD COMMENT

Login before adding your answer.

Traffic: 329 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6