Assessing regional accent barriers for L2 spoken word recognition
Holger Mitterer & Gil Verbeke
2026-03-04
0. Overview
1. Participants
2. Stimulus Materials
2.1. (Lexical)
Properties
2.2. Acoustic Properties
3. Eye-Tracking Data
3.1. Visualization
3.2. Analysis of the Time Course
3.2.1.
Overall Effect of Contrast (FLEECE-KIT vs. DRESS-TRAP)
3.2.2.
FLEECE-KIT Contrast
3.2.3. DRESS-TRAP Contrast
1. Participants
Load datasets and remove data from three participants: (i)
CAT14: born in Dakar, Senegal; (ii) CAT17: .asc-file
with gaze patterns corrupted; (iii) CAT43: L1 Urdu speaker.
Questionnaire <- read.csv("WP3_Main_Questionnaire.csv", row.names = NULL, stringsAsFactors = TRUE)
Familiarity <- read.csv("WP3_Main_Familiarity.csv", row.names = NULL, stringsAsFactors = TRUE)
# Merge datasets by Participant
Questionnaire <- Questionnaire[!(Questionnaire$Part_ID %in% c("CAT14", "CAT17", "CAT43")),]
Questionnaire$Part_ID <- droplevels(Questionnaire$Part_ID)
Familiarity <- Familiarity[!(Familiarity$Part_ID %in% c("CAT14", "CAT17", "CAT43")),]
Familiarity$Part_ID <- droplevels(Familiarity$Part_ID)# Sample Size
str(Questionnaire$Part_ID)## Factor w/ 45 levels " CAT21"," CAT36",..: 3 4 5 6 7 8 9 10 11 12 ...# Age
sapply(list(mean = mean, sd = sd, min = min, max = max), function(f) f(Questionnaire$Age))## mean sd min max
## 19.11111 1.41778 18.00000 24.00000# Gender
table(Questionnaire$Gender)##
## Female Male Non-binary / third gender
## 35 9 1# Language Dominance
table(Questionnaire$LangDomParticipant)##
## Both equally Catalan Spanish
## 0 18 18 9# Programme
table(Questionnaire$FieldStudy)##
## Communication Sciences Engineering
## 1 5
## Linguistics Sea Sciences
## 21 1
## Translation and interpretation
## 17# Age of Onset of Acquisition
sapply(list(mean = mean, sd = sd, min = min, max = max), function(f) f(Questionnaire$AgeOnsetEnglish))## mean sd min max
## 4.788889 1.447917 2.000000 8.000000# Self-reported proficiency - speaking
sapply(list(mean = mean, sd = sd, min = min, max = max), function(f) f(Questionnaire$Prof_Speaking))## mean sd min max
## 5.1777778 0.9837056 3.0000000 7.0000000# Self-reported proficiency - listening
sapply(list(mean = mean, sd = sd, min = min, max = max), function(f) f(Questionnaire$Prof_Listening))## mean sd min max
## 6.0000000 0.7687061 4.0000000 7.0000000# Self-reported proficiency - reading
sapply(list(mean = mean, sd = sd, min = min, max = max), function(f) f(Questionnaire$Prof_Reading))## mean sd min max
## 5.977778 1.097288 3.000000 7.000000# Self-reported proficiency - writing
sapply(list(mean = mean, sd = sd, min = min, max = max), function(f) f(Questionnaire$Prof_Writing))## mean sd min max
## 5.3111111 0.9728641 3.0000000 7.0000000# Self-reported proficiency - overall
sapply(list(mean = mean, sd = sd, min = min, max = max), function(f) f(Questionnaire$Prof_Overall))## mean sd min max
## 5.4666667 0.8146388 4.0000000 7.0000000# LexTALE
sapply(list(mean = mean, sd = sd, min = min, max = max), function(f) f(Questionnaire$LexTALE))## mean sd min max
## 72.111111 9.511589 55.000000 92.500000# Familiarity
library(dplyr)
library(knitr)
Table1 <- Familiarity %>%
group_by(Accent) %>%
summarise(
Mean_Familiarity = mean(Familiarity, na.rm = TRUE),
SD_Familiarity = sd(Familiarity, na.rm = TRUE)
)
kable(Table1)| Accent | Mean_Familiarity | SD_Familiarity |
|---|---|---|
| A.Eng | 51.333333 | 19.522714 |
| AusE | 3.911111 | 5.455698 |
| CanE | 2.622222 | 4.608731 |
| N.Eng | 3.711111 | 5.092875 |
| Other | 4.955556 | 15.265760 |
| S.Eng | 29.577778 | 15.980797 |
| SA.Eng | 3.888889 | 5.296749 |
2. Stimulus Materials
2.1 (Lexical) Properties
library(dplyr)
library(ggplot2)
library(ggrepel)
library(knitr)
Stimuli <- read.csv("WP3_Stimuli.csv", row.names = NULL, stringsAsFactors = TRUE)
contrast_labels <- c(
"/i/-/I/" = "FLEECE-KIT",
"/e/-/a/" = "DRESS-TRAP",
"/u/-/o/" = "GOOSE-LOT"
)
fmt_m_sd <- function(x, digits = 2) {
sprintf(paste0("%.", digits, "f<br>(%.", digits, "f)"),
mean(x, na.rm = TRUE),
sd(x, na.rm = TRUE))
}Table2 <- Stimuli %>%
mutate(Contrast = recode(as.character(Contrast), !!!contrast_labels, .default = as.character(Contrast))) %>%
group_by(Contrast) %>%
summarise(
Freq1 = fmt_m_sd(LexFreq_1),
Freq2 = fmt_m_sd(LexFreq_2),
`ΔFreq` = fmt_m_sd(Delta_LexFreq),
PND1 = fmt_m_sd(PND_1),
PND2 = fmt_m_sd(PND_2),
`ΔPND` = fmt_m_sd(Delta_PND),
Fam1 = fmt_m_sd(Familiarity_1),
Fam2 = fmt_m_sd(Familiarity_2),
`ΔFam` = fmt_m_sd(Delta_Familiarity),
.groups = "drop"
) %>%
arrange(factor(Contrast, levels = c("FLEECE-KIT", "DRESS-TRAP", "GOOSE-LOT")))
kable(Table2)| Contrast | Freq1 | Freq2 | ΔFreq | PND1 | PND2 | ΔPND | Fam1 | Fam2 | ΔFam |
|---|---|---|---|---|---|---|---|---|---|
| FLEECE-KIT | 3.62 (0.70) |
3.39 (0.56) |
0.68 (0.44) |
2.50 (3.27) |
1.80 (3.12) |
1.70 (1.95) |
27.70 (3.13) |
26.50 (2.64) |
1.40 (1.51) |
| DRESS-TRAP | 3.85 (0.65) |
3.90 (0.50) |
0.56 (0.49) |
2.80 (3.91) |
4.20 (3.97) |
3.80 (3.49) |
28.30 (2.50) |
28.30 (2.11) |
2.40 (2.41) |
| GOOSE-LOT | 3.28 (0.44) |
3.93 (0.45) |
0.78 (0.46) |
2.50 (3.87) |
2.60 (3.06) |
2.50 (2.95) |
28.00 (2.00) |
29.00 (1.05) |
1.80 (1.69) |
Table 1 shows the mean values and SDS for lexical frequency (Freq), phonological neighborhood density (PND) and familiarity ratings (Fam) for the item pairs used across the three vowel contrasts. Subscripts 1-2 denote the first and second member of each pair, respectively. Difference scores (ΔFreq, ΔPND, ΔFam) represent the absolute difference between the two words within each pair, averaged across all pairs in the contrast.
2.2 Acoustic Properties
Create an F1-F2 vowel plot showing FLEECE, KIT, DRESS and TRAP for Southern British English (S.Eng) and Australian English (AusE).
library(dplyr)
library(ggplot2)
library(ggrepel)
Stimuli_Acoustics <- read.csv("WP3_Stimuli_Acoustics.csv", row.names = NULL, stringsAsFactors = TRUE)
stimuli_F1_F2 <- Stimuli_Acoustics %>%
mutate(
Accent = factor(Accent, levels = c("S.Eng", "AusE")),
Vowel = factor(toupper(Vowel))
)
#stimuli_F1_F2 <- stimuli_F1_F2 %>%
# filter(stimuli_F1_F2$Contrast != "GOOSE-LOT") %>%
# droplevels()
means_va <- stimuli_F1_F2 %>%
group_by(Accent, Vowel) %>%
summarise(
F1_mean = mean(F1, na.rm = TRUE),
F2_mean = mean(F2, na.rm = TRUE),
.groups = "drop"
)Fig.Stimuli.F1_F2 <- ggplot(stimuli_F1_F2, 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.F1_F2)
Create a bar plot showing vowel duration by accent.
library(ggplot2)
library(plyr)
stimuli_duration <- Stimuli_Acoustics %>%
mutate(
Accent = factor(Accent, levels = c("S.Eng", "AusE")),
Vowel = factor(toupper(Vowel))
)
stimuli_duration <- stimuli_duration %>%
# filter(stimuli_duration$Contrast != "GOOSE-LOT") %>%
mutate(Vowel_Duration = Vowel_Duration * 1000) %>%
droplevels()
stimuli_duration$Vowel <- factor(
stimuli_duration$Vowel,
levels = c("FLEECE", "KIT", "DRESS", "TRAP", "GOOSE", "LOT")
)
data_summary <- function(data, varname, groupnames){
require(plyr)
summary_func <- function(x, col){
c(mean = mean(x[[col]], na.rm=TRUE),
sd = sd(x[[col]], na.rm=TRUE))
}
data_sum<-ddply(data, groupnames, .fun=summary_func,
varname)
data_sum <- rename(data_sum, c("mean" = varname))
return(data_sum)
}
stimuli_duration <- data_summary(stimuli_duration, varname="Vowel_Duration",
groupnames=c("Vowel", "Accent"))Fig.Stimuli.Duration <- ggplot(stimuli_duration, aes(x = Vowel, y = Vowel_Duration, fill = Accent)) +
geom_bar(stat = "identity", color = "black",
position = position_dodge()) +
geom_errorbar(aes(ymin = Vowel_Duration - sd,
ymax = Vowel_Duration + sd),
width = .2,
position = position_dodge(.9)) +
scale_fill_manual(values = c(
"S.Eng" = "#E41A1C",
"AusE" = "#377EB8"
)) +
labs(y = "Duration (ms)") +
theme_grey() +
theme(
strip.background = element_blank(),
strip.text = element_blank(),
legend.position = "bottom",
legend.title = element_text(size = 9, face = "bold"))
plot(Fig.Stimuli.Duration)
library(ggpubr)##
## Attaching package: 'ggpubr'## The following object is masked from 'package:plyr':
##
## mutateFig.Stimuli <- ggarrange(Fig.Stimuli.F1_F2, Fig.Stimuli.Duration,
ncol = 2,
labels = c("A", "B"),
widths = c(3, 2),
common.legend = TRUE, legend = "bottom")
plot(Fig.Stimuli)
ggsave("Stimuli.png",
width = 8.5, height = 3.5, dpi = 600,
bg = "white")detach("package:plyr", unload=TRUE)3. Eye-Tracking Data
3.1 Visualization
rm(list = ls())
library(tidyverse)
library(mgcv)
library(itsadug)
library(stringr)
options(warn = 0)
# Load data
gazeTarget <- read.table("dataCatalan/gazeTarget.txt",
header = TRUE, sep = "\t", quote = "")
gazeComp <- read.table("dataCatalan/gazeCompet.txt",
header = TRUE, sep = "\t", quote = "")
gazeDistr1 <- read.table("dataCatalan/gazeDistract1.txt",
header = TRUE, sep = "\t", quote = "")
gazeDistr2 <- read.table("dataCatalan/gazeDistract2.txt",
header = TRUE, sep = "\t", quote = "")
trialData <- read.table("dataCatalan/trialInfo.txt",
header = TRUE, sep = "\t", quote = "\"")
# Remove fillers
expFilter <- grepl("_", trialData$condition)
gazeTargetE <- gazeTarget %>% filter(expFilter)
gazeCompE <- gazeComp %>% filter(expFilter)
gazeDistr1E <- gazeDistr1 %>% filter(expFilter)
gazeDistr2E <- gazeDistr2 %>% filter(expFilter)
trialDataE <- trialData %>% filter(expFilter)
# Timeline
timeline <- -200 + (0:120) * 10
colnames(gazeTargetE) <- paste0("Tt_", timeline)
colnames(gazeCompE) <- paste0("Tc_", timeline)
colnames(gazeDistr1E) <- paste0("Td1_", timeline)
colnames(gazeDistr2E) <- paste0("Td2_", timeline)
# Accent - Contrast
temp <- str_split_fixed(trialDataE$condition, "_", 2)
trialDataE$accent <- temp[,1]
trialDataE$contrast <- temp[,2]
# Combine datasets
useData <- cbind(trialDataE, gazeTargetE, gazeCompE, gazeDistr1E, gazeDistr2E) %>% as.data.frame()
# Remove distractor responses
useData <- useData %>%
filter(useData$clickedOn != "3") %>%
droplevels()library(dplyr)
library(tidyr)
gamData <- useData %>%
pivot_longer(
cols = matches("^T(t|c|d1|d2)_-?\\d+$"),
names_to = c("which", "Time"),
names_pattern = "^T(t|c|d1|d2)_(-?\\d+)$",
values_to = "fix"
) %>%
pivot_wider(
names_from = which,
values_from = fix
) %>%
rename(
target = t,
competitor = c,
distractor1 = d1,
distractor2 = d2
) %>%
mutate(
target = as.numeric(target),
competitor = as.numeric(competitor),
distractor1 = as.numeric(distractor1),
distractor2 = as.numeric(distractor2),
distractor = (distractor1 + distractor2) / 2,
diffFix = target - competitor,
Time = as.numeric(Time)
) %>%
as.data.frame()
# Sort
gamData <- gamData[order(gamData$pp, gamData$trial, gamData$Time), ]
# Factors
gamData$pp <- as.factor(gamData$pp)
gamData$item <- as.factor(gamData$item)
gamData$contrast2 <- ifelse(gamData$contrast %in% c("FLEECE", "KIT"),
"FLEECE-KIT",
"DRESS-TRAP")
# Time window
gamData <- gamData %>% filter(Time > -150 & Time < 900)
gamData <- start_event(gamData, event = c("pp","trial"))# Average looks by Time x Accent x Contrast Set
plot_data <- gamData %>%
mutate(
accent = factor(accent, levels = c("SEng", "AusE")),
contrast2 = factor(contrast2, levels = c("FLEECE-KIT", "DRESS-TRAP"))
) %>%
group_by(Time, accent, contrast2) %>%
summarise(
target = mean(target, na.rm = TRUE),
competitor = mean(competitor, na.rm = TRUE),
distractor = mean(distractor, na.rm = TRUE),
.groups = "drop"
) %>%
pivot_longer(
cols = c(target, competitor, distractor),
names_to = "Object",
values_to = "Fixation"
) %>%
mutate(
Object = factor(
Object,
levels = c("target", "competitor", "distractor"),
labels = c("Target", "Competitor", "Distractor")
)
)
Fig.Fixations_Overall <- ggplot(
plot_data,
aes(x = Time, y = Fixation, colour = contrast2, linetype = Object)
) +
geom_line(linewidth = 0.6) +
geom_vline(xintercept = 0, colour = "darkgrey", linewidth = 0.3) +
facet_wrap(~ accent, nrow = 1) +
scale_colour_manual(values = c(
"FLEECE-KIT" = "black",
"DRESS-TRAP" = "darkgrey"
)) +
labs(
x = "Time (ms)",
y = "Mean fixation proportion",
colour = "Contrast",
linetype = "Looks to"
) +
theme_light() +
theme(
legend.text = element_text(size = 9),
legend.title = element_text(size = 9, face = "bold"),
legend.title.position = "top",
legend.spacing.y = unit(0.2, "cm"),
legend.box.spacing = unit(0.3, "cm"),
legend.position = "bottom",
axis.title.x = element_text(size = 10),
axis.title.y = element_text(size = 10)
)
print(Fig.Fixations_Overall)
library(ggpubr)
ggsave("Fixations_Overall.png",
width = 8, height = 4.5, dpi = 600,
bg = "white")plot_data_vowels <- gamData %>%
mutate(
accent = factor(accent, levels = c("SEng", "AusE")),
contrast2 = factor(contrast2, levels = c("FLEECE-KIT", "DRESS-TRAP")),
contrast = factor(contrast, levels = c("FLEECE", "KIT", "DRESS", "TRAP"))
) %>%
filter(contrast %in% c("FLEECE", "KIT", "DRESS", "TRAP")) %>%
group_by(Time, accent, contrast2, contrast) %>%
summarise(
target = mean(target, na.rm = TRUE),
competitor = mean(competitor, na.rm = TRUE),
.groups = "drop"
) %>%
pivot_longer(
cols = c(target, competitor),
names_to = "Object",
values_to = "Fixation"
) %>%
mutate(
Object = factor(
Object,
levels = c("target", "competitor"),
labels = c("Target", "Competitor")
)
)
Fig.Vowels <- ggplot(
plot_data_vowels,
aes(x = Time, y = Fixation, colour = contrast, linetype = Object)
) +
geom_line(linewidth = 0.6) +
geom_vline(xintercept = 0, colour = "darkgrey", linewidth = 0.3) +
scale_colour_manual(values = c(
"FLEECE" = "#377EB8",
"KIT" = "#4DAF4A",
"DRESS" = "#E41A1C",
"TRAP" = "#984EA3"
)) +
facet_grid(accent ~ contrast2) +
labs(
x = "Time (ms)",
y = "Mean fixation proportion",
colour = "Vowel",
linetype = "Looks to"
) +
theme_light() +
theme(
legend.text = element_text(size = 8),
legend.title = element_text(size = 8, face = "bold"),
legend.title.position = "top",
legend.spacing.y = unit(0.2, "cm"),
legend.box.spacing = unit(0.3, "cm"),
legend.position = "bottom",
axis.title.x = element_text(size = 10),
axis.title.y = element_text(size = 10)
) +
guides(
colour = guide_legend(title.position = "top", nrow = 1),
linetype = guide_legend(title.position = "top", nrow = 1)
)
print(Fig.Vowels)
library(ggpubr)
ggsave("Fixations_Contrasts.png",
width = 8, height = 6, dpi = 600,
bg = "white")rm(list = ls())
library(tidyverse)
library(mgcv)
library(itsadug)
library(stringr)
library(grid)
options(warn = 0)
# Load data
gazeTarget <- read.table("dataCatalan/gazeTarget.txt",
header = TRUE, sep = "\t", quote = "")
gazeComp <- read.table("dataCatalan/gazeCompet.txt",
header = TRUE, sep = "\t", quote = "")
trialData <- read.table("dataCatalan/trialInfo.txt",
header = TRUE, sep = "\t", quote = "\"")
# Remove fillers
expFilter <- grepl("_", trialData$condition)
gazeTargetE <- gazeTarget %>% filter(expFilter)
gazeCompE <- gazeComp %>% filter(expFilter)
trialDataE <- trialData %>% filter(expFilter)
# Timeline
timeline <- -200 + (0:120) * 10
colnames(gazeTargetE) <- paste0("Tt_", timeline)
colnames(gazeCompE) <- paste0("Tc_", timeline)
# Accent + contrast
temp <- str_split_fixed(trialDataE$condition, "_", 2)
trialDataE$accent <- temp[, 1]
trialDataE$contrast <- temp[, 2]
# Combine
useData <- cbind(trialDataE, gazeTargetE, gazeCompE) %>% as.data.frame()
useData <- useData %>%
filter(useData$clickedOn != "3") %>%
droplevels()
# Reshape
gamData <- useData %>%
pivot_longer(
cols = matches("^T[tc]_-?\\d+$"),
names_to = c("which", "Time"),
names_pattern = "^T([tc])_(-?\\d+)$",
values_to = "fix"
) %>%
pivot_wider(
names_from = which,
values_from = fix
) %>%
rename(target = t, competitor = c) %>%
mutate(
target = as.numeric(target),
competitor = as.numeric(competitor),
diffFix = target - competitor,
Time = as.numeric(Time)
) %>%
as.data.frame()
# Sort
gamData <- gamData[order(gamData$pp, gamData$trial, gamData$Time), ]
# Factors
gamData$pp <- as.factor(gamData$pp)
gamData$item <- as.factor(gamData$item)
# Restrict time window for plotting
gamData <- gamData %>%
filter(Time > -100 & Time < 900)
gamData <- start_event(gamData, event = c("pp", "trial"))
# Clean labels
gamData <- gamData %>%
mutate(
accent = factor(accent, levels = c("SEng", "AusE")),
contrast = toupper(contrast),
contrast2 = case_when(
contrast %in% c("DRESS", "TRAP") ~ "DRESS-TRAP",
contrast %in% c("FLEECE", "KIT") ~ "FLEECE-KIT",
TRUE ~ NA_character_
),
contrast2 = factor(contrast2, levels = c("FLEECE-KIT", "DRESS-TRAP"))
)
overall_data <- gamData %>%
group_by(Time, accent, contrast2) %>%
summarise(
target = mean(target, na.rm = TRUE),
competitor = mean(competitor, na.rm = TRUE),
.groups = "drop"
) %>%
mutate(
diff = target - competitor,
row = "Overall contrasts",
condition = as.character(contrast2)
)
dress_trap_data <- gamData %>%
filter(contrast %in% c("DRESS", "TRAP")) %>%
group_by(Time, accent, contrast) %>%
summarise(
target = mean(target, na.rm = TRUE),
competitor = mean(competitor, na.rm = TRUE),
.groups = "drop"
) %>%
mutate(
diff = target - competitor,
row = "DRESS vs TRAP",
condition = as.character(contrast)
)
fleece_kit_data <- gamData %>%
filter(contrast %in% c("FLEECE", "KIT")) %>%
group_by(Time, accent, contrast) %>%
summarise(
target = mean(target, na.rm = TRUE),
competitor = mean(competitor, na.rm = TRUE),
.groups = "drop"
) %>%
mutate(
diff = target - competitor,
row = "FLEECE vs KIT",
condition = as.character(contrast)
)
# Combine all rows
plot_data <- bind_rows(overall_data, dress_trap_data, fleece_kit_data)
# Long format for plotting
plot_data_long <- plot_data %>%
pivot_longer(
cols = c(target, competitor, diff),
names_to = "Measure",
values_to = "Value"
) %>%
mutate(
row = factor(
row,
levels = c("Overall contrasts", "DRESS vs TRAP", "FLEECE vs KIT")
),
Measure = factor(
Measure,
levels = c("target", "competitor", "diff"),
labels = c("Looks to Target", "Looks to Competitor", "Target - Competitor")
),
condition = factor(
condition,
levels = c("FLEECE-KIT", "DRESS-TRAP", "DRESS", "TRAP", "FLEECE", "KIT")
),
accent = factor(
accent,
levels = c("SEng", "AusE"),
labels = c("S.Eng", "AusE")
)
)
Fig.Fixations <- ggplot(
plot_data_long,
aes(x = Time, y = Value, colour = condition, linetype = accent)
) +
geom_line(linewidth = 0.6) +
geom_vline(xintercept = 0, colour = "darkgrey", linewidth = 0.3) +
facet_grid(row ~ Measure, scales = "free_y") +
labs(
x = "Time (ms)",
y = "Proportion (or difference)",
colour = "Condition",
linetype = "Accent"
) +
theme_light() +
theme(
strip.text = element_text(size = 10),
legend.text = element_text(size = 8),
legend.title = element_text(size = 8, face = "bold"),
legend.title.position = "top",
legend.key.height = unit(0.4, "cm"),
legend.spacing.y = unit(0.1, "cm"),
legend.box.spacing = unit(0.2, "cm"),
legend.position = "bottom",
axis.title.x = element_text(size = 10),
axis.title.y = element_text(size = 10)
) +
guides(
linetype = guide_legend(order = 1, title.position = "top", nrow = 1),
colour = guide_legend(order = 2, title.position = "top", nrow = 2)
)
print(Fig.Fixations)
rm(list = ls())
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)ggplot(useData, aes(x = vowel, y = rt)) +
geom_boxplot() + facet_wrap(.~ accent)
Reaction times (i.e. time interval between vowel onset and mouse click) appear broadly similar across vowels and accents, with comparable medians and variability, though there are several high outliers, particularly for LOT in AusE.
useData_trimmed <- useData %>%
group_by(pp, condition) %>%
mutate(
mean_rt = mean(rt, na.rm = TRUE),
sd_rt = sd(rt, na.rm = TRUE)
) %>%
filter(rt >= mean_rt - 2.5 * sd_rt,
rt <= mean_rt + 2.5 * sd_rt) %>%
ungroup()
ggplot(useData_trimmed, aes(x = rt)) +
geom_histogram() + facet_wrap(.~ condition)
useData_trimmed2 <- useData %>%
filter(rt <= 3000) %>%
droplevels()
ggplot(useData_trimmed2, aes(x = rt)) +
geom_histogram() + facet_wrap(.~ condition)
3.2 Analysis of the time course
The analysis is organised into three main blocks: (i) an overall analysis including both contrasts (FLEECE–KIT and DRESS–TRAP), (ii) a focused analysis of the FLEECE–KIT contrast, and (iii) a focused analysis of the DRESS–TRAP contrast.
Within each block, target and competitor fixations are analysed separately. For each fixation type, we conduct two complementary analyses. First, we apply a time-window analysis to examine differences in fixation proportions within the predefined time window. Second, we fit Generalised Additive Mixed Models (GAMMs) to capture the continuous time course of fixation behaviour.
3.2.1 Overall Effect of Contrast
3.2.1.1. TARGET FIXATIONS
A. TIME-WINDOW ANALYSIS
rm(list = ls())
library(lmerTest)
library(stringr)
library(tidyverse)
library(mgcv)
library(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("KIT|FLEECE|DRESS|TRAP", 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 == "FLEECE" | trialDataE$contrast == "KIT")
table(trialDataE$contrastC)##
## -0.5 0.5
## 3596 3598useData = cbind(trialDataE,gazeDataE)
# Remove distractor responses
useData <- useData %>%
filter(useData$clickedOn != "3") %>%
droplevels()windowData = useData[ ,c(1,2,4:8,10:134)] %>%
pivot_longer(., cols = starts_with("T_"),
names_to = "Time",
values_to = "position") %>% as.data.frame
windowData$position = as.numeric(windowData$position)
windowData$Time = gsub("T_","", windowData$Time) %>% as.numeric()
windowData = windowData[ order(windowData$pp,windowData$trial,windowData$Time), ]
windowData$pp = as.factor(windowData$pp)
windowData$item = as.factor(windowData$item)
windowData$isAustralian = ifelse(windowData$accent == "AusE",1,0)
windowData$isAustralianFleece = windowData$isAustralian *
as.numeric(windowData$contrast == "FLEECE")
windowData$vowelContrast = ifelse(windowData$contrastC < 0, "dressTrap", "kitFleece") %>% as.factor()
#restricting the data to 0-900 after vowel onset
windowData = windowData %>% filter(Time > 0 & Time < 900)
windowData <- start_event(windowData, event=c("pp","trial"))
windowDataPlot = windowData %>% group_by(vowelContrast, Time, accent) %>%
summarize(meanFix = mean(position))
Fig.Target_Overall <- ggplot(windowDataPlot, aes(x = Time, y = meanFix, col = vowelContrast)) + geom_line() + facet_wrap(~accent)
Fig.Target_Overall
The two contrasts show very similar trajectories overall. In AusE, DRESS–TRAP appears to rise slightly more steeply and reaches a marginally higher asymptote than KIT–FLEECE, though the difference seems negligible. In SEng, the curves almost completely overlap across the time course, indicating hardly any differences between the two contrasts.
Next, we restrict the dataset to the time window between 150 and 600 ms.
lmerData = windowData %>% filter(Time > 150 & Time < 600) %>%
group_by(pp, item, trial, isAustralian, contrastC) %>% summarize(meanLooks = mean(position))In addition, proportion looks were logit-transformed prior to analysis. Because the logit function is undefined for proportions of exactly 0 or 1, a small correction was applied to these values before transformation. The transformation maps bounded proportions (0–1) onto the real line, resulting in an approximately linear relationship in the mid-range and greater separation of values near the boundaries.
lmerData$accentC = lmerData$isAustralian - 0.5
outcomes = sort(unique(lmerData$meanLooks))
if (outcomes[1] == 0){
correction = outcomes[2] * 0.5
}
lmerData$logitLooks = ifelse(lmerData$meanLooks == 0, qlogis(correction),
ifelse(lmerData$meanLooks == 1, qlogis(1- correction),
qlogis(lmerData$meanLooks))
)
plot(lmerData$meanLooks, lmerData$logitLooks)
lmerData$trialScaled = scale(lmerData$trial)
target_Overall.lmer = lmer(logitLooks ~ contrastC*accentC*trialScaled +
(1+contrastC+contrastC||pp) +
(1+accentC*trialScaled||item), data = lmerData,
control=lmerControl(optimizer="bobyqa",
optCtrl=list(maxfun=2e5)))summary(target_Overall.lmer)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: logitLooks ~ contrastC * accentC * trialScaled + (1 + contrastC +
## contrastC || pp) + (1 + accentC * trialScaled || item)
## Data: lmerData
## Control: lmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 2e+05))
##
## REML criterion at convergence: 34744.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.28668 -0.60540 -0.00804 0.56125 2.62046
##
## Random effects:
## Groups Name Variance Std.Dev.
## pp (Intercept) 5.503e-01 0.741819
## pp.1 contrastC 0.000e+00 0.000000
## item (Intercept) 1.747e-01 0.417975
## item.1 accentC 1.047e-01 0.323537
## item.2 trialScaled 2.753e-04 0.016591
## item.3 accentC:trialScaled 1.255e-08 0.000112
## Residual 7.114e+00 2.667134
## Number of obs: 7193, groups: pp, 45; item, 40
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -1.445e-01 1.326e-01 6.544e+01 -1.089 0.280
## contrastC -7.488e-03 1.464e-01 3.802e+01 -0.051 0.959
## accentC 4.704e-02 8.110e-02 3.799e+01 0.580 0.565
## trialScaled 4.401e-02 3.168e-02 3.691e+01 1.389 0.173
## contrastC:accentC -3.892e-02 1.622e-01 3.799e+01 -0.240 0.812
## contrastC:trialScaled 8.920e-03 6.360e-02 3.749e+01 0.140 0.889
## accentC:trialScaled -1.898e-02 6.332e-02 7.112e+03 -0.300 0.764
## contrastC:accentC:trialScaled -1.807e-01 1.266e-01 7.112e+03 -1.427 0.154
##
## Correlation of Fixed Effects:
## (Intr) cntrsC accntC trlScl cntC:C cntC:S accC:S
## contrastC 0.000
## accentC 0.000 0.000
## trialScaled 0.000 -0.009 -0.021
## cntrstC:ccC 0.000 0.000 0.000 0.001
## cntrstC:trS -0.005 0.000 0.001 0.000 -0.021
## accntC:trlS -0.006 0.000 0.000 -0.005 -0.016 0.003
## cntrstC:C:S 0.000 -0.012 -0.016 0.004 0.000 -0.005 0.001
## optimizer (bobyqa) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')The linear mixed-effects model revealed no significant main effects or interactions of contrast, accent, or trial, indicating that none of these predictors reliably influenced logit-transformed looks in the analyzed time window.
B. Generalized Additive Mixed Model
gamData <- windowDataBuilding an initial GAMM.
fullModelFile = "target_overall_fullModels.Rdata"
if (file.exists(fullModelFile))
{
load(fullModelFile)
}else{
interaction.target.overall.gamm0 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
s(Time, by = isAustralian) +
s(Time, by = isAustralianFleece) +
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.target.overall.gamm0)
myRhoVal = myAutoCorr[2]
interaction.target.overall.gamm1 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast, k=12) +
s(Time, by = isAustralian, k= 12) +
s(Time, by = isAustralianFleece) +
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.target.overall.gamm0","interaction.target.overall.gamm1"),
file = fullModelFile)
}Model criticism: How are we doing regarding auto-correlation?
myAutoCorr =acf_resid(interaction.target.overall.gamm1)
Is the default for the smooth (k = 10) sufficient?
gam.check(interaction.target.overall.gamm1)
##
## Method: fREML Optimizer: perf chol
## $grad
## [1] 3.296696e-11 1.420641e-12 1.904921e-12 -4.385265e-03 1.573781e-10
## [6] 1.023430e-10 1.928449e-10 6.707523e-11 1.326015e-10 1.773053e-10
## [11] 1.587885e-11 1.391420e-11 1.314238e-11 -2.295053e-11 2.057199e-11
## [16] 4.394440e-11
##
## $hess
## [,1] [,2] [,3] [,4] [,5]
## [1,] 3.5192250527 -4.550710e-02 2.100452e-01 2.626458e-04 2.075795e-01
## [2,] -0.0455070963 2.803017e+00 2.265927e-01 -6.998254e-04 -1.073835e-03
## [3,] 0.2100451791 2.265927e-01 5.258658e-01 -1.852485e-03 5.435939e-03
## [4,] 0.0002626458 -6.998254e-04 -1.852485e-03 4.294458e-03 1.667993e-06
## [5,] 0.2075795041 -1.073835e-03 5.435939e-03 1.667993e-06 2.479997e+01
## [6,] 0.0119472403 -9.586892e-05 9.276785e-05 8.965045e-07 2.650607e-01
## [7,] -0.0409070588 8.909851e-04 -4.240585e-03 -7.496125e-06 -4.773530e-01
## [8,] -0.0029596419 1.449436e-01 1.193502e-02 -6.920561e-05 -1.668679e-04
## [9,] 0.0002745685 4.083678e-02 -1.557307e-04 7.848879e-06 1.640894e-05
## [10,] 0.0006915311 -4.823994e-02 -4.345030e-03 -6.667463e-06 -1.894978e-05
## [11,] 0.0656099373 -6.130468e-04 2.748565e-03 2.455385e-06 -6.641889e-02
## [12,] 0.0021382342 3.209034e-06 -6.600086e-05 7.069456e-08 -5.487546e-02
## [13,] -0.0130772657 2.412666e-04 -1.177022e-03 -2.452385e-06 -4.881541e-02
## [14,] -0.0016689112 6.584013e-02 5.822303e-03 1.240435e-03 -1.407754e-04
## [15,] 0.0004933634 8.026428e-03 1.426500e-03 2.130837e-04 4.342299e-04
## [16,] 0.0002772741 -3.108424e-02 -2.810343e-04 2.825167e-06 7.022105e-04
## [,6] [,7] [,8] [,9] [,10]
## [1,] 1.194724e-02 -4.090706e-02 -2.959642e-03 2.745685e-04 6.915311e-04
## [2,] -9.586892e-05 8.909851e-04 1.449436e-01 4.083678e-02 -4.823994e-02
## [3,] 9.276785e-05 -4.240585e-03 1.193502e-02 -1.557307e-04 -4.345030e-03
## [4,] 8.965045e-07 -7.496125e-06 -6.920561e-05 7.848879e-06 -6.667463e-06
## [5,] 2.650607e-01 -4.773530e-01 -1.668679e-04 1.640894e-05 -1.894978e-05
## [6,] 1.107759e+01 -7.712489e-01 -3.721447e-07 -1.623232e-05 -7.494353e-06
## [7,] -7.712489e-01 1.659373e+01 -3.085594e-06 -7.211144e-06 -5.291425e-05
## [8,] -3.721447e-07 -3.085594e-06 2.081511e+01 -3.072263e-01 -6.163901e-01
## [9,] -1.623232e-05 -7.211144e-06 -3.072263e-01 1.399868e+01 -5.482939e-01
## [10,] -7.494353e-06 -5.291425e-05 -6.163901e-01 -5.482939e-01 1.570422e+01
## [11,] -3.032464e-02 -5.336794e-02 -5.714824e-05 4.117351e-06 1.460817e-06
## [12,] -3.222974e-02 -2.073922e-02 3.376584e-07 -4.467054e-06 -2.989809e-06
## [13,] -1.737475e-02 -5.836679e-02 -7.301666e-06 -3.607158e-06 -1.991309e-05
## [14,] -5.658280e-05 -5.125504e-05 -7.224326e-02 -3.703103e-02 -7.654262e-02
## [15,] 5.535883e-04 3.061314e-04 -4.758906e-02 -4.209796e-02 -2.602548e-02
## [16,] 2.664956e-04 8.400339e-04 -4.253261e-02 -2.311322e-02 -6.469013e-02
## [,11] [,12] [,13] [,14] [,15]
## [1,] 6.560994e-02 2.138234e-03 -1.307727e-02 -1.668911e-03 4.933634e-04
## [2,] -6.130468e-04 3.209034e-06 2.412666e-04 6.584013e-02 8.026428e-03
## [3,] 2.748565e-03 -6.600086e-05 -1.177022e-03 5.822303e-03 1.426500e-03
## [4,] 2.455385e-06 7.069456e-08 -2.452385e-06 1.240435e-03 2.130837e-04
## [5,] -6.641889e-02 -5.487546e-02 -4.881541e-02 -1.407754e-04 4.342299e-04
## [6,] -3.032464e-02 -3.222974e-02 -1.737475e-02 -5.658280e-05 5.535883e-04
## [7,] -5.336794e-02 -2.073922e-02 -5.836679e-02 -5.125504e-05 3.061314e-04
## [8,] -5.714824e-05 3.376584e-07 -7.301666e-06 -7.224326e-02 -4.758906e-02
## [9,] 4.117351e-06 -4.467054e-06 -3.607158e-06 -3.703103e-02 -4.209796e-02
## [10,] 1.460817e-06 -2.989809e-06 -1.991309e-05 -7.654262e-02 -2.602548e-02
## [11,] 6.407871e+00 -1.677157e-01 -2.131213e-01 -4.032319e-05 4.891722e-05
## [12,] -1.677157e-01 3.544196e+00 -6.220722e-02 -1.610628e-05 1.581884e-04
## [13,] -2.131213e-01 -6.220722e-02 4.740807e+00 -2.389891e-05 9.945948e-05
## [14,] -4.032319e-05 -1.610628e-05 -2.389891e-05 8.158307e+00 -1.649808e-01
## [15,] 4.891722e-05 1.581884e-04 9.945948e-05 -1.649808e-01 4.020874e+00
## [16,] 1.870283e-04 8.243050e-05 3.713869e-04 -3.130074e-01 4.104496e-02
## [,16]
## [1,] 2.772741e-04
## [2,] -3.108424e-02
## [3,] -2.810343e-04
## [4,] 2.825167e-06
## [5,] 7.022105e-04
## [6,] 2.664956e-04
## [7,] 8.400339e-04
## [8,] -4.253261e-02
## [9,] -2.311322e-02
## [10,] -6.469013e-02
## [11,] 1.870283e-04
## [12,] 8.243050e-05
## [13,] 3.713869e-04
## [14,] -3.130074e-01
## [15,] 4.104496e-02
## [16,] 5.731466e+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):vowelContrastdressTrap 11.00 8.86 0.99 0.23
## s(Time):vowelContrastkitFleece 11.00 7.76 0.99 0.30
## s(Time):isAustralian 12.00 4.76 0.99 0.22
## s(Time):isAustralianFleece 10.00 2.02 0.99 0.28
## s(Time,pp):vowelContrastdressTrap 450.00 152.00 0.99 0.20
## s(Time,pp):vowelContrastkitFleece 450.00 157.69 0.99 0.24
## s(Time,item):vowelContrastdressTrap 400.00 51.91 0.99 0.23
## s(Time,item):vowelContrastkitFleece 400.00 61.46 0.99 0.27We now examine this model to determine whether there is an overall effect of contrast.
summary(interaction.target.overall.gamm1)##
## Family: binomial
## Link function: logit
##
## Formula:
## position ~ vowelContrast + s(Time, by = vowelContrast, k = 12) +
## s(Time, by = isAustralian, k = 12) + s(Time, by = isAustralianFleece) +
## 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.23207 0.07536 3.079 0.00207 **
## vowelContrastkitFleece 0.01710 0.10355 0.165 0.86882
## ---
## 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):vowelContrastdressTrap 8.858 10.150 692.689 <2e-16 ***
## s(Time):vowelContrastkitFleece 7.765 9.292 457.992 <2e-16 ***
## s(Time):isAustralian 4.757 5.954 8.567 0.196
## s(Time):isAustralianFleece 2.023 2.043 4.509 0.121
## s(Time,pp):vowelContrastdressTrap 151.998 448.000 696.142 <2e-16 ***
## s(Time,pp):vowelContrastkitFleece 157.693 449.000 690.294 <2e-16 ***
## s(Time,item):vowelContrastdressTrap 51.907 198.000 173.950 <2e-16 ***
## s(Time,item):vowelContrastkitFleece 61.460 199.000 235.571 <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.1041e+05 Scale est. = 1 n = 640177There is no evidence for an overall (parametric) effect of contrast, given that the coefficient for vowelContrastkitFleece was not significant (p = .452).
smooth <- plot_smooth(interaction.target.overall.gamm1, view = "Time",
plot_all = "vowelContrast",
rug = FALSE, rm.ranef = TRUE, lwd = 4)## Summary:
## * vowelContrast : factor; set to the value(s): dressTrap, kitFleece.
## * Time : numeric predictor; with 30 values ranging from 10.000000 to 890.000000.
## * isAustralian : numeric predictor; set to the value(s): 0.
## * isAustralianFleece : 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):vowelContrastdressTrap,s(Time,pp):vowelContrastkitFleece,s(Time,item):vowelContrastdressTrap,s(Time,item):vowelContrastkitFleece
## 
We now proceed to model comparison by successively removing effects to evaluate their contribution to model fit.
reducedModelFile = "target_overall_reducedModels.Rdata"
if (file.exists(reducedModelFile))
{
load(reducedModelFile)
}else{
minusInteraction.target.overall.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.target.overall.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.target.overall.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.target.overall.gamm",
"minusAccent.target.overall.gamm",
"timeOnly.target.overall.gamm"),
file = reducedModelFile)
}Model comparisons:
# Full model versus no interaction
compareML(interaction.target.overall.gamm1, minusInteraction.target.overall.gamm)## interaction.target.overall.gamm1: position ~ vowelContrast + s(Time, by = vowelContrast, k = 12) +
## s(Time, by = isAustralian, k = 12) + s(Time, by = isAustralianFleece) +
## s(Time, pp, by = vowelContrast, bs = "fs") + s(Time, item,
## by = vowelContrast, bs = "fs")
##
## minusInteraction.target.overall.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.target.overall.gamm preferred: lower fREML score (200.345), and lower df (1.000).
## -----
## Model Score Edf Difference Df
## 1 interaction.target.overall.gamm1 210409.6 16
## 2 minusInteraction.target.overall.gamm 210209.3 15 -200.345 1.000
##
## AIC difference: 696.16, model minusInteraction.target.overall.gamm has lower AIC.# Both main effects versus minus accent
compareML(minusInteraction.target.overall.gamm, minusAccent.target.overall.gamm)## minusInteraction.target.overall.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.target.overall.gamm: position ~ vowelContrast + s(Time, by = vowelContrast) + s(Time,
## pp, by = vowelContrast, bs = "fs") + s(Time, item, by = vowelContrast,
## bs = "fs")
##
## Model minusAccent.target.overall.gamm preferred: lower fREML score (220.557), and lower df (1.000).
## -----
## Model Score Edf Difference Df
## 1 minusInteraction.target.overall.gamm 210209.3 15
## 2 minusAccent.target.overall.gamm 209988.7 14 -220.557 1.000
##
## AIC difference: 975.34, model minusAccent.target.overall.gamm has lower AIC.# Main effect of Contrast versus time-only
compareML(minusAccent.target.overall.gamm, timeOnly.target.overall.gamm)## minusAccent.target.overall.gamm: position ~ vowelContrast + s(Time, by = vowelContrast) + s(Time,
## pp, by = vowelContrast, bs = "fs") + s(Time, item, by = vowelContrast,
## bs = "fs")
##
## timeOnly.target.overall.gamm: position ~ 1 + s(Time) + s(Time, pp, by = vowelContrast, bs = "fs") +
## s(Time, item, by = vowelContrast, bs = "fs")
##
## Model timeOnly.target.overall.gamm preferred: lower fREML score (139.398), and lower df (1.000).
## -----
## Model Score Edf Difference Df
## 1 minusAccent.target.overall.gamm 209988.7 14
## 2 timeOnly.target.overall.gamm 209849.3 13 -139.398 1.000
##
## AIC difference: 620.82, model timeOnly.target.overall.gamm has lower AIC.Model comparison shows that progressively simpler models fit the data better, with lower fREML scores and substantially lower AIC values at each step. The best-supported model includes only a general smooth of Time (with random smooths for participants and items), providing no evidence that vowel contrast, accent, or their interactions significantly improve model fit.
3.2.1.2. Competitor Fixations
A. TIME-WINDOW ANALYSIS
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("KIT|FLEECE|DRESS|TRAP", 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 == "FLEECE" | trialDataE$contrast == "KIT")
table(trialDataE$contrastC)##
## -0.5 0.5
## 3596 3598useData = cbind(trialDataE,gazeDataE)
# Remove distractor responses
useData <- useData %>%
filter(useData$clickedOn != "3") %>%
droplevels()windowData = useData[ ,c(1,2,4:8,10:134)] %>%
pivot_longer(., cols = starts_with("T_"),
names_to = "Time",
values_to = "position") %>% as.data.frame
windowData$position = as.numeric(windowData$position)
windowData$Time = gsub("T_","", windowData$Time) %>% as.numeric()
windowData = windowData[ order(windowData$pp,windowData$trial,windowData$Time), ]
windowData$pp = as.factor(windowData$pp)
windowData$item = as.factor(windowData$item)
windowData$isAustralian = ifelse(windowData$accent == "AusE",1,0)
windowData$isAustralianFleece = windowData$isAustralian *
as.numeric(windowData$contrast == "FLEECE")
windowData$vowelContrast = ifelse(windowData$contrastC < 0, "dressTrap", "kitFleece") %>% as.factor()
#restricting the data to 0-900 after vowel onset
windowData = windowData %>% filter(Time > 0 & Time < 900)
windowData <- start_event(windowData, event=c("pp","trial"))
windowDataPlot = windowData %>% group_by(vowelContrast, Time, accent) %>%
summarize(meanFix = mean(position))
Fig.Target_Overall <- ggplot(windowDataPlot, aes(x = Time, y = meanFix, col = vowelContrast)) + geom_line() + facet_wrap(~accent)
Fig.Target_Overall
Across both accents, the KIT-FLEECT contrast is consistently slightly higher than the DRESS-TRAP contrast, especially around the peak, but the difference between accents themselves appears minimal.
Next, we restrict the dataset to the time window between 150 and 600 ms.
lmerData = windowData %>% filter(Time > 150 & Time < 600) %>%
group_by(pp, item, trial, isAustralian, contrastC) %>% summarize(meanLooks = mean(position))lmerData$accentC = lmerData$isAustralian - 0.5
outcomes = sort(unique(lmerData$meanLooks))
if (outcomes[1] == 0){
correction = outcomes[2] * 0.5
}
lmerData$logitLooks = ifelse(lmerData$meanLooks == 0, qlogis(correction),
ifelse(lmerData$meanLooks == 1, qlogis(1- correction),
qlogis(lmerData$meanLooks))
)
plot(lmerData$meanLooks, lmerData$logitLooks)
lmerData$trialScaled = scale(lmerData$trial)
compet_overall.lmer = lmer(logitLooks ~ contrastC*accentC*trialScaled +
(1+contrastC+contrastC||pp) +
(1+accentC*trialScaled||item), data = lmerData,
control=lmerControl(optimizer="bobyqa",
optCtrl=list(maxfun=2e5)))summary(compet_overall.lmer)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: logitLooks ~ contrastC * accentC * trialScaled + (1 + contrastC +
## contrastC || pp) + (1 + accentC * trialScaled || item)
## Data: lmerData
## Control: lmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 2e+05))
##
## REML criterion at convergence: 32273.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.46458 -0.95633 0.01079 0.88141 3.14678
##
## Random effects:
## Groups Name Variance Std.Dev.
## pp (Intercept) 0.04883 0.2210
## pp.1 contrastC 0.01842 0.1357
## item (Intercept) 0.05630 0.2373
## item.1 accentC 0.06773 0.2603
## item.2 trialScaled 0.01130 0.1063
## item.3 accentC:trialScaled 0.03034 0.1742
## Residual 5.09251 2.2567
## Number of obs: 7193, groups: pp, 45; item, 40
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -2.20802 0.05660 53.61852 -39.013 < 2e-16 ***
## contrastC 0.26574 0.09424 37.03371 2.820 0.00767 **
## accentC -0.06738 0.06735 38.03649 -1.000 0.32343
## trialScaled -0.02628 0.03159 37.46773 -0.832 0.41072
## contrastC:accentC 0.06442 0.13470 38.03676 0.478 0.63519
## contrastC:trialScaled -0.06007 0.06327 37.68477 -0.949 0.34845
## accentC:trialScaled 0.05946 0.06024 35.86896 0.987 0.33023
## contrastC:accentC:trialScaled 0.03956 0.12047 35.85333 0.328 0.74452
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) cntrsC accntC trlScl cntC:C cntC:S accC:S
## contrastC 0.000
## accentC 0.000 0.001
## trialScaled 0.000 -0.009 -0.018
## cntrstC:ccC 0.001 0.000 -0.001 0.000
## cntrstC:trS -0.008 0.000 0.001 0.001 -0.018
## accntC:trlS -0.011 0.000 0.000 -0.003 -0.014 0.002
## cntrstC:C:S 0.000 -0.014 -0.014 0.003 0.000 -0.003 0.001The analysis revealed a significant main effect of contrast, indicating that fixation proportions differed reliably between the two contrasts (FLEECE–KIT vs. DRESS–TRAP). Specifically, FLEECE-KIT showed higher fixation proportions than DRESS-TRAP.
There was no significant main effect of accent, suggesting that fixation behaviour did not differ between accents. Likewise, there was no effect of trial, indicating no systematic change in fixation proportions over the course of the experiment.
B. Generalized Additive Mixed Model
gamData <- windowDataBuilding an initial GAMM.
fullModelFile = "compet_overall_fullModels.Rdata"
if (file.exists(fullModelFile))
{
load(fullModelFile)
}else{
interaction.compet.overall.gamm0 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
s(Time, by = isAustralian) +
s(Time, by = isAustralianFleece) +
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.compet.overall.gamm0)
myRhoVal = myAutoCorr[2]
interaction.compet.overall.gamm1 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast, k=12) +
s(Time, by = isAustralian, k= 12) +
s(Time, by = isAustralianFleece) +
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.compet.overall.gamm0","interaction.compet.overall.gamm1"),
file = fullModelFile)
}Model criticism: How are we doing regarding auto-correlation?
myAutoCorr =acf_resid(interaction.compet.overall.gamm1)
Is the default for the smooth (k = 10) sufficient?
gam.check(interaction.compet.overall.gamm1)
##
## Method: fREML Optimizer: perf chol
## $grad
## [1] 4.584700e-04 5.426445e-05 -3.329283e-03 -2.421520e-03 1.096973e-03
## [6] 2.767227e-03 1.089927e-03 1.244240e-03 3.661732e-03 1.707183e-03
## [11] -8.276378e-05 1.441878e-03 7.827286e-04 7.409982e-04 1.386816e-03
## [16] 8.104545e-04
##
## $hess
## [,1] [,2] [,3] [,4] [,5]
## [1,] 3.637159e+00 -1.263332e-04 7.565670e-04 -1.730945e-06 1.137561e-01
## [2,] -1.263332e-04 3.880018e+00 6.556908e-04 3.245492e-04 -1.242488e-05
## [3,] 7.565670e-04 6.556908e-04 3.311890e-03 2.702761e-06 2.790554e-05
## [4,] -1.730945e-06 3.245492e-04 2.702761e-06 2.426764e-03 -1.389280e-07
## [5,] 1.137561e-01 -1.242488e-05 2.790554e-05 -1.389280e-07 1.292776e+01
## [6,] -3.698181e-02 -4.732486e-06 -1.505042e-05 2.586537e-08 -9.823530e-01
## [7,] -7.972251e-03 -4.401023e-06 -1.897938e-06 -1.387638e-07 1.250086e-01
## [8,] -7.441482e-06 3.793348e-02 5.543809e-06 -1.952623e-05 -1.791440e-05
## [9,] -4.128245e-06 -1.697266e-02 -3.539104e-06 1.666544e-06 -1.409282e-05
## [10,] -1.954966e-06 -1.271634e-02 -1.853912e-06 -3.341665e-06 -5.273185e-06
## [11,] 2.608905e-02 -1.638988e-06 2.847226e-06 -2.384674e-08 -8.506136e-03
## [12,] -5.589710e-03 -3.570804e-06 -6.941290e-06 -1.002481e-07 -1.865793e-02
## [13,] -7.220000e-03 -3.813719e-06 -1.702436e-06 -1.078706e-07 -1.809707e-02
## [14,] 5.045230e-06 5.560974e-02 3.356883e-05 2.592312e-03 5.273078e-07
## [15,] 4.391103e-05 -6.293344e-03 -1.482192e-06 -2.647517e-04 2.970133e-05
## [16,] 2.822869e-04 -1.760659e-02 -3.829450e-05 -2.018685e-05 8.204527e-04
## [,6] [,7] [,8] [,9] [,10]
## [1,] -3.698181e-02 -7.972251e-03 -7.441482e-06 -4.128245e-06 -1.954966e-06
## [2,] -4.732486e-06 -4.401023e-06 3.793348e-02 -1.697266e-02 -1.271634e-02
## [3,] -1.505042e-05 -1.897938e-06 5.543809e-06 -3.539104e-06 -1.853912e-06
## [4,] 2.586537e-08 -1.387638e-07 -1.952623e-05 1.666544e-06 -3.341665e-06
## [5,] -9.823530e-01 1.250086e-01 -1.791440e-05 -1.409282e-05 -5.273185e-06
## [6,] 1.261291e+01 -9.350702e-01 -1.409617e-05 -1.118615e-05 -4.560996e-06
## [7,] -9.350702e-01 2.657650e+00 -6.896266e-06 -5.290515e-06 -3.378799e-06
## [8,] -1.409617e-05 -6.896266e-06 1.016573e+01 -1.012752e+00 -8.765233e-02
## [9,] -1.118615e-05 -5.290515e-06 -1.012752e+00 1.367874e+01 -1.157978e+00
## [10,] -4.560996e-06 -3.378799e-06 -8.765233e-02 -1.157978e+00 5.059995e+00
## [11,] -9.269353e-03 -2.451294e-03 -2.082300e-06 -1.614729e-06 -8.511872e-07
## [12,] -3.436969e-02 -4.798545e-03 -6.712247e-06 -5.165409e-06 -3.381579e-06
## [13,] -1.633373e-02 -6.628478e-03 -6.327431e-06 -4.880131e-06 -2.902661e-06
## [14,] 9.829808e-06 3.200538e-05 -4.229939e-02 -2.372240e-02 -1.870480e-02
## [15,] 4.969269e-05 1.074157e-04 -1.614626e-02 -3.850110e-02 -1.291256e-02
## [16,] 6.426852e-04 3.007483e-04 -3.683446e-02 -2.477074e-02 -1.604642e-02
## [,11] [,12] [,13] [,14] [,15]
## [1,] 2.608905e-02 -5.589710e-03 -7.220000e-03 5.045230e-06 4.391103e-05
## [2,] -1.638988e-06 -3.570804e-06 -3.813719e-06 5.560974e-02 -6.293344e-03
## [3,] 2.847226e-06 -6.941290e-06 -1.702436e-06 3.356883e-05 -1.482192e-06
## [4,] -2.384674e-08 -1.002481e-07 -1.078706e-07 2.592312e-03 -2.647517e-04
## [5,] -8.506136e-03 -1.865793e-02 -1.809707e-02 5.273078e-07 2.970133e-05
## [6,] -9.269353e-03 -3.436969e-02 -1.633373e-02 9.829808e-06 4.969269e-05
## [7,] -2.451294e-03 -4.798545e-03 -6.628478e-03 3.200538e-05 1.074157e-04
## [8,] -2.082300e-06 -6.712247e-06 -6.327431e-06 -4.229939e-02 -1.614626e-02
## [9,] -1.614729e-06 -5.165409e-06 -4.880131e-06 -2.372240e-02 -3.850110e-02
## [10,] -8.511872e-07 -3.381579e-06 -2.902661e-06 -1.870480e-02 -1.291256e-02
## [11,] 1.096839e+00 -2.704024e-01 5.249934e-02 5.754825e-06 2.052520e-05
## [12,] -2.704024e-01 5.152048e+00 -6.601598e-01 3.340481e-05 1.108091e-04
## [13,] 5.249934e-02 -6.601598e-01 3.291936e+00 2.469080e-05 8.436700e-05
## [14,] 5.754825e-06 3.340481e-05 2.469080e-05 6.554823e+00 -1.042930e-02
## [15,] 2.052520e-05 1.108091e-04 8.436700e-05 -1.042930e-02 4.719261e+00
## [16,] 9.224719e-05 2.913367e-04 2.782832e-04 1.716633e-01 -5.130361e-01
## [,16]
## [1,] 2.822869e-04
## [2,] -1.760659e-02
## [3,] -3.829450e-05
## [4,] -2.018685e-05
## [5,] 8.204527e-04
## [6,] 6.426852e-04
## [7,] 3.007483e-04
## [8,] -3.683446e-02
## [9,] -2.477074e-02
## [10,] -1.604642e-02
## [11,] 9.224719e-05
## [12,] 2.913367e-04
## [13,] 2.782832e-04
## [14,] 1.716633e-01
## [15,] -5.130361e-01
## [16,] 4.889492e+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):vowelContrastdressTrap 11.00 8.12 1 0.99
## s(Time):vowelContrastkitFleece 11.00 8.79 1 0.97
## s(Time):isAustralian 12.00 2.01 1 0.99
## s(Time):isAustralianFleece 10.00 2.01 1 0.98
## s(Time,pp):vowelContrastdressTrap 450.00 121.51 1 0.98
## s(Time,pp):vowelContrastkitFleece 450.00 120.16 1 0.98
## s(Time,item):vowelContrastdressTrap 400.00 35.44 1 0.98
## s(Time,item):vowelContrastkitFleece 400.00 60.80 1 0.99We now examine this model to determine whether there is an overall effect of contrast.
summary(interaction.compet.overall.gamm1)##
## Family: binomial
## Link function: logit
##
## Formula:
## position ~ vowelContrast + s(Time, by = vowelContrast, k = 12) +
## s(Time, by = isAustralian, k = 12) + s(Time, by = isAustralianFleece) +
## 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.60764 0.04977 -32.303 <2e-16 ***
## vowelContrastkitFleece 0.05514 0.07394 0.746 0.456
## ---
## 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):vowelContrastdressTrap 8.117 9.675 172.540 <2e-16 ***
## s(Time):vowelContrastkitFleece 8.790 10.120 232.446 <2e-16 ***
## s(Time):isAustralian 2.012 2.023 4.800 0.101
## s(Time):isAustralianFleece 2.011 2.021 2.791 0.230
## s(Time,pp):vowelContrastdressTrap 121.512 448.000 350.724 <2e-16 ***
## s(Time,pp):vowelContrastkitFleece 120.162 449.000 352.678 <2e-16 ***
## s(Time,item):vowelContrastdressTrap 35.445 198.000 141.809 <2e-16 ***
## s(Time,item):vowelContrastkitFleece 60.799 199.000 196.658 <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.6329e+05 Scale est. = 1 n = 640177The binomial GAMM shows significant nonlinear time effects for both vowel contrasts and substantial participant and item variability over time, but no significant main effect of vowel contrast nor any significant accent-related smooths, indicating that while fixation trajectories change over time, they do not reliably differ by contrast or accent.
smooth <- plot_smooth(interaction.compet.overall.gamm1, view = "Time",
plot_all = "vowelContrast",
rug = FALSE, rm.ranef = TRUE, lwd = 4)## Summary:
## * vowelContrast : factor; set to the value(s): dressTrap, kitFleece.
## * Time : numeric predictor; with 30 values ranging from 10.000000 to 890.000000.
## * isAustralian : numeric predictor; set to the value(s): 0.
## * isAustralianFleece : 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):vowelContrastdressTrap,s(Time,pp):vowelContrastkitFleece,s(Time,item):vowelContrastdressTrap,s(Time,item):vowelContrastkitFleece
## 
We now proceed to model comparison by successively removing effects to evaluate their contribution to model fit.
reducedModelFile = "compet_overall_reducedModels.Rdata"
if (file.exists(reducedModelFile))
{
load(reducedModelFile)
}else{
minusInteraction.compet.overall.gamm <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
s(Time, by = isAustralian) +
#s(Time, by = isAustralianFleece) +
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.compet.overall.gamm <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
#s(Time, by = isAustralian) +
#s(Time, by = isAustralianFleece) +
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.compet.overall.gamm <- bam( position ~ 1 +
s(Time) +
#s(Time, by = isAustralian) +
#s(Time, by = isAustralianFleece) +
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.compet.overall.gamm",
"minusAccent.compet.overall.gamm",
"timeOnly.compet.overall.gamm"),
file = reducedModelFile)
}Model comparisons:
# Full model versus no interaction
compareML(interaction.compet.overall.gamm1, minusInteraction.compet.overall.gamm)## interaction.compet.overall.gamm1: position ~ vowelContrast + s(Time, by = vowelContrast, k = 12) +
## s(Time, by = isAustralian, k = 12) + s(Time, by = isAustralianFleece) +
## s(Time, pp, by = vowelContrast, bs = "fs") + s(Time, item,
## by = vowelContrast, bs = "fs")
##
## minusInteraction.compet.overall.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.compet.overall.gamm preferred: lower fREML score (22.588), and lower df (1.000).
## -----
## Model Score Edf Difference Df
## 1 interaction.compet.overall.gamm1 163293.1 16
## 2 minusInteraction.compet.overall.gamm 163270.5 15 -22.588 1.000
##
## AIC difference: 117.07, model minusInteraction.compet.overall.gamm has lower AIC.# Both main effects versus minus accent
compareML(minusInteraction.compet.overall.gamm, minusAccent.compet.overall.gamm)## minusInteraction.compet.overall.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.compet.overall.gamm: position ~ vowelContrast + s(Time, by = vowelContrast) + s(Time,
## pp, by = vowelContrast, bs = "fs") + s(Time, item, by = vowelContrast,
## bs = "fs")
##
## Model minusAccent.compet.overall.gamm preferred: lower fREML score (99.600), and lower df (1.000).
## -----
## Model Score Edf Difference Df
## 1 minusInteraction.compet.overall.gamm 163270.5 15
## 2 minusAccent.compet.overall.gamm 163170.9 14 -99.600 1.000
##
## AIC difference: 488.39, model minusAccent.compet.overall.gamm has lower AIC.# Main effect of Contrast versus time-only
compareML(minusAccent.compet.overall.gamm, timeOnly.compet.overall.gamm)## minusAccent.compet.overall.gamm: position ~ vowelContrast + s(Time, by = vowelContrast) + s(Time,
## pp, by = vowelContrast, bs = "fs") + s(Time, item, by = vowelContrast,
## bs = "fs")
##
## timeOnly.compet.overall.gamm: position ~ 1 + s(Time) + s(Time, pp, by = vowelContrast, bs = "fs") +
## s(Time, item, by = vowelContrast, bs = "fs")
##
## Model timeOnly.compet.overall.gamm preferred: lower fREML score (435.395), and lower df (1.000).
## -----
## Model Score Edf Difference Df
## 1 minusAccent.compet.overall.gamm 163170.9 14
## 2 timeOnly.compet.overall.gamm 162735.5 13 -435.395 1.000
##
## AIC difference: 1701.89, model timeOnly.compet.overall.gamm has lower AIC.Model comparisons again show that progressively simpler models provide a better fit, with lower fREML scores and substantially lower AIC values at each step. The best-supported model is the time-only model (with random smooths for participants and items), indicating no reliable effects of vowel contrast, accent, or their interaction on the trajectory.
3.2.1.3. TARGET-COMPETITOR FIXATIONS
A. TIME-WINDOW ANALYSIS
rm(list = ls())
library(lmerTest)
library(stringr)
library(tidyverse)
library(mgcv)
library(itsadug)
options(warn = 0)
# Load data
gazeTarget <- read.table("dataCatalan/gazeTarget.txt",
header = TRUE, sep = "\t", quote = "")
gazeComp <- read.table("dataCatalan/gazeCompet.txt",
header = TRUE, sep = "\t", quote = "")
trialData <- read.table("dataCatalan/trialInfo.txt",
header = TRUE, sep = "\t", quote = "\"")
# Keep experimental trials only
expFilter <- grepl("KIT|FLEECE|DRESS|TRAP", trialData$condition)
gazeTargetE <- gazeTarget %>% filter(expFilter)
gazeCompE <- gazeComp %>% filter(expFilter)
trialDataE <- trialData %>% filter(expFilter)
# Timeline
timeline <- -200 + (0:120) * 10
colnames(gazeTargetE) <- paste0("Tt_", timeline)
colnames(gazeCompE) <- paste0("Tc_", timeline)
# Contrast coding variables from condition
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 %in% c("FLEECE", "KIT"))
table(trialDataE$contrastC)##
## -0.5 0.5
## 3596 3598# Combine
useData <- cbind(trialDataE, gazeTargetE, gazeCompE)
# Remove distractor responses
useData <- useData %>%
filter(useData$clickedOn != "3") %>%
droplevels()# Reshape target and competitor into one time-series dataset
windowData <- useData %>%
pivot_longer(
cols = matches("^T(t|c)_-?\\d+$"),
names_to = c("which", "Time"),
names_pattern = "^T(t|c)_(-?\\d+)$",
values_to = "position"
) %>%
pivot_wider(
names_from = which,
values_from = position
) %>%
rename(
target = t,
competitor = c
) %>%
mutate(
target = as.numeric(target),
competitor = as.numeric(competitor),
diffLooks = target - competitor,
Time = as.numeric(Time)
) %>%
arrange(pp, trial, Time) %>%
as.data.frame()
windowData$pp <- as.factor(windowData$pp)
windowData$item <- as.factor(windowData$item)
windowData$isAustralian <- ifelse(windowData$accent == "AusE", 1, 0)
windowData$isAustralianFleece <- windowData$isAustralian *
as.numeric(windowData$contrast == "FLEECE")
windowData$vowelContrast <- ifelse(windowData$contrastC < 0, "dressTrap", "kitFleece") %>%
as.factor()
# Restrict to 0-900 ms after vowel onset
windowData <- windowData %>%
filter(Time > 0 & Time < 900)
windowData <- start_event(windowData, event = c("pp", "trial"))# Quick plot of mean target-competitor difference over time
windowDataPlot <- windowData %>%
group_by(vowelContrast, Time, accent) %>%
summarise(meanDiff = mean(diffLooks, na.rm = TRUE), .groups = "drop")
Fig.Diff_Overall <- ggplot(
windowDataPlot,
aes(x = Time, y = meanDiff, col = vowelContrast)
) +
geom_line() +
facet_wrap(~ accent)
Fig.Diff_Overall
lmerData <- windowData %>%
filter(Time > 150 & Time < 600) %>%
group_by(pp, item, trial, isAustralian, accentC, contrastC) %>%
summarise(
meanDiff = mean(diffLooks, na.rm = TRUE),
.groups = "drop"
)# Optional scaling of trial number
lmerData$trialScaled <- scale(lmerData$trial)
# Linear mixed model on continuous difference score
diff_Overall.lmer <- lmer(
meanDiff ~ contrastC * accentC * trialScaled +
(1+contrastC||pp) +
(1 + accentC * trialScaled || item),
data = lmerData,
control = lmerControl(
optimizer = "bobyqa",
optCtrl = list(maxfun = 2e5)
)
)## boundary (singular) fit: see help('isSingular')summary(diff_Overall.lmer)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: meanDiff ~ contrastC * accentC * trialScaled + (1 + contrastC ||
## pp) + (1 + accentC * trialScaled || item)
## Data: lmerData
## Control: lmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 2e+05))
##
## REML criterion at convergence: 11799.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.67908 -0.78766 0.04043 0.86769 1.99274
##
## Random effects:
## Groups Name Variance Std.Dev.
## pp (Intercept) 0.0093858 0.09688
## pp.1 contrastC 0.0000000 0.00000
## item (Intercept) 0.0061928 0.07869
## item.1 accentC 0.0047420 0.06886
## item.2 trialScaled 0.0002774 0.01665
## item.3 accentC:trialScaled 0.0002296 0.01515
## Residual 0.2930778 0.54137
## Number of obs: 7193, groups: pp, 45; item, 40
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 0.243134 0.020105 67.486703 12.093 <2e-16
## contrastC -0.039095 0.027975 38.005216 -1.398 0.170
## accentC -0.013497 0.016789 38.013475 -0.804 0.426
## trialScaled 0.007454 0.006930 37.097024 1.076 0.289
## contrastC:accentC 0.007677 0.033578 38.013825 0.229 0.820
## contrastC:trialScaled 0.012082 0.013902 37.551140 0.869 0.390
## accentC:trialScaled 0.007706 0.013076 36.672485 0.589 0.559
## contrastC:accentC:trialScaled 0.025948 0.026145 36.638639 0.992 0.327
##
## (Intercept) ***
## contrastC
## accentC
## trialScaled
## contrastC:accentC
## contrastC:trialScaled
## accentC:trialScaled
## contrastC:accentC:trialScaled
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) cntrsC accntC trlScl cntC:C cntC:S accC:S
## contrastC 0.000
## accentC 0.000 0.000
## trialScaled 0.000 -0.008 0.019
## cntrstC:ccC 0.000 0.000 0.000 0.000
## cntrstC:trS -0.006 0.000 0.000 0.000 0.019
## accntC:trlS 0.008 0.000 0.000 0.004 -0.015 -0.003
## cntrstC:C:S 0.000 0.012 -0.015 -0.004 0.000 0.004 0.001
## optimizer (bobyqa) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')anova(diff_Overall.lmer)3.2.2 FLEECE-KIT contrast
3.2.2.1. Target fixations
A. TIME-WINDOW ANALYSIS
rm(list = ls())
library(lmerTest)
library(stringr)
library(tidyverse)
library(mgcv)
library(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("KIT|FLEECE", 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
## 1799 1799trialDataE$contrastC = -0.5 + as.numeric(trialDataE$contrast == "FLEECE")
table(trialDataE$contrastC)##
## -0.5 0.5
## 1799 1799useData = cbind(trialDataE,gazeDataE)
# Remove distractor responses
useData <- useData %>%
filter(useData$clickedOn != "3") %>%
droplevels()We first need to reorganize the data.
windowData = useData[ ,c(1,2,4:8,10:134)] %>%
pivot_longer(., cols = starts_with("T_"),
names_to = "Time",
values_to = "position") %>% as.data.frame
windowData$position = as.numeric(windowData$position)
windowData$Time = gsub("T_","", windowData$Time) %>% as.numeric()
windowData = windowData[ order(windowData$pp,windowData$trial,windowData$Time), ]
windowData$pp = as.factor(windowData$pp)
windowData$item = as.factor(windowData$item)
windowData$isAustralian = ifelse(windowData$accent == "AusE",1,0)
windowData$isAustralianFleece = windowData$isAustralian *
as.numeric(windowData$contrast == "FLEECE")
windowData$vowelContrast = windowData$contrast %>% as.factor()
#restricting the data to 0-900 after vowel onset
windowData = windowData %>% filter(Time > 0 & Time < 900)
windowData <- start_event(windowData, event=c("pp","trial"))
windowDataPlot = windowData %>% group_by(vowelContrast, Time, accent) %>%
summarize(meanFix = mean(position))
Fig.target_kitFleece <- ggplot(windowDataPlot, aes(x = Time, y = meanFix, col = vowelContrast)) + geom_line() + facet_wrap(~accent)
print(Fig.target_kitFleece)
aggregatedData = windowData %>% filter(Time > 200 & Time < 600) %>%
group_by(pp, item, contrast) %>% summarize(meanLooks = mean(position)) %>%
ungroup()
meanLooksValues = sort(unique(aggregatedData$meanLooks))
correction = meanLooksValues[2]/2
aggregatedData$logitLooks = ifelse( aggregatedData$meanLooks == 0,
qlogis(correction),
ifelse(aggregatedData$meanLooks == 1,
qlogis(1 - correction),
qlogis(aggregatedData$meanLooks)
)
)
plot(aggregatedData$meanLooks, aggregatedData$logitLooks)
library(lmerTest)
target_fleeceKit.lmer = lmer(meanLooks ~ contrast +
(1+contrast|pp)+
(1|item), data = aggregatedData)
summary(target_fleeceKit.lmer)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: meanLooks ~ contrast + (1 + contrast | pp) + (1 | item)
## Data: aggregatedData
##
## REML criterion at convergence: -401.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.97962 -0.63052 -0.00039 0.66759 2.62561
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## pp (Intercept) 0.008309 0.09115
## contrastKIT 0.002434 0.04934 -0.36
## item (Intercept) 0.003378 0.05812
## Residual 0.032292 0.17970
## Number of obs: 900, groups: pp, 45; item, 20
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 0.47826 0.02438 32.63326 19.62 <2e-16 ***
## contrastKIT 0.05291 0.02955 19.96635 1.79 0.0886 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## contrastKIT -0.619B. Generalized Additive Mixed Model
gamData <- windowDatafullModelFile = "target_kitFleece_fullModels.Rdata"
if (file.exists(fullModelFile))
{
load(fullModelFile)
}else{
interaction.target_kitFleece.gamm0 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
s(Time, by = isAustralian) +
s(Time, by = isAustralianFleece) +
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.target_kitFleece.gamm0)
myRhoVal = myAutoCorr[2]
interaction.target_kitFleece.gamm1 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast, k=12) +
s(Time, by = isAustralian, k= 12) +
s(Time, by = isAustralianFleece) +
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.target_kitFleece.gamm0","interaction.target_kitFleece.gamm1"), file = fullModelFile)
}Model criticism: How are we doing regarding auto-correlation?
myAutoCorr =acf_resid(interaction.target_kitFleece.gamm1)
Is the default for the smooth (k = 10) sufficient?
gam.check(interaction.target_kitFleece.gamm1)
##
## Method: fREML Optimizer: perf chol
## $grad
## [1] -4.063860e-12 2.150946e-11 -2.189917e-03 -1.377428e-03 -2.582112e-11
## [6] 3.940315e-11 1.943476e-10 -3.286011e-10 1.482423e-10 1.460698e-10
## [11] -1.557332e-11 1.327738e-11 1.055689e-10 -1.530118e-10 5.630341e-11
## [16] -1.641086e-09
##
## $hess
## [,1] [,2] [,3] [,4] [,5]
## [1,] 2.263304e+00 -1.837073e-04 7.009594e-04 3.679420e-04 9.012070e-02
## [2,] -1.837073e-04 3.124999e+00 9.377747e-04 -1.371116e-06 -1.397019e-05
## [3,] 7.009594e-04 9.377747e-04 2.210810e-03 2.028191e-05 6.695191e-05
## [4,] 3.679420e-04 -1.371116e-06 2.028191e-05 1.388854e-03 2.845913e-05
## [5,] 9.012070e-02 -1.397019e-05 6.695191e-05 2.845913e-05 1.111587e+01
## [6,] 5.478076e-02 -4.345524e-06 3.422655e-05 2.682132e-05 -3.568181e-02
## [7,] -6.503518e-02 4.901514e-06 -1.496909e-05 -7.396163e-06 -5.960655e-01
## [8,] -4.758950e-06 1.556852e-01 1.382991e-05 5.263263e-08 -4.360170e-07
## [9,] -1.300686e-06 2.573962e-02 -2.151230e-06 1.725119e-08 3.967844e-08
## [10,] 1.602039e-06 -5.446143e-02 -2.712552e-05 1.061477e-07 1.651173e-07
## [11,] 1.460777e-01 -1.339621e-05 5.475124e-05 2.903803e-05 -4.583996e-02
## [12,] 8.416323e-03 -5.912161e-07 8.215272e-06 7.494480e-06 -2.994154e-02
## [13,] -2.978886e-02 2.181263e-06 -8.809437e-06 -5.027699e-06 -3.515644e-02
## [14,] -2.561218e-06 4.459439e-02 6.223590e-06 2.416267e-08 -1.886760e-07
## [15,] -4.419886e-07 6.805948e-03 4.275629e-06 -1.616029e-08 -1.943994e-08
## [16,] 5.099160e-08 -1.223311e-02 -4.769035e-06 2.563815e-08 8.134387e-09
## [,6] [,7] [,8] [,9] [,10]
## [1,] 5.478076e-02 -6.503518e-02 -4.758950e-06 -1.300686e-06 1.602039e-06
## [2,] -4.345524e-06 4.901514e-06 1.556852e-01 2.573962e-02 -5.446143e-02
## [3,] 3.422655e-05 -1.496909e-05 1.382991e-05 -2.151230e-06 -2.712552e-05
## [4,] 2.682132e-05 -7.396163e-06 5.263263e-08 1.725119e-08 1.061477e-07
## [5,] -3.568181e-02 -5.960655e-01 -4.360170e-07 3.967844e-08 1.651173e-07
## [6,] 8.434833e+00 -3.678398e-01 -6.830323e-08 -2.019510e-08 8.108250e-08
## [7,] -3.678398e-01 1.193362e+01 1.362458e-07 3.586428e-08 -3.054613e-08
## [8,] -6.830323e-08 1.362458e-07 1.300760e+01 -2.998959e-01 -3.091102e-01
## [9,] -2.019510e-08 3.586428e-08 -2.998959e-01 1.099669e+01 -5.484860e-01
## [10,] 8.108250e-08 -3.054613e-08 -3.091102e-01 -5.484860e-01 1.212450e+01
## [11,] -3.082308e-02 -6.105782e-02 -3.424316e-07 -8.235563e-08 1.374512e-07
## [12,] -3.189537e-02 -2.063602e-02 5.045476e-10 -4.580776e-09 1.894908e-08
## [13,] -1.794652e-02 -5.552779e-02 6.073195e-08 1.203298e-08 -1.820422e-08
## [14,] -4.027986e-08 7.321332e-08 -2.475766e-02 -1.637855e-02 -3.362429e-02
## [15,] -1.585963e-08 8.337106e-09 -2.366714e-02 -2.242471e-02 -1.526583e-02
## [16,] 1.412255e-08 -1.289654e-10 -1.068923e-02 -5.424499e-03 -3.289267e-02
## [,11] [,12] [,13] [,14] [,15]
## [1,] 1.460777e-01 8.416323e-03 -2.978886e-02 -2.561218e-06 -4.419886e-07
## [2,] -1.339621e-05 -5.912161e-07 2.181263e-06 4.459439e-02 6.805948e-03
## [3,] 5.475124e-05 8.215272e-06 -8.809437e-06 6.223590e-06 4.275629e-06
## [4,] 2.903803e-05 7.494480e-06 -5.027699e-06 2.416267e-08 -1.616029e-08
## [5,] -4.583996e-02 -2.994154e-02 -3.515644e-02 -1.886760e-07 -1.943994e-08
## [6,] -3.082308e-02 -3.189537e-02 -1.794652e-02 -4.027986e-08 -1.585963e-08
## [7,] -6.105782e-02 -2.063602e-02 -5.552779e-02 7.321332e-08 8.337106e-09
## [8,] -3.424316e-07 5.045476e-10 6.073195e-08 -2.475766e-02 -2.366714e-02
## [9,] -8.235563e-08 -4.580776e-09 1.203298e-08 -1.637855e-02 -2.242471e-02
## [10,] 1.374512e-07 1.894908e-08 -1.820422e-08 -3.362429e-02 -1.526583e-02
## [11,] 4.297629e+00 -2.462764e-02 -1.542217e-01 -1.774305e-07 -3.417733e-08
## [12,] -2.462764e-02 2.060116e+00 -1.788043e-02 -1.789512e-09 -4.008611e-09
## [13,] -1.542217e-01 -1.788043e-02 3.164284e+00 3.205089e-08 4.303616e-09
## [14,] -1.774305e-07 -1.789512e-09 3.205089e-08 2.392398e+00 -7.592085e-02
## [15,] -3.417733e-08 -4.008611e-09 4.303616e-09 -7.592085e-02 1.481738e+00
## [16,] 8.484911e-09 4.086340e-09 -4.811789e-10 -1.819315e-01 8.522086e-02
## [,16]
## [1,] 5.099160e-08
## [2,] -1.223311e-02
## [3,] -4.769035e-06
## [4,] 2.563815e-08
## [5,] 8.134387e-09
## [6,] 1.412255e-08
## [7,] -1.289655e-10
## [8,] -1.068923e-02
## [9,] -5.424499e-03
## [10,] -3.289267e-02
## [11,] 8.484911e-09
## [12,] 4.086340e-09
## [13,] -4.811789e-10
## [14,] -1.819315e-01
## [15,] 8.522086e-02
## [16,] 1.864432e+00
##
## Model rank = 1346 / 1346
##
## 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):vowelContrastFLEECE 11.00 6.65 0.98 0.090 .
## s(Time):vowelContrastKIT 11.00 7.84 0.98 0.030 *
## s(Time):isAustralian 12.00 2.02 0.98 0.080 .
## s(Time):isAustralianFleece 10.00 2.01 0.98 0.055 .
## s(Time,pp):vowelContrastFLEECE 450.00 113.58 0.98 0.050 *
## s(Time,pp):vowelContrastKIT 450.00 122.15 0.98 0.085 .
## s(Time,item):vowelContrastFLEECE 200.00 32.98 0.98 0.095 .
## s(Time,item):vowelContrastKIT 200.00 20.55 0.98 0.050 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1We now examine this model to determine whether there is an overall effect of contrast.
summary(interaction.target_kitFleece.gamm1)##
## Family: binomial
## Link function: logit
##
## Formula:
## position ~ vowelContrast + s(Time, by = vowelContrast, k = 12) +
## s(Time, by = isAustralian, k = 12) + s(Time, by = isAustralianFleece) +
## 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.19097 0.09737 1.961 0.0499 *
## vowelContrastKIT 0.10656 0.12749 0.836 0.4032
## ---
## 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):vowelContrastFLEECE 6.655 8.298 282.646 <2e-16 ***
## s(Time):vowelContrastKIT 7.837 9.474 375.653 <2e-16 ***
## s(Time):isAustralian 2.017 2.033 1.341 0.441
## s(Time):isAustralianFleece 2.010 2.020 5.197 0.099 .
## s(Time,pp):vowelContrastFLEECE 113.582 448.000 374.373 <2e-16 ***
## s(Time,pp):vowelContrastKIT 122.151 448.000 429.997 <2e-16 ***
## s(Time,item):vowelContrastFLEECE 32.983 98.000 139.294 <2e-16 ***
## s(Time,item):vowelContrastKIT 20.545 98.000 79.259 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.25 Deviance explained = 20%
## fREML = 1.0578e+05 Scale est. = 1 n = 320222The model shows nonlinear time effects and substantial participant and item variability, but no significant main effect of vowel contrast and no reliable accent-related effects, indicating that fixation trajectories change over time but do not significantly differ between KIT and FLEECE or across accents.
smooth <- plot_smooth(interaction.target_kitFleece.gamm1, view = "Time",
plot_all = "vowelContrast",
rug = FALSE, rm.ranef = TRUE, lwd = 4)## Summary:
## * vowelContrast : factor; set to the value(s): FLEECE, KIT.
## * Time : numeric predictor; with 30 values ranging from 10.000000 to 890.000000.
## * isAustralian : numeric predictor; set to the value(s): 0.
## * isAustralianFleece : 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): beetle. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(Time,pp):vowelContrastFLEECE,s(Time,pp):vowelContrastKIT,s(Time,item):vowelContrastFLEECE,s(Time,item):vowelContrastKIT
## 
We now proceed to model comparison by successively removing effects to evaluate their contribution to model fit.
reducedModelFile = "target_kitFleece_reducedModel.Rdata"
if (file.exists(reducedModelFile))
{
load(reducedModelFile)
}else{
mainEffect.target_kitFleece.gamm <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
#s(Time, by = isAustralian) +
#s(Time, by = isAustralianFleece) +
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("mainEffect.target_kitFleece.gamm"), file = reducedModelFile)
}summary(mainEffect.target_kitFleece.gamm)##
## Family: binomial
## Link function: logit
##
## Formula:
## position ~ vowelContrast + s(Time, by = vowelContrast) + 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.17262 0.09188 1.879 0.0603 .
## vowelContrastKIT 0.14288 0.11890 1.202 0.2295
## ---
## 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):vowelContrastFLEECE 6.120 7.343 363.1 <2e-16 ***
## s(Time):vowelContrastKIT 6.958 8.121 415.2 <2e-16 ***
## s(Time,pp):vowelContrastFLEECE 112.706 448.000 690.1 <2e-16 ***
## s(Time,pp):vowelContrastKIT 121.372 448.000 708.4 <2e-16 ***
## s(Time,item):vowelContrastFLEECE 32.751 99.000 197.1 <2e-16 ***
## s(Time,item):vowelContrastKIT 20.477 98.000 250.0 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.252 Deviance explained = 20.2%
## fREML = 1.0551e+05 Scale est. = 1 n = 320222There is no significant overall effect of vowel contrast (KIT vs. FLEECE), as the parametric term for vowelContrastKIT is not significant (p = .23). This indicates that, averaged across the time window, fixation proportions do not differ reliably between KIT and FLEECE items.
smooth <- plot_smooth(mainEffect.target_kitFleece.gamm, view = "Time",
plot_all = "vowelContrast",
rug = FALSE, rm.ranef = TRUE, lwd = 4)## Summary:
## * vowelContrast : factor; set to the value(s): FLEECE, KIT.
## * Time : numeric predictor; with 30 values ranging from 10.000000 to 890.000000.
## * pp : factor; set to the value(s): CAT01. (Might be canceled as random effect, check below.)
## * item : factor; set to the value(s): beetle. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(Time,pp):vowelContrastFLEECE,s(Time,pp):vowelContrastKIT,s(Time,item):vowelContrastFLEECE,s(Time,item):vowelContrastKIT
## 
3.2.2.2. Competitor Fixations
A. TIME-WINDOW ANALYSIS
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("KIT|FLEECE", 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
## 1799 1799trialDataE$contrastC = -0.5 + as.numeric(trialDataE$contrast == "FLEECE")
table(trialDataE$contrastC)##
## -0.5 0.5
## 1799 1799useData = cbind(trialDataE,gazeDataE)
# Remove distractor responses
useData <- useData %>%
filter(useData$clickedOn != "3") %>%
droplevels()We first need to reorganize the data.
windowData = useData[ ,c(1,2,4:8,10:134)] %>%
pivot_longer(., cols = starts_with("T_"),
names_to = "Time",
values_to = "position") %>% as.data.frame
windowData$position = as.numeric(windowData$position)
windowData$Time = gsub("T_","", windowData$Time) %>% as.numeric()
windowData = windowData[ order(windowData$pp,windowData$trial,windowData$Time), ]
windowData$pp = as.factor(windowData$pp)
windowData$item = as.factor(windowData$item)
windowData$isAustralian = ifelse(windowData$accent == "AusE",1,0)
windowData$isAustralianFleece = windowData$isAustralian *
as.numeric(windowData$contrast == "FLEECE")
windowData$vowelContrast = windowData$contrast %>% as.factor()
#restricting the data to 0-900 after vowel onset
windowData = windowData %>% filter(Time > 0 & Time < 900)
windowData <- start_event(windowData, event=c("pp","trial"))
windowDataPlot = windowData %>% group_by(vowelContrast, Time, accent) %>%
summarize(meanFix = mean(position))
Fig.compet_kitFleece <- ggplot(windowDataPlot, aes(x = Time, y = meanFix, col = vowelContrast)) + geom_line() + facet_wrap(~accent)
print(Fig.compet_kitFleece)
In both accents (AusE and SEng), fixation proportions to the competitor rise to a peak around 250–300 ms and then decline, with FLEECE showing consistently slightly higher competitor fixation proportions than KIT across most of the time course, and only minimal differences between accents.
lmerData = windowData %>% filter(Time > 150 & Time < 600) %>%
group_by(pp, item, trial, isAustralian, contrastC) %>% summarize(meanLooks = mean(position))
lmerData$accentC = lmerData$isAustralian - 0.5
outcomes = sort(unique(lmerData$meanLooks))
if (outcomes[1] == 0){
correction = outcomes[2] * 0.5
}
lmerData$logitLooks = ifelse(lmerData$meanLooks == 0, qlogis(correction),
ifelse(lmerData$meanLooks == 1, qlogis(1- correction),
qlogis(lmerData$meanLooks))
)
plot(lmerData$meanLooks, lmerData$logitLooks)
lmerData$trialScaled = scale(lmerData$trial)
compet_fleeceKit.lmer = lmer(logitLooks ~ contrastC*accentC*trialScaled +
(1+contrastC+contrastC||pp) +
(1+accentC*trialScaled||item), data = lmerData,
control=lmerControl(optimizer="bobyqa", optCtrl=list(maxfun=2e5)))
summary(compet_fleeceKit.lmer)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: logitLooks ~ contrastC * accentC * trialScaled + (1 + contrastC +
## contrastC || pp) + (1 + accentC * trialScaled || item)
## Data: lmerData
## Control: lmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 2e+05))
##
## REML criterion at convergence: 16295.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.4813 -0.9655 0.1864 0.8378 3.1328
##
## Random effects:
## Groups Name Variance Std.Dev.
## pp (Intercept) 0.05718 0.2391
## pp.1 contrastC 0.05394 0.2322
## item (Intercept) 0.06489 0.2547
## item.1 accentC 0.12455 0.3529
## item.2 trialScaled 0.02393 0.1547
## item.3 accentC:trialScaled 0.03659 0.1913
## Residual 5.27238 2.2962
## Number of obs: 3598, groups: pp, 45; item, 20
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -2.07524 0.07737 25.47510 -26.823 <2e-16 ***
## contrastC 0.39072 0.14163 19.25931 2.759 0.0124 *
## accentC -0.03409 0.11003 18.01479 -0.310 0.7602
## trialScaled -0.05497 0.05181 17.98592 -1.061 0.3028
## contrastC:accentC 0.16490 0.22007 18.01478 0.749 0.4633
## contrastC:trialScaled 0.03384 0.10367 18.01836 0.326 0.7478
## accentC:trialScaled 0.07335 0.08833 16.74850 0.830 0.4180
## contrastC:accentC:trialScaled -0.14668 0.17666 16.75335 -0.830 0.4180
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) cntrsC accntC trlScl cntC:C cntC:S accC:S
## contrastC 0.000
## accentC 0.000 0.000
## trialScaled 0.000 0.001 -0.014
## cntrstC:ccC 0.000 0.000 0.000 -0.013
## cntrstC:trS 0.001 0.000 -0.013 0.002 -0.014
## accntC:trlS -0.011 -0.012 0.001 0.000 0.002 -0.017
## cntrstC:C:S -0.011 -0.012 0.002 -0.017 0.001 0.001 0.003There is a significant main effect of contrast (p = .0137), indicating that logit-transformed looks differ reliably between the two vowel contrasts, but there are no significant effects of accent, trial, or any interactions, meaning the contrast effect does not vary across accents or over trials.
B. Generalized Additive Mixed Model
gamData <- windowDatafullModelFile = "compet_kitFleece_fullModels.Rdata"
if (file.exists(fullModelFile))
{
load(fullModelFile)
}else{
interaction.compet_kitFleece.gamm0 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
s(Time, by = isAustralian) +
s(Time, by = isAustralianFleece) +
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.compet_kitFleece.gamm0)
myRhoVal = myAutoCorr[2]
interaction.compet_kitFleece.gamm1 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast, k=12) +
s(Time, by = isAustralian, k= 12) +
s(Time, by = isAustralianFleece) +
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.compet_kitFleece.gamm0",
"interaction.compet_kitFleece.gamm1"),
file = fullModelFile)
}Model criticism: How are we doing regarding auto-correlation?
myAutoCorr =acf_resid(interaction.compet_kitFleece.gamm1)
Is the default for the smooth (k = 10) sufficient?
gam.check(interaction.compet_kitFleece.gamm1)
##
## Method: fREML Optimizer: perf chol
## $grad
## [1] 3.488698e-10 -5.758505e-12 -2.080648e-03 -1.490578e-03 5.505285e-11
## [6] 1.487330e-08 -3.012283e-07 -1.595239e-09 5.075282e-10 -9.554002e-11
## [11] 8.704504e-11 5.817733e-10 4.767715e-10 6.443379e-11 -3.904450e-09
## [16] -3.568913e-04
##
## $hess
## [,1] [,2] [,3] [,4] [,5]
## [1,] 3.488878e+00 -7.118330e-05 1.883164e-04 2.086634e-04 2.167088e-02
## [2,] -7.118330e-05 3.575813e+00 1.686063e-04 -2.082987e-07 -2.111883e-07
## [3,] 1.883164e-04 1.686063e-04 2.075738e-03 -1.576270e-06 -4.557843e-07
## [4,] 2.086634e-04 -2.082987e-07 -1.576270e-06 1.488017e-03 -1.190370e-06
## [5,] 2.167088e-02 -2.111883e-07 -4.557843e-07 -1.190370e-06 4.348542e+00
## [6,] -2.969726e-02 8.170725e-08 -1.072585e-06 -1.411217e-06 -3.540216e-01
## [7,] -1.232164e-02 2.251190e-07 -6.386259e-07 -7.342654e-07 2.626423e-02
## [8,] -3.876821e-07 4.655849e-02 3.288944e-07 1.321195e-09 -2.374359e-08
## [9,] 1.586598e-06 -1.017152e-02 -4.960540e-06 1.224518e-08 -1.777662e-08
## [10,] 7.072524e-08 -3.039729e-03 -2.175744e-07 4.662594e-10 -9.964996e-10
## [11,] 9.883974e-02 -2.100073e-06 4.650470e-06 4.748376e-06 -2.874770e-02
## [12,] -9.672676e-03 4.173984e-08 -8.184191e-07 -1.253865e-06 -1.359073e-02
## [13,] -1.559061e-02 3.535761e-07 -9.443235e-07 -1.106510e-06 -2.066704e-02
## [14,] -7.443209e-07 4.232443e-02 1.555389e-06 -2.105582e-09 -6.523833e-09
## [15,] 1.731081e-07 -3.814853e-03 -5.098902e-07 1.004885e-09 -1.211248e-09
## [16,] -2.161490e-09 5.364334e-05 4.601878e-09 -1.030844e-11 -4.210237e-11
## [,6] [,7] [,8] [,9] [,10]
## [1,] -2.969726e-02 -1.232164e-02 -3.876821e-07 1.586598e-06 7.072524e-08
## [2,] 8.170725e-08 2.251190e-07 4.655849e-02 -1.017152e-02 -3.039729e-03
## [3,] -1.072585e-06 -6.386259e-07 3.288944e-07 -4.960540e-06 -2.175744e-07
## [4,] -1.411217e-06 -7.342654e-07 1.321195e-09 1.224518e-08 4.662594e-10
## [5,] -3.540216e-01 2.626423e-02 -2.374359e-08 -1.777662e-08 -9.964995e-10
## [6,] 7.666502e+00 -4.441255e-01 -7.843717e-09 -2.837885e-08 -9.067094e-10
## [7,] -4.441255e-01 2.361734e+00 1.572049e-09 -5.117327e-09 -2.365266e-10
## [8,] -7.843717e-09 1.572049e-09 4.904163e+00 -5.875864e-01 1.540887e-02
## [9,] -2.837885e-08 -5.117327e-09 -5.875864e-01 1.047881e+01 -4.271804e-01
## [10,] -9.067095e-10 -2.365266e-10 1.540887e-02 -4.271804e-01 5.181825e-01
## [11,] -1.716472e-02 -1.542121e-02 -3.435465e-08 3.982536e-08 1.254875e-09
## [12,] -2.324477e-02 -7.385966e-03 -6.764359e-09 -2.193674e-08 -7.124700e-10
## [13,] -1.570228e-02 -1.097573e-02 2.274705e-09 -6.669225e-09 -3.537795e-10
## [14,] 2.169655e-10 2.590221e-09 -1.140304e-02 -1.622682e-02 -6.996366e-04
## [15,] -1.961956e-09 -5.337287e-10 -6.512315e-03 -1.355881e-02 6.580466e-05
## [16,] 1.468009e-11 9.562038e-12 -6.481737e-05 -4.328781e-05 -2.642882e-05
## [,11] [,12] [,13] [,14] [,15]
## [1,] 9.883974e-02 -9.672676e-03 -1.559061e-02 -7.443209e-07 1.731081e-07
## [2,] -2.100073e-06 4.173984e-08 3.535761e-07 4.232443e-02 -3.814853e-03
## [3,] 4.650470e-06 -8.184191e-07 -9.443235e-07 1.555389e-06 -5.098902e-07
## [4,] 4.748376e-06 -1.253865e-06 -1.106510e-06 -2.105582e-09 1.004885e-09
## [5,] -2.874770e-02 -1.359073e-02 -2.066704e-02 -6.523833e-09 -1.211248e-09
## [6,] -1.716472e-02 -2.324477e-02 -1.570228e-02 2.169655e-10 -1.961956e-09
## [7,] -1.542121e-02 -7.385966e-03 -1.097573e-02 2.590221e-09 -5.337287e-10
## [8,] -3.435465e-08 -6.764359e-09 2.274705e-09 -1.140304e-02 -6.512315e-03
## [9,] 3.982536e-08 -2.193674e-08 -6.669225e-09 -1.622682e-02 -1.355881e-02
## [10,] 1.254875e-09 -7.124701e-10 -3.537795e-10 -6.996366e-04 6.580466e-05
## [11,] 2.962036e+00 7.366857e-02 9.789900e-02 -3.024778e-08 4.672609e-09
## [12,] 7.366857e-02 2.947388e+00 -2.420712e-01 -4.748003e-10 -1.515864e-09
## [13,] 9.789900e-02 -2.420712e-01 3.108273e+00 4.003006e-09 -7.624866e-10
## [14,] -3.024778e-08 -4.748003e-10 4.003006e-09 1.034150e+00 -3.205520e-02
## [15,] 4.672609e-09 -1.515864e-09 -7.624866e-10 -3.205520e-02 9.420318e-01
## [16,] -1.547068e-10 5.904136e-12 1.255166e-11 6.604327e-04 -2.008086e-03
## [,16]
## [1,] -2.161490e-09
## [2,] 5.364334e-05
## [3,] 4.601878e-09
## [4,] -1.030844e-11
## [5,] -4.210237e-11
## [6,] 1.468009e-11
## [7,] 9.562038e-12
## [8,] -6.481737e-05
## [9,] -4.328781e-05
## [10,] -2.642882e-05
## [11,] -1.547068e-10
## [12,] 5.904135e-12
## [13,] 1.255166e-11
## [14,] 6.604327e-04
## [15,] -2.008086e-03
## [16,] 1.875207e-03
##
## Model rank = 1346 / 1346
##
## 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):vowelContrastFLEECE 11.0 7.9 0.96 0.24
## s(Time):vowelContrastKIT 11.0 7.8 0.96 0.18
## s(Time):isAustralian 12.0 2.0 0.96 0.21
## s(Time):isAustralianFleece 10.0 2.0 0.96 0.27
## s(Time,pp):vowelContrastFLEECE 450.0 76.0 0.96 0.23
## s(Time,pp):vowelContrastKIT 450.0 70.5 0.96 0.24
## s(Time,item):vowelContrastFLEECE 200.0 30.9 0.96 0.23
## s(Time,item):vowelContrastKIT 200.0 13.8 0.96 0.20We now examine this model to determine whether there is an overall effect of contrast.
summary(interaction.compet_kitFleece.gamm1)##
## Family: binomial
## Link function: logit
##
## Formula:
## position ~ vowelContrast + s(Time, by = vowelContrast, k = 12) +
## s(Time, by = isAustralian, k = 12) + s(Time, by = isAustralianFleece) +
## 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.44524 0.09834 -14.70 <2e-16 ***
## vowelContrastKIT -0.13714 0.11237 -1.22 0.222
## ---
## 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):vowelContrastFLEECE 7.901 9.481 119.945 <2e-16 ***
## s(Time):vowelContrastKIT 7.799 9.435 163.014 <2e-16 ***
## s(Time):isAustralian 2.005 2.010 1.451 0.488
## s(Time):isAustralianFleece 2.004 2.007 0.930 0.582
## s(Time,pp):vowelContrastFLEECE 75.973 448.000 346.049 <2e-16 ***
## s(Time,pp):vowelContrastKIT 70.463 449.000 331.189 <2e-16 ***
## s(Time,item):vowelContrastFLEECE 30.876 98.000 141.688 <2e-16 ***
## s(Time,item):vowelContrastKIT 13.791 98.000 84.218 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.0875 Deviance explained = 9.41%
## fREML = 85505 Scale est. = 1 n = 320222The model shows significant nonlinear time effects and substantial participant and item variability, but no significant main effect of vowel contrast and no significant accent-related smooths, indicating that fixation trajectories change over time but do not reliably differ between KIT and FLEECE or across accents.
smooth <- plot_smooth(interaction.compet_kitFleece.gamm1, view = "Time",
plot_all = "vowelContrast",
rug = FALSE, rm.ranef = TRUE, lwd = 4)## Summary:
## * vowelContrast : factor; set to the value(s): FLEECE, KIT.
## * Time : numeric predictor; with 30 values ranging from 10.000000 to 890.000000.
## * isAustralian : numeric predictor; set to the value(s): 0.
## * isAustralianFleece : 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): beetle. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(Time,pp):vowelContrastFLEECE,s(Time,pp):vowelContrastKIT,s(Time,item):vowelContrastFLEECE,s(Time,item):vowelContrastKIT
## 
We now proceed to model comparison by successively removing effects to evaluate their contribution to model fit.
reducedModelFile = "compet_kitFleece_reducedModel.Rdata"
if (file.exists(reducedModelFile))
{
load(reducedModelFile)
}else{
mainEffect.compet_kitFleece.gamm <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
#s(Time, by = isAustralian) +
#s(Time, by = isAustralianFleece) +
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("mainEffect.compet_kitFleece.gamm"), file = reducedModelFile)
}summary(mainEffect.compet_kitFleece.gamm)##
## Family: binomial
## Link function: logit
##
## Formula:
## position ~ vowelContrast + s(Time, by = vowelContrast) + 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.43386 0.09112 -15.735 <2e-16 ***
## vowelContrastKIT -0.17853 0.09986 -1.788 0.0738 .
## ---
## 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):vowelContrastFLEECE 7.172 8.234 125.0 <2e-16 ***
## s(Time):vowelContrastKIT 7.107 8.241 165.5 <2e-16 ***
## s(Time,pp):vowelContrastFLEECE 75.625 448.000 578.8 <2e-16 ***
## s(Time,pp):vowelContrastKIT 70.291 448.000 529.8 <2e-16 ***
## s(Time,item):vowelContrastFLEECE 30.732 99.000 163.3 <2e-16 ***
## s(Time,item):vowelContrastKIT 13.667 98.000 158.2 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.0899 Deviance explained = 9.61%
## fREML = 85259 Scale est. = 1 n = 320222This model shows strong nonlinear time effects and substantial participant and item variability, with only a marginal main effect of vowel contrast (p = .074), suggesting weak evidence that KIT and FLEECE differ overall, but no indication of additional effects beyond contrast-specific time trajectories.
smooth <- plot_smooth(mainEffect.compet_kitFleece.gamm, view = "Time",
plot_all = "vowelContrast",
rug = FALSE, rm.ranef = TRUE, lwd = 4)## Summary:
## * vowelContrast : factor; set to the value(s): FLEECE, KIT.
## * Time : numeric predictor; with 30 values ranging from 10.000000 to 890.000000.
## * pp : factor; set to the value(s): CAT01. (Might be canceled as random effect, check below.)
## * item : factor; set to the value(s): beetle. (Might be canceled as random effect, check below.)
## * NOTE : The following random effects columns are canceled: s(Time,pp):vowelContrastFLEECE,s(Time,pp):vowelContrastKIT,s(Time,item):vowelContrastFLEECE,s(Time,item):vowelContrastKIT
## 
3.2.2.3. TARGET-COMPETITOR FIXATIONS
A. TIME-WINDOW ANALYSIS
rm(list = ls())
library(lmerTest)
library(stringr)
library(tidyverse)
library(mgcv)
library(itsadug)
options(warn = 0)
gazeTarget = read.table("dataCatalan/gazeTarget.txt",
header = TRUE, sep = "\t", quote = "")
gazeComp = read.table("dataCatalan/gazeCompet.txt",
header = TRUE, sep = "\t", quote = "")
trialData = read.table("dataCatalan/trialInfo.txt",
header = TRUE, sep = "\t", quote = "\"")
expFilter = grepl("KIT|FLEECE", trialData$condition)
gazeTargetE = gazeTarget %>% filter(expFilter)
gazeCompE = gazeComp %>% filter(expFilter)
trialDataE = trialData %>% filter(expFilter)
timeline = -200 + (0:120) * 10
colnames(gazeTargetE) = paste0("Tt_", timeline)
colnames(gazeCompE) = paste0("Tc_", 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
## 1799 1799trialDataE$contrastC = -0.5 + as.numeric(trialDataE$contrast == "FLEECE")
table(trialDataE$contrastC)##
## -0.5 0.5
## 1799 1799useData = cbind(trialDataE, gazeTargetE, gazeCompE)
# Remove distractor responses
useData <- useData %>%
filter(useData$clickedOn != "3") %>%
droplevels()windowData = useData %>%
pivot_longer(
cols = matches("^T(t|c)_-?\\d+$"),
names_to = c("which", "Time"),
names_pattern = "^T(t|c)_(-?\\d+)$",
values_to = "position"
) %>%
pivot_wider(
names_from = which,
values_from = position
) %>%
rename(
target = t,
competitor = c
) %>%
mutate(
target = as.numeric(target),
competitor = as.numeric(competitor),
difference = target - competitor,
Time = as.numeric(Time)
) %>%
as.data.frame()
windowData = windowData[order(windowData$pp, windowData$trial, windowData$Time), ]
windowData$pp = as.factor(windowData$pp)
windowData$item = as.factor(windowData$item)
windowData$isAustralian = ifelse(windowData$accent == "AusE", 1, 0)
windowData$isAustralianFleece = windowData$isAustralian *
as.numeric(windowData$contrast == "FLEECE")
windowData$vowelContrast = windowData$contrast %>% as.factor()
# restricting the data to 0-900 after vowel onset
windowData = windowData %>% filter(Time > 0 & Time < 900)
windowData <- start_event(windowData, event = c("pp", "trial"))windowDataPlot = windowData %>%
group_by(vowelContrast, Time, accent) %>%
summarize(meanDiff = mean(difference, na.rm = TRUE), .groups = "drop")
Fig.diff_kitFleece = ggplot(
windowDataPlot,
aes(x = Time, y = meanDiff, col = vowelContrast)
) +
geom_line() +
facet_wrap(~accent)
print(Fig.diff_kitFleece)
aggregatedData = windowData %>%
filter(Time > 0 & Time < 100) %>%
group_by(pp, item, contrast) %>%
summarize(meanDiff = mean(difference, na.rm = TRUE), .groups = "drop") %>%
ungroup()targetCompetitor_fleeceKit.lmer = lmer(
meanDiff ~ contrast +
(1 + contrast | pp) +
(1 | item),
data = aggregatedData
)## boundary (singular) fit: see help('isSingular')summary(targetCompetitor_fleeceKit.lmer)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: meanDiff ~ contrast + (1 + contrast | pp) + (1 | item)
## Data: aggregatedData
##
## REML criterion at convergence: 541.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.1262 -0.7659 0.0046 0.7196 2.9592
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## pp (Intercept) 1.391e-03 3.730e-02
## contrastKIT 1.423e-03 3.772e-02 -1.00
## item (Intercept) 1.702e-10 1.304e-05
## Residual 1.049e-01 3.239e-01
## Number of obs: 900, groups: pp, 45; item, 20
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 0.02358 0.01625 44.00393 1.451 0.154
## contrastKIT -0.02525 0.02231 90.12234 -1.132 0.261
##
## Correlation of Fixed Effects:
## (Intr)
## contrastKIT -0.729
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')3.2.3 DRESS-TRAP contrast
3.2.3.1. Target Fixations
A. TIME-WINDOW ANALYSIS
rm(list = ls())
library(lmerTest)
library(stringr)
library(tidyverse)
library(mgcv)
library(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("DRESS|TRAP", 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
## 1797 1799trialDataE$contrastC = -0.5 + as.numeric(trialDataE$contrast == "TRAP")
table(trialDataE$contrastC)##
## -0.5 0.5
## 1798 1798useData = cbind(trialDataE,gazeDataE)
# Remove distractor responses
useData <- useData %>%
filter(useData$clickedOn != "3") %>%
droplevels()
We first need to reorganize the data.
windowData = useData[ ,c(1,2,4:8,10:134)] %>%
pivot_longer(., cols = starts_with("T_"),
names_to = "Time",
values_to = "position") %>% as.data.frame
windowData$position = as.numeric(windowData$position)
windowData$Time = gsub("T_","", windowData$Time) %>% as.numeric()
windowData = windowData[ order(windowData$pp,windowData$trial,windowData$Time), ]
windowData$pp = as.factor(windowData$pp)
windowData$item = as.factor(windowData$item)
windowData$isAustralian = ifelse(windowData$accent == "AusE",1,0)
windowData$isAustralianTrap = windowData$isAustralian *
as.numeric(windowData$contrast == "TRAP")
windowData$vowelContrast = windowData$contrast %>% as.factor()
#restricting the data to 0-900 after vowel onset
windowData = windowData %>% filter(Time > 0 & Time < 900)
windowData <- start_event(windowData, event=c("pp","trial"))
windowDataPlot = windowData %>% group_by(vowelContrast, Time, accent) %>%
summarize(meanFix = mean(position))
Fig.target_dressTrap <- ggplot(windowDataPlot, aes(x = Time, y = meanFix, col = vowelContrast)) + geom_line() + facet_wrap(~accent)
print(Fig.target_dressTrap)
In both accents (AusE and SEng), fixation proportions increase steadily over time toward the target, reaching high levels by around 700–800 ms, with only minimal differences between accents.
aggregatedData = windowData %>% filter(Time > 150 & Time < 600) %>%
group_by(pp, item, contrast) %>% summarize(meanLooks = mean(position)) %>%
ungroup()
meanLooksValues = sort(unique(aggregatedData$meanLooks))
correction = meanLooksValues[2]/2
aggregatedData$logitLooks = ifelse( aggregatedData$meanLooks == 0,
qlogis(correction),
ifelse(aggregatedData$meanLooks == 1,
qlogis(1 - correction),
qlogis(aggregatedData$meanLooks)
)
)
plot(aggregatedData$meanLooks, aggregatedData$logitLooks)
library(lmerTest)
target.dressTrap.lmer = lmer(meanLooks ~ contrast +
(1+contrast|pp)+
(1|item), data = aggregatedData)summary(target.dressTrap.lmer)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: meanLooks ~ contrast + (1 + contrast | pp) + (1 | item)
## Data: aggregatedData
##
## REML criterion at convergence: -535.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.59744 -0.61093 0.01573 0.69347 2.79841
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## pp (Intercept) 9.349e-03 0.09669
## contrastTRAP 1.109e-05 0.00333 1.00
## item (Intercept) 2.174e-03 0.04663
## Residual 2.789e-02 0.16699
## Number of obs: 900, groups: pp, 45; item, 20
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 0.46624 0.02207 39.84722 21.124 <2e-16 ***
## contrastTRAP 0.03493 0.02364 18.01478 1.477 0.157
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## contrstTRAP -0.522
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')This model shows no significant effect of contrast (p = .144), indicating that mean target looks do not reliably differ between DRESS and TRAP.
B. Generalized Additive Mixed Model
gamData <- windowDatafullModelFile = "target_dressTrap_fullModels.Rdata"
if (file.exists(fullModelFile))
{
load(fullModelFile)
}else{
interaction.target_dressTrap.gamm0 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
s(Time, by = isAustralian) +
s(Time, by = isAustralianTrap) +
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.target_dressTrap.gamm0)
myRhoVal = myAutoCorr[2]
interaction.target_dressTrap.gamm1 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast, k=12) +
s(Time, by = isAustralian, k= 12) +
s(Time, by = isAustralianTrap) +
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.target_dressTrap.gamm0","interaction.target_dressTrap.gamm1"), file = fullModelFile)
}Model criticism: How are we doing regarding auto-correlation?
myAutoCorr =acf_resid(interaction.target_dressTrap.gamm1)
Is the default for the smooth (k = 10) sufficient?
gam.check(interaction.target_dressTrap.gamm1)
##
## Method: fREML Optimizer: perf chol
## $grad
## [1] 1.183631e-10 -7.255085e-12 3.165066e-10 -3.094613e-09 -2.243315e-09
## [6] 9.344348e-11 7.302603e-10 9.935519e-11 -2.609125e-09 6.715339e-10
## [11] -3.522258e-10 2.084111e-11 -1.455384e-03 -2.872671e-10 -4.936589e-10
## [16] 1.972151e-10
##
## $hess
## [,1] [,2] [,3] [,4] [,5]
## [1,] 2.574268e+00 -5.434260e-02 2.649630e-01 -8.803984e-03 1.887122e-01
## [2,] -5.434260e-02 2.450078e+00 1.084269e-01 8.760641e-03 -4.666408e-03
## [3,] 2.649630e-01 1.084269e-01 9.398108e-01 -7.014982e-02 1.006118e-02
## [4,] -8.803984e-03 8.760641e-03 -7.014982e-02 2.701442e-01 -3.473929e-03
## [5,] 1.887122e-01 -4.666408e-03 1.006118e-02 -3.473929e-03 1.271041e+01
## [6,] 5.117804e-02 3.552345e-04 -1.012031e-05 5.945799e-04 1.622991e-01
## [7,] -6.316770e-02 9.099576e-04 -1.145248e-02 -3.272043e-03 -4.124715e-01
## [8,] -3.760602e-03 1.891757e-01 8.092570e-03 1.884310e-03 -3.877303e-04
## [9,] -5.477014e-05 9.016781e-03 -6.444486e-04 -2.487585e-03 4.990290e-05
## [10,] 1.608259e-03 -4.416300e-02 -5.936082e-03 -5.546508e-04 1.097028e-05
## [11,] 6.783641e-02 -1.599474e-03 4.015049e-03 -7.961311e-04 -5.025661e-02
## [12,] 3.726659e-03 2.084986e-04 -6.187737e-04 1.300561e-04 -4.135908e-02
## [13,] -4.593874e-06 2.214217e-07 2.486480e-07 1.330128e-07 -6.432286e-06
## [14,] -1.721535e-03 5.819342e-02 3.637813e-03 1.159901e-03 -1.627980e-04
## [15,] -5.383029e-05 2.447078e-03 2.542820e-05 -4.183228e-04 7.704877e-06
## [16,] 3.664374e-04 -1.101076e-02 -1.171686e-03 9.442004e-05 -7.479329e-06
## [,6] [,7] [,8] [,9] [,10]
## [1,] 5.117804e-02 -6.316770e-02 -3.760602e-03 -5.477014e-05 1.608259e-03
## [2,] 3.552345e-04 9.099576e-04 1.891757e-01 9.016781e-03 -4.416300e-02
## [3,] -1.012031e-05 -1.145248e-02 8.092570e-03 -6.444486e-04 -5.936082e-03
## [4,] 5.945799e-04 -3.272043e-03 1.884310e-03 -2.487585e-03 -5.546508e-04
## [5,] 1.622991e-01 -4.124715e-01 -3.877303e-04 4.990290e-05 1.097028e-05
## [6,] 8.096111e+00 -7.784583e-01 3.348365e-05 -1.208956e-05 1.191196e-05
## [7,] -7.784583e-01 1.308184e+01 2.447178e-05 3.448026e-05 -7.966224e-05
## [8,] 3.348365e-05 2.447178e-05 1.767177e+01 9.969048e-02 -5.843016e-01
## [9,] -1.208956e-05 3.448026e-05 9.969048e-02 5.975782e+00 -6.346378e-01
## [10,] 1.191196e-05 -7.966224e-05 -5.843016e-01 -6.346378e-01 1.190254e+01
## [11,] -2.834727e-02 -4.203984e-02 -9.196062e-05 3.092862e-06 2.871374e-05
## [12,] -2.767100e-02 -2.278798e-02 1.772380e-05 -2.540305e-06 -2.513223e-07
## [13,] -5.124378e-06 -2.332330e-06 1.632504e-08 9.771143e-10 -2.963689e-09
## [14,] 1.311027e-05 1.255238e-05 -7.468771e-02 -2.376021e-02 -5.689310e-02
## [15,] -2.791841e-06 7.956570e-06 -4.272985e-02 -2.081587e-02 -1.568665e-02
## [16,] 4.417882e-06 -2.320583e-05 -6.473264e-02 -1.538814e-02 -6.360023e-02
## [,11] [,12] [,13] [,14] [,15]
## [1,] 6.783641e-02 3.726659e-03 -4.593874e-06 -1.721535e-03 -5.383029e-05
## [2,] -1.599474e-03 2.084986e-04 2.214217e-07 5.819342e-02 2.447078e-03
## [3,] 4.015049e-03 -6.187737e-04 2.486480e-07 3.637813e-03 2.542820e-05
## [4,] -7.961311e-04 1.300561e-04 1.330128e-07 1.159901e-03 -4.183228e-04
## [5,] -5.025661e-02 -4.135908e-02 -6.432286e-06 -1.627980e-04 7.704877e-06
## [6,] -2.834727e-02 -2.767100e-02 -5.124378e-06 1.311027e-05 -2.791841e-06
## [7,] -4.203984e-02 -2.278798e-02 -2.332330e-06 1.255238e-05 7.956570e-06
## [8,] -9.196062e-05 1.772380e-05 1.632504e-08 -7.468771e-02 -4.272985e-02
## [9,] 3.092862e-06 -2.540305e-06 9.771143e-10 -2.376021e-02 -2.081587e-02
## [10,] 2.871374e-05 -2.513223e-07 -2.963689e-09 -5.689310e-02 -1.568665e-02
## [11,] 2.708005e+00 -8.205732e-02 -6.407505e-05 -4.312125e-05 2.201530e-07
## [12,] -8.205732e-02 2.014504e+00 -6.391147e-05 7.179083e-06 -4.001709e-07
## [13,] -6.407505e-05 -6.391147e-05 1.455200e-03 7.100269e-09 4.685097e-10
## [14,] -4.312125e-05 7.179083e-06 7.100269e-09 3.332256e+00 -6.971286e-02
## [15,] 2.201530e-07 -4.001709e-07 4.685097e-10 -6.971286e-02 1.509619e+00
## [16,] 6.603659e-06 3.882613e-07 -6.249999e-10 -1.166019e-01 -3.257821e-02
## [,16]
## [1,] 3.664374e-04
## [2,] -1.101076e-02
## [3,] -1.171686e-03
## [4,] 9.442004e-05
## [5,] -7.479329e-06
## [6,] 4.417882e-06
## [7,] -2.320583e-05
## [8,] -6.473264e-02
## [9,] -1.538814e-02
## [10,] -6.360023e-02
## [11,] 6.603659e-06
## [12,] 3.882613e-07
## [13,] -6.249999e-10
## [14,] -1.166019e-01
## [15,] -3.257821e-02
## [16,] 3.139940e+00
##
## Model rank = 1346 / 1346
##
## 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):vowelContrastDRESS 11.00 7.36 0.97 0.080 .
## s(Time):vowelContrastTRAP 11.00 8.18 0.97 0.050 *
## s(Time):isAustralian 12.00 5.69 0.97 0.060 .
## s(Time):isAustralianTrap 10.00 3.22 0.97 0.055 .
## s(Time,pp):vowelContrastDRESS 450.00 116.62 0.97 0.025 *
## s(Time,pp):vowelContrastTRAP 450.00 124.02 0.97 0.055 .
## s(Time,item):vowelContrastDRESS 200.00 17.26 0.97 0.055 .
## s(Time,item):vowelContrastTRAP 200.00 26.31 0.97 0.060 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1We now examine this model to determine whether there is an overall effect of contrast.
summary(interaction.target_dressTrap.gamm1)##
## Family: binomial
## Link function: logit
##
## Formula:
## position ~ vowelContrast + s(Time, by = vowelContrast, k = 12) +
## s(Time, by = isAustralian, k = 12) + s(Time, by = isAustralianTrap) +
## 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.16880 0.07426 2.273 0.023 *
## vowelContrastTRAP 0.12705 0.12349 1.029 0.304
## ---
## 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):vowelContrastDRESS 7.357 8.984 342.797 <2e-16 ***
## s(Time):vowelContrastTRAP 8.182 9.649 418.874 <2e-16 ***
## s(Time):isAustralian 5.687 7.083 13.663 0.0707 .
## s(Time):isAustralianTrap 3.223 3.875 3.645 0.4976
## s(Time,pp):vowelContrastDRESS 116.621 449.000 670.515 <2e-16 ***
## s(Time,pp):vowelContrastTRAP 124.019 448.000 628.553 <2e-16 ***
## s(Time,item):vowelContrastDRESS 17.261 98.000 344.914 <2e-16 ***
## s(Time,item):vowelContrastTRAP 26.309 98.000 181.079 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.271 Deviance explained = 21.8%
## fREML = 1.0619e+05 Scale est. = 1 n = 319955The model shows strong nonlinear time effects and substantial participant and item variability, but no significant main effect of vowel contrast* and no reliable accent-related effects** (aside from a marginal overall time-by-accent smooth), indicating that while fixation trajectories change robustly over time, they do not significantly differ between DRESS and TRAP or across accents.
smooth <- plot_smooth(interaction.target_dressTrap.gamm1, view = "Time",
plot_all = "vowelContrast",
rug = FALSE, rm.ranef = TRUE, lwd = 4)## Summary:
## * vowelContrast : factor; set to the value(s): DRESS, TRAP.
## * Time : numeric predictor; with 30 values ranging from 10.000000 to 890.000000.
## * isAustralian : numeric predictor; set to the value(s): 0.
## * isAustralianTrap : 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):vowelContrastDRESS,s(Time,pp):vowelContrastTRAP,s(Time,item):vowelContrastDRESS,s(Time,item):vowelContrastTRAP
## 
We now proceed to model comparison by successively removing effects to evaluate their contribution to model fit.
reducedModelFile = "target_dressTrap_reducedModel.Rdata"
if (file.exists(reducedModelFile))
{
load(reducedModelFile)
}else{
mainEffect.target_dressTrap.gamm <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
#s(Time, by = isAustralian) +
#s(Time, by = isAustralianTrap) +
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("mainEffect.target_dressTrap.gamm"), file = reducedModelFile)
}summary(mainEffect.target_dressTrap.gamm)##
## Family: binomial
## Link function: logit
##
## Formula:
## position ~ vowelContrast + s(Time, by = vowelContrast) + 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.20684 0.06785 3.048 0.0023 **
## vowelContrastTRAP 0.11564 0.11565 1.000 0.3173
## ---
## 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):vowelContrastDRESS 7.284 8.341 444.6 <2e-16 ***
## s(Time):vowelContrastTRAP 7.267 8.290 557.0 <2e-16 ***
## s(Time,pp):vowelContrastDRESS 116.395 448.000 671.2 <2e-16 ***
## s(Time,pp):vowelContrastTRAP 124.215 448.000 629.2 <2e-16 ***
## s(Time,item):vowelContrastDRESS 17.155 98.000 345.6 <2e-16 ***
## s(Time,item):vowelContrastTRAP 26.379 98.000 183.5 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.271 Deviance explained = 21.8%
## fREML = 1.0629e+05 Scale est. = 1 n = 319955This model shows strong nonlinear time effects and substantial participant and item variability, but no significant main effect of vowel contrast, indicating that although fixation proportions change robustly over time, they do not reliably differ between DRESS and TRAP overall.
smooth <- plot_smooth(mainEffect.target_dressTrap.gamm, view = "Time",
plot_all = "vowelContrast",
rug = FALSE, rm.ranef = TRUE, lwd = 4)## Summary:
## * vowelContrast : factor; set to the value(s): DRESS, TRAP.
## * Time : numeric predictor; with 30 values ranging from 10.000000 to 890.000000.
## * 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):vowelContrastDRESS,s(Time,pp):vowelContrastTRAP,s(Time,item):vowelContrastDRESS,s(Time,item):vowelContrastTRAP
## 
3.2.3.2. Competitor Fixations
A. TIME-WINDOW ANALYSIS
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("DRESS|TRAP", 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
## 1797 1799trialDataE$contrastC = -0.5 + as.numeric(trialDataE$contrast == "TRAP")
table(trialDataE$contrastC)##
## -0.5 0.5
## 1798 1798useData = cbind(trialDataE,gazeDataE)
# Remove distractor responses
useData <- useData %>%
filter(useData$clickedOn != "3") %>%
droplevels()We first need to reorganize the data.
windowData = useData[ ,c(1,2,4:8,10:134)] %>%
pivot_longer(., cols = starts_with("T_"),
names_to = "Time",
values_to = "position") %>% as.data.frame
windowData$position = as.numeric(windowData$position)
windowData$Time = gsub("T_","", windowData$Time) %>% as.numeric()
windowData = windowData[ order(windowData$pp,windowData$trial,windowData$Time), ]
windowData$pp = as.factor(windowData$pp)
windowData$item = as.factor(windowData$item)
windowData$isAustralian = ifelse(windowData$accent == "AusE",1,0)
windowData$isAustralianTrap = windowData$isAustralian *
as.numeric(windowData$contrast == "TRAP")
windowData$vowelContrast = windowData$contrast %>% as.factor()
#restricting the data to 0-900 after vowel onset
windowData = windowData %>% filter(Time > 0 & Time < 900)
windowData <- start_event(windowData, event=c("pp","trial"))
windowDataPlot = windowData %>% group_by(vowelContrast, Time, accent) %>%
summarize(meanFix = mean(position))
Fig.compet_dressTrap <- ggplot(windowDataPlot, aes(x = Time, y = meanFix, col = vowelContrast)) + geom_line() + facet_wrap(~accent)
print(Fig.compet_dressTrap)
In both accents (AusE and SEng), competitor fixations rise to a peak around 250–300 ms and then decline over time, with DRESS showing slightly higher fixation proportions than TRAP (particularly in SEng), while overall accent differences appear small.
lmerData = windowData %>% filter(Time > 150 & Time < 600) %>%
group_by(pp, item, trial, isAustralian, contrastC) %>% summarize(meanLooks = mean(position))
lmerData$accentC = lmerData$isAustralian - 0.5
outcomes = sort(unique(lmerData$meanLooks))
if (outcomes[1] == 0){
correction = outcomes[2] * 0.5
}
lmerData$logitLooks = ifelse(lmerData$meanLooks == 0, qlogis(correction),
ifelse(lmerData$meanLooks == 1, qlogis(1- correction),
qlogis(lmerData$meanLooks))
)
plot(lmerData$meanLooks, lmerData$logitLooks)
library(lmerTest)
lmerData$trialScaled = scale(lmerData$trial)
compet_dressTrap.lmer = lmer(logitLooks ~ contrastC*accentC*trialScaled +
(1+contrastC+contrastC||pp) +
(1+accentC||item), data = lmerData,
control=lmerControl(optimizer="bobyqa", optCtrl=list(maxfun=2e5)))## boundary (singular) fit: see help('isSingular')summary(compet_dressTrap.lmer)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: logitLooks ~ contrastC * accentC * trialScaled + (1 + contrastC +
## contrastC || pp) + (1 + accentC || item)
## Data: lmerData
## Control: lmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 2e+05))
##
## REML criterion at convergence: 15970.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.2116 -0.9387 -0.2185 0.9229 3.2012
##
## Random effects:
## Groups Name Variance Std.Dev.
## pp (Intercept) 0.049563 0.22263
## pp.1 contrastC 0.011949 0.10931
## item (Intercept) 0.006976 0.08352
## item.1 accentC 0.000000 0.00000
## Residual 4.904899 2.21470
## Number of obs: 3595, groups: pp, 45; item, 20
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -2.341e+00 5.307e-02 2.624e+01 -44.105 <2e-16
## contrastC -1.701e-01 8.441e-02 1.503e+01 -2.016 0.0621
## accentC -9.762e-02 7.392e-02 3.481e+03 -1.321 0.1867
## trialScaled 3.957e-03 3.710e-02 3.567e+03 0.107 0.9151
## contrastC:accentC 2.995e-01 1.478e-01 3.481e+03 2.026 0.0428
## contrastC:trialScaled 1.157e-02 7.421e-02 3.567e+03 0.156 0.8761
## accentC:trialScaled 4.109e-02 7.422e-02 3.571e+03 0.554 0.5799
## contrastC:accentC:trialScaled 1.542e-01 1.485e-01 3.574e+03 1.038 0.2992
##
## (Intercept) ***
## contrastC .
## accentC
## trialScaled
## contrastC:accentC *
## contrastC:trialScaled
## accentC:trialScaled
## contrastC:accentC:trialScaled
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) cntrsC accntC trlScl cntC:C cntC:S accC:S
## contrastC 0.000
## accentC 0.000 0.000
## trialScaled 0.001 -0.015 -0.027
## cntrstC:ccC 0.000 0.000 0.001 -0.008
## cntrstC:trS -0.012 0.001 -0.008 -0.028 -0.027
## accntC:trlS -0.019 -0.007 0.001 -0.005 -0.018 -0.027
## cntrstC:C:S -0.006 -0.023 -0.018 -0.028 0.001 -0.004 -0.028
## optimizer (bobyqa) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')This model shows a significant main effect of contrast (p = .040), indicating that competitor looks differ reliably between DRESS and TRAP, but no significant effects of accent, trial, or higher-order interactions, suggesting that the contrast effect does not vary across accents or over trials (with only a marginal contrast-by-accent interaction, p = .074).
B. Generalized Additive Mixed Model
gamData <- windowDatafullModelFile = "compet_dressTrap_fullModels.Rdata"
if (file.exists(fullModelFile))
{
load(fullModelFile)
}else{
interaction.compet_dressTrap.gamm0 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
s(Time, by = isAustralian) +
s(Time, by = isAustralianTrap) +
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.compet_dressTrap.gamm0)
myRhoVal = myAutoCorr[2]
interaction.compet_dressTrap.gamm1 <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast, k=12) +
s(Time, by = isAustralian, k= 12) +
s(Time, by = isAustralianTrap) +
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.compet_dressTrap.gamm0",
"interaction.compet_dressTrap.gamm1"),
file = fullModelFile)
}Model criticism: How are we doing regarding auto-correlation?
myAutoCorr =acf_resid(interaction.compet_dressTrap.gamm1)
Is the default for the smooth (k = 10) sufficient?
gam.check(interaction.compet_dressTrap.gamm1)
##
## Method: fREML Optimizer: perf chol
## $grad
## [1] -1.178688e-05 -6.989587e-06 -7.486405e-04 -9.609024e-06 -3.310650e-05
## [6] -3.820469e-06 -5.435993e-04 3.309236e-05 5.882892e-06 1.738735e-04
## [11] -1.615667e-06 -2.302186e-06 -1.821877e-06 -7.315706e-04 4.676936e-06
## [16] 1.183560e-06
##
## $hess
## [,1] [,2] [,3] [,4] [,5]
## [1,] 3.133991e+00 -2.823719e-05 5.093472e-05 -9.412601e-05 8.276052e-02
## [2,] -2.823719e-05 1.891323e+00 7.661014e-05 2.487650e-01 -8.345986e-07
## [3,] 5.093472e-05 7.661014e-05 7.484196e-04 -7.946900e-05 -1.568436e-06
## [4,] -9.412601e-05 2.487650e-01 -7.946900e-05 6.345066e-01 -2.682517e-06
## [5,] 8.276052e-02 -8.345986e-07 -1.568436e-06 -2.682517e-06 7.500028e+00
## [6,] -3.876187e-02 2.955259e-07 1.117905e-05 2.377030e-06 -9.652046e-01
## [7,] -5.053342e-06 4.778706e-11 3.747143e-09 6.158777e-10 1.777274e-05
## [8,] -1.613881e-06 7.931302e-02 5.458145e-06 3.482589e-02 -7.785824e-08
## [9,] 6.021735e-08 -3.016254e-02 2.430842e-06 1.928380e-03 5.325279e-09
## [10,] -1.875009e-10 4.448383e-05 -3.336611e-09 4.859196e-06 -4.228263e-11
## [11,] 2.862649e-02 -3.153236e-07 3.901703e-07 -1.130433e-06 -6.475860e-03
## [12,] -2.381963e-03 3.168730e-08 4.393262e-06 4.484872e-07 -1.506467e-02
## [13,] -2.702654e-04 2.921042e-09 3.032800e-07 3.739299e-08 8.070865e-04
## [14,] -1.921152e-09 1.514196e-04 7.371366e-09 2.375318e-05 -6.372019e-11
## [15,] -3.978128e-08 -1.531230e-02 1.284886e-06 8.088122e-04 1.131008e-09
## [16,] 9.546976e-08 -1.015847e-02 3.442315e-07 -1.980713e-03 3.973932e-09
## [,6] [,7] [,8] [,9] [,10]
## [1,] -3.876187e-02 -5.053342e-06 -1.613881e-06 6.021735e-08 -1.875009e-10
## [2,] 2.955259e-07 4.778706e-11 7.931302e-02 -3.016254e-02 4.448383e-05
## [3,] 1.117905e-05 3.747143e-09 5.458145e-06 2.430842e-06 -3.336611e-09
## [4,] 2.377030e-06 6.158777e-10 3.482589e-02 1.928380e-03 4.859196e-06
## [5,] -9.652046e-01 1.777274e-05 -7.785824e-08 5.325279e-09 -4.228263e-11
## [6,] 8.987805e+00 -1.917434e-03 1.837811e-08 -1.571215e-08 2.572062e-11
## [7,] -1.917434e-03 5.566210e-04 3.842400e-12 -4.496739e-12 6.413081e-15
## [8,] 1.837811e-08 3.842400e-12 8.146499e+00 -2.303285e-01 1.032627e-03
## [9,] -1.571215e-08 -4.496739e-12 -2.303285e-01 7.047014e+00 -9.698669e-03
## [10,] 2.572062e-11 6.413081e-15 1.032627e-03 -9.698669e-03 3.713355e-04
## [11,] -5.728497e-03 -1.839729e-06 -2.082894e-08 5.102820e-10 -2.887363e-12
## [12,] -2.662386e-02 5.883356e-06 -5.009482e-09 -6.207114e-09 6.454367e-12
## [13,] -7.275176e-04 1.626211e-06 -1.570550e-10 -4.063459e-10 4.586399e-13
## [14,] 1.219214e-11 2.882146e-15 2.844905e-05 1.761971e-05 -7.462095e-08
## [15,] -7.047874e-09 -2.431001e-12 -1.557562e-02 -1.783414e-02 -7.014699e-07
## [16,] -5.532745e-09 -1.701542e-12 -1.321624e-02 -1.098511e-02 -4.542926e-05
## [,11] [,12] [,13] [,14] [,15]
## [1,] 2.862649e-02 -2.381963e-03 -2.702654e-04 -1.921152e-09 -3.978128e-08
## [2,] -3.153236e-07 3.168730e-08 2.921042e-09 1.514196e-04 -1.531230e-02
## [3,] 3.901703e-07 4.393262e-06 3.032800e-07 7.371366e-09 1.284886e-06
## [4,] -1.130433e-06 4.484872e-07 3.739299e-08 2.375318e-05 8.088122e-04
## [5,] -6.475860e-03 -1.506467e-02 8.070865e-04 -6.372019e-11 1.131008e-09
## [6,] -5.728497e-03 -2.662386e-02 -7.275176e-04 1.219214e-11 -7.047874e-09
## [7,] -1.839729e-06 5.883356e-06 1.626211e-06 2.882146e-15 -2.431001e-12
## [8,] -2.082894e-08 -5.009482e-09 -1.570550e-10 2.844905e-05 -1.557562e-02
## [9,] 5.102820e-10 -6.207114e-09 -4.063459e-10 1.761971e-05 -1.783414e-02
## [10,] -2.887363e-12 6.454367e-12 4.586399e-13 -7.462095e-08 -7.014699e-07
## [11,] 8.150119e-01 -1.597299e-01 1.702626e-02 -2.387114e-11 -1.461755e-10
## [12,] -1.597299e-01 2.413899e+00 -3.713726e-02 -2.669744e-12 -2.750829e-09
## [13,] 1.702626e-02 -3.713726e-02 4.336252e-02 -3.720977e-14 -1.900898e-10
## [14,] -2.387114e-11 -2.669744e-12 -3.720977e-14 7.583464e-04 -1.591605e-03
## [15,] -1.461755e-10 -2.750829e-09 -1.900898e-10 -1.591605e-03 2.553117e+00
## [16,] 1.290874e-09 -1.425838e-09 -1.081010e-10 1.363966e-04 -3.960194e-01
## [,16]
## [1,] 9.546976e-08
## [2,] -1.015847e-02
## [3,] 3.442315e-07
## [4,] -1.980713e-03
## [5,] 3.973932e-09
## [6,] -5.532745e-09
## [7,] -1.701542e-12
## [8,] -1.321624e-02
## [9,] -1.098511e-02
## [10,] -4.542926e-05
## [11,] 1.290874e-09
## [12,] -1.425838e-09
## [13,] -1.081010e-10
## [14,] 1.363966e-04
## [15,] -3.960194e-01
## [16,] 2.467618e+00
##
## Model rank = 1346 / 1346
##
## 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):vowelContrastDRESS 11.00 7.70 0.96 0.28
## s(Time):vowelContrastTRAP 11.00 5.55 0.96 0.22
## s(Time):isAustralian 12.00 2.00 0.96 0.22
## s(Time):isAustralianTrap 10.00 5.19 0.96 0.26
## s(Time,pp):vowelContrastDRESS 450.00 88.75 0.96 0.22
## s(Time,pp):vowelContrastTRAP 450.00 70.01 0.96 0.22
## s(Time,item):vowelContrastDRESS 200.00 13.40 0.96 0.23
## s(Time,item):vowelContrastTRAP 200.00 13.46 0.96 0.32We now examine this model to determine whether there is an overall effect of contrast.
summary(interaction.compet_dressTrap.gamm1)##
## Family: binomial
## Link function: logit
##
## Formula:
## position ~ vowelContrast + s(Time, by = vowelContrast, k = 12) +
## s(Time, by = isAustralian, k = 12) + s(Time, by = isAustralianTrap) +
## 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.54764 0.05172 -29.926 <2e-16 ***
## vowelContrastTRAP -0.12777 0.09708 -1.316 0.188
## ---
## 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):vowelContrastDRESS 7.698 9.347 119.966 < 2e-16 ***
## s(Time):vowelContrastTRAP 5.551 7.134 70.910 < 2e-16 ***
## s(Time):isAustralian 2.004 2.007 2.690 0.261
## s(Time):isAustralianTrap 5.185 6.368 8.971 0.209
## s(Time,pp):vowelContrastDRESS 88.747 448.000 188.714 < 2e-16 ***
## s(Time,pp):vowelContrastTRAP 70.012 449.000 147.311 < 2e-16 ***
## s(Time,item):vowelContrastDRESS 13.396 98.000 34.494 4.93e-05 ***
## s(Time,item):vowelContrastTRAP 13.457 99.000 45.023 3.06e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.0602 Deviance explained = 6.92%
## fREML = 82192 Scale est. = 1 n = 319955This model shows significant nonlinear time effects and substantial participant and item variability, but no significant main effect of vowel contrast (p = .183) and no significant accent-related smooths, indicating that competitor fixation trajectories change over time but do not reliably differ overall between DRESS and TRAP or across accent.
smooth <- plot_smooth(interaction.compet_dressTrap.gamm1, view = "Time",
plot_all = "vowelContrast",
rug = FALSE, rm.ranef = TRUE, lwd = 4)## Summary:
## * vowelContrast : factor; set to the value(s): DRESS, TRAP.
## * Time : numeric predictor; with 30 values ranging from 10.000000 to 890.000000.
## * isAustralian : numeric predictor; set to the value(s): 0.
## * isAustralianTrap : 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):vowelContrastDRESS,s(Time,pp):vowelContrastTRAP,s(Time,item):vowelContrastDRESS,s(Time,item):vowelContrastTRAP
## 
We now proceed to model comparison by successively removing effects to evaluate their contribution to model fit.
reducedModelFile = "compet_dressTrap_reducedModel.Rdata"
if (file.exists(reducedModelFile))
{
load(reducedModelFile)
}else{
mainEffect.compet_dressTrap.gamm <- bam( position ~ vowelContrast+
s(Time, by = vowelContrast) +
#s(Time, by = isAustralian) +
#s(Time, by = isAustralianTrap) +
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("mainEffect.compet_dressTrap.gamm"), file = reducedModelFile)
}summary(mainEffect.compet_dressTrap.gamm)##
## Family: binomial
## Link function: logit
##
## Formula:
## position ~ vowelContrast + s(Time, by = vowelContrast) + 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.59968 0.03859 -41.451 <2e-16 ***
## vowelContrastTRAP -0.10003 0.08329 -1.201 0.23
## ---
## 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):vowelContrastDRESS 6.992 8.148 125.65 < 2e-16 ***
## s(Time):vowelContrastTRAP 6.102 7.457 116.87 < 2e-16 ***
## s(Time,pp):vowelContrastDRESS 88.652 448.000 188.37 < 2e-16 ***
## s(Time,pp):vowelContrastTRAP 70.245 448.000 148.06 < 2e-16 ***
## s(Time,item):vowelContrastDRESS 13.369 98.000 34.38 5.11e-05 ***
## s(Time,item):vowelContrastTRAP 13.496 98.000 45.13 3.01e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.0597 Deviance explained = 6.87%
## fREML = 82265 Scale est. = 1 n = 319955This reduced model shows strong nonlinear time effects and substantial participant and item variability, but no significant main effect of vowel contrast (p = .226), indicating that competitor fixation trajectories change over time but do not reliably differ overall between DRESS and TRAP.
smooth <- plot_smooth(mainEffect.compet_dressTrap.gamm, view = "Time",
plot_all = "vowelContrast",
rug = FALSE, rm.ranef = TRUE, lwd = 4)## Summary:
## * vowelContrast : factor; set to the value(s): DRESS, TRAP.
## * Time : numeric predictor; with 30 values ranging from 10.000000 to 890.000000.
## * 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):vowelContrastDRESS,s(Time,pp):vowelContrastTRAP,s(Time,item):vowelContrastDRESS,s(Time,item):vowelContrastTRAP
## 
3.2.3.3. TARGET-COMPETITOR FIXATIONS
A. TIME-WINDOW ANALYSIS
rm(list = ls())
library(lmerTest)
library(stringr)
library(tidyverse)
library(mgcv)
library(itsadug)
options(warn = 0)
gazeTarget = read.table("dataCatalan/gazeTarget.txt",
header = TRUE, sep = "\t", quote = "")
gazeComp = read.table("dataCatalan/gazeCompet.txt",
header = TRUE, sep = "\t", quote = "")
trialData = read.table("dataCatalan/trialInfo.txt",
header = TRUE, sep = "\t", quote = "\"")
expFilter = grepl("DRESS|TRAP", trialData$condition)
gazeTargetE = gazeTarget %>% filter(expFilter)
gazeCompE = gazeComp %>% filter(expFilter)
trialDataE = trialData %>% filter(expFilter)
timeline = -200 + (0:120) * 10
colnames(gazeTargetE) = paste0("Tt_", timeline)
colnames(gazeCompE) = paste0("Tc_", 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
## 1797 1799trialDataE$contrastC = -0.5 + as.numeric(trialDataE$contrast == "TRAP")
table(trialDataE$contrastC)##
## -0.5 0.5
## 1798 1798useData = cbind(trialDataE, gazeTargetE, gazeCompE)
# Remove distractor responses
useData <- useData %>%
filter(useData$clickedOn != "3") %>%
droplevels()windowData = useData %>%
pivot_longer(
cols = matches("^T(t|c)_-?\\d+$"),
names_to = c("which", "Time"),
names_pattern = "^T(t|c)_(-?\\d+)$",
values_to = "position"
) %>%
pivot_wider(
names_from = which,
values_from = position
) %>%
rename(
target = t,
competitor = c
) %>%
mutate(
target = as.numeric(target),
competitor = as.numeric(competitor),
difference = target - competitor,
Time = as.numeric(Time)
) %>%
as.data.frame()
windowData = windowData[order(windowData$pp, windowData$trial, windowData$Time), ]
windowData$pp = as.factor(windowData$pp)
windowData$item = as.factor(windowData$item)
windowData$isAustralian = ifelse(windowData$accent == "AusE", 1, 0)
windowData$isAustralianFleece = windowData$isAustralian *
as.numeric(windowData$contrast == "TRAP")
windowData$vowelContrast = windowData$contrast %>% as.factor()
# restricting the data to 0-900 after vowel onset
windowData = windowData %>% filter(Time > 0 & Time < 900)
windowData <- start_event(windowData, event = c("pp", "trial"))windowDataPlot = windowData %>%
group_by(vowelContrast, Time, accent) %>%
summarize(meanDiff = mean(difference, na.rm = TRUE), .groups = "drop")
Fig.diff_dressTrap = ggplot(
windowDataPlot,
aes(x = Time, y = meanDiff, col = vowelContrast)
) +
geom_line() +
facet_wrap(~accent)
print(Fig.diff_dressTrap)
aggregatedData = windowData %>%
filter(Time > 150 & Time < 600) %>%
group_by(pp, item, accent, contrast) %>%
summarize(meanDiff = mean(difference, na.rm = TRUE), .groups = "drop") %>%
ungroup()targetCompetitor_fleeceKit.lmer = lmer(
meanDiff ~ contrast + accent +
(1 + contrast | pp) +
(1 | item),
data = aggregatedData
)## boundary (singular) fit: see help('isSingular')summary(targetCompetitor_fleeceKit.lmer)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: meanDiff ~ contrast + accent + (1 + contrast | pp) + (1 | item)
## Data: aggregatedData
##
## REML criterion at convergence: 1668
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.12966 -0.63758 0.01998 0.73536 2.33874
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## pp (Intercept) 0.0090108 0.09493
## contrastTRAP 0.0007775 0.02788 1.00
## item (Intercept) 0.0031585 0.05620
## Residual 0.1396928 0.37376
## Number of obs: 1800, groups: pp, 45; item, 20
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 0.24330 0.02737 37.50773 8.890 9.06e-11 ***
## contrastTRAP 0.05641 0.03097 18.54941 1.821 0.0847 .
## accentSEng -0.01778 0.01762 1735.00036 -1.009 0.3131
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) cnTRAP
## contrstTRAP -0.486
## accentSEng -0.322 0.000
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')