library(tidyverse)
library(tidyr)
library(dplyr)
library(readr)
library(purrr)
library(ggplot2)
library(e1071)
library(emmeans)
library(lme4)
library(lmerTest)
library(patchwork)
library(brms)
library(bayesplot)
library(car)
library(effects)
library(glue)
library(scales)
library(data.table)
library(effects)
# Disable emmeans computation limits for large models
emmeans::emm_options(
lmerTest.limit = Inf,
pbkrtest.limit = Inf
)Group CoM Analysis 3
# --- Load E-Prime file ---
RT <- read.csv("/Users/can/Documents/Uni/Thesis/Data/E-Prime/all_excluded2.csv", sep = ";")
# --- Filter, convert, exclude S12 ---
RTR <- RT %>%
dplyr::filter(procedure == "responsprocedure") %>%
dplyr::mutate(
feedback.ACC = as.numeric(feedback.ACC),
feedback.RT = as.numeric(feedback.RT)
) %>%
dplyr::filter(subject != 12)
# --- Trial numbering within subject × session ---
RTR <- RTR %>%
dplyr::group_by(subject, session) %>%
dplyr::mutate(trial = cumsum(sub.trial.number == 1)) %>%
dplyr::ungroup()
# --- Compute per-trial accuracy & mean RT ---
df <- RTR %>%
dplyr::group_by(subject, session, trial) %>%
dplyr::mutate(
trial.acc = sum(feedback.ACC, na.rm = TRUE) / dplyr::n(),
trial.RT = mean(feedback.RT, na.rm = TRUE)
) %>%
dplyr::ungroup()
# --- keep all trials once, tag correctness ---
df_acc_base <- df %>%
dplyr::mutate(
subject = as.factor(subject),
sub.trial.number = as.factor(sub.trial.number),
# keep session numeric here; convert to labels locally in plots when needed
is_correct = trial.acc == 1,
trial_acc_cat = factor(dplyr::if_else(trial.acc == 1, "correct", "wrong"),
levels = c("correct", "wrong"))
)
# helper to label sessions as "Block 1".."Block 5" on demand
as_block_factor <- function(x) factor(x, levels = 1:5, labels = paste("Block", 1:5))safe_fit <- function(formula, data) {
# Try lmer first; fall back to lm if random effects not identifiable
out <- tryCatch(lme4::lmer(formula, data = data), error = function(e) NULL)
if (is.null(out)) {
message("⚠️ Falling back to lm (random effect not identifiable).")
out <- stats::lm(update(formula, . ~ . - (1 | subject)), data = data)
}
out
}#1 RT per block
# ===== RT per block (correct vs wrong) — one-stop function =====
# - Fits RT ~ trial_acc_cat * Block + (1|subject)
# Helper (safe lmer with fallback)
if (!exists("safe_fit")) {
safe_fit <- function(formula, data) {
m <- tryCatch(lme4::lmer(formula, data = data), error = function(e) NULL)
if (is.null(m)) {
message("⚠️ Falling back to lm (random effect not identifiable).")
m <- stats::lm(update(formula, . ~ . - (1 | subject)), data = data)
}
m
}
}
# Helper (block labels)
if (!exists("as_block_factor")) {
as_block_factor <- function(x) factor(x, levels = 1:5, labels = paste("Block", 1:5))
}
rt_per_block_with_correctness <- function(df_acc_all, make_plot = TRUE) {
# 1) Build trial-level dataset
df_trial_all <- df_acc_all %>%
dplyr::distinct(subject, session, trial, .keep_all = TRUE) %>%
dplyr::transmute(subject, session_f = as_block_factor(session),
rt = trial.RT, trial_acc_cat)
# 2) Fit combined model with correctness factor
cat("\n\n========== RT per block — COMBINED (correctness as factor) ==========\n")
mdl_all <- safe_fit(rt ~ trial_acc_cat * session_f + (1 | subject), data = df_trial_all)
cat("\nModel summary:\n"); print(summary(mdl_all))
cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl_all, type = 2))
# 3) EMMs by block and correctness
emms_all <- emmeans::emmeans(mdl_all, ~ session_f | trial_acc_cat)
cat("\nEMMs (session | correctness):\n"); print(summary(emms_all))
cat("\nPairwise session (within correctness):\n")
print(emmeans::emmeans(mdl_all, pairwise ~ session_f | trial_acc_cat)$contrasts)
# 4) % faster (correct vs wrong) — overall and per block
overall_subtitle <- NULL
em_overall <- tryCatch(emmeans::emmeans(mdl_all, ~ trial_acc_cat) %>% as.data.frame(),
error = function(e) NULL)
if (!is.null(em_overall) && all(c("correct","wrong") %in% em_overall$trial_acc_cat)) {
mean_correct <- em_overall$emmean[em_overall$trial_acc_cat == "correct"]
mean_wrong <- em_overall$emmean[em_overall$trial_acc_cat == "wrong"]
pct_faster_overall <- 100 * (1 - mean_correct / mean_wrong)
cat(sprintf("\nOverall %% faster (correct vs wrong): %.2f%%\n", pct_faster_overall))
overall_subtitle <- sprintf("Overall %% faster (correct vs wrong): %.2f%%", pct_faster_overall)
} else {
cat("\nOverall %% faster: not available (only one correctness level present).\n")
}
em_per_block <- emmeans::emmeans(mdl_all, ~ trial_acc_cat | session_f) %>% as.data.frame()
if (all(c("correct","wrong") %in% em_per_block$trial_acc_cat)) {
pct_tbl <- em_per_block %>%
dplyr::select(session_f, trial_acc_cat, emmean) %>%
tidyr::pivot_wider(names_from = trial_acc_cat, values_from = emmean) %>%
dplyr::mutate(`% faster (correct vs wrong)` = 100 * (1 - correct / wrong)) %>%
dplyr::select(session_f, `% faster (correct vs wrong)`)
cat("\nPer-block %% faster (correct vs wrong):\n"); print(pct_tbl)
} else {
cat("\nPer-block %% faster: not available (only one correctness level present within blocks).\n")
}
# 5) Fitted EMM plot (always printed)
p <- as.data.frame(emms_all) %>%
ggplot2::ggplot(ggplot2::aes(x = session_f, y = emmean,
group = trial_acc_cat, color = trial_acc_cat)) +
ggplot2::geom_point() +
ggplot2::geom_errorbar(ggplot2::aes(ymin = emmean - SE, ymax = emmean + SE), width = .15) +
ggplot2::labs(title = "RT per block — fitted EMMs (correct vs wrong)",
subtitle = overall_subtitle,
x = "Block", y = "Estimated RT", color = "Correctness") +
ggplot2::theme_minimal(base_size = 12) +
ggplot2::theme(panel.grid = ggplot2::element_blank())
if (make_plot) print(p)
# 6) Split (correct-only / wrong-only) models — diagnostics in console (no extra plots)
for (lab in c("correct", "wrong")) {
df_trial <- df_trial_all %>% dplyr::filter(trial_acc_cat == lab)
if (nrow(df_trial) == 0) { cat("\n⚠️ No data for", lab, "\n"); next }
cat("\n\n========== RT per block —", toupper(lab), "==========\n")
mdl <- safe_fit(rt ~ session_f + (1 | subject), data = df_trial)
cat("\nModel summary:\n"); print(summary(mdl))
cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl, type = 2))
emms <- emmeans::emmeans(mdl, ~ session_f)
cat("\nEMMs (session):\n"); print(summary(emms))
cat("\nPairwise (Tukey):\n"); print(emmeans::emmeans(mdl, pairwise ~ session_f)$contrasts)
}
invisible(list(model = mdl_all, emms = emms_all, plot = p))
}
rt_per_block_with_correctness(df_acc_base)
========== RT per block — COMBINED (correctness as factor) ==========
Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ trial_acc_cat * session_f + (1 | subject)
Data: data
REML criterion at convergence: 59183.5
Scaled residuals:
Min 1Q Median 3Q Max
-4.5546 -0.4805 -0.0988 0.3305 17.6628
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 44156 210.1
Residual 51999 228.0
Number of obs: 4320, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 498.49 50.27 9.916
trial_acc_catwrong 362.46 20.05 18.077
session_fBlock 2 -16.31 12.95 -1.260
session_fBlock 3 21.83 13.64 1.600
session_fBlock 4 -36.97 12.41 -2.978
session_fBlock 5 68.17 13.71 4.972
trial_acc_catwrong:session_fBlock 2 -145.39 25.90 -5.613
trial_acc_catwrong:session_fBlock 3 -230.63 25.52 -9.037
trial_acc_catwrong:session_fBlock 4 -210.33 27.09 -7.765
trial_acc_catwrong:session_fBlock 5 -150.08 25.56 -5.872
Correlation of Fixed Effects:
(Intr) trl_c_ sss_B2 sss_B3 sss_B4 sss_B5 t__:_2 t__:_3 t__:_4
trl_cc_ctwr -0.074
sssn_fBlck2 -0.114 0.288
sssn_fBlck3 -0.108 0.271 0.424
sssn_fBlck4 -0.119 0.300 0.462 0.440
sssn_fBlck5 -0.108 0.274 0.419 0.397 0.436
trl_cc_:_B2 0.057 -0.773 -0.506 -0.217 -0.234 -0.213
trl_cc_:_B3 0.058 -0.782 -0.232 -0.541 -0.239 -0.215 0.614
trl_cc_:_B4 0.055 -0.739 -0.215 -0.207 -0.464 -0.202 0.575 0.586
trl_cc_:_B5 0.058 -0.786 -0.227 -0.214 -0.236 -0.545 0.609 0.616 0.580
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
trial_acc_cat 715.960 1 < 2.2e-16 ***
session_f 141.487 4 < 2.2e-16 ***
trial_acc_cat:session_f 92.146 4 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
EMMs (session | correctness):
trial_acc_cat = correct:
session_f emmean SE df lower.CL upper.CL
Block 1 498 50.3 17.9 393 604
Block 2 482 50.5 18.1 376 588
Block 3 520 50.6 18.4 414 627
Block 4 462 50.3 17.9 356 567
Block 5 567 50.7 18.4 460 673
trial_acc_cat = wrong:
session_f emmean SE df lower.CL upper.CL
Block 1 861 52.7 21.6 751 970
Block 2 699 51.3 19.3 592 806
Block 3 652 50.9 18.8 546 759
Block 4 614 52.0 20.4 505 722
Block 5 779 50.8 18.7 672 886
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
Pairwise session (within correctness):
trial_acc_cat = correct:
contrast estimate SE df t.ratio p.value
Block 1 - Block 2 16.3 12.9 4293 1.260 0.7159
Block 1 - Block 3 -21.8 13.6 4293 -1.600 0.4971
Block 1 - Block 4 37.0 12.4 4293 2.978 0.0243
Block 1 - Block 5 -68.2 13.7 4293 -4.972 <.0001
Block 2 - Block 3 -38.1 14.3 4293 -2.672 0.0584
Block 2 - Block 4 20.7 13.2 4293 1.570 0.5172
Block 2 - Block 5 -84.5 14.4 4293 -5.872 <.0001
Block 3 - Block 4 58.8 13.8 4293 4.254 0.0002
Block 3 - Block 5 -46.3 15.0 4294 -3.084 0.0175
Block 4 - Block 5 -105.1 13.9 4293 -7.552 <.0001
trial_acc_cat = wrong:
contrast estimate SE df t.ratio p.value
Block 1 - Block 2 161.7 22.3 4293 7.236 <.0001
Block 1 - Block 3 208.8 21.5 4293 9.731 <.0001
Block 1 - Block 4 247.3 24.0 4293 10.305 <.0001
Block 1 - Block 5 81.9 21.4 4294 3.822 0.0013
Block 2 - Block 3 47.1 17.5 4293 2.697 0.0545
Block 2 - Block 4 85.6 20.5 4293 4.169 0.0003
Block 2 - Block 5 -79.8 17.5 4294 -4.568 <.0001
Block 3 - Block 4 38.5 19.5 4293 1.971 0.2803
Block 3 - Block 5 -126.9 16.4 4294 -7.758 <.0001
Block 4 - Block 5 -165.4 19.6 4294 -8.454 <.0001
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 5 estimates NOTE: Results may be misleading due to involvement in interactions
Overall % faster (correct vs wrong): 29.84%
Per-block %% faster (correct vs wrong):
# A tibble: 5 × 2
session_f `% faster (correct vs wrong)`
<fct> <dbl>
1 Block 1 42.1
2 Block 2 31.0
3 Block 3 20.2
4 Block 4 24.8
5 Block 5 27.3
========== RT per block — CORRECT ==========
Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ session_f + (1 | subject)
Data: data
REML criterion at convergence: 37124.9
Scaled residuals:
Min 1Q Median 3Q Max
-2.8886 -0.5512 -0.1382 0.3723 10.4716
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 37404 193.4
Residual 25812 160.7
Number of obs: 2852, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 499.092 45.988 10.853
session_fBlock 2 -18.479 9.135 -2.023
session_fBlock 3 17.916 9.637 1.859
session_fBlock 4 -38.771 8.756 -4.428
session_fBlock 5 64.359 9.691 6.641
Correlation of Fixed Effects:
(Intr) sss_B2 sss_B3 sss_B4
sssn_fBlck2 -0.088
sssn_fBlck3 -0.083 0.425
sssn_fBlck4 -0.092 0.462 0.441
sssn_fBlck5 -0.083 0.418 0.395 0.435
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
session_f 123.77 4 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
EMMs (session):
session_f emmean SE df lower.CL upper.CL
Block 1 499 46.0 17.4 402 596
Block 2 481 46.1 17.6 384 578
Block 3 517 46.2 17.8 420 614
Block 4 460 46.0 17.5 363 557
Block 5 563 46.2 17.8 466 661
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
Pairwise (Tukey):
contrast estimate SE df t.ratio p.value
Block 1 - Block 2 18.5 9.14 2830 2.023 0.2552
Block 1 - Block 3 -17.9 9.64 2830 -1.859 0.3399
Block 1 - Block 4 38.8 8.76 2830 4.428 0.0001
Block 1 - Block 5 -64.4 9.69 2830 -6.641 <.0001
Block 2 - Block 3 -36.4 10.10 2830 -3.613 0.0028
Block 2 - Block 4 20.3 9.28 2830 2.186 0.1854
Block 2 - Block 5 -82.8 10.20 2830 -8.146 <.0001
Block 3 - Block 4 56.7 9.76 2830 5.811 <.0001
Block 3 - Block 5 -46.4 10.60 2831 -4.368 0.0001
Block 4 - Block 5 -103.1 9.84 2830 -10.480 <.0001
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 5 estimates
========== RT per block — WRONG ==========
Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ session_f + (1 | subject)
Data: data
REML criterion at convergence: 21091.5
Scaled residuals:
Min 1Q Median 3Q Max
-3.3493 -0.4884 -0.0940 0.3321 12.5268
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 58984 242.9
Residual 100214 316.6
Number of obs: 1468, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 862.51 62.58 13.783
session_fBlock 2 -165.20 31.34 -5.271
session_fBlock 3 -214.20 30.07 -7.124
session_fBlock 4 -252.45 33.64 -7.505
session_fBlock 5 -88.35 30.07 -2.938
Correlation of Fixed Effects:
(Intr) sss_B2 sss_B3 sss_B4
sssn_fBlck2 -0.325
sssn_fBlck3 -0.337 0.687
sssn_fBlck4 -0.302 0.611 0.638
sssn_fBlck5 -0.340 0.682 0.706 0.631
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
session_f 86.063 4 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
EMMs (session):
session_f emmean SE df lower.CL upper.CL
Block 1 863 62.6 23.1 733 992
Block 2 697 60.2 19.8 572 823
Block 3 648 59.6 19.1 524 773
Block 4 610 61.4 21.5 482 738
Block 5 774 59.5 18.9 650 899
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
Pairwise (Tukey):
contrast estimate SE df t.ratio p.value
Block 1 - Block 2 165.2 31.3 1448 5.270 <.0001
Block 1 - Block 3 214.2 30.1 1448 7.123 <.0001
Block 1 - Block 4 252.5 33.6 1448 7.504 <.0001
Block 1 - Block 5 88.3 30.1 1448 2.938 0.0277
Block 2 - Block 3 49.0 24.3 1447 2.015 0.2593
Block 2 - Block 4 87.3 28.7 1448 3.036 0.0206
Block 2 - Block 5 -76.8 24.5 1448 -3.134 0.0151
Block 3 - Block 4 38.3 27.3 1447 1.402 0.6263
Block 3 - Block 5 -125.8 23.1 1449 -5.454 <.0001
Block 4 - Block 5 -164.1 27.5 1449 -5.957 <.0001
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 5 estimates rt_block_full_report <- function(df, file = NULL, make_plot = TRUE, open = TRUE, width = 200) {
if (is.null(file)) {
file <- file.path(getwd(), sprintf("rt_block_output_%s.txt",
format(Sys.time(), "%Y-%m-%d_%H-%M-%S")))
}
old_width <- getOption("width")
options(width = width) # wider lines so nothing wraps oddly
con <- file(file, open = "wt")
# capture both stdout and messages
sink(con); sink(con, type = "message")
# ALWAYS restore sinks & width even if an error happens
on.exit({
try(sink(NULL), silent = TRUE)
try(sink(NULL, type = "message"), silent = TRUE)
try(close(con), silent = TRUE)
options(width = old_width)
}, add = TRUE)
# >>> run your existing function (prints go to the file)
res <- rt_per_block_with_correctness(df, make_plot = make_plot)
# restore and open
sink(NULL); sink(NULL, type = "message"); close(con)
options(width = old_width)
if (open) file.show(file)
invisible(list(result = res, path = file))
}
# ===== RUN (this will create & open a scrollable .txt with the full output) =====
rt_block_full_report(df_acc_base, make_plot = TRUE)NOTE: Results may be misleading due to involvement in interactions
#2.1 Correct-(or Wrong-)trials per block (counts, model & plot)
# ===== Combined plot: RT EMMs + % correct on x-axis =====
if (!exists("safe_fit")) {
safe_fit <- function(formula, data) {
m <- tryCatch(lme4::lmer(formula, data = data), error = function(e) NULL)
if (is.null(m)) {
message("⚠️ Falling back to lm (random effect not identifiable).")
m <- stats::lm(update(formula, . ~ . - (1 | subject)), data = data)
}
m
}
}
if (!exists("as_block_factor")) {
as_block_factor <- function(x) factor(x, levels = 1:5, labels = paste("Block", 1:5))
}
plot_rt_with_percent_correct <- function(df_acc_all, total_trials_per_block = 48, show_plot = TRUE, digits = 0) {
# ---------- RT model & EMMs (correct vs wrong by block) ----------
df_trial_all <- df_acc_all %>%
dplyr::distinct(subject, session, trial, .keep_all = TRUE) %>%
dplyr::transmute(subject,
session_f = as_block_factor(session),
rt = trial.RT,
trial_acc_cat)
mdl_rt <- safe_fit(rt ~ trial_acc_cat * session_f + (1 | subject), data = df_trial_all)
em_rt <- emmeans::emmeans(mdl_rt, ~ session_f | trial_acc_cat)
df_rt <- as.data.frame(em_rt)
# Overall % faster (subtitle)
overall_subtitle <- NULL
em_overall <- tryCatch(emmeans::emmeans(mdl_rt, ~ trial_acc_cat) %>% as.data.frame(), error = function(e) NULL)
if (!is.null(em_overall) && all(c("correct","wrong") %in% em_overall$trial_acc_cat)) {
mean_correct <- em_overall$emmean[em_overall$trial_acc_cat == "correct"]
mean_wrong <- em_overall$emmean[em_overall$trial_acc_cat == "wrong"]
pct_faster_overall <- 100 * (1 - mean_correct / mean_wrong)
overall_subtitle <- sprintf("Overall %% faster (correct vs wrong): %.2f%%", pct_faster_overall)
}
# ---------- % correct per block from GLMM (for x-axis labels) ----------
correct_counts <- df_acc_all %>%
dplyr::filter(trial_acc_cat == "correct") %>%
dplyr::distinct(subject, session, trial) %>%
dplyr::count(subject, session, name = "correct_trials") %>%
dplyr::mutate(session_f = as_block_factor(session),
fail_trials = total_trials_per_block - correct_trials)
label_map <- setNames(levels(df_rt$session_f), levels(df_rt$session_f)) # default
base_names <- c("Block 1" = "6 Steps",
"Block 2" = "12 Steps",
"Block 3" = "18 Steps",
"Block 4" = "Familiar",
"Block 5" = "Unfamiliar")
label_map <- base_names[names(label_map)]
if (nrow(correct_counts) > 0) {
mdl_pc <- tryCatch(
lme4::glmer(cbind(correct_trials, fail_trials) ~ session_f + (1 | subject),
data = correct_counts, family = binomial),
error = function(e) NULL
)
if (!is.null(mdl_pc)) {
em_pc <- emmeans::emmeans(mdl_pc, ~ session_f, type = "response")
df_pc <- as.data.frame(em_pc)
val_col <- if ("prob" %in% names(df_pc)) "prob" else if ("response" %in% names(df_pc)) "response" else "emmean"
df_pc <- df_pc %>% dplyr::transmute(session_f, pct = round(.data[[val_col]] * 100, digits))
} else {
# fallback: descriptive mean % correct out of 48
df_pc <- correct_counts %>%
dplyr::mutate(pct = 100 * correct_trials / total_trials_per_block) %>%
dplyr::group_by(session_f) %>%
dplyr::summarise(pct = round(mean(pct, na.rm = TRUE), digits), .groups = "drop")
}
add_map <- setNames(paste0(label_map[df_pc$session_f], " (", df_pc$pct, "%)"),
df_pc$session_f)
label_map[names(add_map)] <- add_map
}
# ---------- Combined plot (clean) ----------
p <- ggplot2::ggplot(df_rt, ggplot2::aes(x = session_f, y = emmean,
group = trial_acc_cat, color = trial_acc_cat)) +
# dotted divider between Block 3 and Block 4 (i.e., after the 3rd discrete position)
ggplot2::geom_vline(xintercept = 3.5, linetype = "dotted") +
ggplot2::geom_point() +
ggplot2::geom_errorbar(ggplot2::aes(ymin = emmean - SE, ymax = emmean + SE), width = 0.15) +
ggplot2::scale_x_discrete(labels = label_map) +
ggplot2::labs(title = "RT across difficulty levels",
subtitle = overall_subtitle,
x = "Difficulty",
y = "Estimated RT (in ms)",
color = "Correctness") +
ggplot2::theme_minimal(base_size = 12) +
ggplot2::theme(panel.grid = ggplot2::element_blank())
if (show_plot) print(p)
invisible(list(m_rt = mdl_rt, em_rt = em_rt,
m_pc = if (exists("mdl_pc")) mdl_pc else NULL,
xaxis_labels = label_map, plot = p))
}
# RUN
plot_rt_with_percent_correct(df_acc_base)NOTE: Results may be misleading due to involvement in interactions
# === Training-only plot: Blocks 1–3 ===
# RT EMMs (correct vs wrong) + % correct on x-axis
if (!exists("safe_fit")) {
safe_fit <- function(formula, data) {
m <- tryCatch(lme4::lmer(formula, data = data), error = function(e) NULL)
if (is.null(m)) {
message("⚠️ Falling back to lm (random effect not identifiable).")
m <- stats::lm(update(formula, . ~ . - (1 | subject)), data = data)
}
m
}
}
if (!exists("as_block_factor")) {
as_block_factor <- function(x) factor(x, levels = 1:5, labels = paste("Block", 1:5))
}
plot_rt_training_blocks <- function(df_acc_all, total_trials_per_block = 48, show_plot = TRUE, digits = 0) {
df_trial_all <- df_acc_all %>%
dplyr::distinct(subject, session, trial, .keep_all = TRUE) %>%
dplyr::transmute(subject,
session_f = as_block_factor(session),
rt = trial.RT,
trial_acc_cat) %>%
dplyr::filter(session_f %in% c("Block 1","Block 2","Block 3")) %>%
dplyr::mutate(session_f = factor(session_f, levels = c("Block 1","Block 2","Block 3")))
# Model & EMMs (training only)
mdl_rt <- safe_fit(rt ~ trial_acc_cat * session_f + (1 | subject), data = df_trial_all)
em_rt <- emmeans::emmeans(mdl_rt, ~ session_f | trial_acc_cat)
df_rt <- as.data.frame(em_rt)
# Overall % faster from training only
overall_subtitle <- NULL
em_overall <- tryCatch(emmeans::emmeans(mdl_rt, ~ trial_acc_cat) %>% as.data.frame(), error = function(e) NULL)
if (!is.null(em_overall) && all(c("correct","wrong") %in% em_overall$trial_acc_cat)) {
mc <- em_overall$emmean[em_overall$trial_acc_cat == "correct"]
mw <- em_overall$emmean[em_overall$trial_acc_cat == "wrong"]
pct_faster_overall <- 100 * (1 - mc / mw)
overall_subtitle <- sprintf("Overall %% faster (correct vs wrong): %.2f%%", pct_faster_overall)
}
# ---------- % correct per training block for x-axis labels ----------
correct_counts <- df_acc_all %>%
dplyr::filter(trial_acc_cat == "correct", session %in% 1:3) %>%
dplyr::distinct(subject, session, trial) %>%
dplyr::count(subject, session, name = "correct_trials") %>%
dplyr::mutate(session_f = factor(as_block_factor(session),
levels = c("Block 1","Block 2","Block 3")),
fail_trials = total_trials_per_block - correct_trials)
label_map <- c("Block 1" = "6 Steps",
"Block 2" = "12 Steps",
"Block 3" = "18 Steps")
if (nrow(correct_counts) > 0) {
mdl_pc <- tryCatch(
lme4::glmer(cbind(correct_trials, fail_trials) ~ session_f + (1 | subject),
data = correct_counts, family = binomial),
error = function(e) NULL
)
if (!is.null(mdl_pc)) {
em_pc <- emmeans::emmeans(mdl_pc, ~ session_f, type = "response")
df_pc <- as.data.frame(em_pc)
val_col <- if ("prob" %in% names(df_pc)) "prob" else if ("response" %in% names(df_pc)) "response" else "emmean"
df_pc <- df_pc %>% dplyr::transmute(session_f, pct = round(.data[[val_col]] * 100, digits))
} else {
df_pc <- correct_counts %>%
dplyr::mutate(pct = 100 * correct_trials / total_trials_per_block) %>%
dplyr::group_by(session_f) %>%
dplyr::summarise(pct = round(mean(pct, na.rm = TRUE), digits), .groups = "drop")
}
label_map[df_pc$session_f] <- paste0(label_map[df_pc$session_f], " (", df_pc$pct, "%)")
} else {
label_map <- paste0(label_map, " (n/a)")
}
pd <- ggplot2::position_dodge(width = 0.35)
p <- ggplot2::ggplot(df_rt, ggplot2::aes(x = session_f, y = emmean,
color = trial_acc_cat, group = trial_acc_cat)) +
ggplot2::geom_point(position = pd) +
ggplot2::geom_errorbar(ggplot2::aes(ymin = emmean - SE, ymax = emmean + SE),
width = 0.15, position = pd) +
ggplot2::scale_x_discrete(labels = label_map) +
ggplot2::labs(title = "RT across difficulty levels",
subtitle = overall_subtitle,
x = "difficulty",
y = "Estimated RT (in ms)",
color = "Correctness") +
ggplot2::theme_minimal(base_size = 12) +
ggplot2::theme(panel.grid = ggplot2::element_blank())
if (show_plot) print(p)
invisible(list(m_rt = mdl_rt, em_rt = em_rt,
xaxis_labels = label_map, plot = p))
}
# RUN (training-only figure)
plot_rt_training_blocks(df_acc_base)NOTE: Results may be misleading due to involvement in interactions
# === Training-only plots: Blocks 1–3 ===
# 1) Correct vs wrong (dodged points) with overall % faster computed on Blocks 1–3
# 2) Correct-only plot (EMMs by block)
# Needs: dplyr, tidyr, ggplot2, lme4, emmeans, car
if (!exists("safe_fit")) {
safe_fit <- function(formula, data) {
m <- tryCatch(lme4::lmer(formula, data = data), error = function(e) NULL)
if (is.null(m)) {
message("⚠️ Falling back to lm (random effect not identifiable).")
m <- stats::lm(update(formula, . ~ . - (1 | subject)), data = data)
}
m
}
}
if (!exists("as_block_factor")) {
as_block_factor <- function(x) factor(x, levels = 1:5, labels = paste("Block", 1:5))
}
plot_rt_training_blocks <- function(df_acc_all, total_trials_per_block = 48, show_plots = TRUE, digits = 0) {
df_trial_all <- df_acc_all %>%
dplyr::distinct(subject, session, trial, .keep_all = TRUE) %>%
dplyr::transmute(subject,
session_f = as_block_factor(session),
rt = trial.RT,
trial_acc_cat) %>%
dplyr::filter(session_f %in% c("Block 1","Block 2","Block 3")) %>%
dplyr::mutate(session_f = factor(session_f, levels = c("Block 1","Block 2","Block 3")))
# Model & EMMs (training only, correct vs wrong)
mdl_rt <- safe_fit(rt ~ trial_acc_cat * session_f + (1 | subject), data = df_trial_all)
em_rt <- emmeans::emmeans(mdl_rt, ~ session_f | trial_acc_cat)
df_rt <- as.data.frame(em_rt)
# Overall % faster from training only
overall_subtitle <- NULL
em_overall <- tryCatch(emmeans::emmeans(mdl_rt, ~ trial_acc_cat) %>% as.data.frame(), error = function(e) NULL)
if (!is.null(em_overall) && all(c("correct","wrong") %in% em_overall$trial_acc_cat)) {
mc <- em_overall$emmean[em_overall$trial_acc_cat == "correct"]
mw <- em_overall$emmean[em_overall$trial_acc_cat == "wrong"]
pct_faster_overall <- 100 * (1 - mc / mw)
overall_subtitle <- sprintf("Overall %% faster (correct vs wrong): %.2f%%", pct_faster_overall)
}
# ---------- % correct per training block for x-axis labels ----------
correct_counts <- df_acc_all %>%
dplyr::filter(trial_acc_cat == "correct", session %in% 1:3) %>%
dplyr::distinct(subject, session, trial) %>%
dplyr::count(subject, session, name = "correct_trials") %>%
dplyr::mutate(session_f = factor(as_block_factor(session),
levels = c("Block 1","Block 2","Block 3")),
fail_trials = total_trials_per_block - correct_trials)
label_map <- c("Block 1" = "6 Steps",
"Block 2" = "12 Steps",
"Block 3" = "18 Steps")
if (nrow(correct_counts) > 0) {
mdl_pc <- tryCatch(
lme4::glmer(cbind(correct_trials, fail_trials) ~ session_f + (1 | subject),
data = correct_counts, family = binomial),
error = function(e) NULL
)
if (!is.null(mdl_pc)) {
em_pc <- emmeans::emmeans(mdl_pc, ~ session_f, type = "response")
df_pc <- as.data.frame(em_pc)
val_col <- if ("prob" %in% names(df_pc)) "prob" else if ("response" %in% names(df_pc)) "response" else "emmean"
df_pc <- df_pc %>% dplyr::transmute(session_f, pct = round(.data[[val_col]] * 100, digits))
} else {
df_pc <- correct_counts %>%
dplyr::mutate(pct = 100 * correct_trials / total_trials_per_block) %>%
dplyr::group_by(session_f) %>%
dplyr::summarise(pct = round(mean(pct, na.rm = TRUE), digits), .groups = "drop")
}
label_map[df_pc$session_f] <- paste0(label_map[df_pc$session_f], " (", df_pc$pct, "%)")
} else {
label_map <- paste0(label_map, " (n/a)")
}
# ---------- Plot 1: correct vs wrong ----------
pd <- ggplot2::position_dodge(width = 0.35)
p_both <- ggplot2::ggplot(df_rt, ggplot2::aes(x = session_f, y = emmean,
color = trial_acc_cat, group = trial_acc_cat)) +
ggplot2::geom_point(position = pd) +
ggplot2::geom_errorbar(ggplot2::aes(ymin = emmean - SE, ymax = emmean + SE),
width = 0.15, position = pd) +
ggplot2::scale_x_discrete(labels = label_map) +
ggplot2::labs(title = "RT across difficulty levels",
subtitle = overall_subtitle,
x = "Difficulty",
y = "Estimated RT (in ms)",
color = "Correctness") +
ggplot2::theme_minimal(base_size = 12) +
ggplot2::theme(panel.grid = ggplot2::element_blank())
if (show_plots) print(p_both)
# ---------- Plot 2: correct-only ----------
df_correct <- df_trial_all %>% dplyr::filter(trial_acc_cat == "correct")
p_correct <- NULL
if (nrow(df_correct) > 0) {
mdl_corr <- safe_fit(rt ~ session_f + (1 | subject), data = df_correct)
em_corr <- emmeans::emmeans(mdl_corr, ~ session_f)
df_corr <- as.data.frame(em_corr)
p_correct <- ggplot2::ggplot(df_corr, ggplot2::aes(x = session_f, y = emmean, group = 1)) +
ggplot2::geom_point() +
ggplot2::geom_errorbar(ggplot2::aes(ymin = emmean - SE, ymax = emmean + SE), width = 0.15) +
ggplot2::scale_x_discrete(labels = label_map) +
ggplot2::labs(title = "RT across difficulty levels — Correct Trials",
x = "Difficulty (correct trials)",
y = "Estimated RT (in ms)") +
ggplot2::theme_minimal(base_size = 12) +
ggplot2::theme(panel.grid = ggplot2::element_blank())
if (show_plots) print(p_correct)
} else {
message("No correct trials found for Blocks 1–3.")
}
invisible(list(
m_rt = mdl_rt, em_rt = em_rt, plot_both = p_both,
plot_correct = p_correct,
xaxis_labels = label_map
))
}
# RUN (training-only figure + correct-only figure)
plot_rt_training_blocks(df_acc_base)NOTE: Results may be misleading due to involvement in interactions
# 2.1 — Statistical output
suppressPackageStartupMessages({
library(dplyr); library(tidyr); library(lme4); library(car); library(emmeans)
})
# Training blocks (1–3): trial-level data
df_train <- df_acc_base %>%
dplyr::distinct(subject, session, trial, .keep_all = TRUE) %>%
dplyr::filter(session %in% 1:3) %>%
dplyr::transmute(
subject,
session_f = factor(session, levels = 1:3, labels = paste("Block", 1:3)),
rt = trial.RT,
trial_acc_cat
)
# --- 2.1A: RT ~ correctness × block (Blocks 1–3) ---
mdl_rt_all <- tryCatch(
lme4::lmer(rt ~ trial_acc_cat * session_f + (1 | subject), data = df_train),
error = function(e) stats::lm(rt ~ trial_acc_cat * session_f, data = df_train)
)
cat("\n## 2.1A RT model (correct vs wrong × block)\n")
## 2.1A RT model (correct vs wrong × block)print(summary(mdl_rt_all))Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ trial_acc_cat * session_f + (1 | subject)
Data: df_train
REML criterion at convergence: 35859.1
Scaled residuals:
Min 1Q Median 3Q Max
-4.8327 -0.4429 -0.0856 0.2722 16.7623
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 57688 240.2
Residual 58826 242.5
Number of obs: 2592, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 497.02 57.35 8.667
trial_acc_catwrong 370.33 21.42 17.290
session_fBlock 2 -11.00 13.80 -0.797
session_fBlock 3 25.18 14.56 1.729
trial_acc_catwrong:session_fBlock 2 -164.19 27.70 -5.928
trial_acc_catwrong:session_fBlock 3 -242.64 27.31 -8.884
Correlation of Fixed Effects:
(Intr) trl_c_ sss_B2 sss_B3 t__:_2
trl_cc_ctwr -0.070
sssn_fBlck2 -0.106 0.288
sssn_fBlck3 -0.100 0.269 0.426
trl_cc_:_B2 0.054 -0.772 -0.508 -0.221
trl_cc_:_B3 0.054 -0.778 -0.235 -0.545 0.616cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl_rt_all, type = 2))
Type II Chi-square ANOVA:Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
trial_acc_cat 402.407 1 < 2.2e-16 ***
session_f 17.823 2 0.0001348 ***
trial_acc_cat:session_f 79.253 2 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1em_rt_all <- emmeans::emmeans(mdl_rt_all, ~ session_f | trial_acc_cat)
cat("\nEMMs by block within correctness:\n"); print(summary(em_rt_all))
EMMs by block within correctness:trial_acc_cat = correct:
session_f emmean SE df lower.CL upper.CL
Block 1 497 57.3 17.6 376 618
Block 2 486 57.5 17.9 365 607
Block 3 522 57.7 18.1 401 643
trial_acc_cat = wrong:
session_f emmean SE df lower.CL upper.CL
Block 1 867 59.8 20.9 743 992
Block 2 692 58.3 18.9 570 814
Block 3 650 58.0 18.4 528 771
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 cat("\nPairwise (Tukey) for blocks within correctness:\n"); print(pairs(em_rt_all, adjust = "tukey"))
Pairwise (Tukey) for blocks within correctness:trial_acc_cat = correct:
contrast estimate SE df t.ratio p.value
Block 1 - Block 2 11.0 13.8 2569 0.797 0.7047
Block 1 - Block 3 -25.2 14.6 2570 -1.729 0.1945
Block 2 - Block 3 -36.2 15.2 2569 -2.379 0.0459
trial_acc_cat = wrong:
contrast estimate SE df t.ratio p.value
Block 1 - Block 2 175.2 23.9 2570 7.343 <.0001
Block 1 - Block 3 217.5 22.9 2570 9.494 <.0001
Block 2 - Block 3 42.3 18.6 2569 2.273 0.0598
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 3 estimates em_acc <- emmeans::emmeans(mdl_rt_all, ~ trial_acc_cat)NOTE: Results may be misleading due to involvement in interactionscat("\nEMMs (correct vs wrong):\n"); print(summary(em_acc))
EMMs (correct vs wrong): trial_acc_cat emmean SE df lower.CL upper.CL
correct 502 56.9 17.1 382 622
wrong 736 57.3 17.6 616 857
Results are averaged over the levels of: session_f
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 cat("\nPairwise (correct vs wrong):\n"); print(pairs(em_acc))
Pairwise (correct vs wrong): contrast estimate SE df t.ratio p.value
correct - wrong -235 11 2572 -21.274 <.0001
Results are averaged over the levels of: session_f
Degrees-of-freedom method: kenward-roger # --- 2.1B: RT ~ block (correct trials only, Blocks 1–3) ---
df_train_corr <- df_train %>% dplyr::filter(trial_acc_cat == "correct")
if (nrow(df_train_corr) > 0 && dplyr::n_distinct(df_train_corr$session_f) > 1) {
mdl_rt_corr <- tryCatch(
lme4::lmer(rt ~ session_f + (1 | subject), data = df_train_corr),
error = function(e) stats::lm(rt ~ session_f, data = df_train_corr)
)
cat("\n## 2.1B RT model (correct-only)\n")
print(summary(mdl_rt_corr))
cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl_rt_corr, type = 2))
em_rt_corr <- emmeans::emmeans(mdl_rt_corr, ~ session_f)
cat("\nEMMs by block (correct-only):\n"); print(summary(em_rt_corr))
cat("\nPairwise (Tukey) between blocks (correct-only):\n"); print(pairs(em_rt_corr, adjust = "tukey"))
} else {
cat("\n[Correct-only model skipped: insufficient data]\n")
}
## 2.1B RT model (correct-only)
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ session_f + (1 | subject)
Data: df_train_corr
REML criterion at convergence: 22560
Scaled residuals:
Min 1Q Median 3Q Max
-3.0473 -0.5630 -0.1433 0.3475 10.7922
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 51782 227.6
Residual 24952 158.0
Number of obs: 1735, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 497.596 53.968 9.220
session_fBlock 2 -13.996 9.010 -1.553
session_fBlock 3 20.260 9.531 2.126
Correlation of Fixed Effects:
(Intr) sss_B2
sssn_fBlck2 -0.073
sssn_fBlck3 -0.069 0.428
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
session_f 11.947 2 0.002546 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
EMMs by block (correct-only):
session_f emmean SE df lower.CL upper.CL
Block 1 498 54.0 17.2 384 611
Block 2 484 54.1 17.4 370 597
Block 3 518 54.1 17.5 404 632
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
Pairwise (Tukey) between blocks (correct-only):
contrast estimate SE df t.ratio p.value
Block 1 - Block 2 14.0 9.01 1715 1.553 0.2664
Block 1 - Block 3 -20.3 9.53 1716 -2.126 0.0850
Block 2 - Block 3 -34.3 9.93 1715 -3.451 0.0017
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 3 estimates # --- 2.1C: Accuracy per block (binomial GLMM, Blocks 1–3) ---
total_trials_per_block <- 48L
counts_train <- df_acc_base %>%
dplyr::distinct(subject, session, trial, .keep_all = TRUE) %>%
dplyr::filter(session %in% 1:3) %>%
dplyr::filter(trial_acc_cat == "correct") %>%
dplyr::count(subject, session, name = "correct_trials") %>%
dplyr::mutate(
session_f = factor(session, levels = 1:3, labels = paste("Block", 1:3)),
fail_trials = total_trials_per_block - correct_trials
)
if (nrow(counts_train) > 0) {
mdl_acc <- lme4::glmer(cbind(correct_trials, fail_trials) ~ session_f + (1 | subject),
data = counts_train, family = binomial)
cat("\n## 2.1C Accuracy model (binomial, Blocks 1–3)\n")
print(summary(mdl_acc))
cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl_acc, type = 2))
em_acc_blk <- emmeans::emmeans(mdl_acc, ~ session_f, type = "response")
cat("\nEMMs (prob correct) by block:\n"); print(summary(em_acc_blk))
cat("\nPairwise (Tukey) between blocks (prob correct):\n"); print(pairs(em_acc_blk, adjust = "tukey"))
} else {
cat("\n[Accuracy model skipped: no counts available]\n")
}
## 2.1C Accuracy model (binomial, Blocks 1–3)
Generalized linear mixed model fit by maximum likelihood (Laplace
Approximation) [glmerMod]
Family: binomial ( logit )
Formula: cbind(correct_trials, fail_trials) ~ session_f + (1 | subject)
Data: counts_train
AIC BIC logLik deviance df.resid
401.8 409.7 -196.9 393.8 50
Scaled residuals:
Min 1Q Median 3Q Max
-2.9722 -1.1243 0.2848 1.0188 3.4491
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 0.4233 0.6506
Number of obs: 54, groups: subject, 18
Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 1.5991 0.1787 8.951 < 2e-16 ***
session_fBlock 2 -0.9251 0.1168 -7.924 2.31e-15 ***
session_fBlock 3 -1.3892 0.1156 -12.014 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Correlation of Fixed Effects:
(Intr) sss_B2
sssn_fBlck2 -0.393
sssn_fBlck3 -0.402 0.606
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: cbind(correct_trials, fail_trials)
Chisq Df Pr(>Chisq)
session_f 144.98 2 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
EMMs (prob correct) by block:
session_f prob SE df asymp.LCL asymp.UCL
Block 1 0.832 0.0250 Inf 0.777 0.875
Block 2 0.662 0.0382 Inf 0.584 0.733
Block 3 0.552 0.0419 Inf 0.470 0.632
Confidence level used: 0.95
Intervals are back-transformed from the logit scale
Pairwise (Tukey) between blocks (prob correct):
contrast odds.ratio SE df null z.ratio p.value
Block 1 / Block 2 2.52 0.294 Inf 1 7.924 <.0001
Block 1 / Block 3 4.01 0.464 Inf 1 12.014 <.0001
Block 2 / Block 3 1.59 0.164 Inf 1 4.499 <.0001
P value adjustment: tukey method for comparing a family of 3 estimates
Tests are performed on the log odds ratio scale # === NEW ANALYSIS (single model): Overall RT per block, ALL trials (Blocks 1–3) ===
# Fixed effects: rt ~ session_f
# Random effects: (1 | subject) + (1 + trial_acc_cat || subject:session_f)
# -> correctness (wrong – correct) varies randomly per participant × block
# -> no fixed correctness effect is included
suppressPackageStartupMessages({
library(dplyr); library(tidyr); library(lme4); library(car); library(emmeans); library(ggplot2)
})
# --- Data (Blocks 1–3; one row per trial) ---
df_overall_train <- df_acc_base %>%
dplyr::distinct(subject, session, trial, .keep_all = TRUE) %>%
dplyr::filter(session %in% 1:3) %>%
dplyr::transmute(
subject,
session_f = factor(session, levels = 1:3, labels = paste("Block", 1:3)),
rt = trial.RT,
trial_acc_cat = factor(trial_acc_cat, levels = c("correct","wrong"))
)
# --- Single, explicit model ---
form_model <- rt ~ session_f + (1 | subject) + (1 + trial_acc_cat | subject:session_f)
mdl_overall <- lme4::lmer(
formula = form_model,
data = df_overall_train,
REML = TRUE,
control = lme4::lmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 1e5))
)
cat("\n=== NEW ANALYSIS: Overall RT per block (ALL trials) ===\n")
=== NEW ANALYSIS: Overall RT per block (ALL trials) ===cat("Formula:\nrt ~ session_f + (1 | subject) + (1 + trial_acc_cat | subject:session_f)\n\n")Formula:
rt ~ session_f + (1 | subject) + (1 + trial_acc_cat | subject:session_f)# (Optional) Note if the fit is singular; we still proceed because this is the chosen model
if (lme4::isSingular(mdl_overall, tol = 1e-6)) {
cat("⚠️ Note: Fit is singular (some variance components ~ 0). Proceeding with the specified model.\n\n")
}
# --- Inference on overall block means (all trials) ---
cat("Model summary:\n"); print(summary(mdl_overall))Model summary:Linear mixed model fit by REML ['lmerMod']
Formula:
rt ~ session_f + (1 | subject) + (1 + trial_acc_cat | subject:session_f)
Data: df_overall_train
Control:
lme4::lmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 1e+05))
REML criterion at convergence: 35686.3
Scaled residuals:
Min 1Q Median 3Q Max
-5.5914 -0.3965 -0.1017 0.2285 15.5099
Random effects:
Groups Name Variance Std.Dev. Corr
subject:session_f (Intercept) 5865 76.58
trial_acc_catwrong 91590 302.64 -0.03
subject (Intercept) 48880 221.09
Residual 50080 223.79
Number of obs: 2592, groups: subject:session_f, 54; subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 506.06 55.78 9.073
session_fBlock 2 -13.14 28.51 -0.461
session_fBlock 3 23.49 29.03 0.809
Correlation of Fixed Effects:
(Intr) sss_B2
sssn_fBlck2 -0.249
sssn_fBlck3 -0.244 0.479cat("\nType II Chi-square ANOVA for session_f:\n"); print(car::Anova(mdl_overall, type = 2))
Type II Chi-square ANOVA for session_f:Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
session_f 1.5901 2 0.4516em_overall <- emmeans::emmeans(mdl_overall, ~ session_f)
cat("\nEstimated marginal means (overall RT by block):\n"); print(summary(em_overall))
Estimated marginal means (overall RT by block): session_f emmean SE df lower.CL upper.CL
Block 1 506 56.0 20.1 389 623
Block 2 493 56.2 20.4 376 610
Block 3 530 56.5 20.8 412 647
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 cat("\nPairwise (Tukey) between blocks:\n"); print(pairs(em_overall, adjust = "tukey"))
Pairwise (Tukey) between blocks: contrast estimate SE df t.ratio p.value
Block 1 - Block 2 13.1 29.3 32.1 0.448 0.8956
Block 1 - Block 3 -23.5 29.9 34.0 -0.787 0.7137
Block 2 - Block 3 -36.6 30.2 35.5 -1.213 0.4535
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 3 estimates # --- Plot: Overall RT per block (all trials) ---
df_em <- as.data.frame(em_overall)
p_overall_train <- ggplot(df_em, aes(x = session_f, y = emmean, group = 1)) +
geom_point() +
geom_errorbar(aes(ymin = emmean - SE, ymax = emmean + SE), width = 0.15) +
labs(title = "Overall RT across difficulty levels (Blocks 1–3, ALL trials)",
subtitle = "No fixed correctness; subject×block random correctness slope ",
x = "Difficulty", y = "Estimated RT (ms)") +
theme_minimal(base_size = 12) +
theme(panel.grid = element_blank())
print(p_overall_train)
#2.2 Correct Trials per block × sequence length — Test phase (Blocks 4–5)
# Test phase (Blocks 4–5)
if (!exists("safe_fit")) {
safe_fit <- function(formula, data) {
m <- tryCatch(lme4::lmer(formula, data = data), error = function(e) NULL)
if (is.null(m)) { message("⚠️ Falling back to lm."); m <- stats::lm(update(formula, . ~ . - (1 | subject)), data = data) }
m
}
}
plot_test_RT_by_length <- function(df_acc_all, total_trials_per_block = 48, digits = 0, show_plots = TRUE) {
trials_per_length <- total_trials_per_block / 3 # 48/3 = 16
df_len <- df_acc_all %>%
dplyr::group_by(subject, session, trial) %>%
dplyr::mutate(seq_length_trial = dplyr::n_distinct(sub.trial.number)) %>%
dplyr::ungroup() %>%
dplyr::filter(session %in% c(4, 5)) %>%
dplyr::mutate(
session_f = factor(session, levels = c(4, 5), labels = c("Familiar", "Unfamiliar")),
seq_length_trial = factor(seq_length_trial, levels = c(6, 12, 18))
)
if (nrow(df_len) == 0) { message("No test-phase data (Blocks 4–5)."); return(invisible(NULL)) }
# % correct per block (denom = 16 per length) -> facet strip labels
corr_len <- df_len %>%
dplyr::filter(trial_acc_cat == "correct") %>%
dplyr::distinct(subject, session_f, trial, seq_length_trial) %>%
dplyr::count(subject, session_f, seq_length_trial, name = "n_correct_len") %>%
dplyr::mutate(n_fail_len = trials_per_length - n_correct_len)
strip_labels <- setNames(levels(df_len$session_f), levels(df_len$session_f))
if (nrow(corr_len) > 0) {
mdl_pc <- tryCatch(
lme4::glmer(cbind(n_correct_len, n_fail_len) ~ session_f + (1 | subject) + (1 | seq_length_trial),
data = corr_len, family = binomial),
error = function(e) NULL
)
pc_tbl <- if (!is.null(mdl_pc)) {
em <- emmeans::emmeans(mdl_pc, ~ session_f, type = "response") %>% as.data.frame()
val <- if ("prob" %in% names(em)) "prob" else if ("response" %in% names(em)) "response" else "emmean"
dplyr::transmute(em, session_f, pct = round(.data[[val]] * 100, digits))
} else {
corr_len %>%
dplyr::mutate(pct_len = 100 * n_correct_len / trials_per_length) %>%
dplyr::group_by(session_f) %>%
dplyr::summarise(pct = round(mean(pct_len, na.rm = TRUE), digits), .groups = "drop")
}
strip_labels[pc_tbl$session_f] <- paste0(pc_tbl$session_f, " (", pc_tbl$pct, "%)")
}
strip_labeller <- ggplot2::as_labeller(strip_labels)
# helper: draw ONE dotted line exactly at the right edge of the left facet
add_separator <- function() {
ggplot2::geom_segment(
data = data.frame(session_f = "Familiar"),
ggplot2::aes(x = Inf, xend = Inf, y = -Inf, yend = Inf),
inherit.aes = FALSE, linetype = "dotted"
)
}
# ===== 1) Combined (correct vs wrong) =====
df_rt_trials <- df_len %>%
dplyr::distinct(subject, session_f, trial, seq_length_trial, .keep_all = TRUE) %>%
dplyr::transmute(subject, session_f, seq_length_trial, trial_acc_cat, rt = trial.RT)
mdl_rt <- safe_fit(rt ~ trial_acc_cat * session_f * seq_length_trial + (1 | subject), data = df_rt_trials)
em_rt <- emmeans::emmeans(mdl_rt, ~ seq_length_trial | session_f * trial_acc_cat)
df_em <- as.data.frame(em_rt)
pd <- ggplot2::position_dodge(width = 0.5)
p_combined <- ggplot2::ggplot(
df_em,
ggplot2::aes(x = seq_length_trial, y = emmean,
color = trial_acc_cat, shape = trial_acc_cat, group = trial_acc_cat)
) +
add_separator() +
ggplot2::geom_point(position = pd) +
ggplot2::geom_errorbar(ggplot2::aes(ymin = emmean - SE, ymax = emmean + SE), width = 0.2, position = pd) +
ggplot2::facet_wrap(~ session_f, nrow = 1, labeller = strip_labeller) +
ggplot2::labs(title = "RT across sequence lengths — test phase",
x = "Sequence length", y = "Estimated RT (in ms)",
color = "Correctness", shape = "Correctness") +
ggplot2::theme_minimal(base_size = 12) +
ggplot2::theme(panel.grid = ggplot2::element_blank())
if (show_plots) print(p_combined)
# ===== 2) Correct-only =====
p_correct <- NULL
df_rt_corr <- df_rt_trials %>% dplyr::filter(trial_acc_cat == "correct")
if (nrow(df_rt_corr) > 0) {
mdl_rt_corr <- safe_fit(rt ~ session_f * seq_length_trial + (1 | subject), data = df_rt_corr)
em_rt_corr <- emmeans::emmeans(mdl_rt_corr, ~ seq_length_trial | session_f)
df_rt_em <- as.data.frame(em_rt_corr)
# % correct per length (denom = 16) -> per-facet tick labels
mdl_pc_len <- tryCatch(
lme4::glmer(cbind(n_correct_len, n_fail_len) ~ session_f * seq_length_trial + (1 | subject),
data = corr_len, family = binomial),
error = function(e) NULL
)
lab_tbl <- if (!is.null(mdl_pc_len)) {
em <- emmeans::emmeans(mdl_pc_len, ~ seq_length_trial | session_f, type = "response") %>% as.data.frame()
val <- if ("prob" %in% names(em)) "prob" else if ("response" %in% names(em)) "response" else "emmean"
dplyr::transmute(em, session_f, seq_length_trial,
tick = paste0(seq_length_trial, " (", round(.data[[val]] * 100, digits), "%)"))
} else {
corr_len %>%
dplyr::mutate(pct = round(100 * n_correct_len / trials_per_length, digits)) %>%
dplyr::group_by(session_f, seq_length_trial) %>%
dplyr::summarise(pct = round(mean(pct, na.rm = TRUE), digits), .groups = "drop") %>%
dplyr::mutate(tick = paste0(seq_length_trial, " (", pct, "%)"))
}
df_rt_em <- df_rt_em %>%
dplyr::left_join(lab_tbl, by = c("session_f", "seq_length_trial")) %>%
dplyr::group_by(session_f) %>%
dplyr::mutate(seq_lab = factor(tick, levels = unique(tick))) %>%
dplyr::ungroup()
p_correct <- ggplot2::ggplot(df_rt_em, ggplot2::aes(x = seq_lab, y = emmean, group = 1)) +
add_separator() +
ggplot2::geom_point() +
ggplot2::geom_errorbar(ggplot2::aes(ymin = emmean - SE, ymax = emmean + SE), width = 0.2) +
ggplot2::facet_wrap(~ session_f, nrow = 1, scales = "free_x", labeller = strip_labeller) +
ggplot2::labs(title = "RT across Sequence Lengths — Correct Trials only (Test Phase)",
x = "Sequence length (Correct Trials)", y = "Estimated RT (in ms)") +
ggplot2::theme_minimal(base_size = 12) +
ggplot2::theme(panel.grid = ggplot2::element_blank())
if (show_plots) print(p_correct)
} else {
message("No correct trials available for the correct-only plot.")
}
invisible(list(m_rt = mdl_rt, em_rt = em_rt,
plot_combined = p_combined, plot_correct = p_correct,
facet_labels = strip_labels))
}
# RUN
plot_test_RT_by_length(df_acc_base)
# 2.2 — Statistical output (inline)
suppressPackageStartupMessages({
library(dplyr); library(tidyr); library(lme4); library(car); library(emmeans)
})
df_test <- df_acc_base %>%
dplyr::group_by(subject, session, trial) %>%
dplyr::mutate(seq_length_trial = dplyr::n_distinct(sub.trial.number)) %>% # TRUE length per trial
dplyr::ungroup() %>%
dplyr::distinct(subject, session, trial, .keep_all = TRUE) %>%
dplyr::filter(session %in% 4:5) %>%
dplyr::mutate(
session_f = factor(session, levels = c(4, 5), labels = paste("Block", 4:5)),
seq_length_trial = factor(seq_length_trial, levels = c(6, 12, 18)),
rt = trial.RT
)
# --- 2.2A: RT ~ block × length × correctness (Blocks 4–5) ---
mdl_rt_test <- tryCatch(
lme4::lmer(rt ~ session_f * seq_length_trial * trial_acc_cat + (1 | subject), data = df_test),
error = function(e) stats::lm(rt ~ session_f * seq_length_trial * trial_acc_cat, data = df_test)
)
cat("\n## 2.2A RT model (block × length × correctness)\n")
## 2.2A RT model (block × length × correctness)print(summary(mdl_rt_test))Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ session_f * seq_length_trial * trial_acc_cat + (1 | subject)
Data: df_test
REML criterion at convergence: 22875.5
Scaled residuals:
Min 1Q Median 3Q Max
-3.9916 -0.5338 -0.0872 0.4120 10.1695
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 32914 181.4
Residual 33321 182.5
Number of obs: 1728, groups: subject, 18
Fixed effects:
Estimate Std. Error
(Intercept) 447.163 44.291
session_fBlock 5 103.205 16.960
seq_length_trial12 20.515 16.793
seq_length_trial18 25.478 17.939
trial_acc_catwrong 246.071 32.480
session_fBlock 5:seq_length_trial12 38.251 26.179
session_fBlock 5:seq_length_trial18 3.666 27.839
session_fBlock 5:trial_acc_catwrong 33.737 41.344
seq_length_trial12:trial_acc_catwrong -109.485 41.456
seq_length_trial18:trial_acc_catwrong -124.187 39.263
session_fBlock 5:seq_length_trial12:trial_acc_catwrong -27.321 53.138
session_fBlock 5:seq_length_trial18:trial_acc_catwrong 26.747 51.517
t value
(Intercept) 10.096
session_fBlock 5 6.085
seq_length_trial12 1.222
seq_length_trial18 1.420
trial_acc_catwrong 7.576
session_fBlock 5:seq_length_trial12 1.461
session_fBlock 5:seq_length_trial18 0.132
session_fBlock 5:trial_acc_catwrong 0.816
seq_length_trial12:trial_acc_catwrong -2.641
seq_length_trial18:trial_acc_catwrong -3.163
session_fBlock 5:seq_length_trial12:trial_acc_catwrong -0.514
session_fBlock 5:seq_length_trial18:trial_acc_catwrong 0.519
Correlation of Fixed Effects:
(Intr) sss_B5 sq__12 sq__18 trl_c_ ss_B5:__12 ss_B5:__18 ss_B5:__
sssn_fBlck5 -0.177
sq_lngth_12 -0.179 0.466
sq_lngth_18 -0.167 0.434 0.442
trl_cc_ctwr -0.094 0.245 0.244 0.226
sss_B5:__12 0.115 -0.643 -0.641 -0.283 -0.158
sss_B5:__18 0.107 -0.602 -0.284 -0.642 -0.142 0.394
sssn_fB5:__ 0.074 -0.421 -0.191 -0.174 -0.783 0.264 0.241
sq_ln_12:__ 0.073 -0.190 -0.408 -0.179 -0.776 0.262 0.114 0.608
sq_ln_18:__ 0.077 -0.198 -0.203 -0.460 -0.819 0.130 0.292 0.638
s_B5:__12:_ -0.057 0.322 0.318 0.138 0.605 -0.495 -0.193 -0.767
s_B5:__18:_ -0.058 0.330 0.154 0.347 0.620 -0.214 -0.540 -0.787
s__12: s__18: s_B5:__12:
sssn_fBlck5
sq_lngth_12
sq_lngth_18
trl_cc_ctwr
sss_B5:__12
sss_B5:__18
sssn_fB5:__
sq_ln_12:__
sq_ln_18:__ 0.641
s_B5:__12:_ -0.779 -0.498
s_B5:__18:_ -0.486 -0.756 0.612 cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl_rt_test, type = 2))
Type II Chi-square ANOVA:Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
session_f 203.5009 1 < 2.2e-16 ***
seq_length_trial 0.5873 2 0.74552
trial_acc_cat 296.8455 1 < 2.2e-16 ***
session_f:seq_length_trial 1.6037 2 0.44850
session_f:trial_acc_cat 2.8424 1 0.09181 .
seq_length_trial:trial_acc_cat 26.4726 2 1.785e-06 ***
session_f:seq_length_trial:trial_acc_cat 1.3761 2 0.50255
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# --- Overall correct vs wrong (averaged over block and length) ---
em_rt_acc_overall <- emmeans::emmeans(mdl_rt_test, ~ trial_acc_cat)NOTE: Results may be misleading due to involvement in interactionscat("\nEMMs (correct vs wrong), averaged over block and length:\n")
EMMs (correct vs wrong), averaged over block and length:print(summary(em_rt_acc_overall)) trial_acc_cat emmean SE df lower.CL upper.CL
correct 521 43.2 17.3 430 612
wrong 706 43.6 18.0 614 798
Results are averaged over the levels of: session_f, seq_length_trial
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 cat("\nPairwise (correct vs wrong):\n")
Pairwise (correct vs wrong):print(pairs(em_rt_acc_overall)) contrast estimate SE df t.ratio p.value
correct - wrong -185 10.6 1701 -17.488 <.0001
Results are averaged over the levels of: session_f, seq_length_trial
Degrees-of-freedom method: kenward-roger em_rt_acc_byblock <- emmeans::emmeans(mdl_rt_test, ~ trial_acc_cat | session_f)NOTE: Results may be misleading due to involvement in interactionscat("\nEMMs (correct vs wrong) within each block:\n")
EMMs (correct vs wrong) within each block:print(summary(em_rt_acc_byblock))session_f = Block 4:
trial_acc_cat emmean SE df lower.CL upper.CL
correct 462 43.4 17.6 371 554
wrong 631 45.0 20.4 537 724
session_f = Block 5:
trial_acc_cat emmean SE df lower.CL upper.CL
correct 580 43.7 18.1 488 671
wrong 781 44.0 18.6 689 874
Results are averaged over the levels of: seq_length_trial
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 cat("\nPairwise (correct vs wrong) within each block:\n")
Pairwise (correct vs wrong) within each block:print(pairs(em_rt_acc_byblock))session_f = Block 4:
contrast estimate SE df t.ratio p.value
correct - wrong -168 16.0 1701 -10.512 <.0001
session_f = Block 5:
contrast estimate SE df t.ratio p.value
correct - wrong -202 13.9 1702 -14.539 <.0001
Results are averaged over the levels of: seq_length_trial
Degrees-of-freedom method: kenward-roger # EMMs and within-block (length) pairwise
em_rt_len <- emmeans::emmeans(mdl_rt_test, ~ seq_length_trial | session_f * trial_acc_cat)
cat("\nEMMs by length within block × correctness:\n"); print(summary(em_rt_len))
EMMs by length within block × correctness:session_f = Block 4, trial_acc_cat = correct:
seq_length_trial emmean SE df lower.CL upper.CL
6 447 44.3 19.1 355 540
12 468 44.5 19.5 375 561
18 473 44.9 20.3 379 566
session_f = Block 5, trial_acc_cat = correct:
seq_length_trial emmean SE df lower.CL upper.CL
6 550 44.5 19.6 457 643
12 609 45.6 21.6 514 704
18 580 46.2 22.6 484 675
session_f = Block 4, trial_acc_cat = wrong:
seq_length_trial emmean SE df lower.CL upper.CL
6 693 52.4 37.4 587 799
12 604 48.6 27.6 505 704
18 595 46.2 22.8 499 690
session_f = Block 5, trial_acc_cat = wrong:
seq_length_trial emmean SE df lower.CL upper.CL
6 830 48.2 26.8 731 929
12 752 45.3 20.9 658 846
18 762 45.0 20.3 668 856
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 cat("\nPairwise (Tukey) across lengths within block × correctness:\n"); print(pairs(em_rt_len, adjust = "tukey"))
Pairwise (Tukey) across lengths within block × correctness:session_f = Block 4, trial_acc_cat = correct:
contrast estimate SE df t.ratio p.value
seq_length_trial6 - seq_length_trial12 -20.52 16.8 1699 -1.222 0.4405
seq_length_trial6 - seq_length_trial18 -25.48 17.9 1699 -1.420 0.3306
seq_length_trial12 - seq_length_trial18 -4.96 18.4 1699 -0.270 0.9606
session_f = Block 5, trial_acc_cat = correct:
contrast estimate SE df t.ratio p.value
seq_length_trial6 - seq_length_trial12 -58.77 20.1 1699 -2.926 0.0097
seq_length_trial6 - seq_length_trial18 -29.14 21.3 1699 -1.366 0.3591
seq_length_trial12 - seq_length_trial18 29.62 23.4 1699 1.263 0.4162
session_f = Block 4, trial_acc_cat = wrong:
contrast estimate SE df t.ratio p.value
seq_length_trial6 - seq_length_trial12 88.97 37.9 1699 2.351 0.0494
seq_length_trial6 - seq_length_trial18 98.71 34.9 1699 2.831 0.0130
seq_length_trial12 - seq_length_trial18 9.74 28.8 1699 0.338 0.9390
session_f = Block 5, trial_acc_cat = wrong:
contrast estimate SE df t.ratio p.value
seq_length_trial6 - seq_length_trial12 78.04 26.5 1700 2.945 0.0092
seq_length_trial6 - seq_length_trial18 68.30 26.0 1700 2.623 0.0239
seq_length_trial12 - seq_length_trial18 -9.74 20.2 1699 -0.483 0.8795
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 3 estimates # Cross-block contrasts within each length × correctness (Block 4 vs Block 5)
em_rt_blk <- emmeans::emmeans(mdl_rt_test, ~ session_f | seq_length_trial * trial_acc_cat)
cat("\nBlock comparisons (RT) within each length × correctness (Block 4 vs Block 5):\n")
Block comparisons (RT) within each length × correctness (Block 4 vs Block 5):print(pairs(em_rt_blk))seq_length_trial = 6, trial_acc_cat = correct:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -103 17.0 1699 -6.085 <.0001
seq_length_trial = 12, trial_acc_cat = correct:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -141 20.1 1699 -7.054 <.0001
seq_length_trial = 18, trial_acc_cat = correct:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -107 22.2 1700 -4.807 <.0001
seq_length_trial = 6, trial_acc_cat = wrong:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -137 37.5 1700 -3.651 0.0003
seq_length_trial = 12, trial_acc_cat = wrong:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -148 27.4 1700 -5.406 <.0001
seq_length_trial = 18, trial_acc_cat = wrong:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -167 22.4 1700 -7.472 <.0001
Degrees-of-freedom method: kenward-roger cat("\nEMMs (RT) by block within length × correctness:\n")
EMMs (RT) by block within length × correctness:print(summary(em_rt_blk))seq_length_trial = 6, trial_acc_cat = correct:
session_f emmean SE df lower.CL upper.CL
Block 4 447 44.3 19.1 355 540
Block 5 550 44.5 19.6 457 643
seq_length_trial = 12, trial_acc_cat = correct:
session_f emmean SE df lower.CL upper.CL
Block 4 468 44.5 19.5 375 561
Block 5 609 45.6 21.6 514 704
seq_length_trial = 18, trial_acc_cat = correct:
session_f emmean SE df lower.CL upper.CL
Block 4 473 44.9 20.3 379 566
Block 5 580 46.2 22.6 484 675
seq_length_trial = 6, trial_acc_cat = wrong:
session_f emmean SE df lower.CL upper.CL
Block 4 693 52.4 37.4 587 799
Block 5 830 48.2 26.8 731 929
seq_length_trial = 12, trial_acc_cat = wrong:
session_f emmean SE df lower.CL upper.CL
Block 4 604 48.6 27.6 505 704
Block 5 752 45.3 20.9 658 846
seq_length_trial = 18, trial_acc_cat = wrong:
session_f emmean SE df lower.CL upper.CL
Block 4 595 46.2 22.8 499 690
Block 5 762 45.0 20.3 668 856
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 # --- 2.2B: RT ~ block × length (correct-only, Blocks 4–5) ---
df_test_corr <- df_test %>% dplyr::filter(trial_acc_cat == "correct")
if (nrow(df_test_corr) > 0) {
mdl_rt_test_corr <- tryCatch(
lme4::lmer(rt ~ session_f * seq_length_trial + (1 | subject), data = df_test_corr),
error = function(e) stats::lm(rt ~ session_f * seq_length_trial, data = df_test_corr)
)
cat("\n## 2.2B RT model (correct-only)\n")
print(summary(mdl_rt_test_corr))
cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl_rt_test_corr, type = 2))
em_rt_len_corr <- emmeans::emmeans(mdl_rt_test_corr, ~ seq_length_trial | session_f)
cat("\nEMMs by length within block (correct-only):\n"); print(summary(em_rt_len_corr))
cat("\nPairwise (Tukey) across lengths within block (correct-only):\n"); print(pairs(em_rt_len_corr, adjust = "tukey"))
# Cross-block contrasts within each length (correct-only)
em_rt_blk_corr <- emmeans::emmeans(mdl_rt_test_corr, ~ session_f | seq_length_trial)
cat("\nBlock comparisons (RT, correct-only) within each length (Block 4 vs Block 5):\n")
print(pairs(em_rt_blk_corr))
cat("\nEMMs (RT, correct-only) by block within length:\n")
print(summary(em_rt_blk_corr))
} else {
cat("\n[Correct-only RT model skipped: insufficient data]\n")
}
## 2.2B RT model (correct-only)
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ session_f * seq_length_trial + (1 | subject)
Data: df_test_corr
REML criterion at convergence: 14254.3
Scaled residuals:
Min 1Q Median 3Q Max
-3.5853 -0.5408 -0.1110 0.4088 7.5602
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 29756 172.5
Residual 19832 140.8
Number of obs: 1117, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 447.2439 41.6235 10.745
session_fBlock 5 105.3359 13.1253 8.025
seq_length_trial12 21.1635 12.9629 1.633
seq_length_trial18 23.7453 13.8623 1.713
session_fBlock 5:seq_length_trial12 35.2801 20.2159 1.745
session_fBlock 5:seq_length_trial18 0.3142 21.5168 0.015
Correlation of Fixed Effects:
(Intr) sss_B5 sq__12 sq__18 s_B5:__12
sssn_fBlck5 -0.145
sq_lngth_12 -0.146 0.464
sq_lngth_18 -0.137 0.429 0.441
sss_B5:__12 0.094 -0.640 -0.641 -0.282
sss_B5:__18 0.088 -0.597 -0.284 -0.641 0.394
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
session_f 174.7094 1 < 2e-16 ***
seq_length_trial 13.5533 2 0.00114 **
session_f:seq_length_trial 3.5825 2 0.16675
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
EMMs by length within block (correct-only):
session_f = Block 4:
seq_length_trial emmean SE df lower.CL upper.CL
6 447 41.6 18.3 360 535
12 468 41.7 18.5 381 556
18 471 42.0 19.0 383 559
session_f = Block 5:
seq_length_trial emmean SE df lower.CL upper.CL
6 553 41.8 18.6 465 640
12 609 42.5 19.8 520 698
18 577 42.9 20.5 487 666
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
Pairwise (Tukey) across lengths within block (correct-only):
session_f = Block 4:
contrast estimate SE df t.ratio p.value
seq_length_trial6 - seq_length_trial12 -21.16 13.0 1094 -1.633 0.2323
seq_length_trial6 - seq_length_trial18 -23.75 13.9 1094 -1.713 0.2008
seq_length_trial12 - seq_length_trial18 -2.58 14.2 1094 -0.182 0.9819
session_f = Block 5:
contrast estimate SE df t.ratio p.value
seq_length_trial6 - seq_length_trial12 -56.44 15.5 1094 -3.638 0.0008
seq_length_trial6 - seq_length_trial18 -24.06 16.5 1094 -1.457 0.3123
seq_length_trial12 - seq_length_trial18 32.38 18.1 1094 1.787 0.1744
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 3 estimates
Block comparisons (RT, correct-only) within each length (Block 4 vs Block 5):
seq_length_trial = 6:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -105 13.1 1094 -8.025 <.0001
seq_length_trial = 12:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -141 15.5 1095 -9.051 <.0001
seq_length_trial = 18:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -106 17.3 1095 -6.117 <.0001
Degrees-of-freedom method: kenward-roger
EMMs (RT, correct-only) by block within length:
seq_length_trial = 6:
session_f emmean SE df lower.CL upper.CL
Block 4 447 41.6 18.3 360 535
Block 5 553 41.8 18.6 465 640
seq_length_trial = 12:
session_f emmean SE df lower.CL upper.CL
Block 4 468 41.7 18.5 381 556
Block 5 609 42.5 19.8 520 698
seq_length_trial = 18:
session_f emmean SE df lower.CL upper.CL
Block 4 471 42.0 19.0 383 559
Block 5 577 42.9 20.5 487 666
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 # --- 2.2C: Accuracy per block × length (binomial GLMM, Blocks 4–5) ---
counts_test <- df_acc_base %>%
dplyr::group_by(subject, session, trial) %>%
dplyr::mutate(seq_length_trial = dplyr::n_distinct(sub.trial.number)) %>% # TRUE length
dplyr::ungroup() %>%
dplyr::distinct(subject, session, trial, .keep_all = TRUE) %>%
dplyr::filter(session %in% 4:5) %>%
dplyr::mutate(
session_f = factor(session, levels = c(4, 5), labels = paste("Block", 4:5)),
seq_length_trial = factor(seq_length_trial, levels = c(6, 12, 18))
) %>%
dplyr::group_by(subject, session_f, seq_length_trial) %>%
dplyr::summarise(
total_trials = dplyr::n(), # actual denominator after any exclusions
correct_trials = sum(trial_acc_cat == "correct", na.rm = TRUE),
.groups = "drop"
) %>%
dplyr::mutate(fail_trials = total_trials - correct_trials) %>%
dplyr::filter(!is.na(seq_length_trial))
if (nrow(counts_test) > 0) {
mdl_acc_test <- lme4::glmer(
cbind(correct_trials, fail_trials) ~ session_f * seq_length_trial + (1 | subject),
data = counts_test, family = binomial
)
cat("\n## 2.2C Accuracy model (binomial, block × length)\n")
print(summary(mdl_acc_test))
cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl_acc_test, type = 2))
# EMMs on probability (response) scale
em_acc_test <- emmeans::emmeans(mdl_acc_test, ~ seq_length_trial | session_f, type = "response")
cat("\nEMMs (prob correct) by length within block:\n"); print(summary(em_acc_test))
cat("\nPairwise (Tukey) across lengths within block (prob correct):\n"); print(pairs(em_acc_test, adjust = "tukey"))
# Cross-block contrasts within each length (prob correct; response scale)
em_acc_blk <- emmeans::emmeans(mdl_acc_test, ~ session_f | seq_length_trial, type = "response")
cat("\nBlock comparisons (prob correct) within each length (Block 4 vs Block 5):\n")
print(pairs(em_acc_blk))
cat("\nEMMs (prob correct) by block within length:\n")
print(summary(em_acc_blk))
} else {
cat("\n[Accuracy model skipped: no counts available]\n")
}
## 2.2C Accuracy model (binomial, block × length)
Generalized linear mixed model fit by maximum likelihood (Laplace
Approximation) [glmerMod]
Family: binomial ( logit )
Formula: cbind(correct_trials, fail_trials) ~ session_f * seq_length_trial +
(1 | subject)
Data: counts_test
AIC BIC logLik deviance df.resid
584.0 602.8 -285.0 570.0 101
Scaled residuals:
Min 1Q Median 3Q Max
-4.9843 -0.8720 0.1909 0.9885 4.0483
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 0.3185 0.5644
Number of obs: 108, groups: subject, 18
Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 2.0307 0.2245 9.044 < 2e-16 ***
session_fBlock 5 -0.8171 0.2284 -3.577 0.000347 ***
seq_length_trial12 -0.6935 0.2310 -3.002 0.002682 **
seq_length_trial18 -1.5158 0.2196 -6.903 5.1e-12 ***
session_fBlock 5:seq_length_trial12 -0.6721 0.2973 -2.261 0.023789 *
session_fBlock 5:seq_length_trial18 -0.1854 0.2890 -0.641 0.521275
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Correlation of Fixed Effects:
(Intr) sss_B5 sq__12 sq__18 s_B5:__12
sssn_fBlck5 -0.631
sq_lngth_12 -0.623 0.610
sq_lngth_18 -0.660 0.644 0.636
sss_B5:__12 0.477 -0.766 -0.775 -0.489
sss_B5:__18 0.493 -0.788 -0.481 -0.753 0.607
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: cbind(correct_trials, fail_trials)
Chisq Df Pr(>Chisq)
session_f 98.2015 1 < 2e-16 ***
seq_length_trial 130.4703 2 < 2e-16 ***
session_f:seq_length_trial 5.9554 2 0.05091 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
EMMs (prob correct) by length within block:
session_f = Block 4:
seq_length_trial prob SE df asymp.LCL asymp.UCL
6 0.884 0.0230 Inf 0.831 0.922
12 0.792 0.0326 Inf 0.721 0.849
18 0.626 0.0429 Inf 0.539 0.706
session_f = Block 5:
seq_length_trial prob SE df asymp.LCL asymp.UCL
6 0.771 0.0344 Inf 0.697 0.831
12 0.462 0.0450 Inf 0.376 0.550
18 0.380 0.0431 Inf 0.300 0.468
Confidence level used: 0.95
Intervals are back-transformed from the logit scale
Pairwise (Tukey) across lengths within block (prob correct):
session_f = Block 4:
contrast odds.ratio SE df null z.ratio
seq_length_trial6 / seq_length_trial12 2.00 0.462 Inf 1 3.002
seq_length_trial6 / seq_length_trial18 4.55 1.000 Inf 1 6.903
seq_length_trial12 / seq_length_trial18 2.28 0.438 Inf 1 4.272
p.value
0.0076
<.0001
0.0001
session_f = Block 5:
contrast odds.ratio SE df null z.ratio
seq_length_trial6 / seq_length_trial12 3.92 0.736 Inf 1 7.268
seq_length_trial6 / seq_length_trial18 5.48 1.040 Inf 1 8.947
seq_length_trial12 / seq_length_trial18 1.40 0.245 Inf 1 1.917
p.value
<.0001
<.0001
0.1339
P value adjustment: tukey method for comparing a family of 3 estimates
Tests are performed on the log odds ratio scale
Block comparisons (prob correct) within each length (Block 4 vs Block 5):
seq_length_trial = 6:
contrast odds.ratio SE df null z.ratio p.value
Block 4 / Block 5 2.26 0.517 Inf 1 3.577 0.0003
seq_length_trial = 12:
contrast odds.ratio SE df null z.ratio p.value
Block 4 / Block 5 4.43 0.848 Inf 1 7.786 <.0001
seq_length_trial = 18:
contrast odds.ratio SE df null z.ratio p.value
Block 4 / Block 5 2.72 0.484 Inf 1 5.639 <.0001
Tests are performed on the log odds ratio scale
EMMs (prob correct) by block within length:
seq_length_trial = 6:
session_f prob SE df asymp.LCL asymp.UCL
Block 4 0.884 0.0230 Inf 0.831 0.922
Block 5 0.771 0.0344 Inf 0.697 0.831
seq_length_trial = 12:
session_f prob SE df asymp.LCL asymp.UCL
Block 4 0.792 0.0326 Inf 0.721 0.849
Block 5 0.462 0.0450 Inf 0.376 0.550
seq_length_trial = 18:
session_f prob SE df asymp.LCL asymp.UCL
Block 4 0.626 0.0429 Inf 0.539 0.706
Block 5 0.380 0.0431 Inf 0.300 0.468
Confidence level used: 0.95
Intervals are back-transformed from the logit scale # 2.2 — Correct Trials per block × sequence length — Test phase (Blocks 4–5)
if (!exists("safe_fit")) {
safe_fit <- function(formula, data) {
m <- tryCatch(lme4::lmer(formula, data = data), error = function(e) NULL)
if (is.null(m)) { message("⚠️ Falling back to lm."); m <- stats::lm(update(formula, . ~ . - (1 | subject)), data = data) }
m
}
}
plot_test_RT_by_length <- function(df_acc_all, total_trials_per_block = 48, digits = 0, show_plots = TRUE) {
df_len <- df_acc_all %>%
dplyr::group_by(subject, session, trial) %>%
dplyr::mutate(seq_length_trial = dplyr::n_distinct(sub.trial.number)) %>%
dplyr::ungroup() %>%
dplyr::filter(session %in% c(4, 5)) %>%
dplyr::mutate(
session_f = factor(session, levels = c(4, 5), labels = paste("Block", 4:5)),
seq_length_trial = factor(seq_length_trial, levels = c(6, 12, 18))
)
if (nrow(df_len) == 0) { message("No test-phase data (Blocks 4–5)."); return(invisible(NULL)) }
# Build per-subject counts with ACTUAL denominators (matches 2.2C stats)
counts_len <- df_len %>%
dplyr::distinct(subject, session_f, trial, seq_length_trial, .keep_all = TRUE) %>%
dplyr::group_by(subject, session_f, seq_length_trial) %>%
dplyr::summarise(
total_trials = dplyr::n(),
correct_trials = sum(trial_acc_cat == "correct", na.rm = TRUE),
.groups = "drop"
) %>%
dplyr::mutate(n_fail_len = total_trials - correct_trials)
# facet strip labels: prob(correct) per block (from GLMM EMMs)
strip_labels <- setNames(levels(df_len$session_f), levels(df_len$session_f))
if (nrow(counts_len) > 0) {
mdl_pc <- lme4::glmer(
cbind(correct_trials, n_fail_len) ~ session_f + (1 | subject) + (1 | seq_length_trial),
data = counts_len, family = binomial
)
em <- emmeans::emmeans(mdl_pc, ~ session_f, type = "response") %>% as.data.frame()
val <- if ("prob" %in% names(em)) "prob" else if ("response" %in% names(em)) "response" else "emmean"
pc_tbl <- dplyr::transmute(em, session_f, pct = round(.data[[val]] * 100, digits))
strip_labels[pc_tbl$session_f] <- paste0(pc_tbl$session_f, " (", pc_tbl$pct, "%)")
}
strip_labeller <- ggplot2::as_labeller(strip_labels)
first_level <- levels(df_len$session_f)[1]
add_separator <- function() {
ggplot2::geom_segment(
data = data.frame(session_f = first_level),
ggplot2::aes(x = Inf, xend = Inf, y = -Inf, yend = Inf),
inherit.aes = FALSE, linetype = "dotted"
)
}
# ===== 1) Combined (correct vs wrong) =====
df_rt_trials <- df_len %>%
dplyr::distinct(subject, session_f, trial, seq_length_trial, .keep_all = TRUE) %>%
dplyr::transmute(subject, session_f, seq_length_trial, trial_acc_cat, rt = trial.RT)
mdl_rt <- safe_fit(rt ~ trial_acc_cat * session_f * seq_length_trial + (1 | subject), data = df_rt_trials)
em_rt <- emmeans::emmeans(mdl_rt, ~ seq_length_trial | session_f * trial_acc_cat)
df_em <- as.data.frame(em_rt)
pd <- ggplot2::position_dodge(width = 0.5)
p_combined <- ggplot2::ggplot(
df_em,
ggplot2::aes(x = seq_length_trial, y = emmean,
color = trial_acc_cat, shape = trial_acc_cat, group = trial_acc_cat)
) +
add_separator() +
ggplot2::geom_point(position = pd) +
ggplot2::geom_errorbar(ggplot2::aes(ymin = emmean - SE, ymax = emmean + SE), width = 0.2, position = pd) +
ggplot2::facet_wrap(~ session_f, nrow = 1, labeller = strip_labeller) +
ggplot2::labs(title = "RT across sequence lengths — Test phase",
x = "Sequence length", y = "Estimated RT (ms)",
color = "Correctness", shape = "Correctness") +
ggplot2::theme_minimal(base_size = 12) +
ggplot2::theme(panel.grid = ggplot2::element_blank())
if (show_plots) print(p_combined)
# ===== 2) Correct-only =====
p_correct <- NULL
df_rt_corr <- df_rt_trials %>% dplyr::filter(trial_acc_cat == "correct")
if (nrow(df_rt_corr) > 0) {
mdl_rt_corr <- safe_fit(rt ~ session_f * seq_length_trial + (1 | subject), data = df_rt_corr)
em_rt_corr <- emmeans::emmeans(mdl_rt_corr, ~ seq_length_trial | session_f)
df_rt_em <- as.data.frame(em_rt_corr)
# tick labels: prob(correct) per length × block (from GLMM EMMs, same counts_len)
mdl_pc_len <- lme4::glmer(
cbind(correct_trials, n_fail_len) ~ session_f * seq_length_trial + (1 | subject),
data = counts_len, family = binomial
)
em <- emmeans::emmeans(mdl_pc_len, ~ seq_length_trial | session_f, type = "response") %>% as.data.frame()
val <- if ("prob" %in% names(em)) "prob" else if ("response" %in% names(em)) "response" else "emmean"
lab_tbl <- dplyr::transmute(
em, session_f, seq_length_trial,
tick = paste0(seq_length_trial, " (", round(.data[[val]] * 100, digits), "%)")
)
df_rt_em <- df_rt_em %>%
dplyr::left_join(lab_tbl, by = c("session_f", "seq_length_trial")) %>%
dplyr::group_by(session_f) %>%
dplyr::mutate(seq_lab = factor(tick, levels = unique(tick))) %>%
dplyr::ungroup()
p_correct <- ggplot2::ggplot(df_rt_em, ggplot2::aes(x = seq_lab, y = emmean, group = 1)) +
add_separator() +
ggplot2::geom_point() +
ggplot2::geom_errorbar(ggplot2::aes(ymin = emmean - SE, ymax = emmean + SE), width = 0.2) +
ggplot2::facet_wrap(~ session_f, nrow = 1, scales = "free_x", labeller = strip_labeller) +
ggplot2::labs(title = "RT across Sequence Lengths — Correct trials only (Test phase)",
x = "Sequence length (Correct trials)", y = "Estimated RT (ms)") +
ggplot2::theme_minimal(base_size = 12) +
ggplot2::theme(panel.grid = ggplot2::element_blank())
if (show_plots) print(p_correct)
} else {
message("No correct trials available for the correct-only plot.")
}
invisible(list(m_rt = mdl_rt, em_rt = em_rt,
plot_combined = p_combined, plot_correct = p_correct,
facet_labels = strip_labels))
}
# RUN
plot_test_RT_by_length(df_acc_base)
# 2.2 — Statistical output
suppressPackageStartupMessages({
library(dplyr); library(tidyr); library(lme4); library(car); library(emmeans)
})
# Derive true sequence length per trial (matches plot code)
df_test <- df_acc_base %>%
dplyr::group_by(subject, session, trial) %>%
dplyr::mutate(seq_length_trial = dplyr::n_distinct(sub.trial.number)) %>% # TRUE length per trial
dplyr::ungroup() %>%
dplyr::distinct(subject, session, trial, .keep_all = TRUE) %>%
dplyr::filter(session %in% 4:5) %>%
dplyr::mutate(
session_f = factor(session, levels = c(4, 5), labels = paste("Block", 4:5)),
seq_length_trial = factor(seq_length_trial, levels = c(6, 12, 18)),
rt = trial.RT
)
# --- 2.2A: RT ~ block × length × correctness (Blocks 4–5) ---
mdl_rt_test <- tryCatch(
lme4::lmer(rt ~ session_f * seq_length_trial * trial_acc_cat + (1 | subject), data = df_test),
error = function(e) stats::lm(rt ~ session_f * seq_length_trial * trial_acc_cat, data = df_test)
)
cat("\n## 2.2A RT model (block × length × correctness)\n")
## 2.2A RT model (block × length × correctness)print(summary(mdl_rt_test))Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ session_f * seq_length_trial * trial_acc_cat + (1 | subject)
Data: df_test
REML criterion at convergence: 22875.5
Scaled residuals:
Min 1Q Median 3Q Max
-3.9916 -0.5338 -0.0872 0.4120 10.1695
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 32914 181.4
Residual 33321 182.5
Number of obs: 1728, groups: subject, 18
Fixed effects:
Estimate Std. Error
(Intercept) 447.163 44.291
session_fBlock 5 103.205 16.960
seq_length_trial12 20.515 16.793
seq_length_trial18 25.478 17.939
trial_acc_catwrong 246.071 32.480
session_fBlock 5:seq_length_trial12 38.251 26.179
session_fBlock 5:seq_length_trial18 3.666 27.839
session_fBlock 5:trial_acc_catwrong 33.737 41.344
seq_length_trial12:trial_acc_catwrong -109.485 41.456
seq_length_trial18:trial_acc_catwrong -124.187 39.263
session_fBlock 5:seq_length_trial12:trial_acc_catwrong -27.321 53.138
session_fBlock 5:seq_length_trial18:trial_acc_catwrong 26.747 51.517
t value
(Intercept) 10.096
session_fBlock 5 6.085
seq_length_trial12 1.222
seq_length_trial18 1.420
trial_acc_catwrong 7.576
session_fBlock 5:seq_length_trial12 1.461
session_fBlock 5:seq_length_trial18 0.132
session_fBlock 5:trial_acc_catwrong 0.816
seq_length_trial12:trial_acc_catwrong -2.641
seq_length_trial18:trial_acc_catwrong -3.163
session_fBlock 5:seq_length_trial12:trial_acc_catwrong -0.514
session_fBlock 5:seq_length_trial18:trial_acc_catwrong 0.519
Correlation of Fixed Effects:
(Intr) sss_B5 sq__12 sq__18 trl_c_ ss_B5:__12 ss_B5:__18 ss_B5:__
sssn_fBlck5 -0.177
sq_lngth_12 -0.179 0.466
sq_lngth_18 -0.167 0.434 0.442
trl_cc_ctwr -0.094 0.245 0.244 0.226
sss_B5:__12 0.115 -0.643 -0.641 -0.283 -0.158
sss_B5:__18 0.107 -0.602 -0.284 -0.642 -0.142 0.394
sssn_fB5:__ 0.074 -0.421 -0.191 -0.174 -0.783 0.264 0.241
sq_ln_12:__ 0.073 -0.190 -0.408 -0.179 -0.776 0.262 0.114 0.608
sq_ln_18:__ 0.077 -0.198 -0.203 -0.460 -0.819 0.130 0.292 0.638
s_B5:__12:_ -0.057 0.322 0.318 0.138 0.605 -0.495 -0.193 -0.767
s_B5:__18:_ -0.058 0.330 0.154 0.347 0.620 -0.214 -0.540 -0.787
s__12: s__18: s_B5:__12:
sssn_fBlck5
sq_lngth_12
sq_lngth_18
trl_cc_ctwr
sss_B5:__12
sss_B5:__18
sssn_fB5:__
sq_ln_12:__
sq_ln_18:__ 0.641
s_B5:__12:_ -0.779 -0.498
s_B5:__18:_ -0.486 -0.756 0.612 cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl_rt_test, type = 2))
Type II Chi-square ANOVA:Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
session_f 203.5009 1 < 2.2e-16 ***
seq_length_trial 0.5873 2 0.74552
trial_acc_cat 296.8455 1 < 2.2e-16 ***
session_f:seq_length_trial 1.6037 2 0.44850
session_f:trial_acc_cat 2.8424 1 0.09181 .
seq_length_trial:trial_acc_cat 26.4726 2 1.785e-06 ***
session_f:seq_length_trial:trial_acc_cat 1.3761 2 0.50255
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# Overall correct vs wrong (averaged over block and length)
em_rt_acc_overall <- emmeans::emmeans(mdl_rt_test, ~ trial_acc_cat)NOTE: Results may be misleading due to involvement in interactionscat("\nEMMs (correct vs wrong), averaged over block and length:\n")
EMMs (correct vs wrong), averaged over block and length:print(summary(em_rt_acc_overall)) trial_acc_cat emmean SE df lower.CL upper.CL
correct 521 43.2 17.3 430 612
wrong 706 43.6 18.0 614 798
Results are averaged over the levels of: session_f, seq_length_trial
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 cat("\nPairwise (correct vs wrong):\n")
Pairwise (correct vs wrong):print(pairs(em_rt_acc_overall)) contrast estimate SE df t.ratio p.value
correct - wrong -185 10.6 1701 -17.488 <.0001
Results are averaged over the levels of: session_f, seq_length_trial
Degrees-of-freedom method: kenward-roger # correct vs wrong within each block
em_rt_acc_byblock <- emmeans::emmeans(mdl_rt_test, ~ trial_acc_cat | session_f)NOTE: Results may be misleading due to involvement in interactionscat("\nEMMs (correct vs wrong) within each block:\n")
EMMs (correct vs wrong) within each block:print(summary(em_rt_acc_byblock))session_f = Block 4:
trial_acc_cat emmean SE df lower.CL upper.CL
correct 462 43.4 17.6 371 554
wrong 631 45.0 20.4 537 724
session_f = Block 5:
trial_acc_cat emmean SE df lower.CL upper.CL
correct 580 43.7 18.1 488 671
wrong 781 44.0 18.6 689 874
Results are averaged over the levels of: seq_length_trial
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 cat("\nPairwise (correct vs wrong) within each block:\n")
Pairwise (correct vs wrong) within each block:print(pairs(em_rt_acc_byblock))session_f = Block 4:
contrast estimate SE df t.ratio p.value
correct - wrong -168 16.0 1701 -10.512 <.0001
session_f = Block 5:
contrast estimate SE df t.ratio p.value
correct - wrong -202 13.9 1702 -14.539 <.0001
Results are averaged over the levels of: seq_length_trial
Degrees-of-freedom method: kenward-roger # EMMs and within-block (length) pairwise
em_rt_len <- emmeans::emmeans(mdl_rt_test, ~ seq_length_trial | session_f * trial_acc_cat)
cat("\nEMMs by length within block × correctness:\n"); print(summary(em_rt_len))
EMMs by length within block × correctness:session_f = Block 4, trial_acc_cat = correct:
seq_length_trial emmean SE df lower.CL upper.CL
6 447 44.3 19.1 355 540
12 468 44.5 19.5 375 561
18 473 44.9 20.3 379 566
session_f = Block 5, trial_acc_cat = correct:
seq_length_trial emmean SE df lower.CL upper.CL
6 550 44.5 19.6 457 643
12 609 45.6 21.6 514 704
18 580 46.2 22.6 484 675
session_f = Block 4, trial_acc_cat = wrong:
seq_length_trial emmean SE df lower.CL upper.CL
6 693 52.4 37.4 587 799
12 604 48.6 27.6 505 704
18 595 46.2 22.8 499 690
session_f = Block 5, trial_acc_cat = wrong:
seq_length_trial emmean SE df lower.CL upper.CL
6 830 48.2 26.8 731 929
12 752 45.3 20.9 658 846
18 762 45.0 20.3 668 856
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 cat("\nPairwise (Tukey) across lengths within block × correctness:\n"); print(pairs(em_rt_len, adjust = "tukey"))
Pairwise (Tukey) across lengths within block × correctness:session_f = Block 4, trial_acc_cat = correct:
contrast estimate SE df t.ratio p.value
seq_length_trial6 - seq_length_trial12 -20.52 16.8 1699 -1.222 0.4405
seq_length_trial6 - seq_length_trial18 -25.48 17.9 1699 -1.420 0.3306
seq_length_trial12 - seq_length_trial18 -4.96 18.4 1699 -0.270 0.9606
session_f = Block 5, trial_acc_cat = correct:
contrast estimate SE df t.ratio p.value
seq_length_trial6 - seq_length_trial12 -58.77 20.1 1699 -2.926 0.0097
seq_length_trial6 - seq_length_trial18 -29.14 21.3 1699 -1.366 0.3591
seq_length_trial12 - seq_length_trial18 29.62 23.4 1699 1.263 0.4162
session_f = Block 4, trial_acc_cat = wrong:
contrast estimate SE df t.ratio p.value
seq_length_trial6 - seq_length_trial12 88.97 37.9 1699 2.351 0.0494
seq_length_trial6 - seq_length_trial18 98.71 34.9 1699 2.831 0.0130
seq_length_trial12 - seq_length_trial18 9.74 28.8 1699 0.338 0.9390
session_f = Block 5, trial_acc_cat = wrong:
contrast estimate SE df t.ratio p.value
seq_length_trial6 - seq_length_trial12 78.04 26.5 1700 2.945 0.0092
seq_length_trial6 - seq_length_trial18 68.30 26.0 1700 2.623 0.0239
seq_length_trial12 - seq_length_trial18 -9.74 20.2 1699 -0.483 0.8795
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 3 estimates # Cross-block contrasts within each length × correctness (Block 4 vs Block 5)
em_rt_blk <- emmeans::emmeans(mdl_rt_test, ~ session_f | seq_length_trial * trial_acc_cat)
cat("\nBlock comparisons (RT) within each length × correctness (Block 4 vs Block 5):\n")
Block comparisons (RT) within each length × correctness (Block 4 vs Block 5):print(pairs(em_rt_blk))seq_length_trial = 6, trial_acc_cat = correct:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -103 17.0 1699 -6.085 <.0001
seq_length_trial = 12, trial_acc_cat = correct:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -141 20.1 1699 -7.054 <.0001
seq_length_trial = 18, trial_acc_cat = correct:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -107 22.2 1700 -4.807 <.0001
seq_length_trial = 6, trial_acc_cat = wrong:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -137 37.5 1700 -3.651 0.0003
seq_length_trial = 12, trial_acc_cat = wrong:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -148 27.4 1700 -5.406 <.0001
seq_length_trial = 18, trial_acc_cat = wrong:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -167 22.4 1700 -7.472 <.0001
Degrees-of-freedom method: kenward-roger cat("\nEMMs (RT) by block within length × correctness:\n")
EMMs (RT) by block within length × correctness:print(summary(em_rt_blk))seq_length_trial = 6, trial_acc_cat = correct:
session_f emmean SE df lower.CL upper.CL
Block 4 447 44.3 19.1 355 540
Block 5 550 44.5 19.6 457 643
seq_length_trial = 12, trial_acc_cat = correct:
session_f emmean SE df lower.CL upper.CL
Block 4 468 44.5 19.5 375 561
Block 5 609 45.6 21.6 514 704
seq_length_trial = 18, trial_acc_cat = correct:
session_f emmean SE df lower.CL upper.CL
Block 4 473 44.9 20.3 379 566
Block 5 580 46.2 22.6 484 675
seq_length_trial = 6, trial_acc_cat = wrong:
session_f emmean SE df lower.CL upper.CL
Block 4 693 52.4 37.4 587 799
Block 5 830 48.2 26.8 731 929
seq_length_trial = 12, trial_acc_cat = wrong:
session_f emmean SE df lower.CL upper.CL
Block 4 604 48.6 27.6 505 704
Block 5 752 45.3 20.9 658 846
seq_length_trial = 18, trial_acc_cat = wrong:
session_f emmean SE df lower.CL upper.CL
Block 4 595 46.2 22.8 499 690
Block 5 762 45.0 20.3 668 856
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 # --- 2.2B: RT ~ block × length (correct-only, Blocks 4–5) ---
df_test_corr <- df_test %>% dplyr::filter(trial_acc_cat == "correct")
if (nrow(df_test_corr) > 0) {
mdl_rt_test_corr <- tryCatch(
lme4::lmer(rt ~ session_f * seq_length_trial + (1 | subject), data = df_test_corr),
error = function(e) stats::lm(rt ~ session_f * seq_length_trial, data = df_test_corr)
)
cat("\n## 2.2B RT model (correct-only)\n")
print(summary(mdl_rt_test_corr))
cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl_rt_test_corr, type = 2))
em_rt_len_corr <- emmeans::emmeans(mdl_rt_test_corr, ~ seq_length_trial | session_f)
cat("\nEMMs by length within block (correct-only):\n"); print(summary(em_rt_len_corr))
cat("\nPairwise (Tukey) across lengths within block (correct-only):\n"); print(pairs(em_rt_len_corr, adjust = "tukey"))
# Cross-block contrasts within each length (correct-only)
em_rt_blk_corr <- emmeans::emmeans(mdl_rt_test_corr, ~ session_f | seq_length_trial)
cat("\nBlock comparisons (RT, correct-only) within each length (Block 4 vs Block 5):\n")
print(pairs(em_rt_blk_corr))
cat("\nEMMs (RT, correct-only) by block within length:\n")
print(summary(em_rt_blk_corr))
} else {
cat("\n[Correct-only RT model skipped: insufficient data]\n")
}
## 2.2B RT model (correct-only)
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ session_f * seq_length_trial + (1 | subject)
Data: df_test_corr
REML criterion at convergence: 14254.3
Scaled residuals:
Min 1Q Median 3Q Max
-3.5853 -0.5408 -0.1110 0.4088 7.5602
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 29756 172.5
Residual 19832 140.8
Number of obs: 1117, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 447.2439 41.6235 10.745
session_fBlock 5 105.3359 13.1253 8.025
seq_length_trial12 21.1635 12.9629 1.633
seq_length_trial18 23.7453 13.8623 1.713
session_fBlock 5:seq_length_trial12 35.2801 20.2159 1.745
session_fBlock 5:seq_length_trial18 0.3142 21.5168 0.015
Correlation of Fixed Effects:
(Intr) sss_B5 sq__12 sq__18 s_B5:__12
sssn_fBlck5 -0.145
sq_lngth_12 -0.146 0.464
sq_lngth_18 -0.137 0.429 0.441
sss_B5:__12 0.094 -0.640 -0.641 -0.282
sss_B5:__18 0.088 -0.597 -0.284 -0.641 0.394
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
session_f 174.7094 1 < 2e-16 ***
seq_length_trial 13.5533 2 0.00114 **
session_f:seq_length_trial 3.5825 2 0.16675
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
EMMs by length within block (correct-only):
session_f = Block 4:
seq_length_trial emmean SE df lower.CL upper.CL
6 447 41.6 18.3 360 535
12 468 41.7 18.5 381 556
18 471 42.0 19.0 383 559
session_f = Block 5:
seq_length_trial emmean SE df lower.CL upper.CL
6 553 41.8 18.6 465 640
12 609 42.5 19.8 520 698
18 577 42.9 20.5 487 666
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
Pairwise (Tukey) across lengths within block (correct-only):
session_f = Block 4:
contrast estimate SE df t.ratio p.value
seq_length_trial6 - seq_length_trial12 -21.16 13.0 1094 -1.633 0.2323
seq_length_trial6 - seq_length_trial18 -23.75 13.9 1094 -1.713 0.2008
seq_length_trial12 - seq_length_trial18 -2.58 14.2 1094 -0.182 0.9819
session_f = Block 5:
contrast estimate SE df t.ratio p.value
seq_length_trial6 - seq_length_trial12 -56.44 15.5 1094 -3.638 0.0008
seq_length_trial6 - seq_length_trial18 -24.06 16.5 1094 -1.457 0.3123
seq_length_trial12 - seq_length_trial18 32.38 18.1 1094 1.787 0.1744
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 3 estimates
Block comparisons (RT, correct-only) within each length (Block 4 vs Block 5):
seq_length_trial = 6:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -105 13.1 1094 -8.025 <.0001
seq_length_trial = 12:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -141 15.5 1095 -9.051 <.0001
seq_length_trial = 18:
contrast estimate SE df t.ratio p.value
Block 4 - Block 5 -106 17.3 1095 -6.117 <.0001
Degrees-of-freedom method: kenward-roger
EMMs (RT, correct-only) by block within length:
seq_length_trial = 6:
session_f emmean SE df lower.CL upper.CL
Block 4 447 41.6 18.3 360 535
Block 5 553 41.8 18.6 465 640
seq_length_trial = 12:
session_f emmean SE df lower.CL upper.CL
Block 4 468 41.7 18.5 381 556
Block 5 609 42.5 19.8 520 698
seq_length_trial = 18:
session_f emmean SE df lower.CL upper.CL
Block 4 471 42.0 19.0 383 559
Block 5 577 42.9 20.5 487 666
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 # --- 2.2C: Accuracy per block × length (binomial GLMM, Blocks 4–5) ---
counts_test <- df_acc_base %>%
dplyr::group_by(subject, session, trial) %>%
dplyr::mutate(seq_length_trial = dplyr::n_distinct(sub.trial.number)) %>%
dplyr::ungroup() %>%
dplyr::distinct(subject, session, trial, .keep_all = TRUE) %>%
dplyr::filter(session %in% 4:5) %>%
dplyr::mutate(
session_f = factor(session, levels = c(4, 5), labels = paste("Block", 4:5)),
seq_length_trial = factor(seq_length_trial, levels = c(6, 12, 18))
) %>%
dplyr::group_by(subject, session_f, seq_length_trial) %>%
dplyr::summarise(
total_trials = dplyr::n(), # actual denominator after any exclusions
correct_trials = sum(trial_acc_cat == "correct", na.rm = TRUE),
.groups = "drop"
) %>%
dplyr::mutate(fail_trials = total_trials - correct_trials) %>%
dplyr::filter(!is.na(seq_length_trial))
if (nrow(counts_test) > 0) {
mdl_acc_test <- lme4::glmer(
cbind(correct_trials, fail_trials) ~ session_f * seq_length_trial + (1 | subject),
data = counts_test, family = binomial
)
cat("\n## 2.2C Accuracy model (binomial, block × length)\n")
print(summary(mdl_acc_test))
cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl_acc_test, type = 2))
# EMMs on probability
em_acc_test <- emmeans::emmeans(mdl_acc_test, ~ seq_length_trial | session_f, type = "response")
cat("\nEMMs (prob correct) by length within block:\n"); print(summary(em_acc_test))
cat("\nPairwise (Tukey) across lengths within block (prob correct):\n"); print(pairs(em_acc_test, adjust = "tukey"))
# Cross-block contrasts within each length
em_acc_blk <- emmeans::emmeans(mdl_acc_test, ~ session_f | seq_length_trial, type = "response")
cat("\nBlock comparisons (prob correct) within each length (Block 4 vs Block 5):\n")
print(pairs(em_acc_blk))
cat("\nEMMs (prob correct) by block within length:\n")
print(summary(em_acc_blk))
} else {
cat("\n[Accuracy model skipped: no counts available]\n")
}
## 2.2C Accuracy model (binomial, block × length)
Generalized linear mixed model fit by maximum likelihood (Laplace
Approximation) [glmerMod]
Family: binomial ( logit )
Formula: cbind(correct_trials, fail_trials) ~ session_f * seq_length_trial +
(1 | subject)
Data: counts_test
AIC BIC logLik deviance df.resid
584.0 602.8 -285.0 570.0 101
Scaled residuals:
Min 1Q Median 3Q Max
-4.9843 -0.8720 0.1909 0.9885 4.0483
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 0.3185 0.5644
Number of obs: 108, groups: subject, 18
Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 2.0307 0.2245 9.044 < 2e-16 ***
session_fBlock 5 -0.8171 0.2284 -3.577 0.000347 ***
seq_length_trial12 -0.6935 0.2310 -3.002 0.002682 **
seq_length_trial18 -1.5158 0.2196 -6.903 5.1e-12 ***
session_fBlock 5:seq_length_trial12 -0.6721 0.2973 -2.261 0.023789 *
session_fBlock 5:seq_length_trial18 -0.1854 0.2890 -0.641 0.521275
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Correlation of Fixed Effects:
(Intr) sss_B5 sq__12 sq__18 s_B5:__12
sssn_fBlck5 -0.631
sq_lngth_12 -0.623 0.610
sq_lngth_18 -0.660 0.644 0.636
sss_B5:__12 0.477 -0.766 -0.775 -0.489
sss_B5:__18 0.493 -0.788 -0.481 -0.753 0.607
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: cbind(correct_trials, fail_trials)
Chisq Df Pr(>Chisq)
session_f 98.2015 1 < 2e-16 ***
seq_length_trial 130.4703 2 < 2e-16 ***
session_f:seq_length_trial 5.9554 2 0.05091 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
EMMs (prob correct) by length within block:
session_f = Block 4:
seq_length_trial prob SE df asymp.LCL asymp.UCL
6 0.884 0.0230 Inf 0.831 0.922
12 0.792 0.0326 Inf 0.721 0.849
18 0.626 0.0429 Inf 0.539 0.706
session_f = Block 5:
seq_length_trial prob SE df asymp.LCL asymp.UCL
6 0.771 0.0344 Inf 0.697 0.831
12 0.462 0.0450 Inf 0.376 0.550
18 0.380 0.0431 Inf 0.300 0.468
Confidence level used: 0.95
Intervals are back-transformed from the logit scale
Pairwise (Tukey) across lengths within block (prob correct):
session_f = Block 4:
contrast odds.ratio SE df null z.ratio
seq_length_trial6 / seq_length_trial12 2.00 0.462 Inf 1 3.002
seq_length_trial6 / seq_length_trial18 4.55 1.000 Inf 1 6.903
seq_length_trial12 / seq_length_trial18 2.28 0.438 Inf 1 4.272
p.value
0.0076
<.0001
0.0001
session_f = Block 5:
contrast odds.ratio SE df null z.ratio
seq_length_trial6 / seq_length_trial12 3.92 0.736 Inf 1 7.268
seq_length_trial6 / seq_length_trial18 5.48 1.040 Inf 1 8.947
seq_length_trial12 / seq_length_trial18 1.40 0.245 Inf 1 1.917
p.value
<.0001
<.0001
0.1339
P value adjustment: tukey method for comparing a family of 3 estimates
Tests are performed on the log odds ratio scale
Block comparisons (prob correct) within each length (Block 4 vs Block 5):
seq_length_trial = 6:
contrast odds.ratio SE df null z.ratio p.value
Block 4 / Block 5 2.26 0.517 Inf 1 3.577 0.0003
seq_length_trial = 12:
contrast odds.ratio SE df null z.ratio p.value
Block 4 / Block 5 4.43 0.848 Inf 1 7.786 <.0001
seq_length_trial = 18:
contrast odds.ratio SE df null z.ratio p.value
Block 4 / Block 5 2.72 0.484 Inf 1 5.639 <.0001
Tests are performed on the log odds ratio scale
EMMs (prob correct) by block within length:
seq_length_trial = 6:
session_f prob SE df asymp.LCL asymp.UCL
Block 4 0.884 0.0230 Inf 0.831 0.922
Block 5 0.771 0.0344 Inf 0.697 0.831
seq_length_trial = 12:
session_f prob SE df asymp.LCL asymp.UCL
Block 4 0.792 0.0326 Inf 0.721 0.849
Block 5 0.462 0.0450 Inf 0.376 0.550
seq_length_trial = 18:
session_f prob SE df asymp.LCL asymp.UCL
Block 4 0.626 0.0429 Inf 0.539 0.706
Block 5 0.380 0.0431 Inf 0.300 0.468
Confidence level used: 0.95
Intervals are back-transformed from the logit scale #3 Stepwise RT — Training blocks (1–3)
suppressPackageStartupMessages({
library(dplyr); library(lme4); library(car); library(emmeans)
})
# fallback in case safe_fit wasn't defined earlier
if (!exists("safe_fit")) {
safe_fit <- function(formula, data) {
tryCatch(lmer(formula, data = data),
error = function(e) { stop("Model failed: ", conditionMessage(e)) })
}
}
stepwise_rt_training_with_correctness <- function(df_acc_all) {
cat("\n\n========== Stepwise RT — Training (1–3) ==========\n")
for (blk in 1:3) {
cat(glue::glue("\n--- Block {blk} ---\n"))
# Per-step dataset (do not collapse steps)
df_blk_all <- df_acc_all %>%
filter(session == blk) %>%
group_by(subject, trial, step = as.integer(as.character(sub.trial.number)), trial_acc_cat) %>%
summarise(rt = mean(feedback.RT, na.rm = TRUE), .groups = "drop")
if (nrow(df_blk_all) == 0) { cat("No data\n"); next }
df_blk_all$step <- factor(df_blk_all$step, levels = sort(unique(df_blk_all$step)), ordered = TRUE)
# Combined model with correctness factor
has_two_steps <- n_distinct(df_blk_all$step) >= 2
has_two_corr <- n_distinct(df_blk_all$trial_acc_cat) >= 2
if (has_two_steps && has_two_corr) {
fmla <- rt ~ step * trial_acc_cat + (1 | subject)
} else if (has_two_steps) {
fmla <- rt ~ step + (1 | subject)
} else if (has_two_corr) {
fmla <- rt ~ trial_acc_cat + (1 | subject)
} else {
cat("⚠️ Not enough levels to model.\n"); next
}
cat("\n[Combined] Model summary:\n")
mdl_all <- safe_fit(fmla, data = df_blk_all)
print(summary(mdl_all))
cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl_all, type = 2))
# EMMs depending on formula
if (has_two_steps && has_two_corr) {
em <- emmeans(mdl_all, ~ step | trial_acc_cat)
cat("\nEMMs (step | correctness):\n"); print(summary(em))
# Pairwise (Tukey) across steps | correctness — print CORRECT only
pw_by <- pairs(em, by = "trial_acc_cat", adjust = "tukey")
if ("correct" %in% names(pw_by)) {
cat("\nPairwise (Tukey) across steps | correctness — correct only:\n")
print(pw_by[["correct"]])
}
# Adjacent-step (“consec”) | correctness — print CORRECT only
adj_by <- contrast(em, method = "consec", by = "trial_acc_cat", adjust = "holm")
if ("correct" %in% names(adj_by)) {
cat("\nAdjacent-step contrasts (consecutive) | correctness — correct only:\n")
print(adj_by[["correct"]])
}
} else if (has_two_steps) {
em <- emmeans(mdl_all, ~ step)
print(summary(em))
} else {
em <- emmeans(mdl_all, ~ trial_acc_cat); print(summary(em))
}
# Split outputs
for (lab in c("correct","wrong")) {
df_blk <- df_blk_all %>% filter(trial_acc_cat == lab)
if (nrow(df_blk) == 0) { cat("\n⚠️ No data for", lab, "\n"); next }
has_two_steps_lab <- n_distinct(df_blk$step) >= 2
if (!has_two_steps_lab) { cat("\n⚠️ Only one step level for", lab, "— skipping.\n"); next }
cat("\n[", toupper(lab), "] Model summary:\n", sep = "")
mdl <- safe_fit(rt ~ step + (1 | subject), data = df_blk)
print(summary(mdl))
cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl, type = 2))
em_lab <- emmeans(mdl, ~ step); print(summary(em_lab))
# >>> Pairwise comparisons ONLY for CORRECT trials <<<
if (lab == "correct") {
cat("\nPairwise (Tukey) across steps — correct:\n")
print(pairs(em_lab, adjust = "tukey"))
cat("\nAdjacent-step contrasts (consecutive) — correct:\n")
print(contrast(em_lab, "consec", adjust = "holm"))
}
}
}
}
stepwise_rt_training_with_correctness(df_acc_base)
========== Stepwise RT — Training (1–3) ==========
--- Block 1 ---
[Combined] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step * trial_acc_cat + (1 | subject)
Data: data
REML criterion at convergence: 83308.2
Scaled residuals:
Min 1Q Median 3Q Max
-3.053 -0.228 -0.068 0.118 33.226
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 72318 268.9
Residual 564412 751.3
Number of obs: 5184, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 495.49 64.44 7.689
step.L -186.72 28.33 -6.590
step.Q 170.62 28.33 6.022
step.C -120.67 28.33 -4.259
step^4 43.37 28.33 1.531
step^5 -12.27 28.33 -0.433
trial_acc_catwrong 378.53 27.67 13.680
step.L:trial_acc_catwrong 701.50 65.64 10.687
step.Q:trial_acc_catwrong 690.47 65.64 10.519
step.C:trial_acc_catwrong 349.61 65.64 5.326
step^4:trial_acc_catwrong 186.16 65.64 2.836
step^5:trial_acc_catwrong 37.65 65.64 0.574
Correlation of Fixed Effects:
(Intr) step.L step.Q step.C step^4 step^5 trl_c_ s.L:__ s.Q:__
step.L 0.000
step.Q 0.000 0.000
step.C 0.000 0.000 0.000
step^4 0.000 0.000 0.000 0.000
step^5 0.000 0.000 0.000 0.000 0.000
trl_cc_ctwr -0.080 0.000 0.000 0.000 0.000 0.000
stp.L:trl__ 0.000 -0.432 0.000 0.000 0.000 0.000 0.000
stp.Q:trl__ 0.000 0.000 -0.432 0.000 0.000 0.000 0.000 0.000
stp.C:trl__ 0.000 0.000 0.000 -0.432 0.000 0.000 0.000 0.000 0.000
stp^4:trl__ 0.000 0.000 0.000 0.000 -0.432 0.000 0.000 0.000 0.000
stp^5:trl__ 0.000 0.000 0.000 0.000 0.000 -0.432 0.000 0.000 0.000
s.C:__ s^4:__
step.L
step.Q
step.C
step^4
step^5
trl_cc_ctwr
stp.L:trl__
stp.Q:trl__
stp.C:trl__
stp^4:trl__ 0.000
stp^5:trl__ 0.000 0.000
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 156.01 5 < 2.2e-16 ***
trial_acc_cat 187.15 1 < 2.2e-16 ***
step:trial_acc_cat 261.61 5 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
EMMs (step | correctness):
trial_acc_cat = correct:
step emmean SE df lower.CL upper.CL
1 754 69.4 23.2 611 898
2 452 69.4 23.2 309 596
3 431 69.4 23.2 288 575
4 443 69.4 23.2 300 587
5 452 69.4 23.2 309 596
6 439 69.4 23.2 296 583
trial_acc_cat = wrong:
step emmean SE df lower.CL upper.CL
1 993 86.9 56.8 819 1167
2 593 86.9 56.8 419 767
3 576 86.9 56.8 402 750
4 594 86.9 56.8 420 768
5 707 86.9 56.8 533 881
6 1782 86.9 56.8 1608 1956
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
[CORRECT] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step + (1 | subject)
Data: data
REML criterion at convergence: 60818.2
Scaled residuals:
Min 1Q Median 3Q Max
-4.4444 -0.3995 -0.1412 0.2213 29.3410
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 66759 258.4
Residual 106031 325.6
Number of obs: 4218, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 498.01 61.11 8.149
step.L -186.72 12.28 -15.204
step.Q 170.62 12.28 13.893
step.C -120.67 12.28 -9.826
step^4 43.37 12.28 3.532
step^5 -12.27 12.28 -0.999
Correlation of Fixed Effects:
(Intr) step.L step.Q step.C step^4
step.L 0.000
step.Q 0.000 0.000
step.C 0.000 0.000 0.000
step^4 0.000 0.000 0.000 0.000
step^5 0.000 0.000 0.000 0.000 0.000
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 534.2 5 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
step emmean SE df lower.CL upper.CL
1 757 62.1 18.2 626 887
2 455 62.1 18.2 324 585
3 434 62.1 18.2 304 564
4 446 62.1 18.2 315 576
5 455 62.1 18.2 324 585
6 442 62.1 18.2 312 572
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
Pairwise (Tukey) across steps — correct:
contrast estimate SE df t.ratio p.value
step1 - step2 301.691 17.4 4195 17.370 <.0001
step1 - step3 322.613 17.4 4195 18.575 <.0001
step1 - step4 310.754 17.4 4195 17.892 <.0001
step1 - step5 301.947 17.4 4195 17.385 <.0001
step1 - step6 314.669 17.4 4195 18.118 <.0001
step2 - step3 20.922 17.4 4195 1.205 0.8348
step2 - step4 9.063 17.4 4195 0.522 0.9953
step2 - step5 0.256 17.4 4195 0.015 1.0000
step2 - step6 12.977 17.4 4195 0.747 0.9759
step3 - step4 -11.859 17.4 4195 -0.683 0.9839
step3 - step5 -20.666 17.4 4195 -1.190 0.8419
step3 - step6 -7.944 17.4 4195 -0.457 0.9975
step4 - step5 -8.806 17.4 4195 -0.507 0.9959
step4 - step6 3.915 17.4 4195 0.225 0.9999
step5 - step6 12.721 17.4 4195 0.732 0.9780
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 6 estimates
Adjacent-step contrasts (consecutive) — correct:
contrast estimate SE df t.ratio p.value
step2 - step1 -301.69 17.4 4195 -17.370 <.0001
step3 - step2 -20.92 17.4 4195 -1.205 0.9137
step4 - step3 11.86 17.4 4195 0.683 1.0000
step5 - step4 8.81 17.4 4195 0.507 1.0000
step6 - step5 -12.72 17.4 4195 -0.732 1.0000
Degrees-of-freedom method: kenward-roger
P value adjustment: holm method for 5 tests
[WRONG] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step + (1 | subject)
Data: data
REML criterion at convergence: 16938.6
Scaled residuals:
Min 1Q Median 3Q Max
-1.4700 -0.2600 -0.0791 0.0668 15.3544
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 129325 359.6
Residual 2548045 1596.3
Number of obs: 966, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 885.50 101.42 8.731
step.L 514.77 125.80 4.092
step.Q 861.09 125.80 6.845
step.C 228.93 125.80 1.820
step^4 229.54 125.80 1.825
step^5 25.38 125.80 0.202
Correlation of Fixed Effects:
(Intr) step.L step.Q step.C step^4
step.L 0.000
step.Q 0.000 0.000
step.C 0.000 0.000 0.000
step^4 0.000 0.000 0.000 0.000
step^5 0.000 0.000 0.000 0.000 0.000
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 70.275 5 8.981e-14 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
step emmean SE df lower.CL upper.CL
1 1004 153 84.2 699 1309
2 604 153 84.2 299 909
3 587 153 84.2 282 892
4 606 153 84.2 301 911
5 719 153 84.2 414 1023
6 1793 153 84.2 1488 2098
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
--- Block 2 ---
[Combined] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step * trial_acc_cat + (1 | subject)
Data: data
REML criterion at convergence: 152133.8
Scaled residuals:
Min 1Q Median 3Q Max
-3.5910 -0.4328 -0.1421 0.1825 13.1480
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 73186 270.5
Residual 139488 373.5
Number of obs: 10368, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 491.831 63.934 7.693
step.L -134.063 15.782 -8.494
step.Q 78.716 15.782 4.988
step.C -183.800 15.782 -11.646
step^4 130.544 15.782 8.271
step^5 -59.893 15.782 -3.795
step^6 31.494 15.782 1.996
step^7 14.042 15.782 0.890
step^8 2.778 15.782 0.176
step^9 -32.254 15.782 -2.044
step^10 22.473 15.782 1.424
step^11 -44.771 15.782 -2.837
trial_acc_catwrong 189.621 8.128 23.329
step.L:trial_acc_catwrong 214.325 26.607 8.055
step.Q:trial_acc_catwrong 155.288 26.607 5.836
step.C:trial_acc_catwrong 51.468 26.607 1.934
step^4:trial_acc_catwrong 157.970 26.607 5.937
step^5:trial_acc_catwrong 24.774 26.607 0.931
step^6:trial_acc_catwrong 40.459 26.607 1.521
step^7:trial_acc_catwrong -20.228 26.607 -0.760
step^8:trial_acc_catwrong -18.810 26.607 -0.707
step^9:trial_acc_catwrong 66.831 26.607 2.512
step^10:trial_acc_catwrong 47.716 26.607 1.793
step^11:trial_acc_catwrong -60.857 26.607 -2.287
Correlation matrix not shown by default, as p = 24 > 12.
Use print(summary(mdl_all), correlation=TRUE) or
vcov(summary(mdl_all)) if you need it
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 582.81 11 < 2.2e-16 ***
trial_acc_cat 544.25 1 < 2.2e-16 ***
step:trial_acc_cat 156.96 11 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
EMMs (step | correctness):
trial_acc_cat = correct:
step emmean SE df lower.CL upper.CL
1 746 65.7 19.0 609 884
2 470 65.7 19.0 333 608
3 439 65.7 19.0 302 577
4 414 65.7 19.0 277 552
5 517 65.7 19.0 380 655
6 480 65.7 19.0 343 618
7 507 65.7 19.0 370 645
8 484 65.7 19.0 347 622
9 511 65.7 19.0 373 648
10 471 65.7 19.0 334 609
11 438 65.7 19.0 300 575
12 424 65.7 19.0 286 561
trial_acc_cat = wrong:
step emmean SE df lower.CL upper.CL
1 950 67.3 20.9 810 1090
2 564 67.3 20.9 424 704
3 524 67.3 20.9 384 664
4 547 67.3 20.9 407 687
5 703 67.3 20.9 563 843
6 573 67.3 20.9 433 713
7 737 67.3 20.9 597 877
8 701 67.3 20.9 561 841
9 664 67.3 20.9 524 804
10 666 67.3 20.9 526 806
11 666 67.3 20.9 526 806
12 884 67.3 20.9 744 1024
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
[CORRECT] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step + (1 | subject)
Data: data
REML criterion at convergence: 92865.7
Scaled residuals:
Min 1Q Median 3Q Max
-4.1593 -0.4831 -0.1817 0.2355 13.2133
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 60609 246.2
Residual 58656 242.2
Number of obs: 6720, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 489.472 58.108 8.423
step.L -134.063 10.234 -13.099
step.Q 78.716 10.234 7.691
step.C -183.800 10.234 -17.959
step^4 130.544 10.234 12.755
step^5 -59.893 10.234 -5.852
step^6 31.494 10.234 3.077
step^7 14.042 10.234 1.372
step^8 2.778 10.234 0.271
step^9 -32.254 10.234 -3.152
step^10 22.473 10.234 2.196
step^11 -44.771 10.234 -4.375
Correlation of Fixed Effects:
(Intr) step.L step.Q step.C step^4 step^5 step^6 step^7 step^8 step^9
step.L 0.000
step.Q 0.000 0.000
step.C 0.000 0.000 0.000
step^4 0.000 0.000 0.000 0.000
step^5 0.000 0.000 0.000 0.000 0.000
step^6 0.000 0.000 0.000 0.000 0.000 0.000
step^7 0.000 0.000 0.000 0.000 0.000 0.000 0.000
step^8 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
step^9 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
step^10 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
step^11 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
stp^10
step.L
step.Q
step.C
step^4
step^5
step^6
step^7
step^8
step^9
step^10
step^11 0.000
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 795.54 11 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
step emmean SE df lower.CL upper.CL
1 744 58.9 18 620 868
2 468 58.9 18 344 592
3 437 58.9 18 313 561
4 412 58.9 18 288 536
5 515 58.9 18 391 639
6 478 58.9 18 354 602
7 505 58.9 18 381 629
8 482 58.9 18 358 606
9 508 58.9 18 385 632
10 469 58.9 18 345 593
11 435 58.9 18 312 559
12 421 58.9 18 298 545
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
Pairwise (Tukey) across steps — correct:
contrast estimate SE df t.ratio p.value
step1 - step2 275.707 14.5 6691 19.049 <.0001
step1 - step3 306.861 14.5 6691 21.201 <.0001
step1 - step4 331.834 14.5 6691 22.927 <.0001
step1 - step5 228.984 14.5 6691 15.821 <.0001
step1 - step6 265.964 14.5 6691 18.376 <.0001
step1 - step7 238.716 14.5 6691 16.493 <.0001
step1 - step8 261.957 14.5 6691 18.099 <.0001
step1 - step9 235.332 14.5 6691 16.259 <.0001
step1 - step10 274.986 14.5 6691 18.999 <.0001
step1 - step11 308.355 14.5 6691 21.305 <.0001
step1 - step12 322.405 14.5 6691 22.275 <.0001
step2 - step3 31.154 14.5 6691 2.152 0.5846
step2 - step4 56.127 14.5 6691 3.878 0.0060
step2 - step5 -46.723 14.5 6691 -3.228 0.0567
step2 - step6 -9.743 14.5 6691 -0.673 0.9999
step2 - step7 -36.991 14.5 6691 -2.556 0.3055
step2 - step8 -13.750 14.5 6691 -0.950 0.9986
step2 - step9 -40.375 14.5 6691 -2.790 0.1845
step2 - step10 -0.721 14.5 6691 -0.050 1.0000
step2 - step11 32.648 14.5 6691 2.256 0.5088
step2 - step12 46.698 14.5 6691 3.226 0.0570
step3 - step4 24.973 14.5 6691 1.725 0.8567
step3 - step5 -77.877 14.5 6691 -5.381 <.0001
step3 - step6 -40.896 14.5 6691 -2.826 0.1693
step3 - step7 -68.145 14.5 6691 -4.708 0.0002
step3 - step8 -44.904 14.5 6691 -3.102 0.0819
step3 - step9 -71.529 14.5 6691 -4.942 0.0001
step3 - step10 -31.875 14.5 6691 -2.202 0.5480
step3 - step11 1.495 14.5 6691 0.103 1.0000
step3 - step12 15.545 14.5 6691 1.074 0.9957
step4 - step5 -102.850 14.5 6691 -7.106 <.0001
step4 - step6 -65.870 14.5 6691 -4.551 0.0003
step4 - step7 -93.118 14.5 6691 -6.434 <.0001
step4 - step8 -69.877 14.5 6691 -4.828 0.0001
step4 - step9 -96.502 14.5 6691 -6.667 <.0001
step4 - step10 -56.848 14.5 6691 -3.928 0.0049
step4 - step11 -23.479 14.5 6691 -1.622 0.9013
step4 - step12 -9.429 14.5 6691 -0.651 1.0000
step5 - step6 36.980 14.5 6691 2.555 0.3059
step5 - step7 9.732 14.5 6691 0.672 0.9999
step5 - step8 32.973 14.5 6691 2.278 0.4924
step5 - step9 6.348 14.5 6691 0.439 1.0000
step5 - step10 46.002 14.5 6691 3.178 0.0657
step5 - step11 79.371 14.5 6691 5.484 <.0001
step5 - step12 93.421 14.5 6691 6.455 <.0001
step6 - step7 -27.248 14.5 6691 -1.883 0.7703
step6 - step8 -4.007 14.5 6691 -0.277 1.0000
step6 - step9 -30.632 14.5 6691 -2.116 0.6109
step6 - step10 9.021 14.5 6691 0.623 1.0000
step6 - step11 42.391 14.5 6691 2.929 0.1309
step6 - step12 56.441 14.5 6691 3.900 0.0055
step7 - step8 23.241 14.5 6691 1.606 0.9074
step7 - step9 -3.384 14.5 6691 -0.234 1.0000
step7 - step10 36.270 14.5 6691 2.506 0.3361
step7 - step11 69.639 14.5 6691 4.811 0.0001
step7 - step12 83.689 14.5 6691 5.782 <.0001
step8 - step9 -26.625 14.5 6691 -1.840 0.7961
step8 - step10 13.029 14.5 6691 0.900 0.9991
step8 - step11 46.398 14.5 6691 3.206 0.0606
step8 - step12 60.448 14.5 6691 4.176 0.0018
step9 - step10 39.654 14.5 6691 2.740 0.2070
step9 - step11 73.023 14.5 6691 5.045 <.0001
step9 - step12 87.073 14.5 6691 6.016 <.0001
step10 - step11 33.370 14.5 6691 2.306 0.4726
step10 - step12 47.420 14.5 6691 3.276 0.0489
step11 - step12 14.050 14.5 6691 0.971 0.9983
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 12 estimates
Adjacent-step contrasts (consecutive) — correct:
contrast estimate SE df t.ratio p.value
step2 - step1 -275.7 14.5 6691 -19.049 <.0001
step3 - step2 -31.2 14.5 6691 -2.152 0.1884
step4 - step3 -25.0 14.5 6691 -1.725 0.2990
step5 - step4 102.8 14.5 6691 7.106 <.0001
step6 - step5 -37.0 14.5 6691 -2.555 0.0851
step7 - step6 27.2 14.5 6691 1.883 0.2990
step8 - step7 -23.2 14.5 6691 -1.606 0.2990
step9 - step8 26.6 14.5 6691 1.840 0.2990
step10 - step9 -39.7 14.5 6691 -2.740 0.0555
step11 - step10 -33.4 14.5 6691 -2.306 0.1482
step12 - step11 -14.1 14.5 6691 -0.971 0.3317
Degrees-of-freedom method: kenward-roger
P value adjustment: holm method for 11 tests
[WRONG] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step + (1 | subject)
Data: data
REML criterion at convergence: 56105.2
Scaled residuals:
Min 1Q Median 3Q Max
-2.6066 -0.4986 -0.1921 0.1556 9.3032
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 90484 300.8
Residual 283253 532.2
Number of obs: 3648, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 683.972 71.647 9.546
step.L 80.262 30.525 2.629
step.Q 234.004 30.525 7.666
step.C -132.333 30.525 -4.335
step^4 288.514 30.525 9.452
step^5 -35.120 30.525 -1.151
step^6 71.953 30.525 2.357
step^7 -6.187 30.525 -0.203
step^8 -16.032 30.525 -0.525
step^9 34.576 30.525 1.133
step^10 70.189 30.525 2.299
step^11 -105.628 30.525 -3.460
Correlation of Fixed Effects:
(Intr) step.L step.Q step.C step^4 step^5 step^6 step^7 step^8 step^9
step.L 0.000
step.Q 0.000 0.000
step.C 0.000 0.000 0.000
step^4 0.000 0.000 0.000 0.000
step^5 0.000 0.000 0.000 0.000 0.000
step^6 0.000 0.000 0.000 0.000 0.000 0.000
step^7 0.000 0.000 0.000 0.000 0.000 0.000 0.000
step^8 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
step^9 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
step^10 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
step^11 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
stp^10
step.L
step.Q
step.C
step^4
step^5
step^6
step^7
step^8
step^9
step^10
step^11 0.000
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 199.56 11 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
step emmean SE df lower.CL upper.CL
1 953 77.4 23.1 793 1113
2 566 77.4 23.1 406 726
3 527 77.4 23.1 367 687
4 550 77.4 23.1 390 710
5 705 77.4 23.1 545 865
6 576 77.4 23.1 416 736
7 739 77.4 23.1 579 899
8 703 77.4 23.1 543 863
9 666 77.4 23.1 506 826
10 668 77.4 23.1 508 828
11 668 77.4 23.1 508 828
12 886 77.4 23.1 726 1046
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
--- Block 3 ---
[Combined] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step * trial_acc_cat + (1 | subject)
Data: data
REML criterion at convergence: 233065.3
Scaled residuals:
Min 1Q Median 3Q Max
-2.6317 -0.4232 -0.1656 0.1450 23.1355
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 45858 214.1
Residual 191532 437.6
Number of obs: 15552, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 534.1761 50.7214 10.532
step.L -74.1687 20.1442 -3.682
step.Q 52.8735 20.1442 2.625
step.C -212.6525 20.1442 -10.557
step^4 148.2090 20.1442 7.357
step^5 -72.2957 20.1442 -3.589
step^6 143.7340 20.1442 7.135
step^7 -39.4603 20.1442 -1.959
step^8 -30.7602 20.1442 -1.527
step^9 -68.5111 20.1442 -3.401
step^10 -65.2846 20.1442 -3.241
step^11 5.9738 20.1442 0.297
step^12 97.7373 20.1442 4.852
step^13 -3.5160 20.1442 -0.175
step^14 -61.4153 20.1442 -3.049
step^15 -91.9082 20.1442 -4.563
step^16 -45.7208 20.1442 -2.270
step^17 -24.5454 20.1442 -1.218
trial_acc_catwrong 101.2860 7.8532 12.897
step.L:trial_acc_catwrong 151.2608 29.9064 5.058
step.Q:trial_acc_catwrong 71.8196 29.9064 2.401
step.C:trial_acc_catwrong 6.5254 29.9064 0.218
step^4:trial_acc_catwrong 161.4595 29.9064 5.399
step^5:trial_acc_catwrong -53.6791 29.9064 -1.795
step^6:trial_acc_catwrong 69.2605 29.9064 2.316
step^7:trial_acc_catwrong 4.7810 29.9064 0.160
step^8:trial_acc_catwrong 116.4322 29.9064 3.893
step^9:trial_acc_catwrong 63.9203 29.9064 2.137
step^10:trial_acc_catwrong 29.9286 29.9064 1.001
step^11:trial_acc_catwrong -27.3549 29.9064 -0.915
step^12:trial_acc_catwrong -18.6902 29.9064 -0.625
step^13:trial_acc_catwrong -25.3792 29.9064 -0.849
step^14:trial_acc_catwrong -1.3827 29.9064 -0.046
step^15:trial_acc_catwrong 66.9411 29.9064 2.238
step^16:trial_acc_catwrong 11.6012 29.9064 0.388
step^17:trial_acc_catwrong 0.2499 29.9064 0.008
Correlation matrix not shown by default, as p = 36 > 12.
Use print(summary(mdl_all), correlation=TRUE) or
vcov(summary(mdl_all)) if you need it
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 742.457 17 < 2.2e-16 ***
trial_acc_cat 166.344 1 < 2.2e-16 ***
step:trial_acc_cat 96.992 17 3.201e-13 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
EMMs (step | correctness):
trial_acc_cat = correct:
step emmean SE df lower.CL upper.CL
1 809 54.4 22.7 696 922
2 501 54.4 22.7 388 613
3 463 54.4 22.7 350 575
4 451 54.4 22.7 339 564
5 529 54.4 22.7 416 641
6 525 54.4 22.7 413 638
7 505 54.4 22.7 392 618
8 493 54.4 22.7 381 606
9 544 54.4 22.7 432 657
10 570 54.4 22.7 457 682
11 514 54.4 22.7 401 626
12 480 54.4 22.7 367 592
13 707 54.4 22.7 595 820
14 565 54.4 22.7 453 678
15 513 54.4 22.7 400 626
16 449 54.4 22.7 337 562
17 508 54.4 22.7 396 621
18 489 54.4 22.7 376 602
trial_acc_cat = wrong:
step emmean SE df lower.CL upper.CL
1 992 55.1 24.0 878 1106
2 486 55.1 24.0 372 600
3 474 55.1 24.0 360 587
4 488 55.1 24.0 374 601
5 574 55.1 24.0 460 688
6 625 55.1 24.0 511 739
7 568 55.1 24.0 454 682
8 588 55.1 24.0 474 702
9 686 55.1 24.0 572 799
10 682 55.1 24.0 568 796
11 642 55.1 24.0 528 756
12 610 55.1 24.0 496 723
13 727 55.1 24.0 613 841
14 655 55.1 24.0 542 769
15 677 55.1 24.0 563 791
16 607 55.1 24.0 494 721
17 590 55.1 24.0 476 704
18 769 55.1 24.0 655 883
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
[CORRECT] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step + (1 | subject)
Data: data
REML criterion at convergence: 121086.8
Scaled residuals:
Min 1Q Median 3Q Max
-2.579 -0.432 -0.161 0.192 33.518
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 40081 200.2
Residual 91189 302.0
Number of obs: 8496, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 530.575 47.335 11.209
step.L -74.169 13.900 -5.336
step.Q 52.873 13.900 3.804
step.C -212.652 13.900 -15.299
step^4 148.209 13.900 10.663
step^5 -72.296 13.900 -5.201
step^6 143.734 13.900 10.341
step^7 -39.460 13.900 -2.839
step^8 -30.760 13.900 -2.213
step^9 -68.511 13.900 -4.929
step^10 -65.285 13.900 -4.697
step^11 5.974 13.900 0.430
step^12 97.737 13.900 7.032
step^13 -3.516 13.900 -0.253
step^14 -61.415 13.900 -4.419
step^15 -91.908 13.900 -6.612
step^16 -45.721 13.900 -3.289
step^17 -24.545 13.900 -1.766
Correlation matrix not shown by default, as p = 18 > 12.
Use print(summary(mdl), correlation=TRUE) or
vcov(summary(mdl)) if you need it
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 710.88 17 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
step emmean SE df lower.CL upper.CL
1 805 49.2 19.9 703 908
2 497 49.2 19.9 395 600
3 459 49.2 19.9 356 562
4 448 49.2 19.9 345 551
5 525 49.2 19.9 422 628
6 522 49.2 19.9 419 624
7 501 49.2 19.9 399 604
8 490 49.2 19.9 387 592
9 541 49.2 19.9 438 643
10 566 49.2 19.9 463 669
11 510 49.2 19.9 407 613
12 476 49.2 19.9 373 579
13 704 49.2 19.9 601 807
14 562 49.2 19.9 459 664
15 509 49.2 19.9 407 612
16 446 49.2 19.9 343 548
17 505 49.2 19.9 402 607
18 485 49.2 19.9 383 588
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
Pairwise (Tukey) across steps — correct:
contrast estimate SE df t.ratio p.value
step1 - step2 308.066 19.7 8461 15.672 <.0001
step1 - step3 346.254 19.7 8461 17.615 <.0001
step1 - step4 357.566 19.7 8461 18.190 <.0001
step1 - step5 280.229 19.7 8461 14.256 <.0001
step1 - step6 283.729 19.7 8461 14.434 <.0001
step1 - step7 303.913 19.7 8461 15.461 <.0001
step1 - step8 315.773 19.7 8461 16.064 <.0001
step1 - step9 264.585 19.7 8461 13.460 <.0001
step1 - step10 239.358 19.7 8461 12.177 <.0001
step1 - step11 295.263 19.7 8461 15.021 <.0001
step1 - step12 329.360 19.7 8461 16.755 <.0001
step1 - step13 101.496 19.7 8461 5.163 <.0001
step1 - step14 243.826 19.7 8461 12.404 <.0001
step1 - step15 296.002 19.7 8461 15.058 <.0001
step1 - step16 359.763 19.7 8461 18.302 <.0001
step1 - step17 300.820 19.7 8461 15.303 <.0001
step1 - step18 319.892 19.7 8461 16.274 <.0001
step2 - step3 38.189 19.7 8461 1.943 0.8954
step2 - step4 49.500 19.7 8461 2.518 0.5162
step2 - step5 -27.837 19.7 8461 -1.416 0.9951
step2 - step6 -24.337 19.7 8461 -1.238 0.9990
step2 - step7 -4.152 19.7 8461 -0.211 1.0000
step2 - step8 7.708 19.7 8461 0.392 1.0000
step2 - step9 -43.481 19.7 8461 -2.212 0.7459
step2 - step10 -68.708 19.7 8461 -3.495 0.0491
step2 - step11 -12.803 19.7 8461 -0.651 1.0000
step2 - step12 21.294 19.7 8461 1.083 0.9998
step2 - step13 -206.570 19.7 8461 -10.509 <.0001
step2 - step14 -64.239 19.7 8461 -3.268 0.0984
step2 - step15 -12.064 19.7 8461 -0.614 1.0000
step2 - step16 51.697 19.7 8461 2.630 0.4311
step2 - step17 -7.246 19.7 8461 -0.369 1.0000
step2 - step18 11.826 19.7 8461 0.602 1.0000
step3 - step4 11.311 19.7 8461 0.575 1.0000
step3 - step5 -66.025 19.7 8461 -3.359 0.0752
step3 - step6 -62.525 19.7 8461 -3.181 0.1256
step3 - step7 -42.341 19.7 8461 -2.154 0.7839
step3 - step8 -30.481 19.7 8461 -1.551 0.9868
step3 - step9 -81.669 19.7 8461 -4.155 0.0043
step3 - step10 -106.896 19.7 8461 -5.438 <.0001
step3 - step11 -50.992 19.7 8461 -2.594 0.4580
step3 - step12 -16.894 19.7 8461 -0.859 1.0000
step3 - step13 -244.758 19.7 8461 -12.451 <.0001
step3 - step14 -102.428 19.7 8461 -5.211 <.0001
step3 - step15 -50.252 19.7 8461 -2.556 0.4867
step3 - step16 13.508 19.7 8461 0.687 1.0000
step3 - step17 -45.434 19.7 8461 -2.311 0.6751
step3 - step18 -26.362 19.7 8461 -1.341 0.9974
step4 - step5 -77.337 19.7 8461 -3.934 0.0103
step4 - step6 -73.837 19.7 8461 -3.756 0.0201
step4 - step7 -53.653 19.7 8461 -2.729 0.3598
step4 - step8 -41.792 19.7 8461 -2.126 0.8012
step4 - step9 -92.981 19.7 8461 -4.730 0.0003
step4 - step10 -118.208 19.7 8461 -6.014 <.0001
step4 - step11 -62.303 19.7 8461 -3.170 0.1295
step4 - step12 -28.206 19.7 8461 -1.435 0.9944
step4 - step13 -256.070 19.7 8461 -13.027 <.0001
step4 - step14 -113.739 19.7 8461 -5.786 <.0001
step4 - step15 -61.564 19.7 8461 -3.132 0.1433
step4 - step16 2.197 19.7 8461 0.112 1.0000
step4 - step17 -56.746 19.7 8461 -2.887 0.2600
step4 - step18 -37.674 19.7 8461 -1.917 0.9059
step5 - step6 3.500 19.7 8461 0.178 1.0000
step5 - step7 23.684 19.7 8461 1.205 0.9993
step5 - step8 35.544 19.7 8461 1.808 0.9420
step5 - step9 -15.644 19.7 8461 -0.796 1.0000
step5 - step10 -40.871 19.7 8461 -2.079 0.8285
step5 - step11 15.034 19.7 8461 0.765 1.0000
step5 - step12 49.131 19.7 8461 2.499 0.5307
step5 - step13 -178.733 19.7 8461 -9.093 <.0001
step5 - step14 -36.403 19.7 8461 -1.852 0.9289
step5 - step15 15.773 19.7 8461 0.802 1.0000
step5 - step16 79.534 19.7 8461 4.046 0.0067
step5 - step17 20.591 19.7 8461 1.048 0.9999
step5 - step18 39.663 19.7 8461 2.018 0.8611
step6 - step7 20.184 19.7 8461 1.027 0.9999
step6 - step8 32.044 19.7 8461 1.630 0.9780
step6 - step9 -19.144 19.7 8461 -0.974 1.0000
step6 - step10 -44.371 19.7 8461 -2.257 0.7144
step6 - step11 11.534 19.7 8461 0.587 1.0000
step6 - step12 45.631 19.7 8461 2.321 0.6677
step6 - step13 -182.233 19.7 8461 -9.271 <.0001
step6 - step14 -39.903 19.7 8461 -2.030 0.8549
step6 - step15 12.273 19.7 8461 0.624 1.0000
step6 - step16 76.034 19.7 8461 3.868 0.0133
step6 - step17 17.091 19.7 8461 0.869 1.0000
step6 - step18 36.163 19.7 8461 1.840 0.9328
step7 - step8 11.860 19.7 8461 0.603 1.0000
step7 - step9 -39.328 19.7 8461 -2.001 0.8694
step7 - step10 -64.555 19.7 8461 -3.284 0.0939
step7 - step11 -8.650 19.7 8461 -0.440 1.0000
step7 - step12 25.447 19.7 8461 1.295 0.9983
step7 - step13 -202.417 19.7 8461 -10.297 <.0001
step7 - step14 -60.087 19.7 8461 -3.057 0.1739
step7 - step15 -7.911 19.7 8461 -0.402 1.0000
step7 - step16 55.850 19.7 8461 2.841 0.2871
step7 - step17 -3.093 19.7 8461 -0.157 1.0000
step7 - step18 15.979 19.7 8461 0.813 1.0000
step8 - step9 -51.189 19.7 8461 -2.604 0.4504
step8 - step10 -76.415 19.7 8461 -3.887 0.0124
step8 - step11 -20.511 19.7 8461 -1.043 0.9999
step8 - step12 13.587 19.7 8461 0.691 1.0000
step8 - step13 -214.278 19.7 8461 -10.901 <.0001
step8 - step14 -71.947 19.7 8461 -3.660 0.0283
step8 - step15 -19.771 19.7 8461 -1.006 0.9999
step8 - step16 43.989 19.7 8461 2.238 0.7281
step8 - step17 -14.953 19.7 8461 -0.761 1.0000
step8 - step18 4.119 19.7 8461 0.210 1.0000
step9 - step10 -25.227 19.7 8461 -1.283 0.9985
step9 - step11 30.678 19.7 8461 1.561 0.9859
step9 - step12 64.775 19.7 8461 3.295 0.0909
step9 - step13 -163.089 19.7 8461 -8.297 <.0001
step9 - step14 -20.759 19.7 8461 -1.056 0.9999
step9 - step15 31.417 19.7 8461 1.598 0.9820
step9 - step16 95.178 19.7 8461 4.842 0.0002
step9 - step17 36.235 19.7 8461 1.843 0.9316
step9 - step18 55.307 19.7 8461 2.814 0.3042
step10 - step11 55.905 19.7 8461 2.844 0.2854
step10 - step12 90.002 19.7 8461 4.579 0.0007
step10 - step13 -137.862 19.7 8461 -7.013 <.0001
step10 - step14 4.468 19.7 8461 0.227 1.0000
step10 - step15 56.644 19.7 8461 2.882 0.2630
step10 - step16 120.405 19.7 8461 6.125 <.0001
step10 - step17 61.462 19.7 8461 3.127 0.1452
step10 - step18 80.534 19.7 8461 4.097 0.0054
step11 - step12 34.097 19.7 8461 1.735 0.9601
step11 - step13 -193.767 19.7 8461 -9.857 <.0001
step11 - step14 -51.436 19.7 8461 -2.617 0.4410
step11 - step15 0.739 19.7 8461 0.038 1.0000
step11 - step16 64.500 19.7 8461 3.281 0.0947
step11 - step17 5.557 19.7 8461 0.283 1.0000
step11 - step18 24.629 19.7 8461 1.253 0.9989
step12 - step13 -227.864 19.7 8461 -11.592 <.0001
step12 - step14 -85.534 19.7 8461 -4.351 0.0019
step12 - step15 -33.358 19.7 8461 -1.697 0.9675
step12 - step16 30.402 19.7 8461 1.547 0.9872
step12 - step17 -28.540 19.7 8461 -1.452 0.9936
step12 - step18 -9.468 19.7 8461 -0.482 1.0000
step13 - step14 142.331 19.7 8461 7.241 <.0001
step13 - step15 194.506 19.7 8461 9.895 <.0001
step13 - step16 258.267 19.7 8461 13.139 <.0001
step13 - step17 199.324 19.7 8461 10.140 <.0001
step13 - step18 218.396 19.7 8461 11.110 <.0001
step14 - step15 52.176 19.7 8461 2.654 0.4132
step14 - step16 115.936 19.7 8461 5.898 <.0001
step14 - step17 56.994 19.7 8461 2.899 0.2529
step14 - step18 76.066 19.7 8461 3.870 0.0132
step15 - step16 63.761 19.7 8461 3.244 0.1054
step15 - step17 4.818 19.7 8461 0.245 1.0000
step15 - step18 23.890 19.7 8461 1.215 0.9992
step16 - step17 -58.943 19.7 8461 -2.999 0.2008
step16 - step18 -39.871 19.7 8461 -2.028 0.8558
step17 - step18 19.072 19.7 8461 0.970 1.0000
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 18 estimates
Adjacent-step contrasts (consecutive) — correct:
contrast estimate SE df t.ratio p.value
step2 - step1 -308.1 19.7 8461 -15.672 <.0001
step3 - step2 -38.2 19.7 8461 -1.943 0.4166
step4 - step3 -11.3 19.7 8461 -0.575 1.0000
step5 - step4 77.3 19.7 8461 3.934 0.0012
step6 - step5 -3.5 19.7 8461 -0.178 1.0000
step7 - step6 -20.2 19.7 8461 -1.027 1.0000
step8 - step7 -11.9 19.7 8461 -0.603 1.0000
step9 - step8 51.2 19.7 8461 2.604 0.0831
step10 - step9 25.2 19.7 8461 1.283 1.0000
step11 - step10 -55.9 19.7 8461 -2.844 0.0491
step12 - step11 -34.1 19.7 8461 -1.735 0.5799
step13 - step12 227.9 19.7 8461 11.592 <.0001
step14 - step13 -142.3 19.7 8461 -7.241 <.0001
step15 - step14 -52.2 19.7 8461 -2.654 0.0796
step16 - step15 -63.8 19.7 8461 -3.244 0.0154
step17 - step16 58.9 19.7 8461 2.999 0.0326
step18 - step17 -19.1 19.7 8461 -0.970 1.0000
Degrees-of-freedom method: kenward-roger
P value adjustment: holm method for 17 tests
[WRONG] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step + (1 | subject)
Data: data
REML criterion at convergence: 109180.4
Scaled residuals:
Min 1Q Median 3Q Max
-2.1545 -0.4740 -0.2235 0.1305 18.0003
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 54101 232.6
Residual 311631 558.2
Number of obs: 7056, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 628.117 55.394 11.339
step.L 77.092 28.195 2.734
step.Q 124.693 28.195 4.422
step.C -206.127 28.195 -7.311
step^4 309.668 28.195 10.983
step^5 -125.975 28.195 -4.468
step^6 212.995 28.195 7.554
step^7 -34.679 28.195 -1.230
step^8 85.672 28.195 3.039
step^9 -4.591 28.195 -0.163
step^10 -35.356 28.195 -1.254
step^11 -21.381 28.195 -0.758
step^12 79.047 28.195 2.804
step^13 -28.895 28.195 -1.025
step^14 -62.798 28.195 -2.227
step^15 -24.967 28.195 -0.886
step^16 -34.120 28.195 -1.210
step^17 -24.296 28.195 -0.862
Correlation matrix not shown by default, as p = 18 > 12.
Use print(summary(mdl), correlation=TRUE) or
vcov(summary(mdl)) if you need it
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 307.92 17 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
step emmean SE df lower.CL upper.CL
1 985 61.8 26.3 858 1112
2 478 61.8 26.3 352 605
3 466 61.8 26.3 339 593
4 480 61.8 26.3 353 607
5 567 61.8 26.3 440 694
6 618 61.8 26.3 491 744
7 561 61.8 26.3 434 688
8 581 61.8 26.3 454 708
9 678 61.8 26.3 551 805
10 675 61.8 26.3 548 802
11 634 61.8 26.3 507 761
12 602 61.8 26.3 475 729
13 719 61.8 26.3 592 846
14 648 61.8 26.3 521 775
15 670 61.8 26.3 543 797
16 600 61.8 26.3 473 727
17 583 61.8 26.3 456 709
18 761 61.8 26.3 634 888
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 # === Stepwise RT — Training (Blocks 1–3), CORRECT trials only ===
suppressPackageStartupMessages({
library(dplyr); library(ggplot2); library(emmeans); library(lme4); library(stringr)
})
if (!exists("safe_fit")) {
safe_fit <- function(formula, data) {
m <- tryCatch(lmer(formula, data = data), error = function(e) NULL)
if (is.null(m)) { message("⚠️ Falling back to lm."); m <- stats::lm(update(formula, . ~ . - (1 | subject)), data = data) }
m
}
}
plot_stepwise_rt_training_correct_only <- function(df_acc_all, show_plot = TRUE, alpha = 0.05) {
emms_list <- list()
sig_list <- list()
for (blk in 1:3) {
df_blk <- df_acc_all %>%
filter(session == blk, trial_acc_cat == "correct") %>%
group_by(subject, trial, step = as.integer(as.character(sub.trial.number))) %>%
summarise(rt = mean(feedback.RT, na.rm = TRUE), .groups = "drop")
if (nrow(df_blk) == 0) next
df_blk <- df_blk %>% mutate(step = factor(step, levels = sort(unique(step)), ordered = TRUE))
mdl <- safe_fit(rt ~ step + (1 | subject), data = df_blk)
em_grid <- emmeans(mdl, ~ step)
# Store EMMs for plotting
em <- as.data.frame(em_grid)
em$block <- blk
em$step_num <- as.integer(as.character(em$step))
emms_list[[blk]] <- em
# Adjacent (consecutive) contrasts (Holm-adjusted)
consec_df <- as.data.frame(contrast(em_grid, method = "consec", adjust = "holm"))
nums <- stringr::str_extract_all(consec_df$contrast, "\\d+")
n1 <- vapply(nums, function(v) as.integer(v[1]), integer(1)) # first number in label
n2 <- vapply(nums, function(v) as.integer(v[2]), integer(1)) # second number in label
# Ensure estimate refers to (next - prev); if label is lower - higher
prev <- pmin(n1, n2)
nextv <- pmax(n1, n2)
sign_adj <- ifelse(nextv == n1, 1, -1)
est_adj <- consec_df$estimate * sign_adj
sig_df <- consec_df %>%
mutate(step_prev = prev,
step_next = nextv,
estimate_adj = est_adj,
sig = !is.na(p.value) & p.value < alpha) %>%
filter(sig) %>%
transmute(
block = blk,
step_num = step_prev, # put star at the previous step (the one followed by a change)
direction = ifelse(estimate_adj > 0, "up", "down"),
p.value = p.value
)
if (nrow(sig_df)) sig_list[[blk]] <- sig_df
}
if (!length(emms_list)) { message("No correct-trial data for Blocks 1–3."); return(invisible(NULL)) }
df_em <- bind_rows(emms_list) %>%
mutate(
seq_label = case_when(block == 1 ~ "6 Steps",
block == 2 ~ "12 Steps",
TRUE ~ "18 Steps"),
seq_label = factor(seq_label, levels = c("6 Steps","12 Steps","18 Steps"), ordered = TRUE)
)
# Star positions
df_stars <- if (length(sig_list)) bind_rows(sig_list) else
tibble(block = integer(), step_num = integer(), direction = character(), p.value = numeric())
if (nrow(df_stars)) {
df_stars <- df_stars %>%
left_join(df_em %>% select(block, step_num, emmean, SE, seq_label), by = c("block","step_num")) %>%
mutate(
y_star = emmean + pmax(SE * 1.2, 35)
)
}
df_stars_up <- df_stars %>% filter(direction == "up")
df_stars_down <- df_stars %>% filter(direction == "down")
# separators at right edge of the first two panels
sep_df <- data.frame(seq_label = factor(c("6 Steps","12 Steps"),
levels = levels(df_em$seq_label), ordered = TRUE))
pal <- c("6 Steps" = "#B22222", "12 Steps" = "#2E7D32", "18 Steps" = "#1E3A8A")
p <- ggplot(df_em, aes(x = step_num, y = emmean, color = seq_label, fill = seq_label, group = 1)) +
geom_segment(data = sep_df,
aes(x = Inf, xend = Inf, y = -Inf, yend = Inf),
inherit.aes = FALSE, linetype = "dotted") +
geom_ribbon(aes(ymin = emmean - SE, ymax = emmean + SE), alpha = 0.15, color = NA) +
geom_line(linewidth = 0.9) +
geom_point(size = 2) +
# star markers: red = following step ↑; blue = following step ↓
{ if (nrow(df_stars_up)) geom_text(data = df_stars_up,
aes(x = step_num, y = y_star, label = "*"),
inherit.aes = FALSE, size = 4.2, fontface = "bold",
color = "#D32F2F") } +
{ if (nrow(df_stars_down)) geom_text(data = df_stars_down,
aes(x = step_num, y = y_star, label = "*"),
inherit.aes = FALSE, size = 4.2, fontface = "bold",
color = "#1F6FEB") } +
facet_grid(cols = vars(seq_label), scales = "free_x", space = "free_x") +
scale_x_continuous(
breaks = function(lims) seq(floor(lims[1]), ceiling(lims[2]), by = 1),
expand = expansion(mult = c(0.02, 0.06))
) +
scale_y_continuous(expand = expansion(mult = c(0.05, 0.10))) +
scale_color_manual(values = pal, breaks = names(pal)) +
scale_fill_manual(values = pal, breaks = names(pal)) +
labs(
title = "RT for each Step across Difficulty levels (Training phase)",
x = "Step", y = "Estimated RT (ms)",
caption = "* Following step is significantly different (red = \u2191, blue = \u2193)"
) +
theme_minimal(base_size = 12) +
theme(panel.grid = element_blank(),
legend.position = "none",
plot.caption = element_text(hjust = 0))
if (show_plot) print(p)
invisible(list(emms = df_em, stars = df_stars, plot = p, palette = pal))
}
# RUN
plot_stepwise_rt_training_correct_only(df_acc_base)Warning in geom_segment(data = sep_df, aes(x = Inf, xend = Inf, y = -Inf, : All aesthetics have length 1, but the data has 2 rows.
ℹ Please consider using `annotate()` or provide this layer with data containing
a single row.
#4 Stepwise RT — Test blocks (4–5), split by sequence length
#4 Stepwise RT — Test blocks (4–5), split by sequence length
stepwise_rt_test_by_length_with_correctness <- function(df_acc_all) {
cat("\n\n====== Stepwise RT — Test (4–5) by length ======\n")
df_len <- df_acc_all %>%
dplyr::group_by(subject, session, trial) %>%
dplyr::mutate(seq_length_trial = dplyr::n_distinct(sub.trial.number)) %>%
dplyr::ungroup()
for (L in c(6,12,18)) {
cat(glue::glue("\n--- Sequence length {L} ---\n"))
df_seq_all <- df_len %>%
dplyr::filter(session %in% c(4,5), seq_length_trial == L) %>%
dplyr::group_by(subject, session, trial, step = as.integer(as.character(sub.trial.number)), trial_acc_cat) %>%
dplyr::summarise(rt = mean(feedback.RT, na.rm = TRUE), .groups = "drop")
if (nrow(df_seq_all) == 0) { cat("No data\n"); next }
df_seq_all$session <- factor(df_seq_all$session, levels = c(4,5), labels = c("Block 4","Block 5"))
df_seq_all$step <- factor(df_seq_all$step, levels = sort(unique(df_seq_all$step)))
has_two_sessions <- dplyr::n_distinct(df_seq_all$session) >= 2
has_two_steps <- dplyr::n_distinct(df_seq_all$step) >= 2
has_two_corr <- dplyr::n_distinct(df_seq_all$trial_acc_cat) >= 2
# Combined model with correctness
if (has_two_sessions && has_two_steps && has_two_corr) {
fmla <- rt ~ step * session * trial_acc_cat + (1 | subject)
} else if (has_two_sessions && has_two_steps) {
fmla <- rt ~ step * session + (1 | subject)
} else if (has_two_steps && has_two_corr) {
fmla <- rt ~ step * trial_acc_cat + (1 | subject)
} else if (has_two_sessions && has_two_corr) {
fmla <- rt ~ session * trial_acc_cat + (1 | subject)
} else if (has_two_steps) {
fmla <- rt ~ step + (1 | subject)
} else if (has_two_sessions) {
fmla <- rt ~ session + (1 | subject)
} else if (has_two_corr) {
fmla <- rt ~ trial_acc_cat + (1 | subject)
} else {
cat("⚠️ Not enough levels to model.\n"); next
}
cat("\n[Combined] Model summary:\n")
mdl_all <- safe_fit(fmla, data = df_seq_all)
print(summary(mdl_all))
cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl_all, type = 2))
# Sensible EMMs depending on fixed effects present
if (grepl("step.*session.*trial_acc_cat", deparse(fmla))) {
em <- emmeans::emmeans(mdl_all, ~ step | session * trial_acc_cat)
print(summary(em))
} else if (grepl("step.*session", deparse(fmla))) {
em <- emmeans::emmeans(mdl_all, ~ step | session); print(summary(em))
} else if (grepl("step.*trial_acc_cat", deparse(fmla))) {
em <- emmeans::emmeans(mdl_all, ~ step | trial_acc_cat); print(summary(em))
} else if (grepl("session.*trial_acc_cat", deparse(fmla))) {
em <- emmeans::emmeans(mdl_all, ~ session | trial_acc_cat); print(summary(em))
} else if (grepl("^rt ~ step", deparse(fmla))) {
em <- emmeans::emmeans(mdl_all, ~ step); print(summary(em))
} else if (grepl("^rt ~ session", deparse(fmla))) {
em <- emmeans::emmeans(mdl_all, ~ session); print(summary(em))
} else {
em <- emmeans::emmeans(mdl_all, ~ trial_acc_cat); print(summary(em))
}
# --- Split outputs: add pairwise ONLY for correct
for (lab in c("correct","wrong")) {
df_seq <- df_seq_all %>% dplyr::filter(trial_acc_cat == lab)
if (nrow(df_seq) == 0) { cat("\n⚠️ No data for", lab, "\n"); next }
has_two_sessions <- dplyr::n_distinct(df_seq$session) >= 2
has_two_steps <- dplyr::n_distinct(df_seq$step) >= 2
if (!has_two_steps) { cat("\n⚠️ Not enough step levels for", lab, "\n"); next }
fmla_s <- if (has_two_sessions) rt ~ step * session + (1 | subject) else rt ~ step + (1 | subject)
cat("\n[", toupper(lab), "] Model summary:\n", sep = "")
mdl <- safe_fit(fmla_s, data = df_seq)
print(summary(mdl))
cat("\nType II Chi-square ANOVA:\n"); print(car::Anova(mdl, type = 2))
# EMMs (always by step)
if (has_two_sessions) {
em_s <- emmeans::emmeans(mdl, ~ step | session)
} else {
em_s <- emmeans::emmeans(mdl, ~ step)
}
print(summary(em_s))
# Pairwise for CORRECT only
if (lab == "correct") {
if (has_two_sessions) {
cat("\nPairwise (Tukey) across steps — correct | session:\n")
} else {
cat("\nPairwise (Tukey) across steps — correct:\n")
}
print(pairs(em_s, adjust = "tukey"))
if (has_two_sessions) {
cat("\nAdjacent-step contrasts (consecutive) — correct | session:\n")
} else {
cat("\nAdjacent-step contrasts (consecutive) — correct:\n")
}
print(emmeans::contrast(em_s, method = "consec", adjust = "holm"))
}
}
}
}
stepwise_rt_test_by_length_with_correctness(df_acc_base)
====== Stepwise RT — Test (4–5) by length ======
--- Sequence length 6 ---
[Combined] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step * session * trial_acc_cat + (1 | subject)
Data: data
REML criterion at convergence: 51434.3
Scaled residuals:
Min 1Q Median 3Q Max
-2.5669 -0.3601 -0.0877 0.1859 27.0864
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 33464 182.9
Residual 179548 423.7
Number of obs: 3456, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 764.911 50.749 15.073
step2 -348.940 37.824 -9.225
step3 -390.669 37.824 -10.329
step4 -402.183 37.824 -10.633
step5 -382.510 37.824 -10.113
step6 -376.163 37.824 -9.945
sessionBlock 5 71.983 39.302 1.832
trial_acc_catwrong 99.424 75.011 1.325
step2:sessionBlock 5 13.248 55.479 0.239
step3:sessionBlock 5 48.610 55.479 0.876
step4:sessionBlock 5 56.422 55.479 1.017
step5:sessionBlock 5 36.051 55.479 0.650
step6:sessionBlock 5 29.439 55.479 0.531
step2:trial_acc_catwrong 50.751 105.527 0.481
step3:trial_acc_catwrong 94.048 105.527 0.891
step4:trial_acc_catwrong -75.519 105.527 -0.716
step5:trial_acc_catwrong 63.591 105.527 0.603
step6:trial_acc_catwrong 700.163 105.527 6.635
sessionBlock 5:trial_acc_catwrong 330.744 95.319 3.470
step2:sessionBlock 5:trial_acc_catwrong -430.244 133.840 -3.215
step3:sessionBlock 5:trial_acc_catwrong -366.302 133.840 -2.737
step4:sessionBlock 5:trial_acc_catwrong -274.605 133.840 -2.052
step5:sessionBlock 5:trial_acc_catwrong 8.282 133.840 0.062
step6:sessionBlock 5:trial_acc_catwrong -682.424 133.840 -5.099
Correlation matrix not shown by default, as p = 24 > 12.
Use print(summary(mdl_all), correlation=TRUE) or
vcov(summary(mdl_all)) if you need it
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 394.1445 5 < 2.2e-16 ***
session 56.3310 1 6.124e-14 ***
trial_acc_cat 172.6516 1 < 2.2e-16 ***
step:session 8.5065 5 0.1304
step:trial_acc_cat 92.9202 5 < 2.2e-16 ***
session:trial_acc_cat 0.9797 1 0.3223
step:session:trial_acc_cat 39.3891 5 1.983e-07 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
session = Block 4, trial_acc_cat = correct:
step emmean SE df lower.CL upper.CL
1 765 50.7 30.8 661 868
2 416 50.7 30.8 312 519
3 374 50.7 30.8 271 478
4 363 50.7 30.8 259 466
5 382 50.7 30.8 279 486
6 389 50.7 30.8 285 492
session = Block 5, trial_acc_cat = correct:
step emmean SE df lower.CL upper.CL
1 837 51.8 33.6 731 942
2 501 51.8 33.6 396 607
3 495 51.8 33.6 389 600
4 491 51.8 33.6 386 597
5 490 51.8 33.6 385 596
6 490 51.8 33.6 385 596
session = Block 4, trial_acc_cat = wrong:
step emmean SE df lower.CL upper.CL
1 864 82.2 205.5 702 1026
2 566 82.2 205.5 404 728
3 568 82.2 205.5 406 730
4 387 82.2 205.5 225 549
5 545 82.2 205.5 383 707
6 1188 82.2 205.5 1026 1350
session = Block 5, trial_acc_cat = wrong:
step emmean SE df lower.CL upper.CL
1 1267 66.9 92.2 1134 1400
2 552 66.9 92.2 419 685
3 653 66.9 92.2 520 786
4 571 66.9 92.2 438 704
5 992 66.9 92.2 860 1125
6 938 66.9 92.2 805 1071
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
[CORRECT] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step * session + (1 | subject)
Data: data
REML criterion at convergence: 39822.1
Scaled residuals:
Min 1Q Median 3Q Max
-2.7780 -0.4318 -0.1148 0.2646 20.8747
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 28464 168.7
Residual 83002 288.1
Number of obs: 2814, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 764.27 43.73 17.475
step2 -348.94 25.72 -13.568
step3 -390.67 25.72 -15.191
step4 -402.18 25.72 -15.639
step5 -382.51 25.72 -14.874
step6 -376.16 25.72 -14.627
sessionBlock 5 74.40 26.75 2.781
step2:sessionBlock 5 13.25 37.72 0.351
step3:sessionBlock 5 48.61 37.72 1.289
step4:sessionBlock 5 56.42 37.72 1.496
step5:sessionBlock 5 36.05 37.72 0.956
step6:sessionBlock 5 29.44 37.72 0.780
Correlation of Fixed Effects:
(Intr) step2 step3 step4 step5 step6 sssnB5 st2:B5 st3:B5
step2 -0.294
step3 -0.294 0.500
step4 -0.294 0.500 0.500
step5 -0.294 0.500 0.500 0.500
step6 -0.294 0.500 0.500 0.500 0.500
sessinBlck5 -0.283 0.481 0.481 0.481 0.481 0.481
stp2:sssnB5 0.200 -0.682 -0.341 -0.341 -0.341 -0.341 -0.705
stp3:sssnB5 0.200 -0.341 -0.682 -0.341 -0.341 -0.341 -0.705 0.500
stp4:sssnB5 0.200 -0.341 -0.341 -0.682 -0.341 -0.341 -0.705 0.500 0.500
stp5:sssnB5 0.200 -0.341 -0.341 -0.341 -0.682 -0.341 -0.705 0.500 0.500
stp6:sssnB5 0.200 -0.341 -0.341 -0.341 -0.341 -0.682 -0.705 0.500 0.500
st4:B5 st5:B5
step2
step3
step4
step5
step6
sessinBlck5
stp2:sssnB5
stp3:sssnB5
stp4:sssnB5
stp5:sssnB5 0.500
stp6:sssnB5 0.500 0.500
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 623.9439 5 <2e-16 ***
session 89.7160 1 <2e-16 ***
step:session 3.1762 5 0.6728
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
session = Block 4:
step emmean SE df lower.CL upper.CL
1 764 43.7 23.9 674 855
2 415 43.7 23.9 325 506
3 374 43.7 23.9 283 464
4 362 43.7 23.9 272 452
5 382 43.7 23.9 291 472
6 388 43.7 23.9 298 478
session = Block 5:
step emmean SE df lower.CL upper.CL
1 839 44.3 25.3 747 930
2 503 44.3 25.3 412 594
3 497 44.3 25.3 405 588
4 493 44.3 25.3 402 584
5 492 44.3 25.3 401 583
6 492 44.3 25.3 401 583
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
Pairwise (Tukey) across steps — correct | session:
session = Block 4:
contrast estimate SE df t.ratio p.value
step1 - step2 348.940 25.7 2785 13.568 <.0001
step1 - step3 390.669 25.7 2785 15.191 <.0001
step1 - step4 402.183 25.7 2785 15.639 <.0001
step1 - step5 382.510 25.7 2785 14.874 <.0001
step1 - step6 376.163 25.7 2785 14.627 <.0001
step2 - step3 41.729 25.7 2785 1.623 0.5835
step2 - step4 53.243 25.7 2785 2.070 0.3032
step2 - step5 33.570 25.7 2785 1.305 0.7822
step2 - step6 27.223 25.7 2785 1.059 0.8976
step3 - step4 11.514 25.7 2785 0.448 0.9977
step3 - step5 -8.159 25.7 2785 -0.317 0.9996
step3 - step6 -14.506 25.7 2785 -0.564 0.9933
step4 - step5 -19.673 25.7 2785 -0.765 0.9733
step4 - step6 -26.020 25.7 2785 -1.012 0.9142
step5 - step6 -6.347 25.7 2785 -0.247 0.9999
session = Block 5:
contrast estimate SE df t.ratio p.value
step1 - step2 335.693 27.6 2785 12.165 <.0001
step1 - step3 342.060 27.6 2785 12.396 <.0001
step1 - step4 345.762 27.6 2785 12.530 <.0001
step1 - step5 346.459 27.6 2785 12.555 <.0001
step1 - step6 346.725 27.6 2785 12.565 <.0001
step2 - step3 6.367 27.6 2785 0.231 0.9999
step2 - step4 10.069 27.6 2785 0.365 0.9992
step2 - step5 10.766 27.6 2785 0.390 0.9988
step2 - step6 11.032 27.6 2785 0.400 0.9987
step3 - step4 3.702 27.6 2785 0.134 1.0000
step3 - step5 4.399 27.6 2785 0.159 1.0000
step3 - step6 4.665 27.6 2785 0.169 1.0000
step4 - step5 0.697 27.6 2785 0.025 1.0000
step4 - step6 0.963 27.6 2785 0.035 1.0000
step5 - step6 0.266 27.6 2785 0.010 1.0000
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 6 estimates
Adjacent-step contrasts (consecutive) — correct | session:
session = Block 4:
contrast estimate SE df t.ratio p.value
step2 - step1 -348.940 25.7 2785 -13.568 <.0001
step3 - step2 -41.729 25.7 2785 -1.623 0.4191
step4 - step3 -11.514 25.7 2785 -0.448 1.0000
step5 - step4 19.673 25.7 2785 0.765 1.0000
step6 - step5 6.347 25.7 2785 0.247 1.0000
session = Block 5:
contrast estimate SE df t.ratio p.value
step2 - step1 -335.693 27.6 2785 -12.165 <.0001
step3 - step2 -6.367 27.6 2785 -0.231 1.0000
step4 - step3 -3.702 27.6 2785 -0.134 1.0000
step5 - step4 -0.697 27.6 2785 -0.025 1.0000
step6 - step5 -0.266 27.6 2785 -0.010 1.0000
Degrees-of-freedom method: kenward-roger
P value adjustment: holm method for 5 tests
[WRONG] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step * session + (1 | subject)
Data: data
REML criterion at convergence: 10242.4
Scaled residuals:
Min 1Q Median 3Q Max
-1.5788 -0.4499 -0.1617 0.1451 14.6358
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 68405 261.5
Residual 600745 775.1
Number of obs: 642, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 863.22 143.94 5.997
step2 -298.19 180.20 -1.655
step3 -296.62 180.20 -1.646
step4 -477.70 180.20 -2.651
step5 -318.92 180.20 -1.770
step6 324.00 180.20 1.798
sessionBlock 5 395.00 160.56 2.460
step2:sessionBlock 5 -417.00 222.79 -1.872
step3:sessionBlock 5 -317.69 222.79 -1.426
step4:sessionBlock 5 -218.18 222.79 -0.979
step5:sessionBlock 5 44.33 222.79 0.199
step6:sessionBlock 5 -652.99 222.79 -2.931
Correlation of Fixed Effects:
(Intr) step2 step3 step4 step5 step6 sssnB5 st2:B5 st3:B5
step2 -0.626
step3 -0.626 0.500
step4 -0.626 0.500 0.500
step5 -0.626 0.500 0.500 0.500
step6 -0.626 0.500 0.500 0.500 0.500
sessinBlck5 -0.726 0.561 0.561 0.561 0.561 0.561
stp2:sssnB5 0.506 -0.809 -0.404 -0.404 -0.404 -0.404 -0.694
stp3:sssnB5 0.506 -0.404 -0.809 -0.404 -0.404 -0.404 -0.694 0.500
stp4:sssnB5 0.506 -0.404 -0.404 -0.809 -0.404 -0.404 -0.694 0.500 0.500
stp5:sssnB5 0.506 -0.404 -0.404 -0.404 -0.809 -0.404 -0.694 0.500 0.500
stp6:sssnB5 0.506 -0.404 -0.404 -0.404 -0.404 -0.809 -0.694 0.500 0.500
st4:B5 st5:B5
step2
step3
step4
step5
step6
sessinBlck5
stp2:sssnB5
stp3:sssnB5
stp4:sssnB5
stp5:sssnB5 0.500
stp6:sssnB5 0.500 0.500
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 59.2825 5 1.71e-11 ***
session 3.5623 1 0.05911 .
step:session 13.8760 5 0.01642 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
session = Block 4:
step emmean SE df lower.CL upper.CL
1 863 144 206 579 1147
2 565 144 206 281 849
3 567 144 206 282 851
4 386 144 206 101 670
5 544 144 206 260 828
6 1187 144 206 903 1471
session = Block 5:
step emmean SE df lower.CL upper.CL
1 1258 114 97 1032 1485
2 543 114 97 317 769
3 644 114 97 418 870
4 562 114 97 336 789
5 984 114 97 757 1210
6 929 114 97 703 1156
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
--- Sequence length 12 ---
[Combined] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step * session * trial_acc_cat + (1 | subject)
Data: data
REML criterion at convergence: 104004.8
Scaled residuals:
Min 1Q Median 3Q Max
-2.2353 -0.4181 -0.1554 0.1457 19.2708
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 36507 191.1
Residual 213173 461.7
Number of obs: 6912, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 739.049 54.599 13.536
step2 -303.513 43.627 -6.957
step3 -343.759 43.627 -7.879
step4 -341.710 43.627 -7.833
step5 -300.915 43.627 -6.897
step6 -266.804 43.627 -6.116
step7 -154.647 43.627 -3.545
step8 -333.433 43.627 -7.643
step9 -247.138 43.627 -5.665
step10 -286.696 43.627 -6.572
step11 -317.045 43.627 -7.267
step12 -334.058 43.627 -7.657
sessionBlock 5 229.706 50.498 4.549
trial_acc_catwrong 21.904 65.621 0.334
step2:sessionBlock 5 -130.300 71.309 -1.827
step3:sessionBlock 5 -110.547 71.309 -1.550
step4:sessionBlock 5 -90.380 71.309 -1.267
step5:sessionBlock 5 -127.450 71.309 -1.787
step6:sessionBlock 5 -106.868 71.309 -1.499
step7:sessionBlock 5 -98.166 71.309 -1.377
step8:sessionBlock 5 -21.619 71.309 -0.303
step9:sessionBlock 5 45.623 71.309 0.640
step10:sessionBlock 5 -85.281 71.309 -1.196
step11:sessionBlock 5 -202.784 71.309 -2.844
step12:sessionBlock 5 -127.218 71.309 -1.784
step2:trial_acc_catwrong -14.940 92.547 -0.161
step3:trial_acc_catwrong 8.118 92.547 0.088
step4:trial_acc_catwrong 22.663 92.547 0.245
step5:trial_acc_catwrong 69.400 92.547 0.750
step6:trial_acc_catwrong -15.196 92.547 -0.164
step7:trial_acc_catwrong 43.491 92.547 0.470
step8:trial_acc_catwrong 273.339 92.547 2.954
step9:trial_acc_catwrong 53.717 92.547 0.580
step10:trial_acc_catwrong 90.493 92.547 0.978
step11:trial_acc_catwrong 316.967 92.547 3.425
step12:trial_acc_catwrong 407.855 92.547 4.407
sessionBlock 5:trial_acc_catwrong 209.009 85.445 2.446
step2:sessionBlock 5:trial_acc_catwrong -133.286 120.479 -1.106
step3:sessionBlock 5:trial_acc_catwrong -153.663 120.479 -1.275
step4:sessionBlock 5:trial_acc_catwrong -111.495 120.479 -0.925
step5:sessionBlock 5:trial_acc_catwrong -134.047 120.479 -1.113
step6:sessionBlock 5:trial_acc_catwrong -165.366 120.479 -1.373
step7:sessionBlock 5:trial_acc_catwrong -111.892 120.479 -0.929
step8:sessionBlock 5:trial_acc_catwrong -399.417 120.479 -3.315
step9:sessionBlock 5:trial_acc_catwrong -216.890 120.479 -1.800
step10:sessionBlock 5:trial_acc_catwrong -265.373 120.479 -2.203
step11:sessionBlock 5:trial_acc_catwrong -256.846 120.479 -2.132
step12:sessionBlock 5:trial_acc_catwrong -420.046 120.479 -3.486
Correlation matrix not shown by default, as p = 48 > 12.
Use print(summary(mdl_all), correlation=TRUE) or
vcov(summary(mdl_all)) if you need it
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 424.9073 11 < 2.2e-16 ***
session 153.1544 1 < 2.2e-16 ***
trial_acc_cat 111.2216 1 < 2.2e-16 ***
step:session 45.6047 11 3.797e-06 ***
step:trial_acc_cat 53.1721 11 1.670e-07 ***
session:trial_acc_cat 0.2094 1 0.64725
step:session:trial_acc_cat 22.5471 11 0.02046 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
session = Block 4, trial_acc_cat = correct:
step emmean SE df lower.CL upper.CL
1 739 54.6 35.6 628 850
2 436 54.6 35.6 325 546
3 395 54.6 35.6 285 506
4 397 54.6 35.6 287 508
5 438 54.6 35.6 327 549
6 472 54.6 35.6 361 583
7 584 54.6 35.6 474 695
8 406 54.6 35.6 295 516
9 492 54.6 35.6 381 603
10 452 54.6 35.6 342 563
11 422 54.6 35.6 311 533
12 405 54.6 35.6 294 516
session = Block 5, trial_acc_cat = correct:
step emmean SE df lower.CL upper.CL
1 969 60.2 52.6 848 1090
2 535 60.2 52.6 414 656
3 514 60.2 52.6 394 635
4 537 60.2 52.6 416 657
5 540 60.2 52.6 420 661
6 595 60.2 52.6 474 716
7 716 60.2 52.6 595 837
8 614 60.2 52.6 493 735
9 767 60.2 52.6 646 888
10 597 60.2 52.6 476 718
11 449 60.2 52.6 328 570
12 507 60.2 52.6 387 628
session = Block 4, trial_acc_cat = wrong:
step emmean SE df lower.CL upper.CL
1 761 73.3 115.0 616 906
2 442 73.3 115.0 297 588
3 425 73.3 115.0 280 571
4 442 73.3 115.0 297 587
5 529 73.3 115.0 384 675
6 479 73.3 115.0 334 624
7 650 73.3 115.0 505 795
8 701 73.3 115.0 556 846
9 568 73.3 115.0 422 713
10 565 73.3 115.0 420 710
11 761 73.3 115.0 616 906
12 835 73.3 115.0 690 980
session = Block 5, trial_acc_cat = wrong:
step emmean SE df lower.CL upper.CL
1 1200 58.5 46.8 1082 1317
2 618 58.5 46.8 500 735
3 600 58.5 46.8 482 717
4 679 58.5 46.8 561 796
5 707 58.5 46.8 589 824
6 645 58.5 46.8 528 763
7 878 58.5 46.8 761 996
8 719 58.5 46.8 601 836
9 835 58.5 46.8 717 953
10 653 58.5 46.8 535 770
11 740 58.5 46.8 622 858
12 726 58.5 46.8 609 844
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
[CORRECT] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step * session + (1 | subject)
Data: data
REML criterion at convergence: 62279.4
Scaled residuals:
Min 1Q Median 3Q Max
-2.3071 -0.3963 -0.1505 0.1613 25.8645
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 34623 186.1
Residual 119957 346.3
Number of obs: 4296, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 738.71 49.60 14.893
step2 -303.51 32.73 -9.274
step3 -343.76 32.73 -10.504
step4 -341.71 32.73 -10.441
step5 -300.92 32.73 -9.195
step6 -266.80 32.73 -8.152
step7 -154.65 32.73 -4.725
step8 -333.43 32.73 -10.188
step9 -247.14 32.73 -7.552
step10 -286.70 32.73 -8.760
step11 -317.04 32.73 -9.688
step12 -334.06 32.73 -10.207
sessionBlock 5 231.17 37.92 6.096
step2:sessionBlock 5 -130.30 53.49 -2.436
step3:sessionBlock 5 -110.55 53.49 -2.067
step4:sessionBlock 5 -90.38 53.49 -1.690
step5:sessionBlock 5 -127.45 53.49 -2.383
step6:sessionBlock 5 -106.87 53.49 -1.998
step7:sessionBlock 5 -98.17 53.49 -1.835
step8:sessionBlock 5 -21.62 53.49 -0.404
step9:sessionBlock 5 45.62 53.49 0.853
step10:sessionBlock 5 -85.28 53.49 -1.594
step11:sessionBlock 5 -202.78 53.49 -3.791
step12:sessionBlock 5 -127.22 53.49 -2.378
Correlation matrix not shown by default, as p = 24 > 12.
Use print(summary(mdl), correlation=TRUE) or
vcov(summary(mdl)) if you need it
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 424.13 11 < 2.2e-16 ***
session 162.02 1 < 2.2e-16 ***
step:session 34.28 11 0.0003255 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
session = Block 4:
step emmean SE df lower.CL upper.CL
1 739 49.6 27.0 637 840
2 435 49.6 27.0 333 537
3 395 49.6 27.0 293 497
4 397 49.6 27.0 295 499
5 438 49.6 27.0 336 540
6 472 49.6 27.0 370 574
7 584 49.6 27.0 482 686
8 405 49.6 27.0 303 507
9 492 49.6 27.0 390 593
10 452 49.6 27.0 350 554
11 422 49.6 27.0 320 523
12 405 49.6 27.0 303 506
session = Block 5:
step emmean SE df lower.CL upper.CL
1 970 53.2 35.5 862 1078
2 536 53.2 35.5 428 644
3 516 53.2 35.5 408 623
4 538 53.2 35.5 430 646
5 542 53.2 35.5 434 649
6 596 53.2 35.5 488 704
7 717 53.2 35.5 609 825
8 615 53.2 35.5 507 723
9 768 53.2 35.5 661 876
10 598 53.2 35.5 490 706
11 450 53.2 35.5 342 558
12 509 53.2 35.5 401 616
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
Pairwise (Tukey) across steps — correct | session:
session = Block 4:
contrast estimate SE df t.ratio p.value
step1 - step2 303.513 32.7 4255 9.274 <.0001
step1 - step3 343.759 32.7 4255 10.504 <.0001
step1 - step4 341.710 32.7 4255 10.441 <.0001
step1 - step5 300.915 32.7 4255 9.195 <.0001
step1 - step6 266.804 32.7 4255 8.152 <.0001
step1 - step7 154.647 32.7 4255 4.725 0.0001
step1 - step8 333.433 32.7 4255 10.188 <.0001
step1 - step9 247.138 32.7 4255 7.552 <.0001
step1 - step10 286.696 32.7 4255 8.760 <.0001
step1 - step11 317.045 32.7 4255 9.688 <.0001
step1 - step12 334.058 32.7 4255 10.207 <.0001
step2 - step3 40.245 32.7 4255 1.230 0.9867
step2 - step4 38.196 32.7 4255 1.167 0.9913
step2 - step5 -2.598 32.7 4255 -0.079 1.0000
step2 - step6 -36.710 32.7 4255 -1.122 0.9938
step2 - step7 -148.866 32.7 4255 -4.549 0.0003
step2 - step8 29.920 32.7 4255 0.914 0.9990
step2 - step9 -56.375 32.7 4255 -1.723 0.8581
step2 - step10 -16.817 32.7 4255 -0.514 1.0000
step2 - step11 13.531 32.7 4255 0.413 1.0000
step2 - step12 30.545 32.7 4255 0.933 0.9988
step3 - step4 -2.049 32.7 4255 -0.063 1.0000
step3 - step5 -42.844 32.7 4255 -1.309 0.9781
step3 - step6 -76.955 32.7 4255 -2.351 0.4399
step3 - step7 -189.112 32.7 4255 -5.778 <.0001
step3 - step8 -10.326 32.7 4255 -0.316 1.0000
step3 - step9 -96.621 32.7 4255 -2.952 0.1234
step3 - step10 -57.062 32.7 4255 -1.744 0.8478
step3 - step11 -26.714 32.7 4255 -0.816 0.9997
step3 - step12 -9.701 32.7 4255 -0.296 1.0000
step4 - step5 -40.795 32.7 4255 -1.247 0.9851
step4 - step6 -74.906 32.7 4255 -2.289 0.4847
step4 - step7 -187.062 32.7 4255 -5.716 <.0001
step4 - step8 -8.277 32.7 4255 -0.253 1.0000
step4 - step9 -94.571 32.7 4255 -2.890 0.1448
step4 - step10 -55.013 32.7 4255 -1.681 0.8771
step4 - step11 -24.665 32.7 4255 -0.754 0.9998
step4 - step12 -7.652 32.7 4255 -0.234 1.0000
step5 - step6 -34.112 32.7 4255 -1.042 0.9967
step5 - step7 -146.268 32.7 4255 -4.469 0.0005
step5 - step8 32.518 32.7 4255 0.994 0.9979
step5 - step9 -53.777 32.7 4255 -1.643 0.8930
step5 - step10 -14.219 32.7 4255 -0.434 1.0000
step5 - step11 16.130 32.7 4255 0.493 1.0000
step5 - step12 33.143 32.7 4255 1.013 0.9975
step6 - step7 -112.156 32.7 4255 -3.427 0.0303
step6 - step8 66.629 32.7 4255 2.036 0.6686
step6 - step9 -19.665 32.7 4255 -0.601 1.0000
step6 - step10 19.893 32.7 4255 0.608 1.0000
step6 - step11 50.241 32.7 4255 1.535 0.9309
step6 - step12 67.254 32.7 4255 2.055 0.6551
step7 - step8 178.786 32.7 4255 5.463 <.0001
step7 - step9 92.491 32.7 4255 2.826 0.1693
step7 - step10 132.049 32.7 4255 4.035 0.0032
step7 - step11 162.397 32.7 4255 4.962 <.0001
step7 - step12 179.411 32.7 4255 5.482 <.0001
step8 - step9 -86.295 32.7 4255 -2.637 0.2593
step8 - step10 -46.737 32.7 4255 -1.428 0.9582
step8 - step11 -16.388 32.7 4255 -0.501 1.0000
step8 - step12 0.625 32.7 4255 0.019 1.0000
step9 - step10 39.558 32.7 4255 1.209 0.9884
step9 - step11 69.906 32.7 4255 2.136 0.5966
step9 - step12 86.920 32.7 4255 2.656 0.2490
step10 - step11 30.348 32.7 4255 0.927 0.9989
step10 - step12 47.362 32.7 4255 1.447 0.9540
step11 - step12 17.013 32.7 4255 0.520 1.0000
session = Block 5:
contrast estimate SE df t.ratio p.value
step1 - step2 433.813 42.3 4255 10.252 <.0001
step1 - step3 454.306 42.3 4255 10.737 <.0001
step1 - step4 432.090 42.3 4255 10.212 <.0001
step1 - step5 428.366 42.3 4255 10.124 <.0001
step1 - step6 373.672 42.3 4255 8.831 <.0001
step1 - step7 252.813 42.3 4255 5.975 <.0001
step1 - step8 355.052 42.3 4255 8.391 <.0001
step1 - step9 201.515 42.3 4255 4.762 0.0001
step1 - step10 371.978 42.3 4255 8.791 <.0001
step1 - step11 519.828 42.3 4255 12.285 <.0001
step1 - step12 461.276 42.3 4255 10.901 <.0001
step2 - step3 20.492 42.3 4255 0.484 1.0000
step2 - step4 -1.724 42.3 4255 -0.041 1.0000
step2 - step5 -5.448 42.3 4255 -0.129 1.0000
step2 - step6 -60.142 42.3 4255 -1.421 0.9596
step2 - step7 -181.000 42.3 4255 -4.278 0.0012
step2 - step8 -78.761 42.3 4255 -1.861 0.7832
step2 - step9 -232.298 42.3 4255 -5.490 <.0001
step2 - step10 -61.836 42.3 4255 -1.461 0.9507
step2 - step11 86.015 42.3 4255 2.033 0.6708
step2 - step12 27.463 42.3 4255 0.649 1.0000
step3 - step4 -22.216 42.3 4255 -0.525 1.0000
step3 - step5 -25.940 42.3 4255 -0.613 1.0000
step3 - step6 -80.634 42.3 4255 -1.906 0.7559
step3 - step7 -201.493 42.3 4255 -4.762 0.0001
step3 - step8 -99.254 42.3 4255 -2.346 0.4440
step3 - step9 -252.791 42.3 4255 -5.974 <.0001
step3 - step10 -82.328 42.3 4255 -1.946 0.7302
step3 - step11 65.522 42.3 4255 1.549 0.9268
step3 - step12 6.970 42.3 4255 0.165 1.0000
step4 - step5 -3.724 42.3 4255 -0.088 1.0000
step4 - step6 -58.418 42.3 4255 -1.381 0.9673
step4 - step7 -179.276 42.3 4255 -4.237 0.0014
step4 - step8 -77.037 42.3 4255 -1.821 0.8069
step4 - step9 -230.575 42.3 4255 -5.449 <.0001
step4 - step10 -60.112 42.3 4255 -1.421 0.9597
step4 - step11 87.739 42.3 4255 2.074 0.6419
step4 - step12 29.187 42.3 4255 0.690 0.9999
step5 - step6 -54.694 42.3 4255 -1.293 0.9802
step5 - step7 -175.552 42.3 4255 -4.149 0.0020
step5 - step8 -73.313 42.3 4255 -1.733 0.8532
step5 - step9 -226.851 42.3 4255 -5.361 <.0001
step5 - step10 -56.388 42.3 4255 -1.333 0.9749
step5 - step11 91.463 42.3 4255 2.162 0.5779
step5 - step12 32.910 42.3 4255 0.778 0.9998
step6 - step7 -120.858 42.3 4255 -2.856 0.1573
step6 - step8 -18.619 42.3 4255 -0.440 1.0000
step6 - step9 -172.157 42.3 4255 -4.069 0.0028
step6 - step10 -1.694 42.3 4255 -0.040 1.0000
step6 - step11 146.157 42.3 4255 3.454 0.0277
step6 - step12 87.605 42.3 4255 2.070 0.6442
step7 - step8 102.239 42.3 4255 2.416 0.3950
step7 - step9 -51.298 42.3 4255 -1.212 0.9881
step7 - step10 119.164 42.3 4255 2.816 0.1733
step7 - step11 267.015 42.3 4255 6.310 <.0001
step7 - step12 208.463 42.3 4255 4.927 0.0001
step8 - step9 -153.537 42.3 4255 -3.629 0.0152
step8 - step10 16.925 42.3 4255 0.400 1.0000
step8 - step11 164.776 42.3 4255 3.894 0.0056
step8 - step12 106.224 42.3 4255 2.510 0.3334
step9 - step10 170.463 42.3 4255 4.029 0.0033
step9 - step11 318.313 42.3 4255 7.523 <.0001
step9 - step12 259.761 42.3 4255 6.139 <.0001
step10 - step11 147.851 42.3 4255 3.494 0.0242
step10 - step12 89.299 42.3 4255 2.110 0.6153
step11 - step12 -58.552 42.3 4255 -1.384 0.9667
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 12 estimates
Adjacent-step contrasts (consecutive) — correct | session:
session = Block 4:
contrast estimate SE df t.ratio p.value
step2 - step1 -303.51 32.7 4255 -9.274 <.0001
step3 - step2 -40.25 32.7 4255 -1.230 1.0000
step4 - step3 2.05 32.7 4255 0.063 1.0000
step5 - step4 40.79 32.7 4255 1.247 1.0000
step6 - step5 34.11 32.7 4255 1.042 1.0000
step7 - step6 112.16 32.7 4255 3.427 0.0055
step8 - step7 -178.79 32.7 4255 -5.463 <.0001
step9 - step8 86.29 32.7 4255 2.637 0.0672
step10 - step9 -39.56 32.7 4255 -1.209 1.0000
step11 - step10 -30.35 32.7 4255 -0.927 1.0000
step12 - step11 -17.01 32.7 4255 -0.520 1.0000
session = Block 5:
contrast estimate SE df t.ratio p.value
step2 - step1 -433.81 42.3 4255 -10.252 <.0001
step3 - step2 -20.49 42.3 4255 -0.484 1.0000
step4 - step3 22.22 42.3 4255 0.525 1.0000
step5 - step4 3.72 42.3 4255 0.088 1.0000
step6 - step5 54.69 42.3 4255 1.293 0.8325
step7 - step6 120.86 42.3 4255 2.856 0.0302
step8 - step7 -102.24 42.3 4255 -2.416 0.0943
step9 - step8 153.54 42.3 4255 3.629 0.0026
step10 - step9 -170.46 42.3 4255 -4.029 0.0006
step11 - step10 -147.85 42.3 4255 -3.494 0.0038
step12 - step11 58.55 42.3 4255 1.384 0.8325
Degrees-of-freedom method: kenward-roger
P value adjustment: holm method for 11 tests
[WRONG] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step * session + (1 | subject)
Data: data
REML criterion at convergence: 40704.9
Scaled residuals:
Min 1Q Median 3Q Max
-1.7587 -0.4847 -0.2199 0.1523 10.2604
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 44324 210.5
Residual 363999 603.3
Number of obs: 2616, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 759.62495 90.66198 8.379
step2 -318.45313 106.65348 -2.986
step3 -335.64063 106.65348 -3.147
step4 -319.04688 106.65348 -2.991
step5 -231.51563 106.65348 -2.171
step6 -282.00000 106.65348 -2.644
step7 -111.15625 106.65348 -1.042
step8 -60.09375 106.65348 -0.563
step9 -193.42188 106.65348 -1.814
step10 -196.20313 106.65348 -1.840
step11 -0.07813 106.65348 -0.001
step12 73.79687 106.65348 0.692
sessionBlock 5 436.34131 90.12844 4.841
step2:sessionBlock 5 -263.58584 126.89454 -2.077
step3:sessionBlock 5 -264.21002 126.89454 -2.082
step4:sessionBlock 5 -201.87520 126.89454 -1.591
step5:sessionBlock 5 -261.49736 126.89454 -2.061
step6:sessionBlock 5 -272.23377 126.89454 -2.145
step7:sessionBlock 5 -210.05804 126.89454 -1.655
step8:sessionBlock 5 -421.03612 126.89454 -3.318
step9:sessionBlock 5 -171.26644 126.89454 -1.350
step10:sessionBlock 5 -350.65402 126.89454 -2.763
step11:sessionBlock 5 -459.62967 126.89454 -3.622
step12:sessionBlock 5 -547.26441 126.89454 -4.313
Correlation matrix not shown by default, as p = 24 > 12.
Use print(summary(mdl), correlation=TRUE) or
vcov(summary(mdl)) if you need it
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 135.160 11 < 2.2e-16 ***
session 30.716 1 2.986e-08 ***
step:session 28.615 11 0.002603 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
session = Block 4:
step emmean SE df lower.CL upper.CL
1 760 90.7 160.7 581 939
2 441 90.7 160.7 262 620
3 424 90.7 160.7 245 603
4 441 90.7 160.7 262 620
5 528 90.7 160.7 349 707
6 478 90.7 160.7 299 657
7 648 90.7 160.7 469 828
8 700 90.7 160.7 520 879
9 566 90.7 160.7 387 745
10 563 90.7 160.7 384 742
11 760 90.7 160.7 580 939
12 833 90.7 160.7 654 1012
session = Block 5:
step emmean SE df lower.CL upper.CL
1 1196 69.7 57.6 1056 1335
2 614 69.7 57.6 474 753
3 596 69.7 57.6 457 736
4 675 69.7 57.6 536 815
5 703 69.7 57.6 563 842
6 642 69.7 57.6 502 781
7 875 69.7 57.6 735 1014
8 715 69.7 57.6 575 854
9 831 69.7 57.6 692 971
10 649 69.7 57.6 510 789
11 736 69.7 57.6 597 876
12 722 69.7 57.6 583 862
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
--- Sequence length 18 ---
[Combined] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step * session * trial_acc_cat + (1 | subject)
Data: data
REML criterion at convergence: 160429.7
Scaled residuals:
Min 1Q Median 3Q Max
-2.048 -0.373 -0.145 0.102 41.159
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 29887 172.9
Residual 328654 573.3
Number of obs: 10368, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 811.48 59.25 13.695
step2 -365.53 60.77 -6.015
step3 -412.41 60.77 -6.787
step4 -420.02 60.77 -6.912
step5 -371.73 60.77 -6.117
step6 -354.89 60.77 -5.840
step7 -311.91 60.77 -5.133
step8 -365.18 60.77 -6.009
step9 -325.52 60.77 -5.357
step10 -329.51 60.77 -5.422
step11 -359.24 60.77 -5.912
step12 -394.91 60.77 -6.499
step13 -235.35 60.77 -3.873
step14 -261.72 60.77 -4.307
step15 -393.92 60.77 -6.482
step16 -415.43 60.77 -6.836
step17 -424.13 60.77 -6.980
step18 -354.54 60.77 -5.834
sessionBlock 5 37.15 69.26 0.536
trial_acc_catwrong -111.73 69.74 -1.602
step2:sessionBlock 5 9.68 97.78 0.099
step3:sessionBlock 5 70.97 97.78 0.726
step4:sessionBlock 5 109.66 97.78 1.121
step5:sessionBlock 5 78.91 97.78 0.807
step6:sessionBlock 5 64.38 97.78 0.658
step7:sessionBlock 5 88.42 97.78 0.904
step8:sessionBlock 5 29.70 97.78 0.304
step9:sessionBlock 5 222.83 97.78 2.279
step10:sessionBlock 5 192.98 97.78 1.974
step11:sessionBlock 5 -16.81 97.78 -0.172
step12:sessionBlock 5 34.57 97.78 0.354
step13:sessionBlock 5 257.23 97.78 2.631
step14:sessionBlock 5 36.57 97.78 0.374
step15:sessionBlock 5 -19.91 97.78 -0.204
step16:sessionBlock 5 38.46 97.78 0.393
step17:sessionBlock 5 93.07 97.78 0.952
step18:sessionBlock 5 -22.04 97.78 -0.225
step2:trial_acc_catwrong 80.78 98.33 0.822
step3:trial_acc_catwrong 157.29 98.33 1.600
step4:trial_acc_catwrong 241.25 98.33 2.454
step5:trial_acc_catwrong 226.18 98.33 2.300
step6:trial_acc_catwrong 154.77 98.33 1.574
step7:trial_acc_catwrong 247.11 98.33 2.513
step8:trial_acc_catwrong 313.71 98.33 3.190
step9:trial_acc_catwrong 389.47 98.33 3.961
step10:trial_acc_catwrong 184.42 98.33 1.876
step11:trial_acc_catwrong 215.97 98.33 2.196
step12:trial_acc_catwrong 232.25 98.33 2.362
step13:trial_acc_catwrong 109.75 98.33 1.116
step14:trial_acc_catwrong 110.25 98.33 1.121
step15:trial_acc_catwrong 212.77 98.33 2.164
step16:trial_acc_catwrong 305.60 98.33 3.108
step17:trial_acc_catwrong 373.30 98.33 3.796
step18:trial_acc_catwrong 642.13 98.33 6.531
sessionBlock 5:trial_acc_catwrong 477.49 98.47 4.849
step2:sessionBlock 5:trial_acc_catwrong -312.80 138.82 -2.253
step3:sessionBlock 5:trial_acc_catwrong -440.40 138.82 -3.172
step4:sessionBlock 5:trial_acc_catwrong -476.65 138.82 -3.434
step5:sessionBlock 5:trial_acc_catwrong -409.56 138.82 -2.950
step6:sessionBlock 5:trial_acc_catwrong -302.10 138.82 -2.176
step7:sessionBlock 5:trial_acc_catwrong -424.65 138.82 -3.059
step8:sessionBlock 5:trial_acc_catwrong -478.26 138.82 -3.445
step9:sessionBlock 5:trial_acc_catwrong -522.35 138.82 -3.763
step10:sessionBlock 5:trial_acc_catwrong -478.88 138.82 -3.450
step11:sessionBlock 5:trial_acc_catwrong -372.63 138.82 -2.684
step12:sessionBlock 5:trial_acc_catwrong -329.98 138.82 -2.377
step13:sessionBlock 5:trial_acc_catwrong -477.01 138.82 -3.436
step14:sessionBlock 5:trial_acc_catwrong -352.75 138.82 -2.541
step15:sessionBlock 5:trial_acc_catwrong -330.58 138.82 -2.381
step16:sessionBlock 5:trial_acc_catwrong -496.37 138.82 -3.576
step17:sessionBlock 5:trial_acc_catwrong -586.13 138.82 -4.222
step18:sessionBlock 5:trial_acc_catwrong -734.34 138.82 -5.290
Correlation matrix not shown by default, as p = 72 > 12.
Use print(summary(mdl_all), correlation=TRUE) or
vcov(summary(mdl_all)) if you need it
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 333.8643 17 < 2.2e-16 ***
session 139.6673 1 < 2.2e-16 ***
trial_acc_cat 152.2053 1 < 2.2e-16 ***
step:session 61.4012 17 6.159e-07 ***
step:trial_acc_cat 57.9133 17 2.313e-06 ***
session:trial_acc_cat 5.9238 1 0.014937 *
step:session:trial_acc_cat 39.5580 17 0.001493 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
session = Block 4, trial_acc_cat = correct:
step emmean SE df lower.CL upper.CL
1 811 59.3 73.0 693 930
2 446 59.3 73.0 328 564
3 399 59.3 73.0 281 517
4 391 59.3 73.0 273 510
5 440 59.3 73.0 322 558
6 457 59.3 73.0 338 575
7 500 59.3 73.0 381 618
8 446 59.3 73.0 328 564
9 486 59.3 73.0 368 604
10 482 59.3 73.0 364 600
11 452 59.3 73.0 334 570
12 417 59.3 73.0 298 535
13 576 59.3 73.0 458 694
14 550 59.3 73.0 432 668
15 418 59.3 73.0 299 536
16 396 59.3 73.0 278 514
17 387 59.3 73.0 269 505
18 457 59.3 73.0 339 575
session = Block 5, trial_acc_cat = correct:
step emmean SE df lower.CL upper.CL
1 849 67.9 125.3 714 983
2 493 67.9 125.3 358 627
3 507 67.9 125.3 373 642
4 538 67.9 125.3 404 673
5 556 67.9 125.3 421 690
6 558 67.9 125.3 424 692
7 625 67.9 125.3 491 759
8 513 67.9 125.3 379 647
9 746 67.9 125.3 612 880
10 712 67.9 125.3 578 846
11 473 67.9 125.3 338 607
12 488 67.9 125.3 354 623
13 870 67.9 125.3 736 1005
14 623 67.9 125.3 489 758
15 435 67.9 125.3 300 569
16 472 67.9 125.3 337 606
17 518 67.9 125.3 383 652
18 472 67.9 125.3 338 606
session = Block 4, trial_acc_cat = wrong:
step emmean SE df lower.CL upper.CL
1 700 68.3 128.2 565 835
2 415 68.3 128.2 280 550
3 445 68.3 128.2 310 580
4 521 68.3 128.2 386 656
5 554 68.3 128.2 419 689
6 500 68.3 128.2 365 635
7 635 68.3 128.2 500 770
8 648 68.3 128.2 513 783
9 764 68.3 128.2 629 899
10 555 68.3 128.2 420 690
11 556 68.3 128.2 421 692
12 537 68.3 128.2 402 672
13 574 68.3 128.2 439 709
14 548 68.3 128.2 413 683
15 519 68.3 128.2 384 654
16 590 68.3 128.2 455 725
17 649 68.3 128.2 514 784
18 987 68.3 128.2 852 1122
session = Block 5, trial_acc_cat = wrong:
step emmean SE df lower.CL upper.CL
1 1214 59.4 73.9 1096 1333
2 627 59.4 73.9 508 745
3 590 59.4 73.9 471 708
4 669 59.4 73.9 550 787
5 738 59.4 73.9 620 857
6 777 59.4 73.9 658 895
7 813 59.4 73.9 695 932
8 714 59.4 73.9 596 833
9 979 59.4 73.9 860 1097
10 783 59.4 73.9 665 902
11 682 59.4 73.9 563 800
12 756 59.4 73.9 638 875
13 869 59.4 73.9 751 987
14 747 59.4 73.9 628 865
15 683 59.4 73.9 564 801
16 647 59.4 73.9 528 765
17 670 59.4 73.9 552 789
18 746 59.4 73.9 627 864
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
[CORRECT] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step * session + (1 | subject)
Data: data
REML criterion at convergence: 76476
Scaled residuals:
Min 1Q Median 3Q Max
-1.8632 -0.3708 -0.1344 0.1428 26.8933
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 31013 176.1
Residual 142434 377.4
Number of obs: 5220, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 814.17 50.27 16.195
step2 -365.53 40.01 -9.137
step3 -412.41 40.01 -10.309
step4 -420.02 40.01 -10.499
step5 -371.73 40.01 -9.292
step6 -354.89 40.01 -8.871
step7 -311.91 40.01 -7.797
step8 -365.18 40.01 -9.128
step9 -325.52 40.01 -8.137
step10 -329.51 40.01 -8.237
step11 -359.24 40.01 -8.980
step12 -394.91 40.01 -9.872
step13 -235.35 40.01 -5.883
step14 -261.72 40.01 -6.542
step15 -393.92 40.01 -9.847
step16 -415.43 40.01 -10.385
step17 -424.13 40.01 -10.602
step18 -354.54 40.01 -8.862
sessionBlock 5 35.90 45.69 0.786
step2:sessionBlock 5 9.68 64.37 0.150
step3:sessionBlock 5 70.97 64.37 1.103
step4:sessionBlock 5 109.66 64.37 1.703
step5:sessionBlock 5 78.91 64.37 1.226
step6:sessionBlock 5 64.38 64.37 1.000
step7:sessionBlock 5 88.42 64.37 1.374
step8:sessionBlock 5 29.70 64.37 0.461
step9:sessionBlock 5 222.83 64.37 3.461
step10:sessionBlock 5 192.98 64.37 2.998
step11:sessionBlock 5 -16.81 64.37 -0.261
step12:sessionBlock 5 34.57 64.37 0.537
step13:sessionBlock 5 257.23 64.37 3.996
step14:sessionBlock 5 36.57 64.37 0.568
step15:sessionBlock 5 -19.91 64.37 -0.309
step16:sessionBlock 5 38.46 64.37 0.597
step17:sessionBlock 5 93.07 64.37 1.446
step18:sessionBlock 5 -22.04 64.37 -0.342
Correlation matrix not shown by default, as p = 36 > 12.
Use print(summary(mdl), correlation=TRUE) or
vcov(summary(mdl)) if you need it
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 385.211 17 < 2.2e-16 ***
session 86.629 1 < 2.2e-16 ***
step:session 54.881 17 7.161e-06 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
session = Block 4:
step emmean SE df lower.CL upper.CL
1 814 50.3 35.3 712 916
2 449 50.3 35.3 347 551
3 402 50.3 35.3 300 504
4 394 50.3 35.3 292 496
5 442 50.3 35.3 340 544
6 459 50.3 35.3 357 561
7 502 50.3 35.3 400 604
8 449 50.3 35.3 347 551
9 489 50.3 35.3 387 591
10 485 50.3 35.3 383 587
11 455 50.3 35.3 353 557
12 419 50.3 35.3 317 521
13 579 50.3 35.3 477 681
14 552 50.3 35.3 450 654
15 420 50.3 35.3 318 522
16 399 50.3 35.3 297 501
17 390 50.3 35.3 288 492
18 460 50.3 35.3 358 562
session = Block 5:
step emmean SE df lower.CL upper.CL
1 850 54.8 49.9 740 960
2 494 54.8 49.9 384 604
3 509 54.8 49.9 398 619
4 540 54.8 49.9 430 650
5 557 54.8 49.9 447 667
6 560 54.8 49.9 449 670
7 627 54.8 49.9 516 737
8 515 54.8 49.9 404 625
9 747 54.8 49.9 637 858
10 714 54.8 49.9 603 824
11 474 54.8 49.9 364 584
12 490 54.8 49.9 380 600
13 872 54.8 49.9 762 982
14 625 54.8 49.9 515 735
15 436 54.8 49.9 326 546
16 473 54.8 49.9 363 583
17 519 54.8 49.9 409 629
18 473 54.8 49.9 363 584
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95
Pairwise (Tukey) across steps — correct | session:
session = Block 4:
contrast estimate SE df t.ratio p.value
step1 - step2 365.528 40.0 5167 9.137 <.0001
step1 - step3 412.410 40.0 5167 10.309 <.0001
step1 - step4 420.022 40.0 5167 10.499 <.0001
step1 - step5 371.730 40.0 5167 9.292 <.0001
step1 - step6 354.888 40.0 5167 8.871 <.0001
step1 - step7 311.910 40.0 5167 7.797 <.0001
step1 - step8 365.180 40.0 5167 9.128 <.0001
step1 - step9 325.522 40.0 5167 8.137 <.0001
step1 - step10 329.511 40.0 5167 8.237 <.0001
step1 - step11 359.242 40.0 5167 8.980 <.0001
step1 - step12 394.910 40.0 5167 9.872 <.0001
step1 - step13 235.354 40.0 5167 5.883 <.0001
step1 - step14 261.719 40.0 5167 6.542 <.0001
step1 - step15 393.916 40.0 5167 9.847 <.0001
step1 - step16 415.433 40.0 5167 10.385 <.0001
step1 - step17 424.135 40.0 5167 10.602 <.0001
step1 - step18 354.539 40.0 5167 8.862 <.0001
step2 - step3 46.882 40.0 5167 1.172 0.9995
step2 - step4 54.494 40.0 5167 1.362 0.9969
step2 - step5 6.202 40.0 5167 0.155 1.0000
step2 - step6 -10.640 40.0 5167 -0.266 1.0000
step2 - step7 -53.618 40.0 5167 -1.340 0.9974
step2 - step8 -0.348 40.0 5167 -0.009 1.0000
step2 - step9 -40.006 40.0 5167 -1.000 0.9999
step2 - step10 -36.017 40.0 5167 -0.900 1.0000
step2 - step11 -6.287 40.0 5167 -0.157 1.0000
step2 - step12 29.382 40.0 5167 0.734 1.0000
step2 - step13 -130.174 40.0 5167 -3.254 0.1026
step2 - step14 -103.809 40.0 5167 -2.595 0.4574
step2 - step15 28.388 40.0 5167 0.710 1.0000
step2 - step16 49.904 40.0 5167 1.247 0.9989
step2 - step17 58.607 40.0 5167 1.465 0.9929
step2 - step18 -10.989 40.0 5167 -0.275 1.0000
step3 - step4 7.612 40.0 5167 0.190 1.0000
step3 - step5 -40.680 40.0 5167 -1.017 0.9999
step3 - step6 -57.523 40.0 5167 -1.438 0.9942
step3 - step7 -100.500 40.0 5167 -2.512 0.5209
step3 - step8 -47.230 40.0 5167 -1.181 0.9995
step3 - step9 -86.888 40.0 5167 -2.172 0.7724
step3 - step10 -82.899 40.0 5167 -2.072 0.8324
step3 - step11 -53.169 40.0 5167 -1.329 0.9977
step3 - step12 -17.500 40.0 5167 -0.437 1.0000
step3 - step13 -177.056 40.0 5167 -4.426 0.0013
step3 - step14 -150.691 40.0 5167 -3.767 0.0194
step3 - step15 -18.494 40.0 5167 -0.462 1.0000
step3 - step16 3.022 40.0 5167 0.076 1.0000
step3 - step17 11.725 40.0 5167 0.293 1.0000
step3 - step18 -57.871 40.0 5167 -1.447 0.9938
step4 - step5 -48.292 40.0 5167 -1.207 0.9993
step4 - step6 -65.135 40.0 5167 -1.628 0.9782
step4 - step7 -108.112 40.0 5167 -2.702 0.3787
step4 - step8 -54.843 40.0 5167 -1.371 0.9967
step4 - step9 -94.500 40.0 5167 -2.362 0.6369
step4 - step10 -90.511 40.0 5167 -2.263 0.7106
step4 - step11 -60.781 40.0 5167 -1.519 0.9894
step4 - step12 -25.112 40.0 5167 -0.628 1.0000
step4 - step13 -184.668 40.0 5167 -4.616 0.0006
step4 - step14 -158.303 40.0 5167 -3.957 0.0095
step4 - step15 -26.107 40.0 5167 -0.653 1.0000
step4 - step16 -4.590 40.0 5167 -0.115 1.0000
step4 - step17 4.112 40.0 5167 0.103 1.0000
step4 - step18 -65.483 40.0 5167 -1.637 0.9771
step5 - step6 -16.843 40.0 5167 -0.421 1.0000
step5 - step7 -59.820 40.0 5167 -1.495 0.9911
step5 - step8 -6.551 40.0 5167 -0.164 1.0000
step5 - step9 -46.208 40.0 5167 -1.155 0.9996
step5 - step10 -42.219 40.0 5167 -1.055 0.9999
step5 - step11 -12.489 40.0 5167 -0.312 1.0000
step5 - step12 23.180 40.0 5167 0.579 1.0000
step5 - step13 -136.376 40.0 5167 -3.409 0.0647
step5 - step14 -110.011 40.0 5167 -2.750 0.3459
step5 - step15 22.185 40.0 5167 0.555 1.0000
step5 - step16 43.702 40.0 5167 1.092 0.9998
step5 - step17 52.404 40.0 5167 1.310 0.9981
step5 - step18 -17.191 40.0 5167 -0.430 1.0000
step6 - step7 -42.977 40.0 5167 -1.074 0.9999
step6 - step8 10.292 40.0 5167 0.257 1.0000
step6 - step9 -29.365 40.0 5167 -0.734 1.0000
step6 - step10 -25.376 40.0 5167 -0.634 1.0000
step6 - step11 4.354 40.0 5167 0.109 1.0000
step6 - step12 40.023 40.0 5167 1.000 0.9999
step6 - step13 -119.534 40.0 5167 -2.988 0.2061
step6 - step14 -93.168 40.0 5167 -2.329 0.6620
step6 - step15 39.028 40.0 5167 0.976 1.0000
step6 - step16 60.545 40.0 5167 1.513 0.9898
step6 - step17 69.247 40.0 5167 1.731 0.9608
step6 - step18 -0.348 40.0 5167 -0.009 1.0000
step7 - step8 53.270 40.0 5167 1.332 0.9976
step7 - step9 13.612 40.0 5167 0.340 1.0000
step7 - step10 17.601 40.0 5167 0.440 1.0000
step7 - step11 47.331 40.0 5167 1.183 0.9995
step7 - step12 83.000 40.0 5167 2.075 0.8310
step7 - step13 -76.556 40.0 5167 -1.914 0.9070
step7 - step14 -50.191 40.0 5167 -1.255 0.9989
step7 - step15 82.006 40.0 5167 2.050 0.8445
step7 - step16 103.522 40.0 5167 2.588 0.4629
step7 - step17 112.225 40.0 5167 2.805 0.3096
step7 - step18 42.629 40.0 5167 1.066 0.9999
step8 - step9 -39.657 40.0 5167 -0.991 1.0000
step8 - step10 -35.669 40.0 5167 -0.892 1.0000
step8 - step11 -5.938 40.0 5167 -0.148 1.0000
step8 - step12 29.730 40.0 5167 0.743 1.0000
step8 - step13 -129.826 40.0 5167 -3.245 0.1051
step8 - step14 -103.461 40.0 5167 -2.586 0.4640
step8 - step15 28.736 40.0 5167 0.718 1.0000
step8 - step16 50.253 40.0 5167 1.256 0.9988
step8 - step17 58.955 40.0 5167 1.474 0.9924
step8 - step18 -10.640 40.0 5167 -0.266 1.0000
step9 - step10 3.989 40.0 5167 0.100 1.0000
step9 - step11 33.719 40.0 5167 0.843 1.0000
step9 - step12 69.388 40.0 5167 1.734 0.9601
step9 - step13 -90.168 40.0 5167 -2.254 0.7167
step9 - step14 -63.803 40.0 5167 -1.595 0.9823
step9 - step15 68.393 40.0 5167 1.710 0.9651
step9 - step16 89.910 40.0 5167 2.247 0.7213
step9 - step17 98.612 40.0 5167 2.465 0.5575
step9 - step18 29.017 40.0 5167 0.725 1.0000
step10 - step11 29.730 40.0 5167 0.743 1.0000
step10 - step12 65.399 40.0 5167 1.635 0.9773
step10 - step13 -94.157 40.0 5167 -2.354 0.6434
step10 - step14 -67.792 40.0 5167 -1.695 0.9679
step10 - step15 64.404 40.0 5167 1.610 0.9806
step10 - step16 85.921 40.0 5167 2.148 0.7878
step10 - step17 94.624 40.0 5167 2.365 0.6346
step10 - step18 25.028 40.0 5167 0.626 1.0000
step11 - step12 35.669 40.0 5167 0.892 1.0000
step11 - step13 -123.888 40.0 5167 -3.097 0.1572
step11 - step14 -97.522 40.0 5167 -2.438 0.5787
step11 - step15 34.674 40.0 5167 0.867 1.0000
step11 - step16 56.191 40.0 5167 1.405 0.9956
step11 - step17 64.893 40.0 5167 1.622 0.9790
step11 - step18 -4.702 40.0 5167 -0.118 1.0000
step12 - step13 -159.556 40.0 5167 -3.988 0.0084
step12 - step14 -133.191 40.0 5167 -3.329 0.0823
step12 - step15 -0.994 40.0 5167 -0.025 1.0000
step12 - step16 20.523 40.0 5167 0.513 1.0000
step12 - step17 29.225 40.0 5167 0.731 1.0000
step12 - step18 -40.371 40.0 5167 -1.009 0.9999
step13 - step14 26.365 40.0 5167 0.659 1.0000
step13 - step15 158.562 40.0 5167 3.964 0.0093
step13 - step16 180.079 40.0 5167 4.501 0.0010
step13 - step17 188.781 40.0 5167 4.719 0.0003
step13 - step18 119.185 40.0 5167 2.979 0.2105
step14 - step15 132.197 40.0 5167 3.305 0.0886
step14 - step16 153.714 40.0 5167 3.842 0.0147
step14 - step17 162.416 40.0 5167 4.060 0.0063
step14 - step18 92.820 40.0 5167 2.320 0.6685
step15 - step16 21.517 40.0 5167 0.538 1.0000
step15 - step17 30.219 40.0 5167 0.755 1.0000
step15 - step18 -39.376 40.0 5167 -0.984 1.0000
step16 - step17 8.702 40.0 5167 0.218 1.0000
step16 - step18 -60.893 40.0 5167 -1.522 0.9892
step17 - step18 -69.596 40.0 5167 -1.740 0.9589
session = Block 5:
contrast estimate SE df t.ratio p.value
step1 - step2 355.848 50.4 5167 7.056 <.0001
step1 - step3 341.438 50.4 5167 6.770 <.0001
step1 - step4 310.366 50.4 5167 6.154 <.0001
step1 - step5 292.821 50.4 5167 5.806 <.0001
step1 - step6 290.509 50.4 5167 5.760 <.0001
step1 - step7 223.491 50.4 5167 4.431 0.0013
step1 - step8 335.482 50.4 5167 6.652 <.0001
step1 - step9 102.696 50.4 5167 2.036 0.8516
step1 - step10 136.536 50.4 5167 2.707 0.3753
step1 - step11 376.054 50.4 5167 7.457 <.0001
step1 - step12 360.339 50.4 5167 7.145 <.0001
step1 - step13 -21.875 50.4 5167 -0.434 1.0000
step1 - step14 225.152 50.4 5167 4.464 0.0011
step1 - step15 413.830 50.4 5167 8.206 <.0001
step1 - step16 376.973 50.4 5167 7.475 <.0001
step1 - step17 331.062 50.4 5167 6.564 <.0001
step1 - step18 376.580 50.4 5167 7.467 <.0001
step2 - step3 -14.411 50.4 5167 -0.286 1.0000
step2 - step4 -45.482 50.4 5167 -0.902 1.0000
step2 - step5 -63.027 50.4 5167 -1.250 0.9989
step2 - step6 -65.339 50.4 5167 -1.296 0.9983
step2 - step7 -132.357 50.4 5167 -2.624 0.4353
step2 - step8 -20.366 50.4 5167 -0.404 1.0000
step2 - step9 -253.152 50.4 5167 -5.020 0.0001
step2 - step10 -219.312 50.4 5167 -4.349 0.0019
step2 - step11 20.205 50.4 5167 0.401 1.0000
step2 - step12 4.491 50.4 5167 0.089 1.0000
step2 - step13 -377.723 50.4 5167 -7.490 <.0001
step2 - step14 -130.696 50.4 5167 -2.591 0.4600
step2 - step15 57.982 50.4 5167 1.150 0.9996
step2 - step16 21.125 50.4 5167 0.419 1.0000
step2 - step17 -24.786 50.4 5167 -0.491 1.0000
step2 - step18 20.732 50.4 5167 0.411 1.0000
step3 - step4 -31.071 50.4 5167 -0.616 1.0000
step3 - step5 -48.616 50.4 5167 -0.964 1.0000
step3 - step6 -50.929 50.4 5167 -1.010 0.9999
step3 - step7 -117.946 50.4 5167 -2.339 0.6547
step3 - step8 -5.955 50.4 5167 -0.118 1.0000
step3 - step9 -238.741 50.4 5167 -4.734 0.0003
step3 - step10 -204.902 50.4 5167 -4.063 0.0063
step3 - step11 34.616 50.4 5167 0.686 1.0000
step3 - step12 18.902 50.4 5167 0.375 1.0000
step3 - step13 -363.312 50.4 5167 -7.204 <.0001
step3 - step14 -116.286 50.4 5167 -2.306 0.6793
step3 - step15 72.393 50.4 5167 1.435 0.9943
step3 - step16 35.536 50.4 5167 0.705 1.0000
step3 - step17 -10.375 50.4 5167 -0.206 1.0000
step3 - step18 35.143 50.4 5167 0.697 1.0000
step4 - step5 -17.545 50.4 5167 -0.348 1.0000
step4 - step6 -19.857 50.4 5167 -0.394 1.0000
step4 - step7 -86.875 50.4 5167 -1.723 0.9625
step4 - step8 25.116 50.4 5167 0.498 1.0000
step4 - step9 -207.670 50.4 5167 -4.118 0.0050
step4 - step10 -173.830 50.4 5167 -3.447 0.0575
step4 - step11 65.688 50.4 5167 1.302 0.9982
step4 - step12 49.973 50.4 5167 0.991 1.0000
step4 - step13 -332.241 50.4 5167 -6.588 <.0001
step4 - step14 -85.214 50.4 5167 -1.690 0.9688
step4 - step15 103.464 50.4 5167 2.052 0.8436
step4 - step16 66.607 50.4 5167 1.321 0.9979
step4 - step17 20.696 50.4 5167 0.410 1.0000
step4 - step18 66.214 50.4 5167 1.313 0.9980
step5 - step6 -2.312 50.4 5167 -0.046 1.0000
step5 - step7 -69.330 50.4 5167 -1.375 0.9966
step5 - step8 42.661 50.4 5167 0.846 1.0000
step5 - step9 -190.125 50.4 5167 -3.770 0.0192
step5 - step10 -156.286 50.4 5167 -3.099 0.1564
step5 - step11 83.232 50.4 5167 1.650 0.9751
step5 - step12 67.518 50.4 5167 1.339 0.9975
step5 - step13 -314.696 50.4 5167 -6.240 <.0001
step5 - step14 -67.670 50.4 5167 -1.342 0.9974
step5 - step15 121.009 50.4 5167 2.399 0.6084
step5 - step16 84.152 50.4 5167 1.669 0.9723
step5 - step17 38.241 50.4 5167 0.758 1.0000
step5 - step18 83.759 50.4 5167 1.661 0.9735
step6 - step7 -67.018 50.4 5167 -1.329 0.9977
step6 - step8 44.973 50.4 5167 0.892 1.0000
step6 - step9 -187.812 50.4 5167 -3.724 0.0226
step6 - step10 -153.973 50.4 5167 -3.053 0.1758
step6 - step11 85.545 50.4 5167 1.696 0.9676
step6 - step12 69.830 50.4 5167 1.385 0.9962
step6 - step13 -312.384 50.4 5167 -6.194 <.0001
step6 - step14 -65.357 50.4 5167 -1.296 0.9983
step6 - step15 123.321 50.4 5167 2.445 0.5729
step6 - step16 86.464 50.4 5167 1.714 0.9641
step6 - step17 40.554 50.4 5167 0.804 1.0000
step6 - step18 86.071 50.4 5167 1.707 0.9656
step7 - step8 111.991 50.4 5167 2.221 0.7400
step7 - step9 -120.795 50.4 5167 -2.395 0.6117
step7 - step10 -86.955 50.4 5167 -1.724 0.9622
step7 - step11 152.562 50.4 5167 3.025 0.1884
step7 - step12 136.848 50.4 5167 2.713 0.3710
step7 - step13 -245.366 50.4 5167 -4.865 0.0002
step7 - step14 1.661 50.4 5167 0.033 1.0000
step7 - step15 190.339 50.4 5167 3.774 0.0189
step7 - step16 153.482 50.4 5167 3.043 0.1801
step7 - step17 107.571 50.4 5167 2.133 0.7969
step7 - step18 153.089 50.4 5167 3.036 0.1836
step8 - step9 -232.786 50.4 5167 -4.616 0.0006
step8 - step10 -198.946 50.4 5167 -3.945 0.0100
step8 - step11 40.571 50.4 5167 0.804 1.0000
step8 - step12 24.857 50.4 5167 0.493 1.0000
step8 - step13 -357.357 50.4 5167 -7.086 <.0001
step8 - step14 -110.330 50.4 5167 -2.188 0.7621
step8 - step15 78.348 50.4 5167 1.554 0.9865
step8 - step16 41.491 50.4 5167 0.823 1.0000
step8 - step17 -4.420 50.4 5167 -0.088 1.0000
step8 - step18 41.098 50.4 5167 0.815 1.0000
step9 - step10 33.839 50.4 5167 0.671 1.0000
step9 - step11 273.357 50.4 5167 5.420 <.0001
step9 - step12 257.643 50.4 5167 5.109 <.0001
step9 - step13 -124.571 50.4 5167 -2.470 0.5536
step9 - step14 122.455 50.4 5167 2.428 0.5862
step9 - step15 311.134 50.4 5167 6.169 <.0001
step9 - step16 274.277 50.4 5167 5.438 <.0001
step9 - step17 228.366 50.4 5167 4.528 0.0009
step9 - step18 273.884 50.4 5167 5.431 <.0001
step10 - step11 239.518 50.4 5167 4.749 0.0003
step10 - step12 223.804 50.4 5167 4.438 0.0013
step10 - step13 -158.411 50.4 5167 -3.141 0.1400
step10 - step14 88.616 50.4 5167 1.757 0.9550
step10 - step15 277.295 50.4 5167 5.498 <.0001
step10 - step16 240.438 50.4 5167 4.767 0.0003
step10 - step17 194.527 50.4 5167 3.857 0.0139
step10 - step18 240.045 50.4 5167 4.760 0.0003
step11 - step12 -15.714 50.4 5167 -0.312 1.0000
step11 - step13 -397.929 50.4 5167 -7.890 <.0001
step11 - step14 -150.902 50.4 5167 -2.992 0.2041
step11 - step15 37.777 50.4 5167 0.749 1.0000
step11 - step16 0.920 50.4 5167 0.018 1.0000
step11 - step17 -44.991 50.4 5167 -0.892 1.0000
step11 - step18 0.527 50.4 5167 0.010 1.0000
step12 - step13 -382.214 50.4 5167 -7.579 <.0001
step12 - step14 -135.188 50.4 5167 -2.681 0.3943
step12 - step15 53.491 50.4 5167 1.061 0.9999
step12 - step16 16.634 50.4 5167 0.330 1.0000
step12 - step17 -29.277 50.4 5167 -0.581 1.0000
step12 - step18 16.241 50.4 5167 0.322 1.0000
step13 - step14 247.027 50.4 5167 4.898 0.0001
step13 - step15 435.705 50.4 5167 8.639 <.0001
step13 - step16 398.848 50.4 5167 7.909 <.0001
step13 - step17 352.938 50.4 5167 6.998 <.0001
step13 - step18 398.455 50.4 5167 7.901 <.0001
step14 - step15 188.679 50.4 5167 3.741 0.0213
step14 - step16 151.821 50.4 5167 3.010 0.1953
step14 - step17 105.911 50.4 5167 2.100 0.8166
step14 - step18 151.429 50.4 5167 3.003 0.1990
step15 - step16 -36.857 50.4 5167 -0.731 1.0000
step15 - step17 -82.768 50.4 5167 -1.641 0.9765
step15 - step18 -37.250 50.4 5167 -0.739 1.0000
step16 - step17 -45.911 50.4 5167 -0.910 1.0000
step16 - step18 -0.393 50.4 5167 -0.008 1.0000
step17 - step18 45.518 50.4 5167 0.903 1.0000
Degrees-of-freedom method: kenward-roger
P value adjustment: tukey method for comparing a family of 18 estimates
Adjacent-step contrasts (consecutive) — correct | session:
session = Block 4:
contrast estimate SE df t.ratio p.value
step2 - step1 -365.53 40.0 5167 -9.137 <.0001
step3 - step2 -46.88 40.0 5167 -1.172 1.0000
step4 - step3 -7.61 40.0 5167 -0.190 1.0000
step5 - step4 48.29 40.0 5167 1.207 1.0000
step6 - step5 16.84 40.0 5167 0.421 1.0000
step7 - step6 42.98 40.0 5167 1.074 1.0000
step8 - step7 -53.27 40.0 5167 -1.332 1.0000
step9 - step8 39.66 40.0 5167 0.991 1.0000
step10 - step9 -3.99 40.0 5167 -0.100 1.0000
step11 - step10 -29.73 40.0 5167 -0.743 1.0000
step12 - step11 -35.67 40.0 5167 -0.892 1.0000
step13 - step12 159.56 40.0 5167 3.988 0.0011
step14 - step13 -26.37 40.0 5167 -0.659 1.0000
step15 - step14 -132.20 40.0 5167 -3.305 0.0144
step16 - step15 -21.52 40.0 5167 -0.538 1.0000
step17 - step16 -8.70 40.0 5167 -0.218 1.0000
step18 - step17 69.60 40.0 5167 1.740 1.0000
session = Block 5:
contrast estimate SE df t.ratio p.value
step2 - step1 -355.85 50.4 5167 -7.056 <.0001
step3 - step2 14.41 50.4 5167 0.286 1.0000
step4 - step3 31.07 50.4 5167 0.616 1.0000
step5 - step4 17.54 50.4 5167 0.348 1.0000
step6 - step5 2.31 50.4 5167 0.046 1.0000
step7 - step6 67.02 50.4 5167 1.329 1.0000
step8 - step7 -111.99 50.4 5167 -2.221 0.2906
step9 - step8 232.79 50.4 5167 4.616 0.0001
step10 - step9 -33.84 50.4 5167 -0.671 1.0000
step11 - step10 -239.52 50.4 5167 -4.749 <.0001
step12 - step11 15.71 50.4 5167 0.312 1.0000
step13 - step12 382.21 50.4 5167 7.579 <.0001
step14 - step13 -247.03 50.4 5167 -4.898 <.0001
step15 - step14 -188.68 50.4 5167 -3.741 0.0022
step16 - step15 36.86 50.4 5167 0.731 1.0000
step17 - step16 45.91 50.4 5167 0.910 1.0000
step18 - step17 -45.52 50.4 5167 -0.903 1.0000
Degrees-of-freedom method: kenward-roger
P value adjustment: holm method for 17 tests
[WRONG] Model summary:
Linear mixed model fit by REML ['lmerMod']
Formula: rt ~ step * session + (1 | subject)
Data: data
REML criterion at convergence: 81979.1
Scaled residuals:
Min 1Q Median 3Q Max
-1.605 -0.411 -0.196 0.089 32.890
Random effects:
Groups Name Variance Std.Dev.
subject (Intercept) 33915 184.2
Residual 516197 718.5
Number of obs: 5148, groups: subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 705.95 81.32 8.682
step2 -284.75 96.88 -2.939
step3 -255.12 96.88 -2.633
step4 -178.77 96.88 -1.845
step5 -145.55 96.88 -1.502
step6 -200.12 96.88 -2.066
step7 -64.80 96.88 -0.669
step8 -51.47 96.88 -0.531
step9 63.95 96.88 0.660
step10 -145.09 96.88 -1.498
step11 -143.27 96.88 -1.479
step12 -162.66 96.88 -1.679
step13 -125.60 96.88 -1.296
step14 -151.47 96.88 -1.564
step15 -181.15 96.88 -1.870
step16 -109.84 96.88 -1.134
step17 -50.84 96.88 -0.525
step18 287.59 96.88 2.969
sessionBlock 5 507.34 87.60 5.791
step2:sessionBlock 5 -303.12 123.50 -2.455
step3:sessionBlock 5 -369.43 123.50 -2.991
step4:sessionBlock 5 -366.99 123.50 -2.972
step5:sessionBlock 5 -330.66 123.50 -2.677
step6:sessionBlock 5 -237.72 123.50 -1.925
step7:sessionBlock 5 -336.23 123.50 -2.723
step8:sessionBlock 5 -448.56 123.50 -3.632
step9:sessionBlock 5 -299.52 123.50 -2.425
step10:sessionBlock 5 -285.91 123.50 -2.315
step11:sessionBlock 5 -389.44 123.50 -3.153
step12:sessionBlock 5 -295.40 123.50 -2.392
step13:sessionBlock 5 -219.79 123.50 -1.780
step14:sessionBlock 5 -316.18 123.50 -2.560
step15:sessionBlock 5 -350.49 123.50 -2.838
step16:sessionBlock 5 -457.91 123.50 -3.708
step17:sessionBlock 5 -493.06 123.50 -3.993
step18:sessionBlock 5 -756.38 123.50 -6.125
Correlation matrix not shown by default, as p = 36 > 12.
Use print(summary(mdl), correlation=TRUE) or
vcov(summary(mdl)) if you need it
Type II Chi-square ANOVA:
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: rt
Chisq Df Pr(>Chisq)
step 135.053 17 < 2.2e-16 ***
session 54.066 1 1.939e-13 ***
step:session 49.136 17 5.747e-05 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
session = Block 4:
step emmean SE df lower.CL upper.CL
1 706 81.3 182.7 546 866
2 421 81.3 182.7 261 582
3 451 81.3 182.7 290 611
4 527 81.3 182.7 367 688
5 560 81.3 182.7 400 721
6 506 81.3 182.7 345 666
7 641 81.3 182.7 481 802
8 654 81.3 182.7 494 815
9 770 81.3 182.7 609 930
10 561 81.3 182.7 400 721
11 563 81.3 182.7 402 723
12 543 81.3 182.7 383 704
13 580 81.3 182.7 420 741
14 554 81.3 182.7 394 715
15 525 81.3 182.7 364 685
16 596 81.3 182.7 436 757
17 655 81.3 182.7 495 816
18 994 81.3 182.7 833 1154
session = Block 5:
step emmean SE df lower.CL upper.CL
1 1213 69.5 98.8 1075 1351
2 625 69.5 98.8 487 763
3 589 69.5 98.8 451 727
4 668 69.5 98.8 530 805
5 737 69.5 98.8 599 875
6 775 69.5 98.8 637 913
7 812 69.5 98.8 674 950
8 713 69.5 98.8 575 851
9 978 69.5 98.8 840 1116
10 782 69.5 98.8 644 920
11 681 69.5 98.8 543 819
12 755 69.5 98.8 617 893
13 868 69.5 98.8 730 1006
14 746 69.5 98.8 608 884
15 682 69.5 98.8 544 820
16 646 69.5 98.8 508 783
17 669 69.5 98.8 531 807
18 744 69.5 98.8 607 882
Degrees-of-freedom method: kenward-roger
Confidence level used: 0.95 # helper (kept)
if (!exists("safe_fit")) {
safe_fit <- function(formula, data) {
m <- tryCatch(lme4::lmer(formula, data = data), error = function(e) NULL)
if (is.null(m)) { message("⚠️ Falling back to lm."); m <- stats::lm(update(formula, . ~ . - (1 | subject)), data = data) }
m
}
}
# === Test phase (Blocks 4–5): Stepwise RT, CORRECT trials only ===
# Facets: rows = Familiar/Unfamiliar, cols = 6/12/18 (with adjacent-step stars)
plot_stepwise_rt_test_correct_only <- function(df_acc_all, show_plot = TRUE, alpha = 0.05) {
df_len0 <- df_acc_all %>%
dplyr::group_by(subject, session, trial) %>%
dplyr::mutate(seq_length_trial = dplyr::n_distinct(sub.trial.number)) %>%
dplyr::ungroup()
df_len <- df_len0 %>%
dplyr::filter(session %in% c(4, 5), trial_acc_cat == "correct") %>%
dplyr::mutate(
session_f = factor(session, levels = c(4, 5), labels = c("Familiar", "Unfamiliar")),
seq_length_trial = factor(seq_length_trial, levels = c(6, 12, 18))
)
if (nrow(df_len) == 0) { message("No correct-trial data for test phase (4–5)."); return(invisible(NULL)) }
lengths <- sort(unique(as.character(df_len$seq_length_trial)))
emms_list <- list(); stars_list <- list()
for (L in lengths) {
df_seq <- df_len %>%
dplyr::filter(as.character(seq_length_trial) == L) %>%
dplyr::group_by(subject, session_f, trial,
step = as.integer(as.character(sub.trial.number))) %>%
dplyr::summarise(rt = mean(feedback.RT, na.rm = TRUE), .groups = "drop") %>%
dplyr::mutate(step = factor(step, levels = sort(unique(step))))
if (nrow(df_seq) == 0) next
has_two_cond <- dplyr::n_distinct(df_seq$session_f) >= 2
has_two_step <- dplyr::n_distinct(df_seq$step) >= 2
if (!has_two_step) next
fmla <- if (has_two_cond) rt ~ step * session_f + (1 | subject) else rt ~ step + (1 | subject)
mdl <- safe_fit(fmla, data = df_seq)
# EMMs
if (has_two_cond) {
em <- emmeans::emmeans(mdl, ~ step | session_f) %>% as.data.frame()
} else {
em <- emmeans::emmeans(mdl, ~ step) %>% as.data.frame()
one_lab <- as.character(unique(df_seq$session_f))
em$session_f <- factor(rep(one_lab, nrow(em)), levels = c("Familiar","Unfamiliar"))
}
em$len_fac <- factor(paste0(L, "-step"), levels = c("6-step","12-step","18-step"))
em$step_num <- as.integer(as.character(em$step))
emms_list[[L]] <- em
# Adjacent-step tests (Holm) -> stars
em_grid <- if (has_two_cond) emmeans::emmeans(mdl, ~ step | session_f) else emmeans::emmeans(mdl, ~ step)
consec <- as.data.frame(emmeans::contrast(em_grid, method = "consec", adjust = "holm"))
nums <- stringr::str_extract_all(consec$contrast, "\\d+")
n1 <- vapply(nums, function(v) as.integer(v[1]), integer(1))
n2 <- vapply(nums, function(v) as.integer(v[2]), integer(1))
prev <- pmin(n1, n2); nextv <- pmax(n1, n2)
sign_adj <- ifelse(nextv == n1, 1, -1)
est_adj <- consec$estimate * sign_adj
st <- consec %>%
dplyr::mutate(step_prev = prev, step_next = nextv,
estimate_adj = est_adj,
sig = !is.na(p.value) & p.value < alpha) %>%
dplyr::filter(sig) %>%
dplyr::transmute(
len_fac = factor(paste0(L,"-step"), levels = c("6-step","12-step","18-step")),
step_num = step_prev,
session_f = if ("session_f" %in% names(.))
as.character(session_f)
else
rep(as.character(unique(df_seq$session_f)), dplyr::n()),
direction = ifelse(estimate_adj > 0, "up", "down"),
p.value = p.value
) %>%
dplyr::mutate(session_f = factor(session_f, levels = c("Familiar","Unfamiliar")))
if (nrow(st)) stars_list[[length(stars_list)+1]] <- st
}
if (!length(emms_list)) { message("No estimable EMMs for test-phase correct trials."); return(invisible(NULL)) }
df_em <- dplyr::bind_rows(emms_list) %>%
dplyr::mutate(
session_f = factor(session_f, levels = c("Familiar","Unfamiliar")),
len_fac = factor(len_fac, levels = c("6-step","12-step","18-step"))
)
df_stars <- if (length(stars_list)) dplyr::bind_rows(stars_list) else
dplyr::tibble(
len_fac = factor(character(), levels = levels(df_em$len_fac)),
step_num = integer(),
session_f = factor(character(), levels = levels(df_em$session_f)),
direction = character(),
p.value = numeric()
)
# y-positions for stars
if (nrow(df_stars)) {
df_stars <- df_stars %>%
dplyr::left_join(df_em %>% dplyr::select(len_fac, session_f, step_num, emmean, SE),
by = c("len_fac","session_f","step_num")) %>%
dplyr::mutate(y_star = emmean + pmax(SE * 1.2, 35))
}
df_stars_up <- df_stars %>% dplyr::filter(direction == "up")
df_stars_down <- df_stars %>% dplyr::filter(direction == "down")
pal <- c("6-step"="#B22222","12-step"="#2E7D32","18-step"="#1E3A8A")
# vertical separators
sep_df <- expand.grid(
session_f = levels(df_em$session_f),
len_fac = c("6-step","12-step"),
KEEP.OUT.ATTRS = FALSE, stringsAsFactors = FALSE
) %>%
dplyr::mutate(
session_f = factor(session_f, levels = levels(df_em$session_f)),
len_fac = factor(len_fac, levels = levels(df_em$len_fac))
)
p <- ggplot2::ggplot(df_em, ggplot2::aes(x = step_num, y = emmean)) +
ggplot2::geom_segment(
data = sep_df,
ggplot2::aes(x = Inf, xend = Inf, y = -Inf, yend = Inf),
inherit.aes = FALSE, linetype = "dotted"
) +
ggplot2::geom_ribbon(
ggplot2::aes(ymin = emmean - SE, ymax = emmean + SE, fill = len_fac),
alpha = 0.15, color = NA
) +
ggplot2::geom_line(ggplot2::aes(color = len_fac), linewidth = 0.9) +
ggplot2::geom_point(ggplot2::aes(fill = len_fac), shape = 21, color = "black", size = 2, stroke = 0.5) +
ggplot2::geom_text(data = df_stars_up,
ggplot2::aes(x = step_num, y = y_star, label = "*"),
inherit.aes = FALSE, size = 4.2, fontface = "bold", color = "#D32F2F") +
ggplot2::geom_text(data = df_stars_down,
ggplot2::aes(x = step_num, y = y_star, label = "*"),
inherit.aes = FALSE, size = 4.2, fontface = "bold", color = "#1F6FEB") +
ggplot2::facet_grid(rows = ggplot2::vars(session_f), cols = ggplot2::vars(len_fac),
scales = "free_x", space = "free_x", drop = FALSE) +
ggplot2::scale_x_continuous(
breaks = function(lims) seq(floor(lims[1]), ceiling(lims[2]), by = 1),
expand = ggplot2::expansion(mult = c(0.02, 0.06))
) +
ggplot2::scale_y_continuous(expand = ggplot2::expansion(mult = c(0.05, 0.12))) +
ggplot2::scale_color_manual(values = pal) +
ggplot2::scale_fill_manual(values = pal) +
ggplot2::labs(
title = "RT Comparison for each step across Difficulty levels (Test phase)",
x = "Step", y = "Estimated RT (ms)",
caption = "* Following step is significantly different (red = ↑, blue = ↓); Holm-adjusted adjacent tests."
) +
ggplot2::guides(color = "none", fill = "none") +
ggplot2::theme_minimal(base_size = 12) +
ggplot2::theme(panel.grid = ggplot2::element_blank(),
plot.caption = ggplot2::element_text(hjust = 0))
if (show_plot) print(p)
invisible(list(emms = df_em, stars = df_stars, plot = p, palette = pal))
}
# === Test phase (Blocks 4–5): Stepwise RT — Familiar & Unfamiliar overlaid ===
# Facets: cols = 6/12/18 (with adjacent-step stars)
plot_stepwise_rt_test_correct_combined <- function(df_acc_all, show_plot = TRUE, alpha = 0.05) {
df_len0 <- df_acc_all %>%
dplyr::group_by(subject, session, trial) %>%
dplyr::mutate(seq_length_trial = dplyr::n_distinct(sub.trial.number)) %>%
dplyr::ungroup()
df_len <- df_len0 %>%
dplyr::filter(session %in% c(4, 5), trial_acc_cat == "correct") %>%
dplyr::mutate(
session_f = factor(session, levels = c(4, 5), labels = c("Familiar", "Unfamiliar")),
seq_length_trial = factor(seq_length_trial, levels = c(6, 12, 18))
)
if (nrow(df_len) == 0) { message("No correct-trial data for test phase (4–5)."); return(invisible(NULL)) }
lengths <- sort(unique(as.character(df_len$seq_length_trial)))
emms_list <- list(); stars_list <- list()
for (L in lengths) {
df_seq <- df_len %>%
dplyr::filter(as.character(seq_length_trial) == L) %>%
dplyr::group_by(
subject, session_f, trial,
step = as.integer(as.character(sub.trial.number))
) %>%
dplyr::summarise(rt = mean(feedback.RT, na.rm = TRUE), .groups = "drop") %>%
dplyr::mutate(step = factor(step, levels = sort(unique(step))))
if (nrow(df_seq) == 0) next
has_two_cond <- dplyr::n_distinct(df_seq$session_f) >= 2
has_two_step <- dplyr::n_distinct(df_seq$step) >= 2
if (!has_two_step) next
fmla <- if (has_two_cond) rt ~ step * session_f + (1 | subject) else rt ~ step + (1 | subject)
mdl <- safe_fit(fmla, data = df_seq)
if (has_two_cond) {
em <- emmeans::emmeans(mdl, ~ step | session_f) %>% as.data.frame()
} else {
em <- emmeans::emmeans(mdl, ~ step) %>% as.data.frame()
one_lab <- as.character(unique(df_seq$session_f))
em$session_f <- factor(rep(one_lab, nrow(em)), levels = c("Familiar","Unfamiliar"))
}
em$len_fac <- factor(paste0(L, "-step"), levels = c("6-step","12-step","18-step"))
em$step_num <- as.integer(as.character(em$step))
emms_list[[L]] <- em
em_grid <- if (has_two_cond) emmeans::emmeans(mdl, ~ step | session_f) else emmeans::emmeans(mdl, ~ step)
consec <- as.data.frame(emmeans::contrast(em_grid, method = "consec", adjust = "holm"))
nums <- stringr::str_extract_all(consec$contrast, "\\d+")
n1 <- vapply(nums, function(v) as.integer(v[1]), integer(1))
n2 <- vapply(nums, function(v) as.integer(v[2]), integer(1))
prev <- pmin(n1, n2); nextv <- pmax(n1, n2)
sign_adj <- ifelse(nextv == n1, 1, -1)
est_adj <- consec$estimate * sign_adj
st <- consec %>%
dplyr::mutate(step_prev = prev, step_next = nextv,
estimate_adj = est_adj,
sig = !is.na(p.value) & p.value < alpha) %>%
dplyr::filter(sig) %>%
dplyr::transmute(
len_fac = factor(paste0(L,"-step"), levels = c("6-step","12-step","18-step")),
step_num = step_prev,
session_f = if ("session_f" %in% names(.))
as.character(session_f)
else
rep(as.character(unique(df_seq$session_f)), dplyr::n()),
direction = ifelse(estimate_adj > 0, "up", "down"),
p.value = p.value
) %>%
dplyr::mutate(session_f = factor(session_f, levels = c("Familiar","Unfamiliar")))
if (nrow(st)) stars_list[[length(stars_list)+1]] <- st
}
df_em <- dplyr::bind_rows(emms_list) %>%
dplyr::mutate(
len_fac = factor(len_fac, levels = c("6-step","12-step","18-step")),
session_f = factor(session_f, levels = c("Familiar","Unfamiliar"))
)
df_stars <- if (length(stars_list)) dplyr::bind_rows(stars_list) else
dplyr::tibble(
len_fac = factor(character(), levels = levels(df_em$len_fac)),
step_num = integer(),
session_f = factor(character(), levels = levels(df_em$session_f)),
direction = character(),
p.value = numeric()
)
if (nrow(df_stars)) {
df_stars <- df_stars %>%
dplyr::left_join(df_em %>% dplyr::select(len_fac, session_f, step_num, emmean, SE),
by = c("len_fac","session_f","step_num")) %>%
dplyr::mutate(y_star = emmean + pmax(SE * 1.2, 35))
}
df_stars_up <- df_stars %>% dplyr::filter(direction == "up")
df_stars_down <- df_stars %>% dplyr::filter(direction == "down")
pal_cond <- c("Familiar"="#1E3A8A", "Unfamiliar"="#C2410C")
sep_df <- data.frame(len_fac = factor(c("6-step","12-step"),
levels = c("6-step","12-step","18-step")))
p <- ggplot2::ggplot(df_em, ggplot2::aes(x = step_num, y = emmean)) +
ggplot2::geom_segment(
data = sep_df,
ggplot2::aes(x = Inf, xend = Inf, y = -Inf, yend = Inf),
inherit.aes = FALSE, linetype = "dotted"
) +
ggplot2::geom_ribbon(
ggplot2::aes(ymin = emmean - SE, ymax = emmean + SE, fill = session_f, group = session_f),
alpha = 0.15, color = NA
) +
ggplot2::geom_line(ggplot2::aes(color = session_f, group = session_f), linewidth = 0.9) +
ggplot2::geom_point(ggplot2::aes(fill = session_f), shape = 21, color = "black", size = 2, stroke = 0.5) +
ggplot2::geom_text(data = df_stars_up,
ggplot2::aes(x = step_num, y = y_star, label = "*"),
inherit.aes = FALSE, size = 4.2, fontface = "bold", color = "#D32F2F") +
ggplot2::geom_text(data = df_stars_down,
ggplot2::aes(x = step_num, y = y_star, label = "*"),
inherit.aes = FALSE, size = 4.2, fontface = "bold", color = "#1F6FEB") +
ggplot2::facet_grid(cols = ggplot2::vars(len_fac), scales = "free_x", space = "free_x", drop = FALSE) +
ggplot2::scale_x_continuous(
breaks = function(lims) seq(floor(lims[1]), ceiling(lims[2]), by = 1),
expand = ggplot2::expansion(mult = c(0.02, 0.06))
) +
ggplot2::scale_y_continuous(expand = ggplot2::expansion(mult = c(0.05, 0.12))) +
ggplot2::scale_color_manual(values = pal_cond, name = NULL) +
ggplot2::scale_fill_manual(values = pal_cond, name = NULL) +
ggplot2::labs(
title = "RT for each Step across Difficulty levels (Test phase)",
x = "Steps", y = "Estimated RT (ms)",
caption = "* Following step is significantly different (red = ↑, blue = ↓)"
) +
ggplot2::theme_minimal(base_size = 12) +
ggplot2::theme(
panel.grid = ggplot2::element_blank(),
legend.position = "bottom",
legend.box = "horizontal",
plot.caption = ggplot2::element_text(hjust = 0)
)
if (show_plot) print(p)
invisible(list(emms = df_em, stars = df_stars, plot = p, palette = pal_cond))
}
plot_stepwise_rt_test_correct_only(df_acc_base)Warning in ggplot2::geom_segment(data = sep_df, ggplot2::aes(x = Inf, xend = Inf, : All aesthetics have length 1, but the data has 4 rows.
ℹ Please consider using `annotate()` or provide this layer with data containing
a single row.
plot_stepwise_rt_test_correct_combined(df_acc_base)Warning in ggplot2::geom_segment(data = sep_df, ggplot2::aes(x = Inf, xend = Inf, : All aesthetics have length 1, but the data has 2 rows.
ℹ Please consider using `annotate()` or provide this layer with data containing
a single row.