####################################################################################
#This is a file with executable R code for chapter 8 of Natalia Levshina's (2015) 
#How to Do Linguistics with R. Amsterdam/Philadelphia: John Benjamins. 
####################################################################################

###Section 8.2


##Main text


install.packages(c("car", "coin", "nparcomp"))

library(Rling); library(car); library(coin); library(nparcomp)

data(NSL)

head(NSL)

str(NSL)

boxplot(NSL$MannerPath ~ NSL$Cohort, xlab = "Cohort", ylab = "Proportion of separate expressions", main = "Path and motion in NSL")

tapply(NSL$MannerPath, NSL$Cohort, mean)

aggregate(MannerPath ~ Cohort, data = NSL, function(x) shapiro.test(x)$p.value)

leveneTest(MannerPath ~ Cohort, data = NSL)

fligner.test(MannerPath ~ Cohort, data = NSL)

NSL.lm <- lm(MannerPath ~ Cohort, data = NSL)

summary(NSL.lm)

NSL.aov <- aov(MannerPath ~ Cohort, data = NSL)

summary(NSL.aov)

oneway.test(MannerPath ~ Cohort, data = NSL)

kruskal.test(MannerPath ~ Cohort, data = NSL)

oneway_test(MannerPath ~ Cohort, data = NSL, distribution = approximate(B = 9999))

TukeyHSD(NSL.aov)

plot(TukeyHSD(NSL.aov))

npar <- nparcomp(MannerPath ~ Cohort, data = NSL, type = "Tukey")

npar$Analysis


###Section 8.3


##Main text


install.packages("car") # if you haven't installed the package yet

library(Rling); library(car)

data(sharedref)

head(sharedref)

ref <- aggregate(mod ~ age + cohort, data = sharedref, FUN = mean)

ref

interaction.plot(ref$age, ref$cohort, ref$mod)

aggregate(mod ~ age + cohort, data = sharedref, function(x) shapiro.test(x)$p.value)

fligner.test(mod ~ interaction(age, cohort), data = sharedref)

contrasts(sharedref$age) <- cbind(c(-2, 1, 1), c(0, -1, 1))

sharedref$age

contrasts(sharedref$cohort) <- c(-1, 1)

sharedref$cohort

sharedref.aov <- aov(mod ~ age*cohort, data = sharedref)

Anova(sharedref.aov, type = "III")

Anova(sharedref.aov, type = "III", white.adjust = TRUE)

TukeyHSD(sharedref.aov, "age")

model.tables(sharedref.aov, type = "means")



###Section 8.4


##Main text


install.packages(c("gplots", "ggplot2", "nlme"))

library(Rling); library(gplots); library(ggplot2); library(nlme)

data(time_exper)

str(time_exper)

boxplot(rt ~ Subj, data = time_exper, xlab = "Subjects", ylab = "Reaction times, in ms", col = c(rep("grey", 10), rep("white", 10)))

legend("topright", c("Chinese", "English"), fill = c("grey", "white"))

m0 <- lme(rt ~ 1, random = ~ 1|Subj/Prime, data = time_exper, method = "ML")

m1 <- lme(rt ~ Lang, random = ~ 1|Subj/Prime, data = time_exper, method = "ML")

m1 <- update(m0,.~. + Lang)

m2 <- update(m1,.~. + Prime)

m3 <- update(m2,.~. + Lang:Prime)

anova(m0, m1, m2, m3)

summary(m3)

intervals(m3)


##Boxes with additional information

ggplot(time_exper, aes(x = Subj, y = rt, fill = Lang)) + geom_boxplot() + labs(x = "Subjects", y = "Reaction times, in ms") + scale_fill_grey(start = 0.5, end = 1)