```
library("betareg")
options(digits = 4)
## Section 4 from Ferrari and Cribari-Neto (2004)
data("GasolineYield", package = "betareg")
data("FoodExpenditure", package = "betareg")
## Table 1
gy <- betareg(yield ~ batch + temp, data = GasolineYield)
summary(gy)
```

```
Call:
betareg(formula = yield ~ batch + temp, data = GasolineYield)
Quantile residuals:
Min 1Q Median 3Q Max
-2.140 -0.570 0.120 0.704 1.751
Coefficients (mean model with logit link):
Estimate Std. Error z value Pr(>|z|)
(Intercept) -6.159571 0.182325 -33.78 < 2e-16 ***
batch1 1.727729 0.101229 17.07 < 2e-16 ***
batch2 1.322597 0.117902 11.22 < 2e-16 ***
batch3 1.572310 0.116105 13.54 < 2e-16 ***
batch4 1.059714 0.102360 10.35 < 2e-16 ***
batch5 1.133752 0.103523 10.95 < 2e-16 ***
batch6 1.040162 0.106036 9.81 < 2e-16 ***
batch7 0.543692 0.109127 4.98 6.3e-07 ***
batch8 0.495901 0.108926 4.55 5.3e-06 ***
batch9 0.385793 0.118593 3.25 0.0011 **
temp 0.010967 0.000413 26.58 < 2e-16 ***
Phi coefficients (precision model with identity link):
Estimate Std. Error z value Pr(>|z|)
(phi) 440 110 4 6.3e-05 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Type of estimator: ML (maximum likelihood)
Log-likelihood: 84.8 on 12 Df
Pseudo R-squared: 0.962
Number of iterations: 51 (BFGS) + 3 (Fisher scoring)
```

```
## Table 2
fe_lin <- lm(I(food/income) ~ income + persons, data = FoodExpenditure)
library("lmtest")
bptest(fe_lin)
```

```
studentized Breusch-Pagan test
data: fe_lin
BP = 5.9, df = 2, p-value = 0.05
```

```
Call:
betareg(formula = I(food/income) ~ income + persons, data = FoodExpenditure)
Quantile residuals:
Min 1Q Median 3Q Max
-2.533 -0.460 0.170 0.642 1.773
Coefficients (mean model with logit link):
Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.62255 0.22385 -2.78 0.0054 **
income -0.01230 0.00304 -4.05 5.1e-05 ***
persons 0.11846 0.03534 3.35 0.0008 ***
Phi coefficients (precision model with identity link):
Estimate Std. Error z value Pr(>|z|)
(phi) 35.61 8.08 4.41 1e-05 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Type of estimator: ML (maximum likelihood)
Log-likelihood: 45.3 on 4 Df
Pseudo R-squared: 0.388
Number of iterations: 28 (BFGS) + 4 (Fisher scoring)
```

```
## nested model comparisons via Wald and LR tests
fe_beta2 <- betareg(I(food/income) ~ income, data = FoodExpenditure)
lrtest(fe_beta, fe_beta2)
```

```
Likelihood ratio test
Model 1: I(food/income) ~ income + persons
Model 2: I(food/income) ~ income
#Df LogLik Df Chisq Pr(>Chisq)
1 4 45.3
2 3 40.5 -1 9.65 0.0019 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```

`waldtest(fe_beta, fe_beta2)`

```
Wald test
Model 1: I(food/income) ~ income + persons
Model 2: I(food/income) ~ income
Res.Df Df Chisq Pr(>Chisq)
1 34
2 35 -1 11.2 8e-04 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```

```
## Section 3 from online supplements to Simas et al. (2010)
## mean model as in gy above
## precision model with regressor temp
gy2 <- betareg(yield ~ batch + temp | temp, data = GasolineYield)
## MLE column in Table 19
summary(gy2)
```

```
Call:
betareg(formula = yield ~ batch + temp | temp, data = GasolineYield)
Quantile residuals:
Min 1Q Median 3Q Max
-2.104 -0.585 -0.143 0.690 2.520
Coefficients (mean model with logit link):
Estimate Std. Error z value Pr(>|z|)
(Intercept) -5.923236 0.183526 -32.27 < 2e-16 ***
batch1 1.601988 0.063856 25.09 < 2e-16 ***
batch2 1.297266 0.099100 13.09 < 2e-16 ***
batch3 1.565338 0.099739 15.69 < 2e-16 ***
batch4 1.030072 0.063288 16.28 < 2e-16 ***
batch5 1.154163 0.065643 17.58 < 2e-16 ***
batch6 1.019445 0.066351 15.36 < 2e-16 ***
batch7 0.622259 0.065632 9.48 < 2e-16 ***
batch8 0.564583 0.060185 9.38 < 2e-16 ***
batch9 0.359439 0.067141 5.35 8.6e-08 ***
temp 0.010359 0.000436 23.75 < 2e-16 ***
Phi coefficients (precision model with log link):
Estimate Std. Error z value Pr(>|z|)
(Intercept) 1.36409 1.22578 1.11 0.27
temp 0.01457 0.00362 4.03 5.7e-05 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Type of estimator: ML (maximum likelihood)
Log-likelihood: 87 on 13 Df
Pseudo R-squared: 0.952
Number of iterations: 33 (BFGS) + 28 (Fisher scoring)
```

```
## LRT row in Table 18
lrtest(gy, gy2)
```

```
Likelihood ratio test
Model 1: yield ~ batch + temp
Model 2: yield ~ batch + temp | temp
#Df LogLik Df Chisq Pr(>Chisq)
1 12 84.8
2 13 87.0 1 4.36 0.037 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```