methylKit: Error in lm.fit when including covariates

[yaaminiv]

Posted in the methylKit Github, but adding here too in case anyone has some insight!

Any suggestions for successfully running calculateDiffMeth with a covariate matrix without an lm.fit error?

I'm running methylKit with the Manchester samples in this R Markdown script to identify DML. I'm also using a covariate matrix to account for differences in maturation stage, and correcting for over dispersion. When I run this code:

differentialMethylationStatsTreatment <- methylKit::calculateDiffMeth(methylationInformationFilteredCov5, covariates = covariateMetadata, overdispersion = "MN", test = "Chisq", mc.cores = 2). #test = "Chisq" is default

I get this error about model fit:

two groups detected:
 will calculate methylation difference as the difference of
treatment (group: 1) - control (group: 0)
Error in lm.fit(x, y, offset = offset, singular.ok = singular.ok, ...) : 
  0 (non-NA) cases

I dug into the methylKit Github repo, and found a seemingly related issue. The OP may have had issues with too many NAs after unite since they used a min.per.group argument, but I didn't do that! Interestingly, when I run the code without the covariate matrix:

differentialMethylationStatsTreatment <- methylKit::calculateDiffMeth(methylationInformationFilteredCov5, mc.cores = 4)

I don't get an error because R crashes, so I can't figure out if there's an issue with my covariate matrix, or if I don't have a large enough sample size to try these additional model parameters.

Possible solutions:

1. Running the code above with mc.cores = 1
2. Running calculateDiffMeth with overdispersion = "MN" but no covariate matrix
3. Giving up and trying DSS

[kubu4]

The OP may have had issues with too many NAs after unite

Does your data end up with NAs? If yes, any way to filter those out without impacting analysis?

[yaaminiv]

There are no NAs after filtering, normalizing, and uniting! So I don't think there's a problem with the input. I think it's a modeling issue of some sort

[yaaminiv]

Possible solutions: Running the code above with mc.cores = 1

Ran differentialMethylationStatsTreatment <- methylKit::calculateDiffMeth(methylationInformationFilteredCov5) and got an output dataframe! I wonder if there's an issue with specifying mc.cores (default is 1).

Now trying differentialMethylationStatsTreatment <- methylKit::calculateDiffMeth(methylationInformationFilteredCov5, overdispersion = "MN", test = "Chisq")

[yaaminiv]

I have a theory that mc.cores is the issue.

I tried differentialMethylationStatsTreatment <- methylKit::calculateDiffMeth(methylationInformationFilteredCov5, overdispersion = "MN", test = "Chisq", mc.cores = 2) and got the same lm.fit error.

Ran differentialMethylationStatsTreatment <- methylKit::calculateDiffMeth(methylationInformationFilteredCov5, overdispersion = "MN", test = "Chisq") and got output!

Going to try code below to see if I can get output without an lm.fit error:

differentialMethylationStatsTreatment <- methylKit::calculateDiffMeth(methylationInformationFilteredCov5, covariates = covariateMetadata, overdispersion = "MN", test = "Chisq")