FrederickHuangLin / ANCOMBC

Differential abundance (DA) and correlation analyses for microbial absolute abundance data
https://www.nature.com/articles/s41467-020-17041-7
107 stars 29 forks source link

Using rand_formula results in warning from lme4::lmer #166

Closed evetb closed 1 year ago

evetb commented 1 year ago

Hello there,

I keep having an issue where, when I include rand_formula when running ancombc2, I receive the following warning several times:

In lme4::lmer(formula = tformula, data = df, control = lme_control) : restarting interrupted promise evaluation

I have looked at other issues & discussion in ANCOMBC and I have noticed that past issues have been due to incorrect formatting of rand_formula. However, as seen with my code below, I believe I am formatting it correctly. As such, I do not understand why I continue to receive this error.

ancom_output <- ancombc2(phylo_obj, 
                                 assay_name = "counts",
                                 tax_level = "Phylum",
                                 fix_formula = "hs",
                                 rand_formula = "(1|cage)",
                                 p_adj_method = "bonferroni",
                                 pseudo = 0,
                                 pseudo_sens = TRUE,
                                 prv_cut = 0.10, lib_cut = 1000, s0_perc = 0.05,
                                 group = "hs",
                                 struc_zero = TRUE,
                                 alpha = 0.05, n_cl = 2, verbose = TRUE,
                                 global = TRUE,
                                 dunnet = TRUE,
                                 iter_control = list(tol = 1e-2, max_iter = 20, 
                                                     verbose = TRUE),
                                 em_control = list(tol = 1e-5, max_iter = 100),
                                 lme_control = lme4::lmerControl(),
                                 mdfdr_control = list(fwer_ctrl_method = "bonferroni", B = 100))

Please let me know if you need any other information to help my case. I am using R version 4.2.3 and running in RStudio version 2023.03.0.

Many thanks in advance!

FrederickHuangLin commented 1 year ago

Hi @evetb,

Thanks for your feedback!

To clarify, the message you encountered, "restarting interrupted promise evaluation" is a warning message, not an error, right?

To debug, I suggest fitting the mixed effects model to your RAW data using (lmerTest)[https://cran.r-project.org/web/packages/lmerTest/index.html] using the same fixed effects and random effects you used for the ancombc2 function and see if you would encounter the same error message. Usually, this will give you more informative error messages.

Lastly, could you try downloading and implementing the devel version of ANCOMBC (version 2.1.4) from GitHub? The warning message you received could result from the bug in the parallel computing process, which has been fixed in the devel version. The stable version of ANCOMBC 2.1.4 will be released on Bioconductor 3.17.

Best, Huang

evetb commented 1 year ago

Hello Huang,

Thank you for your reply! I tried using n_cl =1 instead of n_cl = 2 and this got rid of the warning (so it likely did have to do with the parallel processing). I also made sure that the values I was using for rand_formula were numeric. However, now I am encountering a different warning, which I will described below. Please let me know if I should create a different issue for this.

When running ANCOMBC2 with rand_formula (as below - note max_iter and B are set to = 1 for speed while testing):

ancom_output <- ancombc2(phylo_obj, 
                                 assay_name = "counts",
                                 tax_level = "Phylum",
                                 fix_formula = "hs",
                                 rand_formula = "(day | sample)",
                                 p_adj_method = "bonferroni",
                                 pseudo = 0,
                                 pseudo_sens = TRUE,
                                 prv_cut = 0.10, lib_cut = 1000, s0_perc = 0.05,
                                 group = "hs",
                                 struc_zero = TRUE,
                                 alpha = 0.05, n_cl = 1, verbose = TRUE,
                                 global = TRUE,
                                 dunnet = TRUE,
                                 iter_control = list(tol = 1e-2, max_iter = 1, 
                                                     verbose = TRUE),
                                 em_control = list(tol = 1e-5, max_iter = 1),
                                 lme_control = lme4::lmerControl(),
                                 mdfdr_control = list(fwer_ctrl_method = "bonferroni", B = 1)) 

I receive the following type of warnings:

Warning messages:
1: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  ... :
  Model failed to converge with max|grad| = 0.01015 (tol = 0.002, component 1)
...
21: Model failed to converge with 1 negative eigenvalue: -5.9e+01

I tried running lmer on the raw data as instructed - I'm not sure I did it correctly, but here is what I did:

phylo_obj_df <- psmelt(phylo_obj)

model <- lmer(abundance ~ hs + (day|sample),
                   data = phylo_obj_df)

It gave me similar warnings:

Warning messages:
  1: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  :
      unable to evaluate scaled gradient
  2: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  :
      Model failed to converge: degenerate  Hessian with 1 negative eigenvalues
  3: Model failed to converge with 1 negative eigenvalue: -1.9e+03

Looking into it a bit, it seems these issuesn could be due to a variety of issues, including singularity (https://rstudio-pubs-static.s3.amazonaws.com/33653_57fc7b8e5d484c909b615d8633c01d51.html)

FrederickHuangLin commented 1 year ago

Hi @evetb,

Thank you for your feedback. We appreciate your concern regarding the warning message you encountered. The warning message you saw is from the lmerTest package, which we have imported into the ANCOMBC package for fitting mixed effects models. This warning message can be triggered by various factors, as you have mentioned, and it may not necessarily indicate a problem with the analysis.

To address this issue, you can try fine-tuning the lme_control argument in the ancombc2 function. By adjusting the control parameters, you may be able to mitigate the warning message or reduce its frequency.

We would also like to inform you that we have submitted a revised version of the ANCOM-BC2 paper, and we expect the preprint of the paper to be available in approximately 2 weeks.

Best regards, Huang