DMRcate Error in if (nsig == 0) { : missing value where TRUE/FALSE needed
1
0
Entering edit mode
Jeni • 0
@jeni-23673
Last seen 7 days ago
Spain

Hi!

I am using DMRcate to find diferencially methylated regions between two groups of samples, from the Beta Values of a set of samples data that I downloaded from TCGA repository.

I have performed the probe wise analysis first, using limma, so I have generated the design and contrast matrix with this package. Now, I am using this command with DMRcate:

myAnnotation <- cpg.annotate(object = as.matrix(bVals_filtered), datatype = "array", arraytype = c("450K"), what = "Beta",
analysis.type = "differential", design = design,
contrasts = TRUE, cont.matrix = contMatrix,
coef = "B - A", fdr = 0.01)


After that, I get this error:

Error in if (nsig == 0) { : missing value where TRUE/FALSE needed


Is any way to handle NA in the Beta values matrix? If I filter rows containing some NA then I keep just a few probes, compared with the total amount of probes in the array.

DMRcate differencial methylation • 884 views
0
Entering edit mode

Why do you have so many NA values in your Beta matrix? I've never actually (so far as I recall) even seen one NA value, not to mention enough that you filter out most of your probes.

0
Entering edit mode

I don't know that, I have downloaded this data from the TCGA. But, if I remove rows with some NA, then I get just 47495 probes, which is a quite low amount.

3
Entering edit mode
@james-w-macdonald-5106
Last seen 4 days ago
United States

The test for nsig == 0 comes from a test for the adjusted p-values being < your significance cutoff. The only way you will have NA values for your adjusted p-values is if you have NA values for all the subjects in a single group. So you don't have to get rid of all the probes with any NA value, just those that have all NA values for one or more of the groups you are comparing.

You should have seen a warning that said

Warning message:
Partial NA coefficients for N probe(s)



When you did the limma analysis. And you could find out which probes those are by doing

naind <- is.na(rowSums(fit$p.value))  where 'fit' is the name of the MArrayLM object you got after calling eBayes. ADD COMMENT 1 Entering edit mode Hi Jeni, Yep Jim's right - run your limma fit manually with his code to see which probes you need to filter out. Alternatively, you could impute the NAs using knn.impute or some such method, taking into account the risks inherent to imputation. Cheers, Tim ADD REPLY 0 Entering edit mode Dear James, I have the same problem that the OP, and I have followed your steps to identify the probes with NA values after fitting, but the output is 0. The error message says: Error in if (nsig == 0) { : missing value where TRUE/FALSE needed In addition: Warning messages: 1: In logit2(assay(object, "Beta")) : NaNs produced 2: In logit2(assay(object, "Beta")) : NaNs produced 3: Partial NA coefficients for 72970 probe(s)  This is what I have done: design2 <- model.matrix(~0+location + sex + age) cm2 <- makeContrasts( ProximalvsDistal = proximal - distal, ProximalvsRectum = proximal - rectum, DistalvsRectum = distal - rectum, levels=design2) fit2 <- lmFit(MVals, design2) fit2 <- contrasts.fit(fit2, cm2) efit2 <- eBayes(fit2, robust = TRUE) naind <- is.na(rowSums(efit2$p.value))
naind[which(naind=='TRUE')]
#named logical(0)


Also I have checked that my M-values matrix has no NAs

apply(is.na(MVals), 2, which)
#integer(0)


Any ideas why I am getting that error?

0
Entering edit mode

Please don't add questions to old posts. Just make a new one.

Anyway, you don't show how you got your MVals object, but this part:

1: In logit2(assay(object, "Beta")) : NaNs produced
2: In logit2(assay(object, "Beta")) : NaNs produced


comes from getM, and indicates that the GenomicRatioSet that you used for cpg.annotate is problematic because it has beta values outside of (0,1)

> cbind(seq(-0.1, 1.1, 0.1), logit2(seq(-0.1, 1.1, 0.1)))
[,1]          [,2]
[1,] -0.1           NaN
[2,]  0.0          -Inf
[3,]  0.1 -3.169925e+00
[4,]  0.2 -2.000000e+00
[5,]  0.3 -1.222392e+00
[6,]  0.4 -5.849625e-01
[7,]  0.5  5.551115e-16
[8,]  0.6  5.849625e-01
[9,]  0.7  1.222392e+00
[10,]  0.8  2.000000e+00
[11,]  0.9  3.169925e+00
[12,]  1.0           Inf
[13,]  1.1           NaN
Warning messages:
1: In logit2(seq(-0.1, 1.1, 0.1)) : NaNs produced
2: In logit2(seq(-0.1, 1.1, 0.1)) : NaNs produced


And if you have a group that is all NA or NaN, you get that error from lmFit

> z <- matrix(rnorm(1000), 100)
> z[3,1:5] <- NaN
> zz <- lmFit(z, model.matrix(~gl(2,5)))
Warning message:
Partial NA coefficients for 1 probe(s)

0
Entering edit mode

Thanks a lot James, I realised the issue about posting in a comment and I asked a new question but forgot to delete the comment. My MVals object is a matrix of M-values not Beta values, so yes, it goes outside of (0,1). I just thought that I could use a matrix of M-values, according to cpg.annotate help

object.          Either: - A matrix of M-values, with unique Illumina probe IDs as rownames and unique sample IDs as column names or, - Output from DSS:::DMLtest().Either: - A matrix of M-values, with unique Illumina probe IDs as rownames and unique sample IDs as column names or, - Output from DSS:::DMLtest().


The command that triggers the error is:

myannotation_proxvsdis <- cpg.annotate(datatype="array",
object=MVals,
arraytype = "EPIC",
annotation=c(array = "IlluminaHumanMethylationEPIC", annotation = "ilm10b4.hg19"),
analysis.type="differential",
design=design2,
contrasts=TRUE,
cont.matrix=cm2,
coef="ProximalvsDistal")


Should I use a matrix of Beta values instead of M-values or matrices cannot be used as input for cpg.annotate?

0
Entering edit mode

I think you missed this part of the help page:

Arguments:

datatype: Character string representing the type of data being
analysed.

object: Either:

- A matrix of _M_-values, with unique Illumina probe IDs as
rownames and unique sample IDs as column names or,

- A GenomicRatioSet, appropriately annotated.

what: Does the data matrix contain Beta or M-values? Not needed if
object is a GenomicRatioSet.


Specifically the 'what' argument, which has 'Beta' as the default.

0
Entering edit mode

And you are right James, I feel very silly right now... Thanks a lot!! and sorry for bothering