DMRcate cpg.annotate error
1
1
Entering edit mode
@giovanni-calice-6415
Last seen 5 weeks ago
Italy

Hi all,

I am using DMRcate to find differentially methylated regions in illumina 450k data between one group of 4 samples and another group of 9 samples.

I followed all the steps of DMRcate user's guide and I build the design matrix in a similar manner,

and this is it:

str(design)
 num [1:13, 1:14] 1 1 1 1 1 1 1 1 1 1 ...
 - attr(*, "dimnames")=List of 2
  ..$ : chr [1:13] "1" "2" "3" "4" ...
  ..$ : chr [1:14] "(Intercept)" "patient020" "patient021" "patient023" ...
 - attr(*, "assign")= int [1:14] 0 1 1 1 1 1 1 1 1 1 ...
 - attr(*, "contrasts")=List of 2
  ..$ patient: chr "contr.treatment"
  ..$ type   : chr "contr.treatment"

 

But when I call the cpg.annotate(), I have this error:

 

myannotation <- cpg.annotate(myMvalue, analysis.type = "differential", design = design, coef = 14)
Coefficients not estimable: typeGroup2
Error in ebayes(fit = fit, proportion = proportion, stdev.coef.lim = stdev.coef.lim,  :
  No residual degrees of freedom in linear model fits
Inoltre: Warning message:
Partial NA coefficients for 425827 probe(s) 

 

In the DMRcate example dataset, that I have tested, the differential analysis is paired 38 vs 38

and I think that the error probably is in the design matrix building of my dataset.

So how to set the design matrix correctly?

 

Thanks in advance, Regards

 

Giovanni

 
 
 
dmrcate design matrix • 2.4k views
ADD COMMENT
1
Entering edit mode
Tim Peters ▴ 170
@tim-peters-7579
Last seen 1 day ago
Australia

Hi Giovanni,

 

You're right, there is definitely a problem with the design matrix, especially one that has 14 columns when you only have 13 samples. This means that limma is interpreting each sample as having its own phenotype, hence none have replicates, and therefore you don't have any residual degrees of freedom in the model, just like the error says.

 

What you need to do is create a vector describing the groups you want to contrast, corresponding to the matching columns in your M-matrix. Say the first 4 columns are group 1 and the remaining 9 columns are Group 2.

So to calculate DMRs that are Group 2 - Group 1, input like so, and make sure the coef argument corresponds to the column of the design matrix which contains the contrasts of interest, i.e. 2. :

> type <- c(rep("1", 4), rep("2", 9))
> type
[1] "1" "1" "1" "1" "2" "2" "2" "2" "2" "2" "2" "2" "2"
> design <- model.matrix(~type)
> design
   (Intercept) type2
1            1     0
2            1     0
3            1     0
4            1     0
5            1     1
6            1     1
7            1     1
8            1     1
9            1     1
10           1     1
11           1     1
12           1     1
13           1     1
attr(,"assign")
[1] 0 1
attr(,"contrasts")
attr(,"contrasts")$type
[1] "contr.treatment"

myannotation <- cpg.annotate(myMvalue, analysis.type = "differential", design = design, coef = 2)

Good luck,

Tim

ADD COMMENT
0
Entering edit mode

Hi Tim,

 

Many Thanks, It works fine!

In particular I had already tried to build the matrix by

myPhenoData$Sample_Group

used to calculate DMPs, but I wrong in the coef argument in cpg.annotate call

(coef = 14)

having again error. Maybe I misinterpreted the cpg.annotate call in the paired analysis of example dataset.

Then using the coef argument

(coef = 2)

that reflects my originary

myPhenoData$Sample_Group

I got DMRs and this is the design matrix:

   (Intercept)      myPhenoData$SampleGroup
1            1                            0
2            1                            0
3            1                            0
4            1                            0
5            1                            0
6            1                            0
7            1                            0
8            1                            0
9            1                            0
10           1                            1
11           1                            1
12           1                            1
13           1                            1
attr(,"assign")
[1] 0 1
attr(,"contrasts")
attr(,"contrasts")$`myPhenoData$SampleGroup`
[1] "contr.treatment"

 

Good luck,

Giovanni

ADD REPLY

Login before adding your answer.

Traffic: 871 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