Question

Methylation Analysis experiment design using Limma

0

Entering edit mode

hortense96 • 0

@f94ed731

Last seen 6 months ago

France

Hello everyone,

I'm trying to analyze DNA methylation data with the Limma package to identify potential differentially methylated CpGs between my conditions using M-values. I would like to take into account the proportion of different blood cell fractions in my samples, age, sex and possible batch effects as covariates in my linear analysis model.

This is the code I'm trying to use with the following error message when using lmFit

head(design) fDisease fControl Age SexM CD8T CD4T NK Bcell Mono Neu 1 1 0 53 0 0.0162 0.1618 0.0484 0.0435 0.0853 0.6471 2 0 1 53 0 0.0305 0.1034 0.0430 0.0409 0.1516 0.6391 3 1 0 37 1 0.0814 0.1255 0.0761 0.0718 0.0970 0.5559 4 0 1 36 1 0.0114 0.1454 0.0457 0.0499 0.1328 0.6277 5 1 0 37 1 0.0912 0.1006 0.0620 0.0948 0.0636 0.5912 6 0 1 36 1 0.0703 0.1295 0.0474 0.0359 0.1037 0.6249

fit<-lmFit(mvalues,design)

Coefficients not estimable: NK Bcell Mono Neu Warning message: Partial NA coefficients for 869576 probe(s) ```

I'm not sure to understand why some of my cell fractions are problematic and not others in my model?
I'm also wondering how to integrate into this model any technical batch effects linked in particular to the use of different arrays?

Thank you in advance for any suggestions and help.

Hortense

DNAMethylation limma methylationArrayAnalysis • 1.3k views

ADD COMMENT • link 8 months ago • updated 7 months ago hortense96 • 0

0

Entering edit mode

I think adding , myLoad$pd in model.matrix() is misleading/incorrect. Remove it and try again.

Try a simple model such as ~0 +f+NK in order to test if the same error oocurs.

ADD REPLY • link 8 months ago SamGG ▴ 360

0

Entering edit mode

Hello Sam,

Thank you for your feedback ! So I have tried your suggestions by removing myLoad$pd and the cell populations that cause problems : NK, Bcell Mono and the Neu.

While this allows me to use the lmFit function without an error message, it causes problems in the rest of the code when I try to apply eBayes :

> fit3<-eBayes(fit2)
Error in .ebayes(fit = fit, proportion = proportion, stdev.coef.lim = stdev.coef.lim,  : 
  No finite residual standard deviations

The strange thing is that I had no problems when I used these lines of code for another data set of the same nature...

ADD REPLY • link 8 months ago hortense96 • 0

0

Entering edit mode

I said that , myLoad$pd was incorrect, because the variables you transformed are not reintroduced in myLoad$pd. It should be myLoad$pd$Neu<- as.numeric(myLoad$pd$Neu). The design shows that those numerical variables are not considered as factors, which is OK.

What are dim(design)? dim(mvalues)?

ADD REPLY • link 8 months ago SamGG ▴ 360

score 0 · Answer 1 · 2024-04-17

0

Entering edit mode

Gordon Smyth 52k

@gordon-smyth

Last seen 1 hour ago

WEHI, Melbourne, Australia

You don't explain how many samples you have (columns of data), but the warning message shown in your post would arise if there are in fact only 6 samples in the analysis. In other words, the matrix of M-values mvalues has only six columns.

If there are only six samples in the analysis, then it is only possible to estimate six coefficients. Hence limma removes columns 7-10 from the design matrix, which is what the warning message is telling you.

If you want to analyse a design matrix with so many covariates then you need much more data.

ADD COMMENT • link 8 months ago Gordon Smyth 52k

0

Entering edit mode

Hi,

Thank you very much for your answer. It's true that I only used 6 samples for this analysis, which is consistent with your explanation. When I include more samples, the problem no longer occurs.

I'd also like to take this opportunity to ask you, when using the eBayes function in this type of methylation analysis (with a sufficient number of samples), how should the TREND argument of this function be considered? I'm not sure how to explore my mean-variance trend with this type of data and whether it's even relevant to consider? In any case, when I try the eBayes function with trend = true or not in a completely empirical way, I don't get the same number of differentially methylated CpGs (900 versus 600). I have to admit that I don't understand the subtleties of what this difference in CpG numbers means in this case ?

ADD REPLY • link 8 months ago hortense96 • 0

0

Entering edit mode

I don't analyse methylation arrays so cannot advise on specific limma settings for that type of data. Please follow the advice of the methylation packages that produce M-values, for example the missMethyl package.

ADD REPLY • link 8 months ago Gordon Smyth 52k

0

Entering edit mode

Thank you very much for your reply and your suggestion, as I'm having a bit of trouble finding information on how to manage this part of the analysis. Hortense

ADD REPLY • link 7 months ago hortense96 • 0