Vowel contrast: GAMM analysis
Holger Mitterer
2026-02-18
Note that this Markdown requires the file “mouse_data.txt”, which is generated by running the preprocessing markdown as well as the mousePreProcess.R script.
TARGET FIXATIONS
rm(list = ls())
library(lmerTest)## Loading required package: lme4## Loading required package: Matrix##
## Attaching package: 'lmerTest'## The following object is masked from 'package:lme4':
##
## lmer## The following object is masked from 'package:stats':
##
## steplibrary(stringr)
library(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.2.0 ✔ purrr 1.2.1
## ✔ forcats 1.0.1 ✔ readr 2.2.0
## ✔ ggplot2 4.0.2 ✔ tibble 3.3.1
## ✔ lubridate 1.9.5 ✔ tidyr 1.3.2## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ tidyr::expand() masks Matrix::expand()
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ✖ tidyr::pack() masks Matrix::pack()
## ✖ tidyr::unpack() masks Matrix::unpack()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorslibrary(mgcv)## Loading required package: nlme
##
## Attaching package: 'nlme'
##
## The following object is masked from 'package:dplyr':
##
## collapse
##
## The following object is masked from 'package:lme4':
##
## lmList
##
## This is mgcv 1.9-3. For overview type 'help("mgcv-package")'.library(itsadug)## Loading required package: plotfunctions
##
## Attaching package: 'plotfunctions'
##
## The following object is masked from 'package:ggplot2':
##
## alpha
##
## Loaded package itsadug 2.5 (see 'help("itsadug")' ).options(warn = 0)
gazeData = read.table("dataCatalan/gazeTarget.txt", T, sep = "\t", quote = "")
trialData = read.table("dataCatalan/trialInfo.txt", T, sep = "\t", quote = "\"")
expFilter = grepl("_", trialData$condition)
gazeDataE = gazeData %>% filter(expFilter)
trialDataE = trialData %>% filter(expFilter)
timeline = -200 + (0:120) * 10
colnames(gazeDataE) = paste0("T_",timeline)
#contrast coding
temp = str_split_fixed(trialDataE$condition,"_",2)
trialDataE$accent = temp[, 1]
trialDataE$contrast = temp[, 2]
trialDataE$accentC = -0.5 + as.numeric(trialDataE$accent == "SEng")
table(trialDataE$accentC)##
## -0.5 0.5
## 3596 3598trialDataE$contrastC = -0.5 + as.numeric(trialDataE$contrast == "DRESS" |trialDataE$contrast == "TRAP" )
table(trialDataE$contrastC)##
## -0.5 0.5
## 3598 3596useData = cbind(trialDataE,gazeDataE)Analysis of the time course using GAMMs
We first need to reorganize the data.
gamData = useData[ ,c(1,2,4:8,10:13,14:134)] %>%
pivot_longer(., cols = starts_with("T_"),
names_to = "Time",
values_to = "position") %>% as.data.frame
gamData$position = as.numeric(gamData$position)
gamData$Time = gsub("T_","", gamData$Time) %>% as.numeric()
gamData = gamData[ order(gamData$pp,gamData$trial,gamData$Time), ]
gamData$pp = as.factor(gamData$pp)
gamData$item = as.factor(gamData$item)
gamData$isAustralian = ifelse(gamData$accent == "AusE",1,0)
gamData$isAustralianDress = gamData$isAustralian *
as.numeric(gamData$contrast == "DRESS" |gamData$contrast == "TRAP" )
gamData$vowelContrast = ifelse(gamData$contrast == "DRESS" |gamData$contrast == "TRAP",
"midFront","tensenessI") %>% as.factor()
#restricting the data to 0-900 after vowel onset
gamData = gamData %>% filter(Time > 0 & Time < 900)
gamData <- start_event(gamData, event=c("pp","trial"))
gamDataPlot = gamData %>% group_by(vowelContrast, Time, accent) %>%
summarize(meanFix = mean(position))## `summarise()` has regrouped the output.
## ℹ Summaries were computed grouped by vowelContrast, Time, and accent.
## ℹ Output is grouped by vowelContrast and Time.
## ℹ Use `summarise(.groups = "drop_last")` to silence this message.
## ℹ Use `summarise(.by = c(vowelContrast, Time, accent))` for per-operation
## grouping (`?dplyr::dplyr_by`) instead.ggplot(gamDataPlot, aes(x = Time, y = meanFix, col = vowelContrast)) +
geom_line() + facet_wrap(~accent)
This is the descriptive data
Building an initial GAMM
fullModelFile = "target_overall_fullModels.Rdata"
if (file.exists(fullModelFile))
{
load(fullModelFile)
}else{
interaction.Targ.gamm0 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
s(Time, by = isAustralian) +
s(Time, by = isAustralianDress) +
s(Time, pp, by = vowelContrast, bs="fs") +
s(Time, item, by = vowelContrast, bs="fs"),
family = "binomial",
data=gamData, discrete = T, nThreads = 3)
#dealing with autocorrelation
myAutoCorr =acf_resid(interaction.Targ.gamm0)
myRhoVal = myAutoCorr[2]
interaction.Targ.gamm1 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast, k=12) +
s(Time, by = isAustralian, k= 12) +
s(Time, by = isAustralianDress) +
s(Time, pp, by = vowelContrast, bs="fs") +
s(Time, item, by = vowelContrast, bs="fs"),
family = "binomial",
data=gamData, rho = myRhoVal,
AR.start = gamData$start.event, nThreads = 3, discrete = T)
save(list= c("interaction.Targ.gamm0","interaction.Targ.gamm1"), file = fullModelFile)
}Model criticism
How are we doing regarding auto-correlation?
myAutoCorr =acf_resid(interaction.Targ.gamm1)
Is the default for the smooth (k = 10) sufficient?
gam.check(interaction.Targ.gamm1)
##
## Method: fREML Optimizer: perf chol
## $grad
## [1] 0.0069406337 0.0024180552 -0.0022415194 -0.0019898415 0.0314626836
## [6] 0.0153926254 0.0354470408 0.0445616939 0.0195573102 0.0377280812
## [11] 0.0048887509 0.0066139624 0.0002983201 0.0079791663 0.0053302942
## [16] 0.0053596874
##
## $hess
## [,1] [,2] [,3] [,4] [,5]
## [1,] 3.340892e+00 -1.798569e-04 9.577181e-04 2.132693e-01 2.147290e-01
## [2,] -1.798569e-04 3.053790e+00 2.571085e-03 -7.880004e-04 -4.850833e-06
## [3,] 9.577181e-04 2.571085e-03 2.311638e-03 2.669913e-03 3.149169e-05
## [4,] 2.132693e-01 -7.880004e-04 2.669913e-03 7.238158e-01 7.384829e-03
## [5,] 2.147290e-01 -4.850833e-06 3.149169e-05 7.384829e-03 2.474965e+01
## [6,] 1.514245e-02 -4.698830e-07 4.643919e-06 -4.979736e-04 2.700870e-01
## [7,] -3.966766e-02 4.387085e-06 -2.615465e-05 -8.140549e-03 -4.895610e-01
## [8,] -1.141793e-05 1.428172e-01 8.473813e-05 -2.807682e-05 -6.192438e-07
## [9,] 7.968583e-07 4.237776e-02 2.217294e-05 -5.403145e-06 1.060016e-07
## [10,] 1.838134e-06 -4.716512e-02 -5.144949e-05 1.067204e-05 3.869345e-08
## [11,] 7.242385e-02 -2.295805e-06 1.297510e-05 1.937084e-03 -7.943797e-02
## [12,] 1.582360e-03 1.666258e-07 3.421469e-07 -4.328364e-04 -6.932169e-02
## [13,] -1.222761e-02 1.182772e-06 -7.226080e-06 -2.399865e-03 -4.764104e-02
## [14,] -5.168744e-06 6.239933e-02 5.907710e-05 -1.778927e-05 -2.573258e-07
## [15,] 3.484781e-07 8.836137e-03 4.127427e-07 -1.037451e-06 3.380736e-08
## [16,] 7.547161e-07 -3.197547e-02 -3.306720e-05 7.285799e-06 -7.595624e-09
## [,6] [,7] [,8] [,9] [,10]
## [1,] 1.514245e-02 -3.966766e-02 -1.141793e-05 7.968583e-07 1.838134e-06
## [2,] -4.698830e-07 4.387085e-06 1.428172e-01 4.237776e-02 -4.716512e-02
## [3,] 4.643919e-06 -2.615465e-05 8.473813e-05 2.217294e-05 -5.144949e-05
## [4,] -4.979736e-04 -8.140549e-03 -2.807682e-05 -5.403145e-06 1.067204e-05
## [5,] 2.700870e-01 -4.895610e-01 -6.192438e-07 1.060016e-07 3.869345e-08
## [6,] 1.108683e+01 -7.695496e-01 5.950836e-08 -2.945867e-08 1.745073e-08
## [7,] -7.695496e-01 1.662183e+01 1.237807e-07 4.088049e-08 -6.623594e-08
## [8,] 5.950836e-08 1.237807e-07 2.069404e+01 -3.315983e-01 -6.395533e-01
## [9,] -2.945867e-08 4.088049e-08 -3.315983e-01 1.399071e+01 -5.609156e-01
## [10,] 1.745073e-08 -6.623594e-08 -6.395533e-01 -5.609156e-01 1.570437e+01
## [11,] -3.558707e-02 -5.819340e-02 -2.064639e-07 2.566554e-08 2.269722e-08
## [12,] -3.876018e-02 -2.676317e-02 3.358399e-08 -8.862064e-09 1.955237e-09
## [13,] -1.524750e-02 -5.407243e-02 2.361057e-08 1.411396e-08 -1.895829e-08
## [14,] 1.108600e-08 8.871286e-08 -8.433821e-02 -4.216662e-02 -8.025604e-02
## [15,] -7.207956e-09 7.932692e-09 -5.520802e-02 -4.532412e-02 -2.846319e-02
## [16,] 1.572271e-08 -4.623023e-08 -4.613904e-02 -2.238309e-02 -6.695430e-02
## [,11] [,12] [,13] [,14] [,15]
## [1,] 7.242385e-02 1.582360e-03 -1.222761e-02 -5.168744e-06 3.484781e-07
## [2,] -2.295805e-06 1.666258e-07 1.182772e-06 6.239933e-02 8.836137e-03
## [3,] 1.297510e-05 3.421469e-07 -7.226080e-06 5.907710e-05 4.127427e-07
## [4,] 1.937084e-03 -4.328364e-04 -2.399865e-03 -1.778927e-05 -1.037451e-06
## [5,] -7.943797e-02 -6.932169e-02 -4.764104e-02 -2.573258e-07 3.380736e-08
## [6,] -3.558707e-02 -3.876018e-02 -1.524750e-02 1.108600e-08 -7.207956e-09
## [7,] -5.819340e-02 -2.676317e-02 -5.407243e-02 8.871286e-08 7.932692e-09
## [8,] -2.064639e-07 3.358399e-08 2.361057e-08 -8.433821e-02 -5.520802e-02
## [9,] 2.566554e-08 -8.862064e-09 1.411396e-08 -4.216662e-02 -4.532412e-02
## [10,] 2.269722e-08 1.955237e-09 -1.895829e-08 -8.025604e-02 -2.846319e-02
## [11,] 6.437702e+00 -1.742088e-01 -2.092983e-01 -8.719682e-08 8.840444e-09
## [12,] -1.742088e-01 3.545384e+00 -6.236342e-02 1.148506e-08 -2.224622e-09
## [13,] -2.092983e-01 -6.236342e-02 4.774143e+00 2.164207e-08 2.697182e-09
## [14,] -8.719682e-08 1.148506e-08 2.164207e-08 8.121517e+00 -1.602243e-01
## [15,] 8.840444e-09 -2.224622e-09 2.697182e-09 -1.602243e-01 3.910068e+00
## [16,] 4.075270e-09 3.163292e-09 -1.405559e-08 -3.110829e-01 4.323348e-02
## [,16]
## [1,] 7.547161e-07
## [2,] -3.197547e-02
## [3,] -3.306720e-05
## [4,] 7.285799e-06
## [5,] -7.595624e-09
## [6,] 1.572271e-08
## [7,] -4.623023e-08
## [8,] -4.613904e-02
## [9,] -2.238309e-02
## [10,] -6.695430e-02
## [11,] 4.075270e-09
## [12,] 3.163292e-09
## [13,] -1.405559e-08
## [14,] -3.110829e-01
## [15,] 4.323348e-02
## [16,] 6.009899e+00
##
## Model rank = 1746 / 1746
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(Time):vowelContrastmidFront 11.00 8.67 1 0.35
## s(Time):vowelContrasttensenessI 11.00 8.15 1 0.43
## s(Time):isAustralian 12.00 2.03 1 0.49
## s(Time):isAustralianDress 10.00 5.16 1 0.45
## s(Time,pp):vowelContrastmidFront 450.00 152.01 1 0.45
## s(Time,pp):vowelContrasttensenessI 450.00 157.48 1 0.46
## s(Time,item):vowelContrastmidFront 400.00 51.94 1 0.48
## s(Time,item):vowelContrasttensenessI 400.00 61.49 1 0.47Let’s have a look at this model no whether there is an overall effect of contrast
summary(interaction.Targ.gamm1)##
## Family: binomial
## Link function: logit
##
## Formula:
## position ~ vowelContrast + s(Time, by = vowelContrast, k = 12) +
## s(Time, by = isAustralian, k = 12) + s(Time, by = isAustralianDress) +
## s(Time, pp, by = vowelContrast, bs = "fs") + s(Time, item,
## by = vowelContrast, bs = "fs")
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.23484 0.07650 3.070 0.00214 **
## vowelContrasttensenessI 0.01743 0.10737 0.162 0.87103
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(Time):vowelContrastmidFront 8.672 10.003 628.280 <2e-16 ***
## s(Time):vowelContrasttensenessI 8.146 9.644 472.640 <2e-16 ***
## s(Time):isAustralian 2.032 2.059 0.389 0.714
## s(Time):isAustralianDress 5.159 6.370 5.516 0.511
## s(Time,pp):vowelContrastmidFront 152.015 449.000 714.637 <2e-16 ***
## s(Time,pp):vowelContrasttensenessI 157.481 449.000 714.361 <2e-16 ***
## s(Time,item):vowelContrastmidFront 51.936 199.000 210.473 <2e-16 ***
## s(Time,item):vowelContrasttensenessI 61.485 199.000 269.413 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.254 Deviance explained = 20.3%
## fREML = 2.1037e+05 Scale est. = 1 n = 640266There is a main effect of competitor, meaning that, overall, the mouse position is closer to target when the competitor is dissimilar to the target. Note that the significant effects of the tensors only indicate changes over time, they do not constitute evidence for a significant difference between conditions.
smooth <- plot_smooth(interaction.Targ.gamm1, view = "Time",
plot_all = "vowelContrast",
rug = FALSE, rm.ranef = TRUE, lwd = 4)## Summary:
## * vowelContrast : factor; set to the value(s): midFront, tensenessI.
## * Time : numeric predictor; with 30 values ranging from 10.000000 to 890.000000.
## * isAustralian : numeric predictor; set to the value(s): 0.
## * isAustralianDress : numeric predictor; set to the value(s): 0.
## * pp : factor; set to the value(s): CAT01. (Might be canceled as random effect, check below.)
## * item : factor; set to the value(s): baggage. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(Time,pp):vowelContrastmidFront,s(Time,pp):vowelContrasttensenessI,s(Time,item):vowelContrastmidFront,s(Time,item):vowelContrasttensenessI
## 
Now model comparison
We succesively take effects out to prepare model comparisons.
reducedModelFile = "target_overall_reducedModels.Rdata"
if (file.exists(reducedModelFile))
{
load(reducedModelFile)
}else{
minusInteraction.Targ.gamm <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
s(Time, by = isAustralian) +
#s(Time, by = isAustralianDress) +
s(Time, pp, by = vowelContrast, bs="fs") +
s(Time, item, by = vowelContrast, bs="fs"),
family = "binomial",
data=gamData, rho = myRhoVal,
AR.start = gamData$start.event, nThreads = 3, discrete = T)
minusAccent.Targ.gamm <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
#s(Time, by = isAustralian) +
#s(Time, by = isAustralianDress) +
s(Time, pp, by = vowelContrast, bs="fs") +
s(Time, item, by = vowelContrast, bs="fs"),
family = "binomial",
data=gamData, rho = myRhoVal,
AR.start = gamData$start.event, nThreads = 3, discrete = T)
timeOnly.Targ.gamm <- bam( position ~ 1 +
s(Time) +
#s(Time, by = isAustralian) +
#s(Time, by = isAustralianDress) +
s(Time, pp, by = vowelContrast, bs="fs") +
s(Time, item, by = vowelContrast, bs="fs"), family = "binomial",
data=gamData, rho = myRhoVal,
AR.start = gamData$start.event, nThreads = 3, discrete = T)
save(list= c("minusInteraction.Targ.gamm","minusAccent.Targ.gamm","timeOnly.Targ.gamm"), file = reducedModelFile)
}Model comparisons
Full model versus no interaction
compareML(interaction.Targ.gamm1, minusInteraction.Targ.gamm)## interaction.Targ.gamm1: position ~ vowelContrast + s(Time, by = vowelContrast, k = 12) +
## s(Time, by = isAustralian, k = 12) + s(Time, by = isAustralianDress) +
## s(Time, pp, by = vowelContrast, bs = "fs") + s(Time, item,
## by = vowelContrast, bs = "fs")
##
## minusInteraction.Targ.gamm: position ~ vowelContrast + s(Time, by = vowelContrast) + s(Time,
## by = isAustralian) + s(Time, pp, by = vowelContrast, bs = "fs") +
## s(Time, item, by = vowelContrast, bs = "fs")
##
## Model minusInteraction.Targ.gamm preferred: lower fREML score (145.423), and lower df (1.000).
## -----
## Model Score Edf Difference Df
## 1 interaction.Targ.gamm1 210367.6 16
## 2 minusInteraction.Targ.gamm 210222.2 15 -145.423 1.000
##
## AIC difference: 492.28, model minusInteraction.Targ.gamm has lower AIC.Both main effects versus minus accent
compareML(minusInteraction.Targ.gamm, minusAccent.Targ.gamm)## minusInteraction.Targ.gamm: position ~ vowelContrast + s(Time, by = vowelContrast) + s(Time,
## by = isAustralian) + s(Time, pp, by = vowelContrast, bs = "fs") +
## s(Time, item, by = vowelContrast, bs = "fs")
##
## minusAccent.Targ.gamm: position ~ vowelContrast + s(Time, by = vowelContrast) + s(Time,
## pp, by = vowelContrast, bs = "fs") + s(Time, item, by = vowelContrast,
## bs = "fs")
##
## Model minusAccent.Targ.gamm preferred: lower fREML score (219.395), and lower df (1.000).
## -----
## Model Score Edf Difference Df
## 1 minusInteraction.Targ.gamm 210222.2 15
## 2 minusAccent.Targ.gamm 210002.8 14 -219.395 1.000
##
## AIC difference: 972.60, model minusAccent.Targ.gamm has lower AIC.Main effect of Contrast versus time-only
compareML(minusAccent.Targ.gamm, timeOnly.Targ.gamm)## minusAccent.Targ.gamm: position ~ vowelContrast + s(Time, by = vowelContrast) + s(Time,
## pp, by = vowelContrast, bs = "fs") + s(Time, item, by = vowelContrast,
## bs = "fs")
##
## timeOnly.Targ.gamm: position ~ 1 + s(Time) + s(Time, pp, by = vowelContrast, bs = "fs") +
## s(Time, item, by = vowelContrast, bs = "fs")
##
## Model timeOnly.Targ.gamm preferred: lower fREML score (137.133), and lower df (1.000).
## -----
## Model Score Edf Difference Df
## 1 minusAccent.Targ.gamm 210002.8 14
## 2 timeOnly.Targ.gamm 209865.7 13 -137.133 1.000
##
## AIC difference: 613.35, model timeOnly.Targ.gamm has lower AIC.COMPETITOR FIXATIONS
Setting up the environment and then Getting all the mouse positions as a matrix
rm(list = ls())
library(lmerTest)
library(stringr)
library(tidyverse)
library(mgcv)
library(itsadug)
options(warn = 0)
gazeData = read.table("dataCatalan/gazeCompet.txt", T, sep = "\t", quote = "")
trialData = read.table("dataCatalan/trialInfo.txt", T, sep = "\t", quote = "\"")
expFilter = grepl("_", trialData$condition)
gazeDataE = gazeData %>% filter(expFilter)
trialDataE = trialData %>% filter(expFilter)
timeline = -200 + (0:120) * 10
colnames(gazeDataE) = paste0("T_",timeline)
#contrast coding
temp = str_split_fixed(trialDataE$condition,"_",2)
trialDataE$accent = temp[, 1]
trialDataE$contrast = temp[, 2]
trialDataE$accentC = -0.5 + as.numeric(trialDataE$accent == "SEng")
table(trialDataE$accentC)##
## -0.5 0.5
## 3596 3598trialDataE$contrastC = -0.5 + as.numeric(trialDataE$contrast == "DRESS" |trialDataE$contrast == "TRAP" )
table(trialDataE$contrastC)##
## -0.5 0.5
## 3598 3596useData = cbind(trialDataE,gazeDataE)Analysis of the time course using GAMMs
We first need to reorganize the data.
gamData = useData[ ,c(1,2,4:8,10:13,14:134)] %>%
pivot_longer(., cols = starts_with("T_"),
names_to = "Time",
values_to = "position") %>% as.data.frame
gamData$position = as.numeric(gamData$position)
gamData$Time = gsub("T_","", gamData$Time) %>% as.numeric()
gamData = gamData[ order(gamData$pp,gamData$trial,gamData$Time), ]
gamData$pp = as.factor(gamData$pp)
gamData$item = as.factor(gamData$item)
gamData$isAustralian = ifelse(gamData$accent == "AusE",1,0)
gamData$isAustralianDress = gamData$isAustralian *
as.numeric(gamData$contrast == "DRESS" |gamData$contrast == "TRAP" )
gamData$vowelContrast = ifelse(gamData$contrast == "DRESS" |gamData$contrast == "TRAP",
"midFront","tensenessI") %>% as.factor()
#restricting the data to 0-900 after vowel onset
gamData = gamData %>% filter(Time > 0 & Time < 900)
gamData <- start_event(gamData, event=c("pp","trial"))
gamDataPlot = gamData %>% group_by(vowelContrast, Time, accent) %>%
summarize(meanFix = mean(position))## `summarise()` has regrouped the output.
## ℹ Summaries were computed grouped by vowelContrast, Time, and accent.
## ℹ Output is grouped by vowelContrast and Time.
## ℹ Use `summarise(.groups = "drop_last")` to silence this message.
## ℹ Use `summarise(.by = c(vowelContrast, Time, accent))` for per-operation
## grouping (`?dplyr::dplyr_by`) instead.ggplot(gamDataPlot, aes(x = Time, y = meanFix, col = vowelContrast)) +
geom_line() + facet_wrap(~accent)
This is the descriptive data
Building an initial GAMM
fullModelFile = "competitor_overall_fullModels.Rdata"
if (file.exists(fullModelFile))
{
load(fullModelFile)
}else{
interaction.Comp.gamm0 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
s(Time, by = isAustralian) +
s(Time, by = isAustralianDress) +
s(Time, pp, by = vowelContrast, bs="fs") +
s(Time, item, by = vowelContrast, bs="fs"),
family = "binomial",
data=gamData, discrete = T, nThreads = 3)
#dealing with autocorrelation
myAutoCorr =acf_resid(interaction.Comp.gamm0)
myRhoVal = myAutoCorr[2]
interaction.Comp.gamm1 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast, k=12) +
s(Time, by = isAustralian, k= 12) +
s(Time, by = isAustralianDress) +
s(Time, pp, by = vowelContrast, bs="fs") +
s(Time, item, by = vowelContrast, bs="fs"),
family = "binomial",
data=gamData, rho = myRhoVal,
AR.start = gamData$start.event, nThreads = 3, discrete = T)
save(list= c("interaction.Comp.gamm0","interaction.Comp.gamm1"), file = fullModelFile)
}Model criticism
How are we doing regarding auto-correlation?
myAutoCorr =acf_resid(interaction.Comp.gamm1)
Is the default for the smooth (k = 10) sufficient?
gam.check(interaction.Comp.gamm1)
##
## Method: fREML Optimizer: perf chol
## $grad
## [1] 4.344936e-09 -6.668310e-11 -3.617088e-03 -2.361695e-03 -1.437751e-08
## [6] 4.397094e-09 -1.403622e-07 -3.270612e-08 1.320329e-08 -1.223585e-07
## [11] -3.114113e-10 -9.537457e-09 -3.080776e-08 -5.261800e-09 -1.131491e-08
## [16] -3.280221e-09
##
## $hess
## [,1] [,2] [,3] [,4] [,5]
## [1,] 3.635382e+00 -1.376290e-04 8.619980e-04 8.586773e-04 1.138355e-01
## [2,] -1.376290e-04 3.880525e+00 7.466937e-04 -1.558598e-06 -4.275009e-06
## [3,] 8.619980e-04 7.466937e-04 3.597892e-03 -4.456656e-06 3.025956e-05
## [4,] 8.586773e-04 -1.558598e-06 -4.456656e-06 2.349919e-03 3.030709e-05
## [5,] 1.138355e-01 -4.275009e-06 3.025956e-05 3.030709e-05 1.292494e+01
## [6,] -3.700456e-02 2.276860e-06 -1.621932e-05 -1.369423e-05 -9.831543e-01
## [7,] -7.906943e-03 2.875409e-07 -2.106202e-06 -2.068759e-06 1.255385e-01
## [8,] -1.602161e-06 3.773216e-02 6.891158e-06 -8.127057e-09 -7.072587e-08
## [9,] 1.045072e-06 -1.702599e-02 -6.433961e-06 3.504716e-08 6.856727e-08
## [10,] 4.089598e-07 -1.273209e-02 -2.499181e-06 7.974350e-09 2.095914e-08
## [11,] 2.629652e-02 -5.265767e-07 3.552015e-06 4.319390e-06 -8.710586e-03
## [12,] -5.603244e-03 1.024203e-06 -7.201901e-06 -5.408272e-06 -1.884878e-02
## [13,] -7.166771e-03 2.895458e-07 -1.989687e-06 -1.841388e-06 -1.822608e-02
## [14,] -2.529303e-06 5.614114e-02 1.268922e-05 -2.578752e-08 -1.389110e-07
## [15,] 6.679473e-07 -7.128001e-03 -4.334080e-06 1.994009e-08 3.573136e-08
## [16,] 2.588168e-07 -1.940266e-02 -1.579194e-06 -1.970647e-09 1.683026e-09
## [,6] [,7] [,8] [,9] [,10]
## [1,] -3.700456e-02 -7.906943e-03 -1.602161e-06 1.045072e-06 4.089598e-07
## [2,] 2.276860e-06 2.875409e-07 3.773216e-02 -1.702599e-02 -1.273209e-02
## [3,] -1.621932e-05 -2.106202e-06 6.891158e-06 -6.433961e-06 -2.499181e-06
## [4,] -1.369423e-05 -2.068759e-06 -8.127057e-09 3.504716e-08 7.974350e-09
## [5,] -9.831543e-01 1.255385e-01 -7.072587e-08 6.856727e-08 2.095914e-08
## [6,] 1.261535e+01 -9.316762e-01 1.032658e-08 -7.923694e-08 -1.505821e-08
## [7,] -9.316762e-01 2.657327e+00 2.262021e-09 -6.058032e-09 -1.609988e-09
## [8,] 1.032658e-08 2.262021e-09 1.015215e+01 -1.013432e+00 -8.810777e-02
## [9,] -7.923694e-08 -6.058032e-09 -1.013432e+00 1.368019e+01 -1.156318e+00
## [10,] -1.505821e-08 -1.609988e-09 -8.810777e-02 -1.156318e+00 5.072175e+00
## [11,] -9.425580e-03 -2.471591e-03 -1.264574e-08 -4.784365e-09 9.058073e-10
## [12,] -3.369710e-02 -4.802659e-03 6.349869e-11 -4.345179e-08 -7.465993e-09
## [13,] -1.632764e-02 -6.561908e-03 1.299653e-09 -5.977144e-09 -1.389967e-09
## [14,] 4.858758e-08 6.524317e-09 -4.271093e-02 -2.416663e-02 -1.883049e-02
## [15,] -4.188182e-08 -3.449844e-09 -1.691265e-02 -3.940526e-02 -1.346172e-02
## [16,] 9.881644e-09 1.251340e-10 -4.067646e-02 -2.646676e-02 -1.665774e-02
## [,11] [,12] [,13] [,14] [,15]
## [1,] 2.629652e-02 -5.603244e-03 -7.166771e-03 -2.529303e-06 6.679473e-07
## [2,] -5.265767e-07 1.024203e-06 2.895458e-07 5.614114e-02 -7.128001e-03
## [3,] 3.552015e-06 -7.201901e-06 -1.989687e-06 1.268922e-05 -4.334080e-06
## [4,] 4.319390e-06 -5.408272e-06 -1.841388e-06 -2.578752e-08 1.994009e-08
## [5,] -8.710586e-03 -1.884878e-02 -1.822608e-02 -1.389110e-07 3.573136e-08
## [6,] -9.425580e-03 -3.369710e-02 -1.632764e-02 4.858758e-08 -4.188182e-08
## [7,] -2.471591e-03 -4.802659e-03 -6.561908e-03 6.524317e-09 -3.449844e-09
## [8,] -1.264574e-08 6.349865e-11 1.299653e-09 -4.271093e-02 -1.691265e-02
## [9,] -4.784365e-09 -4.345179e-08 -5.977144e-09 -2.416663e-02 -3.940526e-02
## [10,] 9.058073e-10 -7.465993e-09 -1.389967e-09 -1.883049e-02 -1.346172e-02
## [11,] 1.097995e+00 -2.705718e-01 5.232447e-02 -1.882327e-08 -2.334746e-09
## [12,] -2.705718e-01 5.151456e+00 -6.577427e-01 1.434919e-08 -2.315007e-08
## [13,] 5.232447e-02 -6.577427e-01 3.283239e+00 3.786390e-09 -3.514842e-09
## [14,] -1.882327e-08 1.434919e-08 3.786390e-09 6.572665e+00 -3.769508e-03
## [15,] -2.334746e-09 -2.315007e-08 -3.514842e-09 -3.769508e-03 4.759903e+00
## [16,] 3.295020e-09 6.427254e-09 2.429660e-10 1.597832e-01 -4.868616e-01
## [,16]
## [1,] 2.588168e-07
## [2,] -1.940266e-02
## [3,] -1.579194e-06
## [4,] -1.970647e-09
## [5,] 1.683026e-09
## [6,] 9.881644e-09
## [7,] 1.251340e-10
## [8,] -4.067646e-02
## [9,] -2.646676e-02
## [10,] -1.665774e-02
## [11,] 3.295020e-09
## [12,] 6.427254e-09
## [13,] 2.429660e-10
## [14,] 1.597832e-01
## [15,] -4.868616e-01
## [16,] 5.361220e+00
##
## Model rank = 1746 / 1746
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(Time):vowelContrastmidFront 11.00 8.11 0.94 0.115
## s(Time):vowelContrasttensenessI 11.00 8.79 0.94 0.190
## s(Time):isAustralian 12.00 2.01 0.94 0.105
## s(Time):isAustralianDress 10.00 2.01 0.94 0.120
## s(Time,pp):vowelContrastmidFront 450.00 121.50 0.94 0.075 .
## s(Time,pp):vowelContrasttensenessI 450.00 120.12 0.94 0.165
## s(Time,item):vowelContrastmidFront 400.00 35.44 0.94 0.110
## s(Time,item):vowelContrasttensenessI 400.00 61.45 0.94 0.095 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Let’s have a look at this model no whether there is an overall effect of contrast
summary(interaction.Comp.gamm1)##
## Family: binomial
## Link function: logit
##
## Formula:
## position ~ vowelContrast + s(Time, by = vowelContrast, k = 12) +
## s(Time, by = isAustralian, k = 12) + s(Time, by = isAustralianDress) +
## s(Time, pp, by = vowelContrast, bs = "fs") + s(Time, item,
## by = vowelContrast, bs = "fs")
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.61502 0.05197 -31.076 <2e-16 ***
## vowelContrasttensenessI 0.06879 0.08114 0.848 0.397
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(Time):vowelContrastmidFront 8.114 9.673 168.256 <2e-16 ***
## s(Time):vowelContrasttensenessI 8.790 10.120 229.128 <2e-16 ***
## s(Time):isAustralian 2.013 2.025 0.704 0.772
## s(Time):isAustralianDress 2.009 2.018 0.327 0.880
## s(Time,pp):vowelContrastmidFront 121.499 449.000 352.710 <2e-16 ***
## s(Time,pp):vowelContrasttensenessI 120.124 448.000 354.939 <2e-16 ***
## s(Time,item):vowelContrastmidFront 35.440 199.000 146.981 <2e-16 ***
## s(Time,item):vowelContrasttensenessI 61.451 198.000 209.584 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.0731 Deviance explained = 8.17%
## fREML = 1.6334e+05 Scale est. = 1 n = 640266There is a main effect of competitor, meaning that, overall, the mouse position is closer to target when the competitor is dissimilar to the target. Note that the significant effects of the tensors only indicate changes over time, they do not constitute evidence for a significant difference between conditions.
smooth <- plot_smooth(interaction.Comp.gamm1, view = "Time",
plot_all = "vowelContrast",
rug = FALSE, rm.ranef = TRUE, lwd = 4)## Summary:
## * vowelContrast : factor; set to the value(s): midFront, tensenessI.
## * Time : numeric predictor; with 30 values ranging from 10.000000 to 890.000000.
## * isAustralian : numeric predictor; set to the value(s): 0.
## * isAustralianDress : numeric predictor; set to the value(s): 0.
## * pp : factor; set to the value(s): CAT01. (Might be canceled as random effect, check below.)
## * item : factor; set to the value(s): baggage. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(Time,pp):vowelContrastmidFront,s(Time,pp):vowelContrasttensenessI,s(Time,item):vowelContrastmidFront,s(Time,item):vowelContrasttensenessI
## 
Now model comparison
We succesively take effects out to prepare model comparisons.
reducedModelFile = "competitor_overall_reducedModels.Rdata"
if (file.exists(reducedModelFile))
{
load(reducedModelFile)
}else{
minusInteraction.Comp.gamm <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
s(Time, by = isAustralian) +
#s(Time, by = isAustralianDress) +
s(Time, pp, by = vowelContrast, bs="fs") +
s(Time, item, by = vowelContrast, bs="fs"),
family = "binomial",
data=gamData, rho = myRhoVal,
AR.start = gamData$start.event, nThreads = 3, discrete = T)
minusAccent.Comp.gamm <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
#s(Time, by = isAustralian) +
#s(Time, by = isAustralianDress) +
s(Time, pp, by = vowelContrast, bs="fs") +
s(Time, item, by = vowelContrast, bs="fs"),
family = "binomial",
data=gamData, rho = myRhoVal,
AR.start = gamData$start.event, nThreads = 3, discrete = T)
timeOnly.Comp.gamm <- bam( position ~ 1 +
s(Time) +
#s(Time, by = isAustralian) +
#s(Time, by = isAustralianDress) +
s(Time, pp, by = vowelContrast, bs="fs") +
s(Time, item, by = vowelContrast, bs="fs"),
family = "binomial",
data=gamData, rho = myRhoVal,
AR.start = gamData$start.event, nThreads = 3, discrete = T)
save(list= c("minusInteraction.Comp.gamm","minusAccent.Comp.gamm","timeOnly.Comp.gamm"), file = reducedModelFile)
}Model comparisons
Full model versus no interaction
compareML(interaction.Comp.gamm1, minusInteraction.Comp.gamm)## interaction.Comp.gamm1: position ~ vowelContrast + s(Time, by = vowelContrast, k = 12) +
## s(Time, by = isAustralian, k = 12) + s(Time, by = isAustralianDress) +
## s(Time, pp, by = vowelContrast, bs = "fs") + s(Time, item,
## by = vowelContrast, bs = "fs")
##
## minusInteraction.Comp.gamm: position ~ vowelContrast + s(Time, by = vowelContrast) + s(Time,
## by = isAustralian) + s(Time, pp, by = vowelContrast, bs = "fs") +
## s(Time, item, by = vowelContrast, bs = "fs")
##
## Model minusInteraction.Comp.gamm preferred: lower fREML score (111.479), and lower df (1.000).
## -----
## Model Score Edf Difference Df
## 1 interaction.Comp.gamm1 163337.1 16
## 2 minusInteraction.Comp.gamm 163225.6 15 -111.479 1.000
##
## AIC difference: 221.24, model minusInteraction.Comp.gamm has lower AIC.Both main effects versus minus accent
compareML(minusInteraction.Comp.gamm, minusAccent.Comp.gamm)## minusInteraction.Comp.gamm: position ~ vowelContrast + s(Time, by = vowelContrast) + s(Time,
## by = isAustralian) + s(Time, pp, by = vowelContrast, bs = "fs") +
## s(Time, item, by = vowelContrast, bs = "fs")
##
## minusAccent.Comp.gamm: position ~ vowelContrast + s(Time, by = vowelContrast) + s(Time,
## pp, by = vowelContrast, bs = "fs") + s(Time, item, by = vowelContrast,
## bs = "fs")
##
## Model minusAccent.Comp.gamm preferred: lower fREML score (14.380), and lower df (1.000).
## -----
## Model Score Edf Difference Df
## 1 minusInteraction.Comp.gamm 163225.6 15
## 2 minusAccent.Comp.gamm 163211.3 14 -14.380 1.000
##
## AIC difference: 374.14, model minusAccent.Comp.gamm has lower AIC.Main effect of Contrast versus time-only
compareML(minusAccent.Comp.gamm, timeOnly.Comp.gamm)## minusAccent.Comp.gamm: position ~ vowelContrast + s(Time, by = vowelContrast) + s(Time,
## pp, by = vowelContrast, bs = "fs") + s(Time, item, by = vowelContrast,
## bs = "fs")
##
## timeOnly.Comp.gamm: position ~ 1 + s(Time) + s(Time, pp, by = vowelContrast, bs = "fs") +
## s(Time, item, by = vowelContrast, bs = "fs")
##
## Model timeOnly.Comp.gamm preferred: lower fREML score (436.489), and lower df (1.000).
## -----
## Model Score Edf Difference Df
## 1 minusAccent.Comp.gamm 163211.3 14
## 2 timeOnly.Comp.gamm 162774.8 13 -436.489 1.000
##
## AIC difference: 1703.05, model timeOnly.Comp.gamm has lower AIC.FILLER ITEMS
Target fixations
library(tidyverse)
library(mgcv)
library(itsadug)
library(stringr)
options(warn = 0)
# Load data
gazeTarget <- read.table("dataCatalan/gazeTarget.txt",
header = TRUE, sep = "\t", quote = "")
trialData <- read.table("dataCatalan/trialInfo.txt",
header = TRUE, sep = "\t", quote = "\"")
# Timeline
timeline <- -200 + (0:120) * 10
colnames(gazeTarget) <- paste0("T_", timeline)
# Combine datasets (same style as original)
useData <- cbind(trialData, gazeTarget) %>% as.data.frame()
# Remove practice trials
useData <- useData %>%
filter(useData$condition != "practice") %>%
droplevels()
# Rename levels of condition
accent <- ifelse(grepl("^SEng", useData$wavfile), "SEng",
ifelse(grepl("^AusE", useData$wavfile), "AusE", NA))
vowel <- ifelse(grepl("oo", useData$item, ignore.case = TRUE),
"GOOSE", "LOT")
useData$condition <- ifelse(
useData$condition == "filler" & !is.na(accent),
paste(accent, vowel, sep = "_"),
useData$condition
)
# Add accent and vowel variables
useData$accent <- sub("_.*", "", useData$condition)
useData$vowel <- sub(".*_", "", useData$condition)gamData = useData[ ,c(1:132)] %>%
pivot_longer(., cols = starts_with("T_"),
names_to = "Time",
values_to = "position") %>% as.data.frame
gamData$position = as.numeric(gamData$position)
gamData$Time = gsub("T_","", gamData$Time) %>% as.numeric()
# Sort (like original)
gamData <- gamData[order(gamData$pp, gamData$trial, gamData$Time), ]
# Factors (like original)
gamData$pp <- as.factor(gamData$pp)
gamData$item <- as.factor(gamData$item)
# Add variable Accent
# Time window
gamData <- gamData %>% filter(Time > 0 & Time < 900)
gamData <- start_event(gamData, event = c("pp","trial"))gamDataPlot = gamData %>% group_by(accent, vowel, Time) %>%
summarize(meanFix = mean(position))## `summarise()` has regrouped the output.
## ℹ Summaries were computed grouped by accent, vowel, and Time.
## ℹ Output is grouped by accent and vowel.
## ℹ Use `summarise(.groups = "drop_last")` to silence this message.
## ℹ Use `summarise(.by = c(accent, vowel, Time))` for per-operation grouping
## (`?dplyr::dplyr_by`) instead.ggplot(gamDataPlot, aes(x = Time, y = meanFix, col = vowel)) +
geom_line() + facet_wrap(. ~ accent)
Across all three conditions, looks to the target increase over time while looks to the competitor decrease in a broadly similar pattern. The two accents (solid vs. dotted lines) also show very similar time courses, with only minor divergence between them within each condition.
# Create variable block
gamData$block <- NA
gamData$block[gamData$trial >= 7] <-
((gamData$trial[gamData$trial >= 7] - 7) %/% 50) + 1
gamDataPlot = gamData %>% group_by(block, accent, vowel, Time) %>%
summarize(meanFix = mean(position))## `summarise()` has regrouped the output.
## ℹ Summaries were computed grouped by block, accent, vowel, and Time.
## ℹ Output is grouped by block, accent, and vowel.
## ℹ Use `summarise(.groups = "drop_last")` to silence this message.
## ℹ Use `summarise(.by = c(block, accent, vowel, Time))` for per-operation
## grouping (`?dplyr::dplyr_by`) instead.ggplot(gamDataPlot, aes(x = Time, y = meanFix, col = vowel)) +
geom_line() + facet_wrap(accent ~ block, nrow = 2)
Competitor fixations
library(tidyverse)
library(mgcv)
library(itsadug)
library(stringr)
options(warn = 0)
# Load data
gazeCompet <- read.table("dataCatalan/gazeCompet.txt",
header = TRUE, sep = "\t", quote = "")
trialData <- read.table("dataCatalan/trialInfo.txt",
header = TRUE, sep = "\t", quote = "\"")
# Timeline
timeline <- -200 + (0:120) * 10
colnames(gazeCompet) <- paste0("T_", timeline)
# Combine datasets (same style as original)
useData <- cbind(trialData, gazeCompet) %>% as.data.frame()
# Remove practice trials
useData <- useData %>%
filter(useData$condition != "practice") %>%
droplevels()
# Rename levels of condition
accent <- ifelse(grepl("^SEng", useData$wavfile), "SEng",
ifelse(grepl("^AusE", useData$wavfile), "AusE", NA))
vowel <- ifelse(grepl("oo", useData$item, ignore.case = TRUE),
"GOOSE", "LOT")
useData$condition <- ifelse(
useData$condition == "filler" & !is.na(accent),
paste(accent, vowel, sep = "_"),
useData$condition
)
# Add accent and vowel variables
useData$accent <- sub("_.*", "", useData$condition)
useData$vowel <- sub(".*_", "", useData$condition)gamData = useData[ ,c(1,2,4:132)] %>%
pivot_longer(., cols = starts_with("T_"),
names_to = "Time",
values_to = "position") %>% as.data.frame
gamData$position = as.numeric(gamData$position)
gamData$Time = gsub("T_","", gamData$Time) %>% as.numeric()
# Sort (like original)
gamData <- gamData[order(gamData$pp, gamData$trial, gamData$Time), ]
# Factors (like original)
gamData$pp <- as.factor(gamData$pp)
gamData$item <- as.factor(gamData$item)
# Time window
gamData <- gamData %>% filter(Time > 0 & Time < 900)
gamData <- start_event(gamData, event = c("pp","trial"))gamDataPlot = gamData %>% group_by(accent, vowel, Time) %>%
summarize(meanFix = mean(position))## `summarise()` has regrouped the output.
## ℹ Summaries were computed grouped by accent, vowel, and Time.
## ℹ Output is grouped by accent and vowel.
## ℹ Use `summarise(.groups = "drop_last")` to silence this message.
## ℹ Use `summarise(.by = c(accent, vowel, Time))` for per-operation grouping
## (`?dplyr::dplyr_by`) instead.ggplot(gamDataPlot, aes(x = Time, y = meanFix, col = vowel)) +
geom_line() + facet_wrap(. ~ accent)
Across all three conditions, looks to the target increase over time while looks to the competitor decrease in a broadly similar pattern. The two accents (solid vs. dotted lines) also show very similar time courses, with only minor divergence between them within each condition.
gamData <- gamData %>%
filter(gamData$vowel != "GOOSE") %>%
droplevels()
gamData <- gamData %>%
filter(gamData$vowel != "LOT") %>%
droplevels()
# Create variable block
gamData$block <- NA
gamData$block[gamData$trial >= 7] <-
((gamData$trial[gamData$trial >= 7] - 7) %/% 50) + 1
gamDataPlot = gamData %>% group_by(block, accent, vowel, Time) %>%
summarize(meanFix = mean(position))## `summarise()` has regrouped the output.
## ℹ Summaries were computed grouped by block, accent, vowel, and Time.
## ℹ Output is grouped by block, accent, and vowel.
## ℹ Use `summarise(.groups = "drop_last")` to silence this message.
## ℹ Use `summarise(.by = c(block, accent, vowel, Time))` for per-operation
## grouping (`?dplyr::dplyr_by`) instead.ggplot(gamDataPlot, aes(x = Time, y = meanFix, col = vowel)) +
geom_line() + facet_wrap(accent ~ block, nrow = 2)
RESPONSE TIMES
library(tidyverse)
library(mgcv)
library(itsadug)
library(stringr)
options(warn = 0)
# Load data
gazeTarget <- read.table("dataCatalan/gazeTarget.txt",
header = TRUE, sep = "\t", quote = "")
trialData <- read.table("dataCatalan/trialInfo.txt",
header = TRUE, sep = "\t", quote = "\"")
# Timeline
timeline <- -200 + (0:120) * 10
colnames(gazeCompet) <- paste0("T_", timeline)
# Combine datasets (same style as original)
useData <- cbind(trialData, gazeCompet) %>% as.data.frame()
# Remove practice trials
useData <- useData %>%
filter(useData$condition != "practice") %>%
droplevels()
# Rename levels of condition
accent <- ifelse(grepl("^SEng", useData$wavfile), "SEng",
ifelse(grepl("^AusE", useData$wavfile), "AusE", NA))
vowel <- ifelse(grepl("oo", useData$item, ignore.case = TRUE),
"GOOSE", "LOT")
useData$condition <- ifelse(
useData$condition == "filler" & !is.na(accent),
paste(accent, vowel, sep = "_"),
useData$condition
)
# Add accent and vowel variables
useData$accent <- sub("_.*", "", useData$condition)
useData$vowel <- sub(".*_", "", useData$condition)ggplot(useData, aes(x = vowel, y = rt)) +
geom_boxplot() + facet_wrap(.~ accent)
STIMULUS MATERIALS
library(dplyr)
library(ggplot2)
library(ggrepel)
Stimuli_Acoustics <- read.csv("WP3_Stimuli_Acoustics.csv", row.names = NULL, stringsAsFactors = TRUE)
stimuli_plot <- Stimuli_Acoustics %>%
mutate(
Accent = factor(Accent, levels = c("S.Eng", "AusE")),
Vowel = factor(toupper(Vowel))
)
stimuli_plot <- stimuli_plot %>%
filter(stimuli_plot$Contrast != "GOOSE-LOT") %>%
droplevels()
means_va <- stimuli_plot %>%
group_by(Accent, Vowel) %>%
summarise(
F1_mean = mean(F1, na.rm = TRUE),
F2_mean = mean(F2, na.rm = TRUE),
.groups = "drop"
)
Fig.Stimuli <- ggplot(stimuli_plot, aes(x = F2, y = F1, colour = Accent)) +
stat_ellipse(aes(group = interaction(Accent, Vowel)),
type = "norm",
level = 0.95,
linewidth = 0.6,
alpha = 0.5) +
geom_point(data = means_va,
aes(x = F2_mean, y = F1_mean),
size = 2) +
geom_text_repel(data = means_va,
aes(x = F2_mean, y = F1_mean, label = Vowel),
show.legend = FALSE,
size = 3,
point.padding = 0.3,
box.padding = 0.4,
min.segment.length = 0,
segment.alpha = 0.35) +
scale_colour_manual(values = c(
"S.Eng" = "#E41A1C",
"AusE" = "#377EB8"
)) +
scale_x_reverse() +
scale_y_reverse() +
labs(x = "F2 (Hz)", y = "F1 (Hz)", colour = "Accent") +
theme_grey() +
facet_wrap(.~ Accent) +
theme(
strip.background = element_blank(),
strip.text = element_blank(),
legend.position = "bottom",
legend.title = element_text(size = 9, face = "bold"))
plot(Fig.Stimuli)