# R code for the Llamas (2010) Gretna (r) study, from Gorman and Johnson 
# in press

source('ky.R')
library(Design)
library(mlogit)
options(contrasts = c('contr.treatment', 'contr.poly'))

gretna <- read.csv("gretna.csv")
gretna$response <- as.factor(gretna$response)
levels(gretna$response) <- c("tap_or_trill", "approximant", "zero")

# mlogit
m.gretna <- mlogit.data(gretna, "response", "wide")
contrasts(m.gretna$agegroup) <- contr.sum(2)
contrasts(m.gretna$gender) <- contr.sum(2)
contrasts(m.gretna$class) <- contr.sum(2)
contrasts(m.gretna$preceding) <- contr.sum(11)
contrasts(m.gretna$stress) <- contr.sum(2)
contrasts(m.gretna$style) <- contr.sum(3)
contrasts(m.gretna$folseg) <- contr.sum(17)
contrasts(m.gretna$folbound) <- contr.sum(2)
summary(mlogit(response ~ 1 | agegroup + gender+ class + preceding + stress + style + folseg + folbound, data=m.gretna, ref.level="0"))

# proportional odds assumption table
summary(response ~ agegroup + gender + class + preceding + stress + style + folseg + folbound, data=gretna, fun=sf)