How to improve the fit of a beta zero-inflated regression model (GAMLSS)?

Question

I'm working with a response variable with values between 0.0 and 1.0. I have a lot of zero. Thus, I'm using beta zero-inflated regression model. Specifically, I'm using the function gamlss from the gamlss package, with family = BEZI, to test the effect of a factor (three levels) on the response variable. With the wp function, I get the worm plot to evaluate the adjustment (following this document: http://www.gamlss.com/wp-content/uploads/2013/01/gamlss-manual.pdf). However, the adjustment is not good, and it is better when I use family = BEINF (which consideres 1 values. But, my data do not show 1s).

The AIC is lower and global deviance is higher when BEZI is used, but the fit (follow worm plot) is better when I consider BEINF.

Do you know what can I do to choose the best model or how to improve the fit? I appreciate any help you can bring me!

RESULT FOR BEZI

m1<-gamlss(resposta ~fase, data = ciliados, family = BEZI(),trace = F);summary(m1)

---

Family: c("BEZI", "Zero Inflated Beta")

Call: gamlss(formula = resposta ~ fase, family = BEZI(), data = ciliados, trace = F)

Fitting method: RS()

---

Mu link function: logit Mu Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) -3.3939 0.1812 -18.735 7.85e-16 *** faseBefore 0.5548 0.2065 2.686 0.0129 * faseFill 0.4391 0.2091 2.100 0.0464 *

Signif. codes: 0 ‘’ 0.001 ‘’ 0.01 ‘’ 0.05 ‘.’ 0.1 ‘ ’ 1

---

Sigma link function: log Sigma Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 6.8640 0.6332 10.84 9.93e-11 ***

Signif. codes: 0 ‘’ 0.001 ‘’ 0.01 ‘’ 0.05 ‘.’ 0.1 ‘ ’ 1

---

Nu link function: logit Nu Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 1.5686 0.4916 3.191 0.00393 **

Signif. codes: 0 ‘’ 0.001 ‘’ 0.01 ‘’ 0.05 ‘.’ 0.1 ‘ ’ 1

---

No. of observations in the fit: 29 Degrees of Freedom for the fit: 5 Residual Deg. of Freedom: 24 at cycle: 8

Global Deviance: -9.016796 AIC: 0.9832042 SBC: 7.819683

---

RESULT FOR BEINF

m1<-gamlss(resposta ~fase, data = ciliados, family = BEINF(),trace = F);summary(m1)

---

Family: c("BEINF", "Beta Inflated")

Call: gamlss(formula = resposta ~ fase, family = BEINF(), data = ciliados, trace = F)

Fitting method: RS()

---

Mu link function: logit Mu Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) -3.3938 0.1812 -18.735 2e-15 *** faseBefore 0.5548 0.2065 2.686 0.0132 * faseFill 0.4391 0.2091 2.100 0.0469 *

Signif. codes: 0 ‘’ 0.001 ‘’ 0.01 ‘’ 0.05 ‘.’ 0.1 ‘ ’ 1

---

Sigma link function: logit Sigma Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) -3.3997 0.3268 -10.4 3.61e-10 ***

Signif. codes: 0 ‘’ 0.001 ‘’ 0.01 ‘’ 0.05 ‘.’ 0.1 ‘ ’ 1

---

Nu link function: log Nu Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 1.5686 0.4916 3.191 0.00407 **

Signif. codes: 0 ‘’ 0.001 ‘’ 0.01 ‘’ 0.05 ‘.’ 0.1 ‘ ’ 1

---

Tau link function: log Tau Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) -15.27 924.26 -0.017 0.987

---

No. of observations in the fit: 29 Degrees of Freedom for the fit: 6 Residual Deg. of Freedom: 23 at cycle: 6

Global Deviance: -9.016793 AIC: 2.983207 SBC: 11.18698

---

Answer 1

For a zero inflated beta you can use BEINF0 in gamlss, [which is a submodel of BEINF].

The reason for your different residual plots, may be that BEZI and BEINF0 use slightly different parameterisations.

In your BEINF fit you will find that the fitted probability that Y=1 is effectively zero, so best to use BEINF0.

Note also that you are only modelling mu in BEINF. In BEINF0(mu,sigma,mu) you could try including models for mu, sigma and nu.