Limma analysis - decomposition problem
1
0
Entering edit mode
7.5 years ago
iside ▴ 20

Dear all,

I am running a differential expression analysis with limma. My dataset includes unrelated individuals as well as twins, and I am aware that limma is not designed to deal with familiar data, but, as explained in the users' guide I may use the family as a blocking variable (in this way limma will estimate a within-family correlation from the expression data and will use it for the differential analysis).

Unfortunately, when running lmFit with the correlation factor I get this error:

Error in chol.default(V) :
the leading minor of order 2 is not positive definite


I understand that this depends on my data structure, and clearly I do not have a positive definite matrix, therefore it doesn't respect the assumptions of cholesky decomposition. But do you have any suggestion on how to overcome the problem?

This is the code I am using:

formula <- ~ condition
design <- model.matrix(formula, data=info)
corfit <- duplicateCorrelation(eset, design, block=info$famID) fit <- lmFit(eset, design, block=info$famID, correlation=corfit.corr\$consensus)


duplicateCorrelation lmFit Cholesky Limma • 2.3k views
0
Entering edit mode

I have the same problem with lmFit! Have you found any solution yet? I have log2 transformed data as well.

Thanks a lot!

1
Entering edit mode
7.5 years ago

Is your normalised eset in log2 space? It'd be useful to see the output of:

range(eset)

edit: It'd also be useful to see a pheno table so we can gauge the experimental design.

0
Entering edit mode

Hi Andrew,

The data had already been normalised, and corrected for batch effect with the package Combat. The output you asked for is the following:

range(expression)
[1]  6.109361 16.553681


My design matrix is composed by two columns: Intercept (136 rows set to 1) and conditioncase, where control is set to 0 (68 rows) and case is set to 1 (68 rows).

0
Entering edit mode

Thanks, how many levels are there to your condition variable? i.e. How many conditions do you have?

0
Entering edit mode

I have only 2 conditions. Either case or control.

0
Entering edit mode

OK, that explains your intercept model. What's the output of corfit?

0
Entering edit mode

I get only a matrix of NAs!

0
Entering edit mode

Ah, no errors when running duplicateCorrelation?