Ajuste de modelos aditivos generalizados para posición, escala y forma (GAMLSS) con datos Dengue-Clima en Cali

Modelo Poisson

library(gamlss)
modelo_PO <- gamlss(
  formula = casos ~ cs(prec1, df = 8)+cs(tmax1, df = 8)+cs(hum3, df = 8)+cs(ao, df = 8)+cs(mes, df = 8)+cs(lncasos1, df = 8),
  family  = PO,
  data    = data,
  trace   = FALSE
)
summary(modelo_PO)
## ******************************************************************
## Family:  c("PO", "Poisson") 
## 
## Call:  gamlss(formula = casos ~ cs(prec1, df = 8) + cs(tmax1,  
##     df = 8) + cs(hum3, df = 8) + cs(ao, df = 8) + cs(mes,  
##     df = 8) + cs(lncasos1, df = 8), family = PO, data = data,  
##     trace = FALSE) 
## 
## Fitting method: RS() 
## 
## ------------------------------------------------------------------
## Mu link function:  log
## Mu Coefficients:
##                        Estimate Std. Error t value Pr(>|t|)    
## (Intercept)          -6.680e+01  2.539e+00 -26.306  < 2e-16 ***
## cs(prec1, df = 8)     1.801e-02  4.109e-03   4.382 3.82e-05 ***
## cs(tmax1, df = 8)     1.820e-01  7.765e-03  23.445  < 2e-16 ***
## cs(hum3, df = 8)     -8.295e-03  1.376e-03  -6.030 5.94e-08 ***
## cs(ao, df = 8)        3.174e-02  6.658e-04  47.675  < 2e-16 ***
## cs(mes, df = 8)      -4.372e-02  1.742e-03 -25.093  < 2e-16 ***
## cs(lncasos1, df = 8)  7.199e-01  7.177e-03 100.311  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## ------------------------------------------------------------------
## NOTE: Additive smoothing terms exist in the formulas: 
##  i) Std. Error for smoothers are for the linear effect only. 
## ii) Std. Error for the linear terms maybe are not accurate. 
## ------------------------------------------------------------------
## No. of observations in the fit:  129 
## Degrees of Freedom for the fit:  55.00002
##       Residual Deg. of Freedom:  73.99998 
##                       at cycle:  20 
##  
## Global Deviance:     2337.892 
##             AIC:     2447.892 
##             SBC:     2605.182 
## ******************************************************************

Modelo Yule

modelo_YU <- gamlss(
  formula = casos ~ cs(prec1, df = 8)+cs(tmax1, df = 8)+cs(hum3, df = 8)+cs(ao, df = 8)+cs(mes, df = 8)+cs(lncasos1, df = 8),
  family  = YULE,
  mu.fix = TRUE,
  method = mixed(1,2),
  data    = data,
  trace   = FALSE
)
summary(modelo_YU)
## ******************************************************************
## Family:  c("YULE", "Yule") 
## 
## Call:  gamlss(formula = casos ~ cs(prec1, df = 8) + cs(tmax1,  
##     df = 8) + cs(hum3, df = 8) + cs(ao, df = 8) + cs(mes,  
##     df = 8) + cs(lncasos1, df = 8), family = YULE,  
##     data = data, method = mixed(1, 2), mu.fix = TRUE,      trace = FALSE) 
## 
## Fitting method: mixed(1, 2) 
## 
## ------------------------------------------------------------------
## Mu parameter is fixed 
## Mu =  251.2713 
## ------------------------------------------------------------------
## No. of observations in the fit:  129 
## Degrees of Freedom for the fit:  0
##       Residual Deg. of Freedom:  129 
##                       at cycle:  1 
##  
## Global Deviance:     2628.616 
##             AIC:     2628.616 
##             SBC:     2628.616 
## ******************************************************************

Modelo Binomial negativo II

modelo_NB <- gamlss(
  formula = casos ~ cs(prec1, df = 8)+cs(tmax1, df = 8)+cs(hum3, df = 8)+cs(ao, df = 8)+cs(mes, df = 8)+cs(lncasos1, df = 8),
  sigma.formula = ~ cs(prec1, df = 8)+cs(tmax1, df = 8)+cs(hum3, df = 8)+cs(ao, df = 8)+cs(mes, df = 8)+cs(lncasos1, df = 8),
  family  = NBII,
  mu.fix = TRUE,
  sigma.fix = TRUE,
  method = mixed(1,2),
  data    = data,
  trace   = FALSE
)
summary(modelo_NB)
## ******************************************************************
## Family:  c("NBII", "Negative Binomial type II") 
## 
## Call:  gamlss(formula = casos ~ cs(prec1, df = 8) + cs(tmax1,  
##     df = 8) + cs(hum3, df = 8) + cs(ao, df = 8) + cs(mes,  
##     df = 8) + cs(lncasos1, df = 8), sigma.formula = ~cs(prec1,  
##     df = 8) + cs(tmax1, df = 8) + cs(hum3, df = 8) +  
##     cs(ao, df = 8) + cs(mes, df = 8) + cs(lncasos1,  
##     df = 8), family = NBII, data = data, method = mixed(1,  
##     2), mu.fix = TRUE, sigma.fix = TRUE, trace = FALSE) 
## 
## Fitting method: mixed(1, 2) 
## 
## ------------------------------------------------------------------
## Mu parameter is fixed 
## Mu is equal with the vector ( 177.6357 , 171.1357 , 189.1357 , 209.6357 , ...) 
## ------------------------------------------------------------------
## Sigma parameter is fixed 
## Sigma =  446.419 
## ------------------------------------------------------------------
## No. of observations in the fit:  129 
## Degrees of Freedom for the fit:  0
##       Residual Deg. of Freedom:  129 
##                       at cycle:  1 
##  
## Global Deviance:     1694.316 
##             AIC:     1694.316 
##             SBC:     1694.316 
## ******************************************************************

Modelo Poisson-inv Gaussiano

modelo_PIG <- gamlss(
  formula = casos ~ cs(prec1, df = 8)+cs(tmax1, df = 8)+cs(hum3, df = 8)+cs(ao, df = 8)+cs(mes, df = 8)+cs(lncasos1, df = 8),
  sigma.formula = ~ cs(prec1, df = 8)+cs(tmax1, df = 8)+cs(hum3, df = 8)+cs(ao, df = 8)+cs(mes, df = 8)+cs(lncasos1, df = 8),
  family  = PIG,
  mu.fix = TRUE,
  sigma.fix = TRUE,
  method = mixed(1,2),
  data    = data,
  trace   = FALSE
)
summary(modelo_PIG)
## ******************************************************************
## Family:  c("PIG", "Poisson.Inverse.Gaussian") 
## 
## Call:  gamlss(formula = casos ~ cs(prec1, df = 8) + cs(tmax1,  
##     df = 8) + cs(hum3, df = 8) + cs(ao, df = 8) + cs(mes,  
##     df = 8) + cs(lncasos1, df = 8), sigma.formula = ~cs(prec1,  
##     df = 8) + cs(tmax1, df = 8) + cs(hum3, df = 8) +  
##     cs(ao, df = 8) + cs(mes, df = 8) + cs(lncasos1,  
##     df = 8), family = PIG, data = data, method = mixed(1,  
##     2), mu.fix = TRUE, sigma.fix = TRUE, trace = FALSE) 
## 
## Fitting method: mixed(1, 2) 
## 
## ------------------------------------------------------------------
## Mu parameter is fixed 
## Mu is equal with the vector ( 177.6357 , 171.1357 , 189.1357 , 209.6357 , ...) 
## ------------------------------------------------------------------
## Sigma parameter is fixed 
## Sigma =  1.776641 
## ------------------------------------------------------------------
## No. of observations in the fit:  129 
## Degrees of Freedom for the fit:  0
##       Residual Deg. of Freedom:  129 
##                       at cycle:  1 
##  
## Global Deviance:     1602.082 
##             AIC:     1602.082 
##             SBC:     1602.082 
## ******************************************************************

Modelo Poisson doble

modelo_DPO <- gamlss(
  formula = casos ~ cs(prec1, df = 8)+cs(tmax1, df = 8)+cs(hum3, df = 8)+cs(ao, df = 8)+cs(mes, df = 8)+cs(lncasos1, df = 8),
  sigma.formula = ~ cs(prec1, df = 8)+cs(tmax1, df = 8)+cs(hum3, df = 8)+cs(ao, df = 8)+cs(mes, df = 8)+cs(lncasos1, df = 8),
  family  = DPO,
  mu.fix = TRUE,
  sigma.fix = TRUE,
  method = mixed(1,2),
  data    = data,
  trace   = FALSE
)
summary(modelo_DPO)
## ******************************************************************
## Family:  c("DPO", "Double Poisson") 
## 
## Call:  gamlss(formula = casos ~ cs(prec1, df = 8) + cs(tmax1,  
##     df = 8) + cs(hum3, df = 8) + cs(ao, df = 8) + cs(mes,  
##     df = 8) + cs(lncasos1, df = 8), sigma.formula = ~cs(prec1,  
##     df = 8) + cs(tmax1, df = 8) + cs(hum3, df = 8) +  
##     cs(ao, df = 8) + cs(mes, df = 8) + cs(lncasos1,  
##     df = 8), family = DPO, data = data, method = mixed(1,  
##     2), mu.fix = TRUE, sigma.fix = TRUE, trace = FALSE) 
## 
## Fitting method: mixed(1, 2) 
## 
## ------------------------------------------------------------------
## Mu parameter is fixed 
## Mu is equal with the vector ( 177.6357 , 171.1357 , 189.1357 , 209.6357 , ...) 
## ------------------------------------------------------------------
## Sigma parameter is fixed 
## Sigma =  1 
## ------------------------------------------------------------------
## No. of observations in the fit:  129 
## Degrees of Freedom for the fit:  0
##       Residual Deg. of Freedom:  129 
##                       at cycle:  1 
##  
## Global Deviance:     7713.224 
##             AIC:     7713.224 
##             SBC:     7713.224 
## ******************************************************************

Modelo Poisson generalizado

modelo_GPO <- gamlss(
  formula = casos ~ cs(prec1, df = 8)+cs(tmax1, df = 8)+cs(hum3, df = 8)+cs(ao, df = 8)+cs(mes, df = 8)+cs(lncasos1, df = 8),
  sigma.formula = ~ cs(prec1, df = 8)+cs(tmax1, df = 8)+cs(hum3, df = 8)+cs(ao, df = 8)+cs(mes, df = 8)+cs(lncasos1, df = 8),
  family  = GPO,
  mu.fix = TRUE,
  sigma.fix = TRUE,
  method = mixed(1,2),
  data    = data,
  trace   = FALSE
)
summary(modelo_GPO)
## ******************************************************************
## Family:  c("GPO", "Generalised Poisson") 
## 
## Call:  gamlss(formula = casos ~ cs(prec1, df = 8) + cs(tmax1,  
##     df = 8) + cs(hum3, df = 8) + cs(ao, df = 8) + cs(mes,  
##     df = 8) + cs(lncasos1, df = 8), sigma.formula = ~cs(prec1,  
##     df = 8) + cs(tmax1, df = 8) + cs(hum3, df = 8) +  
##     cs(ao, df = 8) + cs(mes, df = 8) + cs(lncasos1,  
##     df = 8), family = GPO, data = data, method = mixed(1,  
##     2), mu.fix = TRUE, sigma.fix = TRUE, trace = FALSE) 
## 
## Fitting method: mixed(1, 2) 
## 
## ------------------------------------------------------------------
## Mu parameter is fixed 
## Mu is equal with the vector ( 177.6357 , 171.1357 , 189.1357 , 209.6357 , ...) 
## ------------------------------------------------------------------
## Sigma parameter is fixed 
## Sigma =  0.1 
## ------------------------------------------------------------------
## No. of observations in the fit:  129 
## Degrees of Freedom for the fit:  0
##       Residual Deg. of Freedom:  129 
##                       at cycle:  1 
##  
## Global Deviance:     1648.157 
##             AIC:     1648.157 
##             SBC:     1648.157 
## ******************************************************************

Selección de modelos a través de GAIC

GAIC(modelo_PO, modelo_YU, modelo_NB, modelo_PIG, modelo_DPO, modelo_GPO)
##                  df      AIC
## modelo_PIG  0.00000 1602.082
## modelo_GPO  0.00000 1648.157
## modelo_NB   0.00000 1694.316
## modelo_PO  55.00002 2447.892
## modelo_YU   0.00000 2628.616
## modelo_DPO  0.00000 7713.224