Bitterness of wine
wine.RdThe wine data set is adopted from Randall(1989) and from a
factorial experiment on factors determining the bitterness of
wine. Two treatment factors (temperature and contact) each have two
levels. Temperature and contact between juice and skins can be
controlled when cruching grapes during wine production. Nine judges
each assessed wine from two bottles from each of the four treatment
conditions, hence there are 72 observations in all.
Format
responsescorings of wine bitterness on a 0—100 continuous scale.
ratingordered factor with 5 levels; a grouped version of
response.temptemperature: factor with two levels.
contactfactor with two levels (
"no"and"yes").bottlefactor with eight levels.
judgefactor with nine levels.
References
Randall, J (1989). The analysis of sensory data by generalised linear model. Biometrical journal 7, pp. 781–793.
Tutz, G. and W. Hennevogl (1996). Random effects in ordinal regression models. Computational Statistics & Data Analysis 22, pp. 537–557.
Examples
head(wine)
#> response rating temp contact bottle judge
#> 1 36 2 cold no 1 1
#> 2 48 3 cold no 2 1
#> 3 47 3 cold yes 3 1
#> 4 67 4 cold yes 4 1
#> 5 77 4 warm no 5 1
#> 6 60 4 warm no 6 1
str(wine)
#> 'data.frame': 72 obs. of 6 variables:
#> $ response: num 36 48 47 67 77 60 83 90 17 22 ...
#> $ rating : Ord.factor w/ 5 levels "1"<"2"<"3"<"4"<..: 2 3 3 4 4 4 5 5 1 2 ...
#> $ temp : Factor w/ 2 levels "cold","warm": 1 1 1 1 2 2 2 2 1 1 ...
#> $ contact : Factor w/ 2 levels "no","yes": 1 1 2 2 1 1 2 2 1 1 ...
#> $ bottle : Factor w/ 8 levels "1","2","3","4",..: 1 2 3 4 5 6 7 8 1 2 ...
#> $ judge : Factor w/ 9 levels "1","2","3","4",..: 1 1 1 1 1 1 1 1 2 2 ...
## Variables 'rating' and 'response' are related in the following way:
(intervals <- seq(0,100, by = 20))
#> [1] 0 20 40 60 80 100
all(wine$rating == findInterval(wine$response, intervals)) ## ok
#> [1] TRUE
## A few illustrative tabulations:
## Table matching Table 5 in Randall (1989):
temp.contact.bottle <- with(wine, temp:contact:bottle)[drop=TRUE]
xtabs(response ~ temp.contact.bottle + judge, data = wine)
#> judge
#> temp.contact.bottle 1 2 3 4 5 6 7 8 9
#> cold:no:1 36 17 36 46 26 46 13 25 12
#> cold:no:2 48 22 50 27 45 30 19 32 29
#> cold:yes:3 47 14 42 48 61 54 31 39 47
#> cold:yes:4 67 50 23 32 41 37 29 40 28
#> warm:no:5 77 30 80 57 48 32 22 51 47
#> warm:no:6 60 51 81 37 41 60 43 45 38
#> warm:yes:7 83 90 73 84 58 88 32 42 72
#> warm:yes:8 90 70 62 58 55 73 49 67 65
## Table matching Table 6 in Randall (1989):
with(wine, {
tcb <- temp:contact:bottle
tcb <- tcb[drop=TRUE]
table(tcb, rating)
})
#> rating
#> tcb 1 2 3 4 5
#> cold:no:1 3 4 2 0 0
#> cold:no:2 1 5 3 0 0
#> cold:yes:3 1 2 5 1 0
#> cold:yes:4 0 5 3 1 0
#> warm:no:5 0 3 4 1 1
#> warm:no:6 0 2 4 2 1
#> warm:yes:7 0 1 2 2 4
#> warm:yes:8 0 0 3 5 1
## or simply: with(wine, table(bottle, rating))
## Table matching Table 1 in Tutz & Hennevogl (1996):
tab <- xtabs(as.numeric(rating) ~ judge + temp.contact.bottle,
data = wine)
colnames(tab) <-
paste(rep(c("c","w"), each = 4), rep(c("n", "n", "y", "y"), 2),
1:8, sep=".")
tab
#> temp.contact.bottle
#> judge c.n.1 c.n.2 c.y.3 c.y.4 w.n.5 w.n.6 w.y.7 w.y.8
#> 1 2 3 3 4 4 4 5 5
#> 2 1 2 1 3 2 3 5 4
#> 3 2 3 3 2 5 5 4 4
#> 4 3 2 3 2 3 2 5 3
#> 5 2 3 4 3 3 3 3 3
#> 6 3 2 3 2 2 4 5 4
#> 7 1 1 2 2 2 3 2 3
#> 8 2 2 2 3 3 3 3 4
#> 9 1 2 3 2 3 2 4 4
## A simple model:
m1 <- clm(rating ~ temp * contact, data = wine)
summary(m1)
#> formula: rating ~ temp * contact
#> data: wine
#>
#> link threshold nobs logLik AIC niter max.grad cond.H
#> logit flexible 72 -86.42 186.83 6(0) 5.22e-12 5.1e+01
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> tempwarm 2.3212 0.7009 3.311 0.000928 ***
#> contactyes 1.3475 0.6604 2.041 0.041300 *
#> tempwarm:contactyes 0.3595 0.9238 0.389 0.697129
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#>
#> Threshold coefficients:
#> Estimate Std. Error z value
#> 1|2 -1.4113 0.5454 -2.588
#> 2|3 1.1436 0.5097 2.244
#> 3|4 3.3771 0.6382 5.292
#> 4|5 4.9420 0.7509 6.581