Code
cat('\014') # clean terminalCode
rm(list = ls())
library(tidyverse)
library(GGally)
library(afex)
library(emmeans)
library(performance)Robot 2023
Recordings with Emotiv EPOC Flex (gel)
Code
theme_set(
theme_minimal()
)
a_posteriori <- function(afex_aov, sig_level = .05) {
factors <- as.list(rownames(afex_aov$anova_table))
for (j in 1:length(factors)) {
if (grepl(':', factors[[j]])) {
factors[[j]] <- unlist(strsplit(factors[[j]], ':'))
}
}
p_values <- afex_aov$anova_table$`Pr(>F)`
for (i in 1:length(p_values)) {
if (p_values[i] <= sig_level) {
cat(rep('_', 60), '\n', sep = '')
print(emmeans(afex_aov, factors[[i]], contr = 'pairwise'))
}
}
}Code
robot_files <- c('erplab/average_voltage_100_to_150_robot_2023.txt',
'erplab/average_voltage_290_to_350_robot_2023.txt',
'erplab/average_voltage_600_to_800_robot_2023.txt')
robot_data <- read_tsv(robot_files, show_col_types = FALSE)
robot_data <- robot_data |>
rename(component = mlabel,
amplitude = value,
stimulus = binlabel,
electrode = chlabel) |>
separate(stimulus, c('valence', 'face'), sep = '_', remove = FALSE) |>
mutate(valence = factor(valence, levels = c('neutral', 'happy'))) |>
mutate(face = factor(face, levels = c('robot', 'human'))) |>
mutate(component = factor(component, levels = c('P1', 'P3', 'Late'))) |>
separate(worklat, c('start', 'end'), sep = ' ', remove = FALSE) |>
mutate(start = parse_number(start)) |>
mutate(end = parse_number(end)) |>
mutate(num_id = parse_number(ERPset)) |>
mutate_if(is.character, as.factor) |>
mutate(electrode = factor(electrode, levels = c('P3', 'Pz', 'P4', 'O1', 'Oz', 'O2'))) |>
mutate(antero_posterior = ifelse(electrode %in% c('P3', 'Pz', 'P4'), 'parietal', 'occipital')) |>
mutate(antero_posterior = factor(antero_posterior, levels = c('parietal', 'occipital'))) |>
mutate(medio_lateral = ifelse(electrode %in% c('P3', 'O1'), 'left', ifelse(electrode %in% c('Pz', 'Oz'), 'central', 'right'))) |>
mutate(medio_lateral = factor(medio_lateral, levels = c('left', 'central', 'right'))) |>
mutate_at(c('num_id', 'chindex', 'bini'), as.factor)
write_csv(robot_data, 'mrs_robot_2023_data_clean.csv')
n170_file <- 'erplab/average_voltage_150_to_200_robot_2023.txt'
n170_data <- read_tsv(n170_file, show_col_types = FALSE) |>
rename(component = mlabel,
amplitude = value,
stimulus = binlabel,
electrode = chlabel) |>
separate(stimulus, c('valence', 'face'), sep = '_', remove = FALSE) |>
mutate(valence = factor(valence, levels = c('neutral', 'happy'))) |>
mutate(face = factor(face, levels = c('robot', 'human'))) |>
separate(worklat, c('start', 'end'), sep = ' ', remove = FALSE) |>
mutate(start = parse_number(start)) |>
mutate(end = parse_number(end)) |>
mutate(num_id = parse_number(ERPset)) |>
mutate_if(is.character, as.factor) |>
mutate(electrode = factor(electrode, levels = c('P7', 'P8', 'PO9', 'PO10'))) |>
mutate(antero_posterior = ifelse(electrode %in% c('P7', 'P8'), 'parietal', 'parieto_occipital')) |>
mutate(antero_posterior = factor(antero_posterior, levels = c('parietal', 'parieto_occipital'))) |>
mutate(medio_lateral = ifelse(electrode %in% c('P7', 'PO9'), 'left', 'right')) |>
mutate(medio_lateral = factor(medio_lateral, levels = c('left', 'right'))) |>
mutate_at(c('num_id', 'chindex', 'bini'), as.factor)
write_csv(n170_data, 'mrs_robot_2023_n170_data_clean.csv')
connectivity_file <- 'erplab/median_connectivity.csv'
connectivity_data <- read_csv(connectivity_file, show_col_types = FALSE) |>
mutate(log_median_connectivity = log(median_connectivity)) |>
mutate(face = factor(face, levels = c('robot', 'human'))) |>
mutate(num_id = parse_number(ERPset)) |>
mutate_if(is.character, as.factor) |>
mutate(region_b = factor(region_b, levels = c('parietal_ipsi', 'parietal_contra')))
write_csv(connectivity_data, 'mrs_robot_2023_connectivity_data_clean.csv')ERP plots
EDF files were pre-processed in Matlab using the EEGLAB 2024.0 toolbox (Delorme and Makeig 2004), the ERPLAB 10.0 toolbox (Lopez-Calderon and Luck 2014), and automated with in-house scripts. EEG data was high-pass filtered at 0.5 Hz with a 12 dB/oct roll-off. A working time window was selected from 2 seconds before the first stimulus to 2 seconds after the last one. Defective channels were identified by eye inspection and omitted from the preprocessing steps. Artifacts originating from eye blinks or movements, heart beats, muscle contraction, channel noise, or electrical interference were identified and removed by means of independent component analysis (ICA) and using the ICLabel 1.4 classifier (Pion-Tonachini, Kreutz-Delgado, and Makeig 2019), components with an score between 0.8 and 1 in the aforementioned artifactual categories were removed. At this point defective channels were spherically interpolated. Data was re-referenced to infinity with the REST 1.2 toolbox (Dong et al. 2017), low-pass filtered at 30Hz with a 12 dB/oct roll-off, cut in [-200ms 1000ms] epochs around stimuli, and baseline-corrected relative to the mean voltage in the pre-stimulus time.
Face & Valence ERPs:
Topographic layout:
General description
Code
summary(robot_data) worklat start end component amplitude
[100.0 150.0]:552 Min. :100 Min. :150.0 P1 :552 Min. :-5.184
[290.0 350.0]:552 1st Qu.:100 1st Qu.:150.0 P3 :552 1st Qu.: 1.218
[600.0 800.0]:552 Median :290 Median :350.0 Late:552 Median : 3.093
Mean :330 Mean :433.3 Mean : 3.540
3rd Qu.:600 3rd Qu.:800.0 3rd Qu.: 5.561
Max. :600 Max. :800.0 Max. :20.933
chindex electrode bini stimulus valence face
13:368 P3:368 1:414 happy_human :414 neutral:828 robot:828
16:184 Pz:368 2:414 happy_robot :414 happy :828 human:828
17:368 P4:368 3:414 neutral_human:414
18:184 O1:184 4:414 neutral_robot:414
19:184 Oz:184
22:368 O2:184
ERPset num_id antero_posterior medio_lateral
MRS_04_ct_FT9: 36 4 : 36 parietal :1104 left :552
MRS_07_ct_FT9: 36 7 : 36 occipital: 552 central:552
MRS_08_ct_FT9: 36 8 : 36 right :552
MRS_09_ct_FT9: 36 9 : 36
MRS_10_ct_FT9: 36 10 : 36
MRS_11_ct_FT9: 36 11 : 36
(Other) :1440 (Other):1440 Code
ggplot(
robot_data, aes(x = amplitude, fill = component, color = component)) +
geom_histogram(alpha = .4) +
facet_wrap(~component)
robot_data_wide <- robot_data[c('ERPset', 'face', 'valence', 'medio_lateral', 'component', 'amplitude')] |>
pivot_wider(names_from = component, values_from = amplitude)
components <- c('P1', 'P3', 'Late')
components_pairs <- ggpairs(robot_data_wide,
aes(colour = face, alpha = .5, linewidth = NA),
columns = components,
lower = list(continuous = wrap('points', alpha = .4)),
progress = FALSE
)
suppressWarnings(print(components_pairs))
P1 by Face & Valence, average 100 to 150 ms
P1 Scalp voltage distribution:
P1 ANOVA:
Code
robot_p1_rep_anova <- aov_ez('num_id', 'amplitude',
subset(robot_data, component == 'P1',),
within = c('face', 'valence', 'electrode'))
robot_p1_afex_plot <- afex_plot(robot_p1_rep_anova,
x = 'electrode',
trace = 'face',
panel = 'valence',
error = 'within',
error_arg = list(width = .3, lwd = .75),
dodge = .5,
data_arg = list(
position =
position_jitterdodge(
jitter.width = .2,
dodge.width = .5 ## needs to be same as dodge
)),
mapping = c('color'),, data_alpha = .3,
point_arg = list(size = 3)
)
suppressWarnings(print(robot_p1_afex_plot))
Code
mytable <- xtabs(~ face + valence + electrode, data = robot_p1_rep_anova$data$long)
ftable(mytable) electrode O1 Oz O2
face valence
robot neutral 46 46 46
happy 46 46 46
human neutral 46 46 46
happy 46 46 46Code
print(robot_p1_rep_anova)Anova Table (Type 3 tests)
Response: amplitude
Effect df MSE F ges p.value
1 face 1, 45 11.89 13.50 *** .021 <.001
2 valence 1, 45 2.96 0.78 <.001 .383
3 electrode 1.96, 88.09 5.03 0.26 <.001 .764
4 face:valence 1, 45 6.94 0.41 <.001 .525
5 face:electrode 1.53, 68.94 1.42 2.32 <.001 .118
6 valence:electrode 1.40, 63.16 0.67 3.41 + <.001 .055
7 face:valence:electrode 1.67, 75.35 0.68 1.97 <.001 .154
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
Sphericity correction method: GG Code
a_posteriori(robot_p1_rep_anova)____________________________________________________________
$emmeans
face emmean SE df lower.CL upper.CL
robot 5.59 0.572 45 4.44 6.75
human 4.52 0.425 45 3.66 5.37
Results are averaged over the levels of: electrode, valence
Confidence level used: 0.95
$contrasts
contrast estimate SE df t.ratio p.value
robot - human 1.08 0.293 45 3.674 0.0006
Results are averaged over the levels of: electrode, valence P1 Assumptions:
Code
norm_p1 <- check_normality(robot_p1_rep_anova)
hetero_p1 <- check_heteroscedasticity(robot_p1_rep_anova)
print(norm_p1)
print('', quote = FALSE)
print(hetero_p1)
print('', quote = FALSE)
check_sphericity(robot_p1_rep_anova)
plot(norm_p1, type = 'density')
plot(norm_p1, type = 'qq')
plot(hetero_p1)Warning: Non-normality of residuals detected (p < .001).[1] Warning: Heteroscedasticity (non-constant error variance) detected (p < .001).[1] Warning: Sphericity violated for:
- face:electrode (p < .001)
- valence:electrode (p < .001)
- face:valence:electrode (p = 0.009).
P3 by Face & Valence, average 290 to 350 ms
P3 Scalp voltage distribution:
P3 ANOVA:
Code
robot_p3_rep_anova = aov_ez('num_id', 'amplitude',
subset(robot_data, component == 'P3'),
within = c('face', 'valence', 'electrode'))
robot_p3_afex_plot <- afex_plot(robot_p3_rep_anova,
x = 'electrode',
trace = 'face',
panel = 'valence',
error = 'within',
error_arg = list(width = .3, lwd = .75),
dodge = .5,
data_arg = list(
position =
position_jitterdodge(
jitter.width = .2,
dodge.width = .5 ## needs to be same as dodge
)),
mapping = c('color'),, data_alpha = .3,
point_arg = list(size = 3.5)
)
suppressWarnings(print(robot_p3_afex_plot))
Code
mytable <- xtabs(~ face + valence + electrode, data = robot_p3_rep_anova$data$long)
ftable(mytable) electrode P3 Pz P4
face valence
robot neutral 46 46 46
happy 46 46 46
human neutral 46 46 46
happy 46 46 46Code
print(robot_p3_rep_anova)Anova Table (Type 3 tests)
Response: amplitude
Effect df MSE F ges p.value
1 face 1, 45 9.56 10.59 ** .028 .002
2 valence 1, 45 3.90 7.80 ** .008 .008
3 electrode 1.86, 83.63 5.36 11.13 *** .030 <.001
4 face:valence 1, 45 4.91 0.00 <.001 .956
5 face:electrode 1.96, 88.26 1.55 0.30 <.001 .740
6 valence:electrode 2.00, 89.83 0.42 1.00 <.001 .370
7 face:valence:electrode 1.98, 88.94 0.55 0.51 <.001 .599
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
Sphericity correction method: GG Code
a_posteriori(robot_p3_rep_anova)____________________________________________________________
$emmeans
face emmean SE df lower.CL upper.CL
robot 4.55 0.319 45 3.90 5.19
human 3.69 0.316 45 3.05 4.33
Results are averaged over the levels of: electrode, valence
Confidence level used: 0.95
$contrasts
contrast estimate SE df t.ratio p.value
robot - human 0.857 0.263 45 3.255 0.0022
Results are averaged over the levels of: electrode, valence
____________________________________________________________
$emmeans
valence emmean SE df lower.CL upper.CL
neutral 4.35 0.308 45 3.73 4.98
happy 3.88 0.294 45 3.29 4.48
Results are averaged over the levels of: electrode, face
Confidence level used: 0.95
$contrasts
contrast estimate SE df t.ratio p.value
neutral - happy 0.47 0.168 45 2.792 0.0077
Results are averaged over the levels of: electrode, face
____________________________________________________________
$emmeans
electrode emmean SE df lower.CL upper.CL
P3 3.69 0.294 45 3.10 4.28
Pz 3.93 0.363 45 3.20 4.66
P4 4.74 0.296 45 4.14 5.33
Results are averaged over the levels of: valence, face
Confidence level used: 0.95
$contrasts
contrast estimate SE df t.ratio p.value
P3 - Pz -0.242 0.203 45 -1.192 0.4642
P3 - P4 -1.048 0.232 45 -4.509 0.0001
Pz - P4 -0.806 0.259 45 -3.113 0.0089
Results are averaged over the levels of: valence, face
P value adjustment: tukey method for comparing a family of 3 estimates P3 Assumptions:
Code
norm_p3 <- check_normality(robot_p3_rep_anova)
hetero_p3 <- check_heteroscedasticity(robot_p3_rep_anova)
print(norm_p3)
print('', quote = FALSE)
print(hetero_p3)
print('', quote = FALSE)
check_sphericity(robot_p3_rep_anova)
plot(norm_p3, type = 'density')
plot(norm_p3, type = 'qq')
plot(hetero_p3)OK: residuals appear as normally distributed (p = 0.338).[1] OK: Error variance appears to be homoscedastic (p = 0.862).[1] OK: Data seems to be spherical (p > 0.175).
Late by Face & Valence, average 600 to 800 ms
Late Scalp voltage distribution:
Late ANOVA:
Code
robot_late_rep_anova = aov_ez('num_id', 'amplitude',
subset(robot_data, component == 'Late'),
within = c('face', 'valence', 'electrode'))
robot_late_afex_plot <- afex_plot(robot_late_rep_anova,
x = 'electrode',
trace = 'face',
panel = 'valence',
error = 'within',
error_arg = list(width = .3, lwd = .75),
dodge = .5,
data_arg = list(
position =
position_jitterdodge(
jitter.width = .2,
dodge.width = .5 ## needs to be same as dodge
)),
mapping = c('color'),, data_alpha = .3,
point_arg = list(size = 3.5)
)
suppressWarnings(print(robot_late_afex_plot))
Code
mytable <- xtabs(~ face + valence + electrode, data = robot_late_rep_anova$data$long)
ftable(mytable) electrode P3 Pz P4
face valence
robot neutral 46 46 46
happy 46 46 46
human neutral 46 46 46
happy 46 46 46Code
print(robot_late_rep_anova)Anova Table (Type 3 tests)
Response: amplitude
Effect df MSE F ges p.value
1 face 1, 45 4.28 54.35 *** .128 <.001
2 valence 1, 45 2.69 3.14 + .005 .083
3 electrode 1.99, 89.36 1.08 2.44 + .003 .094
4 face:valence 1, 45 3.76 0.02 <.001 .897
5 face:electrode 1.99, 89.60 0.57 3.18 * .002 .046
6 valence:electrode 1.86, 83.87 0.41 3.91 * .002 .026
7 face:valence:electrode 1.78, 80.30 0.68 1.11 <.001 .330
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
Sphericity correction method: GG Code
cat(rep('_', 60), '\n', sep = '')____________________________________________________________Code
emmeans(robot_late_rep_anova, pairwise ~ face)$emmeans
face emmean SE df lower.CL upper.CL
robot 0.796 0.205 45 0.383 1.21
human 2.095 0.208 45 1.677 2.51
Results are averaged over the levels of: electrode, valence
Confidence level used: 0.95
$contrasts
contrast estimate SE df t.ratio p.value
robot - human -1.3 0.176 45 -7.373 <.0001
Results are averaged over the levels of: electrode, valence Code
cat(rep('_', 60), '\n', sep = '')____________________________________________________________Code
emmeans(robot_late_rep_anova, pairwise ~ electrode | face)$emmeans
face = robot:
electrode emmean SE df lower.CL upper.CL
P3 0.912 0.201 45 0.507 1.32
Pz 0.797 0.232 45 0.330 1.26
P4 0.679 0.220 45 0.237 1.12
face = human:
electrode emmean SE df lower.CL upper.CL
P3 2.027 0.203 45 1.618 2.44
Pz 2.304 0.253 45 1.794 2.81
P4 1.953 0.207 45 1.536 2.37
Results are averaged over the levels of: valence
Confidence level used: 0.95
$contrasts
face = robot:
contrast estimate SE df t.ratio p.value
P3 - Pz 0.1153 0.119 45 0.967 0.6012
P3 - P4 0.2325 0.133 45 1.745 0.1999
Pz - P4 0.1172 0.132 45 0.889 0.6500
face = human:
contrast estimate SE df t.ratio p.value
P3 - Pz -0.2778 0.140 45 -1.984 0.1278
P3 - P4 0.0735 0.136 45 0.541 0.8517
Pz - P4 0.3513 0.138 45 2.541 0.0380
Results are averaged over the levels of: valence
P value adjustment: tukey method for comparing a family of 3 estimates Code
cat(rep('_', 60), '\n', sep = '')____________________________________________________________Code
emmeans(robot_late_rep_anova, pairwise ~ electrode | valence)$emmeans
valence = neutral:
electrode emmean SE df lower.CL upper.CL
P3 1.58 0.204 45 1.173 2.00
Pz 1.59 0.247 45 1.093 2.09
P4 1.53 0.202 45 1.127 1.94
valence = happy:
electrode emmean SE df lower.CL upper.CL
P3 1.35 0.208 45 0.936 1.77
Pz 1.51 0.215 45 1.078 1.95
P4 1.10 0.192 45 0.712 1.49
Results are averaged over the levels of: face
Confidence level used: 0.95
$contrasts
valence = neutral:
contrast estimate SE df t.ratio p.value
P3 - Pz -0.00475 0.118 45 -0.040 0.9991
P3 - P4 0.05088 0.120 45 0.424 0.9059
Pz - P4 0.05564 0.127 45 0.438 0.8998
valence = happy:
contrast estimate SE df t.ratio p.value
P3 - Pz -0.15783 0.133 45 -1.188 0.4667
P3 - P4 0.25511 0.130 45 1.957 0.1349
Pz - P4 0.41294 0.123 45 3.357 0.0045
Results are averaged over the levels of: face
P value adjustment: tukey method for comparing a family of 3 estimates Late Assumptions:
Code
norm_late <- check_normality(robot_late_rep_anova)
hetero_late <- check_heteroscedasticity(robot_late_rep_anova)
print(norm_late)
print('', quote = FALSE)
print(hetero_late)
print('', quote = FALSE)
check_sphericity(robot_late_rep_anova)
plot(norm_late, type = 'density')
plot(norm_late, type = 'qq')
plot(hetero_late)Warning: Non-normality of residuals detected (p < .001).[1] OK: Error variance appears to be homoscedastic (p = 0.178).[1] OK: Data seems to be spherical (p > 0.059).
Factorial Mass Univariate analysis (FMUT, cluster-mass)
Groppe, Urbach, and Kutas (2011), Fields and Kuperberg (2020)
Interaction
Face
Valence
Multivariate pattern analysis (decoding)
Done in the ADAM 1.14 toolbox (Fahrenfort 2020) with FieldTrip 20211209 (Oostenveld et al. 2011).
Exemplar ERPs on O2
Decoding over time
Individual results
Activation pattern
Temporal generalization plot
Temporal generalization over time
Coda: N170 by Face & Valence, average 150 to 200 ms
Code
summary(n170_data) worklat start end component
[150.0 200.0]:736 Min. :150 Min. :200 N170:736
1st Qu.:150 1st Qu.:200
Median :150 Median :200
Mean :150 Mean :200
3rd Qu.:150 3rd Qu.:200
Max. :150 Max. :200
amplitude chindex electrode bini stimulus
Min. :-14.282380 14:184 P7 :184 1:184 happy_human :184
1st Qu.: -2.164432 15:184 P8 :184 2:184 happy_robot :184
Median : -0.003614 20:184 PO9 :184 3:184 neutral_human:184
Mean : -0.224474 21:184 PO10:184 4:184 neutral_robot:184
3rd Qu.: 1.938663
Max. : 10.430920
valence face ERPset num_id
neutral:368 robot:368 MRS_04_ct_FT9: 16 4 : 16
happy :368 human:368 MRS_07_ct_FT9: 16 7 : 16
MRS_08_ct_FT9: 16 8 : 16
MRS_09_ct_FT9: 16 9 : 16
MRS_10_ct_FT9: 16 10 : 16
MRS_11_ct_FT9: 16 11 : 16
(Other) :640 (Other):640
antero_posterior medio_lateral
parietal :368 left :368
parieto_occipital:368 right:368
Code
ggplot(
n170_data, aes(x = amplitude, fill = face, color = face)) +
geom_histogram(alpha = .4) +
facet_wrap(~face)
Face & Valence ERPs:
N170 Scalp voltage distribution:
N170 ANOVA:
Code
robot_n170_rep_anova <- aov_ez('num_id', 'amplitude',
n170_data,
within = c('face', 'valence', 'medio_lateral', 'antero_posterior'))
robot_n170_afex_plot <- afex_plot(robot_n170_rep_anova,
x = 'medio_lateral',
trace = 'face',
panel = 'valence',
error = 'within',
error_arg = list(width = .3, lwd = .75),
dodge = .5,
data_arg = list(
position =
position_jitterdodge(
jitter.width = .2,
dodge.width = .5 ## needs to be same as dodge
)),
mapping = c('color'),, data_alpha = .3,
point_arg = list(size = 3)
)
suppressWarnings(print(robot_n170_afex_plot))
Code
mytable <- xtabs(~ face + valence + medio_lateral + antero_posterior, data = robot_n170_rep_anova$data$long)
ftable(mytable) antero_posterior parietal parieto_occipital
face valence medio_lateral
robot neutral left 46 46
right 46 46
happy left 46 46
right 46 46
human neutral left 46 46
right 46 46
happy left 46 46
right 46 46Code
print(robot_n170_rep_anova)Anova Table (Type 3 tests)
Response: amplitude
Effect df MSE F ges p.value
1 face 1, 45 5.18 1.97 .001 .168
2 valence 1, 45 4.43 0.88 <.001 .354
3 medio_lateral 1, 45 17.31 0.64 .002 .427
4 antero_posterior 1, 45 12.30 1.28 .002 .263
5 face:valence 1, 45 5.29 0.01 <.001 .905
6 face:medio_lateral 1, 45 2.64 6.34 * .002 .015
7 valence:medio_lateral 1, 45 2.04 0.35 <.001 .557
8 face:antero_posterior 1, 45 3.01 0.00 <.001 .959
9 valence:antero_posterior 1, 45 4.19 0.90 <.001 .348
10 medio_lateral:antero_posterior 1, 45 3.78 1.62 <.001 .209
11 face:valence:medio_lateral 1, 45 2.74 0.13 <.001 .716
12 face:valence:antero_posterior 1, 45 1.34 0.47 <.001 .497
13 face:medio_lateral:antero_posterior 1, 45 0.80 0.01 <.001 .909
14 valence:medio_lateral:antero_posterior 1, 45 0.96 1.75 <.001 .193
15 face:valence:medio_lateral:antero_posterior 1, 45 1.37 0.06 <.001 .804
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1Code
a_posteriori(robot_n170_rep_anova)____________________________________________________________
$emmeans
face medio_lateral emmean SE df lower.CL upper.CL
robot left 0.167 0.368 45 -0.573 0.908
human left -0.370 0.366 45 -1.108 0.368
robot right -0.381 0.436 45 -1.258 0.496
human right -0.314 0.434 45 -1.188 0.559
Results are averaged over the levels of: antero_posterior, valence
Confidence level used: 0.95
$contrasts
contrast estimate SE df t.ratio p.value
robot left - human left 0.5372 0.201 45 2.673 0.0493
robot left - robot right 0.5478 0.313 45 1.752 0.3096
robot left - human right 0.4812 0.354 45 1.359 0.5312
human left - robot right 0.0106 0.345 45 0.031 1.0000
human left - human right -0.0560 0.345 45 -0.162 0.9985
robot right - human right -0.0666 0.211 45 -0.315 0.9890
Results are averaged over the levels of: antero_posterior, valence
P value adjustment: tukey method for comparing a family of 4 estimates N170 Assumptions:
Code
norm_n170 <- check_normality(robot_n170_rep_anova)
hetero_n170 <- check_heteroscedasticity(robot_n170_rep_anova)
print(norm_n170)
print('', quote = FALSE)
print(hetero_n170)
print('', quote = FALSE)
check_sphericity(robot_n170_rep_anova)
plot(norm_n170, type = 'density')
plot(norm_n170, type = 'qq')
plot(hetero_n170)Warning: Non-normality of residuals detected (p < .001).[1] Warning: Heteroscedasticity (non-constant error variance) detected (p = 0.002).[1] OK: Data seems to be spherical (p > .999).
Encore: Occipital-parietal theta-alpha (4-12 Hz) Functional Connectivity (ROIconnect)
Multivariate Interaction Measure, MIM (Pellegrini et al. 2023).
1 sec post-stimulus.
Code
summary(connectivity_data) ERPset face region_a region_b
MRS_04_ct_FT9: 8 robot:184 left_occipital :184 parietal_ipsi :184
MRS_07_ct_FT9: 8 human:184 right_occipital:184 parietal_contra:184
MRS_08_ct_FT9: 8
MRS_09_ct_FT9: 8
MRS_10_ct_FT9: 8
MRS_11_ct_FT9: 8
(Other) :320
median_connectivity log_median_connectivity num_id
Min. :0.04429 Min. :-3.117 Min. : 4.00
1st Qu.:0.07753 1st Qu.:-2.557 1st Qu.:17.00
Median :0.08960 Median :-2.412 Median :30.50
Mean :0.09370 Mean :-2.397 Mean :29.65
3rd Qu.:0.10655 3rd Qu.:-2.239 3rd Qu.:42.00
Max. :0.17459 Max. :-1.745 Max. :53.00
Code
connectivity_data_wide_occipital <- connectivity_data[c('ERPset', 'region_a', 'region_b', 'face', 'log_median_connectivity')] |>
pivot_wider(names_from = region_a, values_from = log_median_connectivity)
connectivity_data_wide_temporal <- connectivity_data[c('ERPset', 'region_a', 'region_b', 'face', 'log_median_connectivity')] |>
pivot_wider(names_from = region_b, values_from = log_median_connectivity)
connectivity_data_wide_face <- connectivity_data[c('ERPset', 'region_a', 'region_b', 'face', 'log_median_connectivity')] |>
pivot_wider(names_from = face, values_from = log_median_connectivity)
face <- c('robot' , 'human')
parietal <- c('parietal_ipsi' , 'parietal_contra')
occipital <- c('left_occipital', 'right_occipital')
occipital_pairs <- ggpairs(connectivity_data_wide_occipital,
aes(colour = face, alpha = .5, linewidth = NA),
columns = occipital,
lower = list(continuous = wrap('points', alpha = .4)),
progress = FALSE
)
parietal_pairs <- ggpairs(connectivity_data_wide_temporal,
aes(colour = face, alpha = .5, linewidth = NA),
columns = parietal,
lower = list(continuous = wrap('points', alpha = .4)),
progress = FALSE
)
face_pairs <- ggpairs(connectivity_data_wide_face ,
aes(colour = region_b, alpha = .5, linewidth = NA),
columns = face,
lower = list(continuous = wrap('points', alpha = .4)),
progress = FALSE
)
suppressWarnings(print(occipital_pairs))
suppressWarnings(print(parietal_pairs))
suppressWarnings(print(face_pairs))
Connectivity ANOVA:
Code
connectivity_rep_anova <- aov_ez('num_id', 'log_median_connectivity',
connectivity_data,
within = c('face', 'region_a', 'region_b'))
connectivity_afex_plot <- afex_plot(connectivity_rep_anova,
x = 'region_b',
trace = 'face',
panel = 'region_a',
error = 'within',
error_arg = list(width = .3, lwd = .75),
dodge = .5,
data_arg = list(
position =
position_jitterdodge(
jitter.width = .2,
dodge.width = .5 ## needs to be same as dodge
)),
mapping = c('color'),, data_alpha = .3,
point_arg = list(size = 3)
)
suppressWarnings(print(connectivity_afex_plot))
Code
mytable <- xtabs(~ face + region_a + region_b, data = connectivity_rep_anova$data$long)
ftable(mytable) region_b parietal_ipsi parietal_contra
face region_a
robot left_occipital 46 46
right_occipital 46 46
human left_occipital 46 46
right_occipital 46 46Code
print(connectivity_rep_anova)Anova Table (Type 3 tests)
Response: log_median_connectivity
Effect df MSE F ges p.value
1 face 1, 45 0.10 8.74 ** .040 .005
2 region_a 1, 45 0.01 0.12 <.001 .731
3 region_b 1, 45 0.01 24.96 *** .008 <.001
4 face:region_a 1, 45 0.01 0.01 <.001 .905
5 face:region_b 1, 45 0.00 0.87 <.001 .357
6 region_a:region_b 1, 45 0.04 1.01 .002 .319
7 face:region_a:region_b 1, 45 0.03 0.72 .001 .400
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1Code
a_posteriori(connectivity_rep_anova)____________________________________________________________
$emmeans
face emmean SE df lower.CL upper.CL
robot -2.45 0.0328 45 -2.51 -2.38
human -2.35 0.0286 45 -2.41 -2.29
Results are averaged over the levels of: region_b, region_a
Confidence level used: 0.95
$contrasts
contrast estimate SE df t.ratio p.value
robot - human -0.0969 0.0328 45 -2.956 0.0050
Results are averaged over the levels of: region_b, region_a
____________________________________________________________
$emmeans
region_b emmean SE df lower.CL upper.CL
parietal_ipsi -2.42 0.0274 45 -2.47 -2.36
parietal_contra -2.38 0.0254 45 -2.43 -2.32
Results are averaged over the levels of: region_a, face
Confidence level used: 0.95
$contrasts
contrast estimate SE df t.ratio p.value
parietal_ipsi - parietal_contra -0.0433 0.00866 45 -4.996 <.0001
Results are averaged over the levels of: region_a, face Connectivity Assumptions:
Code
norm_connectivity <- check_normality(connectivity_rep_anova)
hetero_connectivity <- check_heteroscedasticity(connectivity_rep_anova)
print(norm_connectivity)
print('', quote = FALSE)
print(hetero_connectivity)
print('', quote = FALSE)
check_sphericity(connectivity_rep_anova)
plot(norm_connectivity, type = 'density')
plot(norm_connectivity, type = 'qq')
plot(hetero_connectivity)OK: residuals appear as normally distributed (p = 0.128).[1] Warning: Heteroscedasticity (non-constant error variance) detected (p = 0.005).[1] OK: Data seems to be spherical (p > .999).
References
Delorme, Arnaud, and Scott Makeig. 2004. “EEGLAB: An Open Source Toolbox for Analysis of Single-Trial EEG Dynamics Including Independent Component Analysis.” Journal of Neuroscience Methods 134 (March): 9–21. https://doi.org/10.1016/J.JNEUMETH.2003.10.009.
Dong, Li, Fali Li, Qiang Liu, Xin Wen, Yongxiu Lai, Peng Xu, and Dezhong Yao. 2017. “MATLAB Toolboxes for Reference Electrode Standardization Technique (REST) of Scalp EEG.” Frontiers in Neuroscience 11 (October): 601. https://doi.org/10.3389/fnins.2017.00601.
Fahrenfort, Johannes Jacobus. 2020. “Multivariate Methods to Track the Spatiotemporal Profile of Feature-Based Attentional Selection Using EEG.” Neuromethods 151: 129–56. https://doi.org/10.1007/7657_2019_26/COVER.
Fields, Eric C., and Gina R. Kuperberg. 2020. “Having Your Cake and Eating It Too: Flexibility and Power with Mass Univariate Statistics for ERP Data.” Psychophysiology 57 (February): e13468. https://doi.org/10.1111/PSYP.13468.
Groppe, David M., Thomas P. Urbach, and Marta Kutas. 2011. “Mass Univariate Analysis of Event-Related Brain Potentials/Fields i: A Critical Tutorial Review.” Psychophysiology 48 (December): 1711–25. https://doi.org/10.1111/J.1469-8986.2011.01273.X.
Lopez-Calderon, Javier, and Steven J. Luck. 2014. “ERPLAB: An Open-Source Toolbox for the Analysis of Event-Related Potentials.” Frontiers in Human Neuroscience 8 (April): 75729. https://doi.org/10.3389/FNHUM.2014.00213/BIBTEX.
Oostenveld, Robert, Pascal Fries, Eric Maris, and Jan Mathijs Schoffelen. 2011. “FieldTrip: Open Source Software for Advanced Analysis of MEG, EEG, and Invasive Electrophysiological Data.” Computational Intelligence and Neuroscience 2011. https://doi.org/10.1155/2011/156869.
Pellegrini, Franziska, Arnaud Delorme, Vadim Nikulin, and Stefan Haufe. 2023. “Identifying Good Practices for Detecting Inter-Regional Linear Functional Connectivity from EEG.” NeuroImage 277 (August): 120218. https://doi.org/10.1016/J.NEUROIMAGE.2023.120218.
Pion-Tonachini, Luca, Ken Kreutz-Delgado, and Scott Makeig. 2019. “ICLabel: An Automated Electroencephalographic Independent Component Classifier, Dataset, and Website.” NeuroImage 198: 181–97. https://doi.org/10.1016/j.neuroimage.2019.05.026.