Error de cálculo al ejecutar el modelo de regresión

Question

Este problema me ha retenido durante tres días, así que espero de verdad que alguien aquí tenga una solución para el problema.

Tengo un modelo con un número excesivo de ceros, por lo que utilizo un modelo de regresión de Poisson cero-inflado con el siguiente código y resumen.

cr_f1 = formula(cr ~ depth + habtype2 + month + year + lightregime + depth*month + depth*lightregime + depth*habtype2 + habtype2*year | depth + habtype2 + month + year + lightregime + depth*month + depth*lightregime + depth*habtype2)
summary(zeroinfl(cr_f1, dist = "poisson", link = "logit", data = allUVCdata))

Call:
zeroinfl(formula = cr_f1, data = allUVCdata, dist = "poisson", link = "logit")

Pearson residuals:
    Min      1Q  Median      3Q     Max 
-1.6430 -0.5680 -0.2893  0.1426 16.8090 

Count model coefficients (poisson with log link):
                            Estimate Std. Error z value Pr(>|z|)   
(Intercept)                -4.515522   2.182503  -2.069  0.03855 * 
depth                       0.108941   0.072278   1.507  0.13175   
habtype2Pinnacles           0.879765   0.791166   1.112  0.26614   
habtype2Unexposed          -0.604246   0.786129  -0.769  0.44211   
month2                      0.628468   0.380450   1.652  0.09855 . 
month3                      0.309282   0.367690   0.841  0.40026   
month4                      0.649411   0.371667   1.747  0.08059 . 
month5                      0.758717   0.364079   2.084  0.03717 * 
month6                      0.467611   0.341024   1.371  0.17031   
month7                      0.523043   0.343363   1.523  0.12768   
month8                      0.563272   0.356843   1.578  0.11445   
month9                      0.204509   0.400398   0.511  0.60952   
month10                     0.662415   0.341616   1.939  0.05249 . 
month11                     0.934844   0.335077   2.790  0.00527 **
month12                     0.252216   0.360512   0.700  0.48417   
year2013                   -1.271010   1.282158  -0.991  0.32154   
year2014                    1.221887   0.753644   1.621  0.10495   
year2015                   -0.463176   0.771131  -0.601  0.54808   
lightregimeLight            2.754925   1.948779   1.414  0.15746   
depth:month2               -0.019864   0.008906  -2.230  0.02572 * 
depth:month3               -0.014157   0.008106  -1.747  0.08071 . 
depth:month4               -0.020553   0.008332  -2.467  0.01364 * 
depth:month5               -0.021213   0.008373  -2.533  0.01129 * 
depth:month6               -0.013561   0.007393  -1.834  0.06663 . 
depth:month7               -0.015043   0.007544  -1.994  0.04615 * 
depth:month8               -0.017383   0.008011  -2.170  0.03003 * 
depth:month9               -0.012340   0.008990  -1.373  0.16988   
depth:month10              -0.019631   0.007629  -2.573  0.01008 * 
depth:month11              -0.024101   0.007611  -3.167  0.00154 **
depth:month12              -0.014319   0.007952  -1.801  0.07174 . 
depth:lightregimeLight     -0.079860   0.071024  -1.124  0.26084   
depth:habtype2Pinnacles    -0.006819   0.011178  -0.610  0.54182   
depth:habtype2Unexposed     0.014857   0.011103   1.338  0.18086   
habtype2Pinnacles:year2013  1.351509   1.277930   1.058  0.29025   
habtype2Unexposed:year2013  1.538282   1.256047   1.225  0.22069   
habtype2Pinnacles:year2014 -1.213233   0.754305  -1.608  0.10775   
habtype2Unexposed:year2014 -0.495275   0.726863  -0.681  0.49563   
habtype2Pinnacles:year2015  0.389117   0.775476   0.502  0.61582   
habtype2Unexposed:year2015  0.659117   0.750396   0.878  0.37975   

Zero-inflation model coefficients (binomial with logit link):
                        Estimate Std. Error z value Pr(>|z|)    
(Intercept)             -4.61555    7.04621  -0.655 0.512442    
depth                    0.28728    0.28211   1.018 0.308524    
habtype2Pinnacles        9.41037    3.82210   2.462 0.013813 *  
habtype2Unexposed        2.11213    1.46465   1.442 0.149282    
month2                   8.67847    3.91193   2.218 0.026523 *  
month3                   7.12210    3.86428   1.843 0.065320 .  
month4                   4.10296    2.41285   1.700 0.089044 .  
month5                  12.76919    4.28035   2.983 0.002852 ** 
month6                   3.57695    2.49820   1.432 0.152198    
month7                   5.85534    3.27394   1.788 0.073700 .  
month8                   5.59503    3.33054   1.680 0.092974 .  
month9                   4.22953    3.76919   1.122 0.261807    
month10                  6.35022    3.59424   1.767 0.077265 .  
month11                  5.92079    3.36405   1.760 0.078404 .  
month12                  4.36214    3.17233   1.375 0.169113    
year2013                -0.18722    0.42651  -0.439 0.660688    
year2014                -1.50194    0.45263  -3.318 0.000906 ***
year2015                -9.79773    4.87536  -2.010 0.044469 *  
lightregimeLight         0.79826    5.62419   0.142 0.887133    
depth:month2            -0.39212    0.16795  -2.335 0.019557 *  
depth:month3            -0.36363    0.16695  -2.178 0.029397 *  
depth:month4            -0.21521    0.10211  -2.108 0.035059 *  
depth:month5            -0.57543    0.16933  -3.398 0.000678 ***
depth:month6            -0.24336    0.10398  -2.341 0.019256 *  
depth:month7            -0.33704    0.13975  -2.412 0.015874 *  
depth:month8            -0.35343    0.14683  -2.407 0.016082 *  
depth:month9            -0.31787    0.16903  -1.881 0.060026 .  
depth:month10           -0.37550    0.16021  -2.344 0.019087 *  
depth:month11           -0.34650    0.14821  -2.338 0.019397 *  
depth:month12           -0.29639    0.14221  -2.084 0.037142 *  
depth:lightregimeLight   0.08117    0.21795   0.372 0.709571    
depth:habtype2Pinnacles -0.57765    0.17049  -3.388 0.000704 ***
depth:habtype2Unexposed -0.17897    0.06252  -2.863 0.004200 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Number of iterations in BFGS optimization: 146 
Log-likelihood: -3977 on 72 Df

Incluí la interacción "habtype2*year" en la parte del recuento de la fórmula, pero ahora quiero incluirla también en el segundo modelo (el binomial), pero si lo hago obtengo el siguiente error:

cr_f1 = formula(cr ~ depth + habtype2 + month + year + lightregime + depth*month + depth*lightregime + depth*habtype2 + habtype2*year | depth + habtype2 + month + year + lightregime + depth*month + depth*lightregime + depth*habtype2 + habtype2*year)

summary(zeroinfl(cr_f1, dist = "poisson", link = "logit", data = allUVCdata))

Error in solve.default(as.matrix(fit$hessian)) : 
  system is computationally singular: reciprocal condition number = 2.08629e-37

Lo mismo ocurre si intento incluir alguno de los otros términos de interacción que todavía quiero incluir en el modelo ("mes año", "mes lightregime" y "month*habtype2").

He buscado aquí en el foro y en google, parece que más gente se ha encontrado con este error (también en otras funciones que haciendo un zeroinfl), pero no he encontrado ninguna solución adecuada.

Datos: avistamientos de as especies en 29 lugares diferentes, >5200 observaciones (incluidos ceros).

¿Qué podría solucionar esto para que pueda ejecutar el modelo con los términos de interacción que quiero?

EDITAR : ha añadido algunos resultados nuevos para dar una idea del problema.

allUVCdata$year = as.numeric(as.character(allUVCdata$year))
cr_f1 = formula(cr ~ depth + lightregime + month + year + habtype2 + month*habtype2 + month*year + habtype2*year + depth*month)
summary(hurdle(cr_f1, dist = "poisson", link = "logit", data = allUVCdata))

Error in solve.default(as.matrix(fit_count$hessian)) : 
  system is computationally singular: reciprocal condition number = 6.23277e-26
> allUVCdata$year = as.factor(as.character(allUVCdata$year))
> table(allUVCdata$year, allUVCdata$cr)

          0    1    2    3    4    5    6    7
  2012  750  149   25   12    3    0    0    0
  2013 1133  209   69   16    4    1    1    0
  2014  844  387  142   42   11    7    0    1
  2015  833  401  125   31    5    3    1    2
> table(allUVCdata$month, allUVCdata$cr)

       0   1   2   3   4   5   6   7
  1  299  53  18   7   1   1   1   2
  10 346 104  40   9   4   4   0   0
  11 328 114  43  17   5   0   0   0
  12 350 112  29  10   2   0   0   0
  2  248  80  16   1   1   0   0   1
  3  303  82  24   6   0   0   1   0
  4  329  93  32   4   0   0   0   0
  5  277 105  28   9   1   1   0   0
  6  312 111  36  12   3   2   0   0
  7  362 113  46  14   5   2   0   0
  8  213 100  25   8   1   1   0   0
  9  193  79  24   4   0   0   0   0
> table(allUVCdata$month, allUVCdata$year)

     2012 2013 2014 2015
  1    41  129   72  140
  10   91  149  152  115
  11  112  121  150  124
  12  112  124  154  113
  2    35  108   79  125
  3    33  149  101  133
  4    88  105  150  115
  5   101  108   95  117
  6    94  115  142  125
  7   118  133  153  138
  8    61  114  100   73
  9    53   78   86   83

table(allUVCdata$habtype2, allUVCdata$year)

            2012 2013 2014 2015
  Exposed     93  138  120  144
  Pinnacles  274  386  339  338
  Unexposed  572  909  975  919

Answer 1

Es probable que el problema se deba a que la matriz modelo de la parte de recuento o de inflado cero (aunque por tu código parece que son la misma) está cerca del singular.

Una forma de avanzar es inspeccionar las matrices del modelo y ver si se puede observar si existe dependencia lineal entre alguna de las filas o columnas, o algo más causando singularidad .

Eso podría resultar complicado, así que antes de intentarlo, prueba a simplificar el modelo. ¿Tiene month ¿realmente tiene que ser un factor? A mí no me lo parece. En la parte de recuento del modelo, las estimaciones de los efectos principales y las interacciones con depth son todos de una magnitud similar, por lo que el mes podría ser numérico. Con la parte inflada a cero, las estimaciones de los efectos principales de month indican cierta no linealidad, pero los de la interacción son similares, por lo que volvería a utilizar month como numérico, e incluir un término cuadrático (quizás después de centrar para evitar la colinealidad entre los términos lineales y cuadráticos).