Search
Question: LIMMA - problem with block analysis
0
10.4 years ago by
Milena Gongora80 wrote:
Dear LIMMA users, Hello. I am analysing single colour microarrays in LIMMA and am facing some obstacles, which I am not sure have a solution at all. So wondering if someone can confirm this or illuminate me on why these errors. I have a really bad experimental design. It is a simple pair-wise comparison between mutant and a control conditions. I have 7 samples, 3 control samples (comming from independent individuals) and 4 mutant samples where 3 come from the same individual and 1 from its mother. So in essence, I have a high level of correlation between the mutant samples, but independence in the control samples. What I thought I could do to get the most realistic information out of this data, was to do a "blocked" analysis where I would look at the differences between mutant and control, but blocking the samples by patient first. I tried doing this by following the example in the LIMMA guide for section 8.3 "paired samples" as below: > samples Sample Mutation Patient Related 1 Mu_1A mu 1 R1 2 Mu_1B mu 1 R1 3 Mu_1C mu 1 R1 4 C_3 ctrl 3 R2 5 C_4 ctrl 4 R3 6 C_5 ctrl 5 R4 7 Mu_2 mu 2 R1 # extract the variables of interest as factors > mutation <- factor(samples$Mutation) > patient <- factor(samples$Patient) > design_c1 <- model.matrix(~patient+mutation) > fit <- lmFit(data_Qnorm_log2, design_c1) Coefficients not estimable: mutationmu My question here is: Why can the coefficient not be estimated? So I tried a second approach. Instead of putting my blocking vairable in the design matrix. I tried fitting the model specifying the correlation using duplicateCorrelation() like in section 8.5 of LIMMA guide "Technical Replication". In this approach, the problem is that my correlation is negative and can't fit the model. See below: > block_p <- as.vector(patient) > corr_patient <- duplicateCorrelation(exprs(data_Qnorm_log2),design_1, block=block_p) > corr_patient$consensus.correlation [1] -1 > fit_corr_pat <- lmFit(data_Qnorm_log2, design_1, block=block_p, correlation=corr_patient$consensus) Error in chol.default(V) : the leading minor of order 2 is not positive definite Why is it producing a negative correlation? if I do cor() on the expression values, I get positive correlations as seen here > cor(exprs(data_Qnorm_log2)) Mu_1A Mu_1B Mu_1C C_3 C_4 C_5 Mu_2 Mu_1A 1.0000000 0.9814710 0.9768791 0.9753479 0.9783894 0.9744542 0.9710319 Mu_1B 0.9814710 1.0000000 0.9785653 0.9765288 0.9783644 0.9755863 0.9706648 Mu_1C 0.9768791 0.9785653 1.0000000 0.9727115 0.9750562 0.9756139 0.9673867 C_3 0.9753479 0.9765288 0.9727115 1.0000000 0.9762357 0.9765884 0.9673692 C_4 0.9783894 0.9783644 0.9750562 0.9762357 1.0000000 0.9761812 0.9707479 C_5 0.9744542 0.9755863 0.9756139 0.9765884 0.9761812 1.0000000 0.9710843 Mu_2 0.9710319 0.9706648 0.9673867 0.9673692 0.9707479 0.9710843 1.0000000 > In summary is it possible to run this analysis or is there just too many constraints in this data? If possible, how could I do so? I appreciate every insight and thanks very much for reading this long post! Milena