statdivlab / corncob

Count Regression for Correlated Observations with the Beta-binomial
102 stars 22 forks source link

Model not converging #139

Closed Auerilas closed 1 year ago

Auerilas commented 2 years ago

Hi, I'm trying to run difftest and I can't get my model to converge (even though all the examples work):

library(phyloseq)
library(corncob)
library(ggplot2)

# ASV data
seq_data = readRDS('ps_object.rds')
# filter out dates after the experiment
seq_data = seq_data %>% phyloseq::subset_samples(Date %in% c("2020.09.01", "2020.06.11", "2020.08.06", "2020.08.19", "2020.07.08", "2020.07.07", "2020.06.25", "2020.07.22", "2020.05.11", "2020.05.29"))
# collapse to phylum data
phy_data = seq_data %>% tax_glom('phylum')

test = bbdml(formula=ASV1~Date,                                                                                                                                            
+              phi.formula=~Date,                                                                                                                                            
+              data=phy_data,                                                                                                                                                
+              taxa_are_rows=T)                                                                                                                                              
> test                                                                                                                                                                       

Call:                                                                                                                                                                        
bbdml(formula = ASV1 ~ Date, phi.formula = ~Date, data = phy_data,                                                                                                           
    ... = pairlist(taxa_are_rows = T))
Coefficients associated with abundance:
               Estimate Std. Error t value Pr(>|t|)    
(Intercept)    -3.48050    0.06874 -50.635  < 2e-16 ***
Date2020.05.29  0.14888    0.11188   1.331    0.187    
Date2020.06.11  0.49810    0.09036   5.512  4.2e-07 ***
Date2020.06.25 -0.05654    0.08979  -0.630    0.531    
Date2020.07.07  0.04556    0.09530   0.478    0.634    
Date2020.07.08  0.07992    0.09455   0.845    0.401    
Date2020.07.22  0.09872    0.09886   0.999    0.321    
Date2020.08.06  0.00730    0.09463   0.077    0.939    
Date2020.08.19  0.08539    0.08293   1.030    0.306    
Date2020.09.01  0.15254    0.09270   1.646    0.104    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Coefficients associated with dispersion:
               Estimate Std. Error t value Pr(>|t|)    
(Intercept)    -6.62462    0.46342 -14.295   <2e-16 ***
Date2020.05.29  0.65342    0.65156   1.003    0.319    
Date2020.06.11  0.11408    0.66311   0.172    0.864    
Date2020.06.25 -0.40465    0.65670  -0.616    0.540    
Date2020.07.07 -0.03213    0.65362  -0.049    0.961    
Date2020.07.08 -0.03100    0.65301  -0.047    0.962    
Date2020.07.22  0.16200    0.65495   0.247    0.805    
Date2020.08.06 -0.09950    0.65437  -0.152    0.880    
Date2020.08.19 -0.73318    0.66298  -1.106    0.272    
Date2020.09.01 -0.05732    0.65613  -0.087    0.931    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Log-likelihood: -648.38

> da_analysis = differentialTest(formula = ~ Date,
+                                phi.formula = ~ Date,
+                                formula_null = ~ 1,
+                                phi.formula.null = ~ Date,
+                                test='Wald',
+                                boot=F,
+                                data=common_phyla)
Error in differentialTest(formula = ~Date, phi.formula = ~Date, formula_null = ~1,  : 
  All models failed to converge! 

           If you are seeing this, it is likely that your model is overspecified. This occurs when your sample size is not large enough to estimate all the parameters of your model. This is most commonly due to categorical variables that include many categories. 

           Alternatively, double-check your values for the arguments `link`, `phi.link`, and `method` to makes sure that they follow the specified options. 

           To confirm you have fixed the issue, try running a model for a single taxon with bbdml.

It doesn't matter how I specify the model, they all fail:

> da_analysis = differentialTest(formula = ~ 1,
+                                phi.formula = ~ 1,
+                                formula_null = ~ 1,
+                                phi.formula.null = ~ 1,
+                                test='Wald',
+                                boot=F,
+                                data=common_phyla)
Error in differentialTest(formula = ~1, phi.formula = ~1, formula_null = ~1,  : 
  All models failed to converge! 

           If you are seeing this, it is likely that your model is overspecified. This occurs when your sample size is not large enough to estimate all the parameters of your model. This is most commonly due to categorical variables that include many categories. 

           Alternatively, double-check your values for the arguments `link`, `phi.link`, and `method` to makes sure that they follow the specified options. 

           To confirm you have fixed the issue, try running a model for a single taxon with bbdml.

Thoughts?

bryandmartin commented 1 year ago

Hi @Auerilas ,

I believe this should be fixed in 365add6. Let me know if you are still encountering issues!