## [1] "Data Loaded. Checking Survey Labels:"
##
## Baseline Post-Training Early Learning Late Learning
## 48 48 48 48
# ==============================================================================
# DEMOGRAPHICS: Footedness Check
# ==============================================================================
# Load demographics
library(readxl)
demo_df <- read_excel("Demographics.xlsx")
View(demo_df)
# Filter for the 3 main groups
demo_df <- demo_df %>%
filter(Group %in% c("CON", "FAM", "OMM")) %>%
mutate(Group = factor(Group, levels = c("CON", "FAM", "OMM")))
# Create a Contingency Table
foot_table <- table(demo_df$Group, demo_df$Footedness)
print("--- Footedness Counts by Group ---")
## [1] "--- Footedness Counts by Group ---"
print(foot_table)
##
## Comfortable with both feet Left-foot Right-foot
## CON 1 3 12
## FAM 2 1 13
## OMM 0 3 13
# Run Chi-Square Test
print("--- Chi-Square Test for Footedness ---")
## [1] "--- Chi-Square Test for Footedness ---"
foot_chi <- chisq.test(foot_table)
## Warning in chisq.test(foot_table): Chi-squared approximation may be incorrect
print(foot_chi)
##
## Pearson's Chi-squared test
##
## data: foot_table
## X-squared = 3.1955, df = 4, p-value = 0.5257
# 1. Gaming Hours Per Week (ANOVA)
print("--- ANOVA: Gaming Hours Per Week ---")
## [1] "--- ANOVA: Gaming Hours Per Week ---"
m_gaming_hours <- lm(GamingHoursPerWeek ~ Group, data = demo_df)
print(anova(m_gaming_hours))
## Analysis of Variance Table
##
## Response: GamingHoursPerWeek
## Df Sum Sq Mean Sq F value Pr(>F)
## Group 2 100.5 50.250 0.9325 0.401
## Residuals 45 2424.8 53.885
# 2. Gaming Level Self-Rating (Chi-Square)
print("--- Chi-Square Test: Gaming Level Self-Rating ---")
## [1] "--- Chi-Square Test: Gaming Level Self-Rating ---"
gaming_table <- table(demo_df$Group, demo_df$GamingLevelSelfRating)
print(gaming_table)
##
## Advanced Beginner Complete beginner/I do not game Expert Intermediate
## CON 4 1 4 0 7
## FAM 3 4 6 1 2
## OMM 5 2 6 0 3
gaming_chi <- chisq.test(gaming_table)
## Warning in chisq.test(gaming_table): Chi-squared approximation may be incorrect
print(gaming_chi)
##
## Pearson's Chi-squared test
##
## data: gaming_table
## X-squared = 8.5, df = 8, p-value = 0.3862
# Optional: Kruskal-Wallis for Hours (Robust to Skew/Outliers)
kruskal.test(GamingHoursPerWeek ~ Group, data = demo_df)
##
## Kruskal-Wallis rank sum test
##
## data: GamingHoursPerWeek by Group
## Kruskal-Wallis chi-squared = 1.0407, df = 2, p-value = 0.5943
# ==============================================================================
# 2. ACCURACY ANALYSIS (Linear "Lukas" Method - Consistent)
# ==============================================================================
# --- A. Learning Phase (Blocks 1-6) ---
ACC_Learn <- trial_df %>% filter(Session <= 6)
ACC_Learn$Subject <- factor(ACC_Learn$Subject)
ACC_Learn$Group <- factor(ACC_Learn$Group, levels = c("CON", "FAM", "OMM"))
ACC_Learn$Session <- factor(ACC_Learn$Session, levels=c('1', '2', '3', '4', '5', '6'))
ACC_Learn$feedback.ACC.trial <- as.numeric(ACC_Learn$feedback.ACC.trial)
print("--- MODEL: ACCURACY LEARNING (Linear LMER) ---")
## [1] "--- MODEL: ACCURACY LEARNING (Linear LMER) ---"
accuracymodel <- lmer(feedback.ACC.trial ~ Session * Group + (1|Subject),
data=ACC_Learn, REML = FALSE)
print(Anova(accuracymodel, type="II"))
## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: feedback.ACC.trial
## Chisq Df Pr(>Chisq)
## Session 551.3683 5 < 2.2e-16 ***
## Group 0.4717 2 0.789895
## Session:Group 24.7158 10 0.005911 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# --- B. Test Phase (Retention & Transfer) ---
ACC_Test <- trial_df %>% filter(Session > 6)
ACC_Test$Subject <- factor(ACC_Test$Subject)
ACC_Test$Group <- factor(ACC_Test$Group, levels = c("CON", "FAM", "OMM"))
ACC_Test$Condition <- factor(ACC_Test$Condition, levels=c("Retention", "Transfer"))
ACC_Test$feedback.ACC.trial <- as.numeric(ACC_Test$feedback.ACC.trial)
print("--- MODEL: ACCURACY TEST (Linear LMER) ---")
## [1] "--- MODEL: ACCURACY TEST (Linear LMER) ---"
accuracyTestModel <- lmer(feedback.ACC.trial ~ Condition * Group + (1|Subject),
data=ACC_Test, REML = FALSE)
print(Anova(accuracyTestModel, type="II"))
## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: feedback.ACC.trial
## Chisq Df Pr(>Chisq)
## Condition 71.5140 1 <2e-16 ***
## Group 0.5147 2 0.7731
## Condition:Group 0.8815 2 0.6435
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# --- C. Pairwise T-Tests (Learning Blocks) ---
run_pairwise_ttests_acc <- function(data_subset, value_col) {
groups <- c("OMM", "FAM", "CON")
pairs <- combn(groups, 2, simplify = FALSE)
results <- list()
for(pair in pairs) {
g1 <- pair[1]; g2 <- pair[2]
d1 <- data_subset %>% filter(Group == g1)
d2 <- data_subset %>% filter(Group == g2)
if(nrow(d1) > 1 & nrow(d2) > 1) {
t_res <- t.test(d1[[value_col]], d2[[value_col]])
results[[paste(g1, "vs", g2)]] <- data.frame(
Comparison = paste(g1, "vs", g2),
p_value = t_res$p.value,
Mean_G1 = mean(d1[[value_col]]),
Mean_G2 = mean(d2[[value_col]]),
Diff_Raw = mean(d1[[value_col]]) - mean(d2[[value_col]])
)
}
}
bind_rows(results)
}
print("--- T-TESTS: Accuracy Learning (Block 1-6) ---")
## [1] "--- T-TESTS: Accuracy Learning (Block 1-6) ---"
for(b in 1:6) {
block_res <- run_pairwise_ttests_acc(ACC_Learn %>% filter(Session == b), "feedback.ACC.trial")
print(paste(">> BLOCK", b))
print(block_res %>% mutate(p_value = round(p_value, 4), Diff_Raw = round(Diff_Raw, 3)))
}
## [1] ">> BLOCK 1"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.0910 0.6562500 0.6966146 -0.040
## 2 OMM vs CON 0.0049 0.6562500 0.7226562 -0.066
## 3 FAM vs CON 0.2612 0.6966146 0.7226562 -0.026
## [1] ">> BLOCK 2"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.1718 0.8906250 0.9114583 -0.021
## 2 OMM vs CON 0.1030 0.8906250 0.8632812 0.027
## 3 FAM vs CON 0.0028 0.9114583 0.8632812 0.048
## [1] ">> BLOCK 3"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.7500 0.8867188 0.8815104 0.005
## 2 OMM vs CON 0.3107 0.8867188 0.8697917 0.017
## 3 FAM vs CON 0.4868 0.8815104 0.8697917 0.012
## [1] ">> BLOCK 4"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.8785 0.8736979 0.8710938 0.003
## 2 OMM vs CON 0.2670 0.8736979 0.8919271 -0.018
## 3 FAM vs CON 0.2068 0.8710938 0.8919271 -0.021
## [1] ">> BLOCK 5"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.0554 0.8567708 0.8893229 -0.033
## 2 OMM vs CON 0.5543 0.8567708 0.8671875 -0.010
## 3 FAM vs CON 0.1849 0.8893229 0.8671875 0.022
## [1] ">> BLOCK 6"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.0323 0.8424479 0.8802083 -0.038
## 2 OMM vs CON 0.5224 0.8424479 0.8541667 -0.012
## 3 FAM vs CON 0.1328 0.8802083 0.8541667 0.026
print("--- T-TESTS: Accuracy Test Phase (Retention & Transfer) ---")
## [1] "--- T-TESTS: Accuracy Test Phase (Retention & Transfer) ---"
for(cond in c("Retention", "Transfer")) {
test_res <- run_pairwise_ttests_acc(ACC_Test %>% filter(Condition == cond), "feedback.ACC.trial")
print(paste(">> CONDITION:", cond))
print(test_res %>% mutate(p_value = round(p_value, 4), Diff_Raw = round(Diff_Raw, 3)))
}
## [1] ">> CONDITION: Retention"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.9431 0.8489583 0.8502604 -0.001
## 2 OMM vs CON 0.1609 0.8489583 0.8736979 -0.025
## 3 FAM vs CON 0.1832 0.8502604 0.8736979 -0.023
## [1] ">> CONDITION: Transfer"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.2849 0.7682292 0.7447917 0.023
## 2 OMM vs CON 0.4630 0.7682292 0.7838542 -0.016
## 3 FAM vs CON 0.0714 0.7447917 0.7838542 -0.039
# --- D. Step-Level Accuracy (Checking for Speed-Accuracy Trade-offs) ---
print("--- MODEL: STEP-LEVEL ACCURACY (Aggregated LMM) ---")
## [1] "--- MODEL: STEP-LEVEL ACCURACY (Aggregated LMM) ---"
# 1. Calculate % Accuracy per Subject per Step (collapsing across blocks)
acc_step_summary <- foot_df %>%
filter(Session <= 6) %>%
group_by(Subject, Group, foot.position) %>%
summarise(Step_Accuracy = mean(feedback.ACC, na.rm = TRUE), .groups = "drop")
# 2. Run a standard Linear Mixed Model (lmer)
m_acc_step_linear <- lmer(Step_Accuracy ~ Group * as.factor(foot.position) + (1|Subject),
data = acc_step_summary)
print(Anova(m_acc_step_linear, type="II"))
## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Step_Accuracy
## Chisq Df Pr(>Chisq)
## Group 0.2330 2 0.8900
## as.factor(foot.position) 288.4034 5 <2e-16 ***
## Group:as.factor(foot.position) 6.1003 10 0.8068
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# ==============================================================================
# 3. FILTERING LOGIC (For RT Models)
# ==============================================================================
# --- A. ACCURACY ONLY (For Stats - LogRT) ---
apply_filtering_acc_only <- function(df, target_col) {
valid_df <- df %>% filter(feedback.ACC.trial == 1, .data[[target_col]] > 0)
valid_df$Log_RT <- log(valid_df[[target_col]])
return(valid_df)
}
df_learn <- apply_filtering_acc_only(trial_df %>% filter(Condition == "Learning"), "feedback.RT.mean")
df_test <- apply_filtering_acc_only(trial_df %>% filter(Condition %in% c("Retention", "Transfer")), "feedback.RT.mean")
df_foot_clean <- apply_filtering_acc_only(foot_df, "feedback.RT")
# --- B. 2.5 SD CLEANING (For Plots & % Improvement) ---
apply_filtering_plot <- function(df, target_col) {
valid_df <- df %>% filter(feedback.ACC.trial == 1, .data[[target_col]] > 0)
valid_df$Log_RT <- log(valid_df[[target_col]])
mean_log <- mean(valid_df$Log_RT)
sd_log <- sd(valid_df$Log_RT)
final_df <- valid_df %>%
filter(Log_RT >= (mean_log - 2.5 * sd_log) & Log_RT <= (mean_log + 2.5 * sd_log))
return(final_df)
}
df_learn_plot <- apply_filtering_plot(trial_df %>% filter(Condition == "Learning"), "feedback.RT.mean")
df_test_plot <- apply_filtering_plot(trial_df %>% filter(Condition %in% c("Retention", "Transfer")), "feedback.RT.mean")
df_foot_plot <- apply_filtering_plot(foot_df, "feedback.RT")
# ==============================================================================
# 4. STATISTICAL MODELS (RT ONLY)
# ==============================================================================
# --- A. RT Learning (LMM) ---
print("--- MODEL: RT LEARNING (Slope) ---")
## [1] "--- MODEL: RT LEARNING (Slope) ---"
m_rt_learn <- lmer(Log_RT ~ Group * as.numeric(Session) + (1|Subject), data = df_learn)
print(Anova(m_rt_learn, type="II"))
## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Log_RT
## Chisq Df Pr(>Chisq)
## Group 0.2357 2 0.8888
## as.numeric(Session) 3785.6317 1 < 2.2e-16 ***
## Group:as.numeric(Session) 32.2258 2 1.005e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# --- B. RT Test (LMM) ---
print("--- MODEL: RT TEST (Retention vs Transfer) ---")
## [1] "--- MODEL: RT TEST (Retention vs Transfer) ---"
m_rt_test <- lmer(Log_RT ~ Group * Condition + (1|Subject), data = df_test)
print(Anova(m_rt_test, type="II"))
## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Log_RT
## Chisq Df Pr(>Chisq)
## Group 0.2312 2 0.8908276
## Condition 746.9451 1 < 2.2e-16 ***
## Group:Condition 17.7158 2 0.0001423 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# --- C. Concatenation Learning (LMM - Step Level) ---
print("--- MODEL: CONCATENATION LEARNING (Group x Block x Step) ---")
## [1] "--- MODEL: CONCATENATION LEARNING (Group x Block x Step) ---"
m_foot_learn <- lmer(Log_RT ~ Group * as.factor(Session) * as.factor(foot.position) + (1|Subject),
data = df_foot_clean %>% filter(Session <= 6))
print(Anova(m_foot_learn, type="II"))
## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Log_RT
## Chisq Df Pr(>Chisq)
## Group 0.4362 2 0.8041
## as.factor(Session) 13721.1997 5 <2e-16 ***
## as.factor(foot.position) 17995.4985 5 <2e-16 ***
## Group:as.factor(Session) 194.3912 10 <2e-16 ***
## Group:as.factor(foot.position) 302.4145 10 <2e-16 ***
## as.factor(Session):as.factor(foot.position) 192.4877 25 <2e-16 ***
## Group:as.factor(Session):as.factor(foot.position) 228.4567 50 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# --- D. Concatenation Test (LMM - Step Level) ---
print("--- MODEL: CONCATENATION TEST (Group x Condition x Step) ---")
## [1] "--- MODEL: CONCATENATION TEST (Group x Condition x Step) ---"
m_foot_test <- lmer(Log_RT ~ Group * Panel_Label * as.factor(foot.position) + (1|Subject),
data = df_foot_clean %>% filter(Session > 6))
print(Anova(m_foot_test, type="II"))
## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Log_RT
## Chisq Df Pr(>Chisq)
## Group 0.2763 2 0.8709
## Panel_Label 1672.2449 1 < 2.2e-16 ***
## as.factor(foot.position) 6979.9902 5 < 2.2e-16 ***
## Group:Panel_Label 39.3200 2 2.896e-09 ***
## Group:as.factor(foot.position) 97.5077 10 < 2.2e-16 ***
## Panel_Label:as.factor(foot.position) 78.6621 5 1.598e-15 ***
## Group:Panel_Label:as.factor(foot.position) 44.1915 10 3.041e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# ==============================================================================
# 5. EXPLORATORY T-TESTS (RT)
# ==============================================================================
run_pairwise_ttests <- function(data_subset, value_col, log_col) {
groups <- c("OMM", "FAM", "CON")
pairs <- combn(groups, 2, simplify = FALSE)
results <- list()
for(pair in pairs) {
g1 <- pair[1]; g2 <- pair[2]
d1 <- data_subset %>% filter(Group == g1)
d2 <- data_subset %>% filter(Group == g2)
if(nrow(d1) > 1 & nrow(d2) > 1) {
t_res <- t.test(d1[[log_col]], d2[[log_col]])
results[[paste(g1, "vs", g2)]] <- data.frame(
Comparison = paste(g1, "vs", g2),
p_value = t_res$p.value,
Mean_G1 = mean(d1[[value_col]]), Mean_G2 = mean(d2[[value_col]]),
Diff_Raw = mean(d1[[value_col]]) - mean(d2[[value_col]])
)
}
}
bind_rows(results)
}
# --- A1. Learning RT ---
print("--- T-TESTS: RT Learning (Block 1-6) ---")
## [1] "--- T-TESTS: RT Learning (Block 1-6) ---"
for(b in 1:6) {
block_res <- run_pairwise_ttests(df_learn %>% filter(Session == b), "feedback.RT.mean", "Log_RT")
print(paste(">> BLOCK", b))
print(block_res %>% mutate(p_value = round(p_value, 4), Diff_Raw = round(Diff_Raw)))
}
## [1] ">> BLOCK 1"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.0031 687.8442 646.7439 41
## 2 OMM vs CON 0.0000 687.8442 637.8715 50
## 3 FAM vs CON 0.1789 646.7439 637.8715 9
## [1] ">> BLOCK 2"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.0000 482.5470 427.8462 55
## 2 OMM vs CON 0.0000 482.5470 409.9186 73
## 3 FAM vs CON 0.4202 427.8462 409.9186 18
## [1] ">> BLOCK 3"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.0002 411.6951 375.8831 36
## 2 OMM vs CON 0.0003 411.6951 381.0337 31
## 3 FAM vs CON 0.9739 375.8831 381.0337 -5
## [1] ">> BLOCK 4"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.0008 389.6975 354.3286 35
## 2 OMM vs CON 0.0071 389.6975 358.7377 31
## 3 FAM vs CON 0.4312 354.3286 358.7377 -4
## [1] ">> BLOCK 5"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.0116 377.8642 357.0178 21
## 2 OMM vs CON 0.0007 377.8642 342.5105 35
## 3 FAM vs CON 0.4705 357.0178 342.5105 15
## [1] ">> BLOCK 6"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.7336 355.2022 353.8190 1
## 2 OMM vs CON 0.1474 355.2022 334.7414 20
## 3 FAM vs CON 0.0627 353.8190 334.7414 19
# --- A2. Test Phase RT ---
print("--- T-TESTS: Test Phase RT ---")
## [1] "--- T-TESTS: Test Phase RT ---"
for(cond in c("Retention", "Transfer")) {
test_res <- run_pairwise_ttests(df_test %>% filter(Condition == cond), "feedback.RT.mean", "Log_RT")
print(paste(">> CONDITION:", cond))
print(test_res %>% mutate(p_value = round(p_value, 4), Diff_Raw = round(Diff_Raw)))
}
## [1] ">> CONDITION: Retention"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.3561 364.0999 349.1460 15
## 2 OMM vs CON 0.0082 364.0999 324.2983 40
## 3 FAM vs CON 0.0801 349.1460 324.2983 25
## [1] ">> CONDITION: Transfer"
## Comparison p_value Mean_G1 Mean_G2 Diff_Raw
## 1 OMM vs FAM 0.0008 468.8932 411.0204 58
## 2 OMM vs CON 0.0001 468.8932 401.4023 67
## 3 FAM vs CON 0.7250 411.0204 401.4023 10
# --- B. Concatenation RT ---
concat_results <- list()
for(b in 1:6) { # Learning
for(p in 1:6) {
res <- run_pairwise_ttests(df_foot_clean %>% filter(Session == b, foot.position == p), "feedback.RT", "Log_RT")
if(nrow(res) > 0) { res$Phase <- paste("Block", b); res$Step <- p; concat_results[[paste("B", b, p)]] <- res }
}
}
for(cond in c("Retention", "Transfer")) { # Test
for(p in 1:6) {
res <- run_pairwise_ttests(df_foot_clean %>% filter(Panel_Label == cond, foot.position == p), "feedback.RT", "Log_RT")
if(nrow(res) > 0) { res$Phase <- cond; res$Step <- p; concat_results[[paste(cond, p)]] <- res }
}
}
print("--- T-TESTS: Concatenation (Significant Only) ---")
## [1] "--- T-TESTS: Concatenation (Significant Only) ---"
final_concat <- bind_rows(concat_results) %>% filter(p_value < 0.05) %>% select(Phase, Step, Comparison, p_value, Diff_Raw)
print(head(final_concat, 100))
## Phase Step Comparison p_value Diff_Raw
## 1 Block 1 1 FAM vs CON 2.926226e-02 -209.7592321
## 2 Block 1 2 OMM vs FAM 2.620306e-05 68.6304703
## 3 Block 1 2 OMM vs CON 4.670809e-11 107.0418704
## 4 Block 1 2 FAM vs CON 2.978048e-02 38.4114002
## 5 Block 1 3 OMM vs FAM 5.146015e-03 45.0401276
## 6 Block 1 3 OMM vs CON 2.880810e-06 118.0499785
## 7 Block 1 4 OMM vs FAM 1.102689e-02 18.5755526
## 8 Block 1 4 OMM vs CON 4.294532e-08 106.7172887
## 9 Block 1 4 FAM vs CON 8.955553e-03 88.1417361
## 10 Block 1 5 OMM vs FAM 6.230306e-08 103.3303108
## 11 Block 1 5 OMM vs CON 1.652877e-11 155.3145324
## 12 Block 1 6 OMM vs FAM 1.156872e-05 70.9329625
## 13 Block 1 6 OMM vs CON 1.780815e-06 82.3797083
## 14 Block 2 1 OMM vs FAM 7.780598e-03 -5.2963409
## 15 Block 2 1 OMM vs CON 1.494462e-06 80.2488886
## 16 Block 2 2 OMM vs FAM 7.274168e-05 45.0920217
## 17 Block 2 2 OMM vs CON 5.705961e-11 68.0049284
## 18 Block 2 2 FAM vs CON 1.602813e-03 22.9129067
## 19 Block 2 3 OMM vs FAM 3.355918e-08 67.3701838
## 20 Block 2 3 OMM vs CON 9.433374e-12 86.9931465
## 21 Block 2 4 OMM vs FAM 1.437731e-09 72.9670927
## 22 Block 2 4 OMM vs CON 1.546658e-13 89.8249186
## 23 Block 2 5 OMM vs FAM 1.704973e-10 81.9947285
## 24 Block 2 5 OMM vs CON 5.350715e-11 87.4419769
## 25 Block 2 6 OMM vs FAM 1.067780e-09 66.0773350
## 26 Block 2 6 OMM vs CON 3.182774e-03 23.2569924
## 27 Block 2 6 FAM vs CON 3.131489e-03 -42.8203426
## 28 Block 3 1 OMM vs FAM 2.520562e-06 33.8301893
## 29 Block 3 1 OMM vs CON 1.684813e-03 20.3063235
## 30 Block 3 1 FAM vs CON 4.855950e-02 -13.5238659
## 31 Block 3 2 OMM vs CON 7.666155e-05 17.1280171
## 32 Block 3 2 FAM vs CON 6.346276e-08 27.3582112
## 33 Block 3 3 OMM vs FAM 5.226842e-03 29.6381570
## 34 Block 3 3 OMM vs CON 4.021453e-08 49.8556961
## 35 Block 3 3 FAM vs CON 2.646057e-03 20.2175391
## 36 Block 3 4 OMM vs FAM 1.009982e-03 27.5799860
## 37 Block 3 4 OMM vs CON 8.244858e-07 28.9661140
## 38 Block 3 5 OMM vs FAM 4.168149e-07 74.3059841
## 39 Block 3 5 OMM vs CON 6.465741e-08 67.5205712
## 40 Block 3 6 OMM vs FAM 2.163221e-05 59.7478402
## 41 Block 3 6 FAM vs CON 4.476545e-05 -59.5563202
## 42 Block 4 1 OMM vs FAM 3.503340e-03 54.8205365
## 43 Block 4 1 OMM vs CON 2.444126e-03 61.7071938
## 44 Block 4 2 FAM vs CON 2.795740e-02 4.1968119
## 45 Block 4 3 OMM vs FAM 1.796622e-03 29.5190789
## 46 Block 4 3 OMM vs CON 1.289094e-06 47.1184723
## 47 Block 4 4 OMM vs FAM 2.435674e-03 29.0621387
## 48 Block 4 4 OMM vs CON 2.611729e-03 26.1094956
## 49 Block 4 5 OMM vs FAM 1.828391e-07 73.2048902
## 50 Block 4 5 OMM vs CON 1.504925e-05 61.9291960
## 51 Block 4 6 OMM vs FAM 1.103089e-02 35.4796580
## 52 Block 4 6 FAM vs CON 5.220804e-05 -40.9092032
## 53 Block 5 1 OMM vs FAM 6.197775e-04 49.7937781
## 54 Block 5 1 OMM vs CON 6.235256e-13 128.4866188
## 55 Block 5 1 FAM vs CON 2.688475e-04 78.6928407
## 56 Block 5 3 OMM vs FAM 3.365879e-02 12.1547704
## 57 Block 5 3 OMM vs CON 1.012254e-04 36.0594394
## 58 Block 5 4 OMM vs CON 1.406548e-02 18.4980604
## 59 Block 5 5 OMM vs FAM 9.837092e-07 50.9223722
## 60 Block 5 5 OMM vs CON 2.939599e-07 54.6118824
## 61 Block 5 6 FAM vs CON 1.696904e-02 -24.2855755
## 62 Block 6 1 OMM vs FAM 1.072453e-02 38.9089654
## 63 Block 6 1 OMM vs CON 2.706324e-08 86.0333340
## 64 Block 6 1 FAM vs CON 2.263673e-03 47.1243686
## 65 Block 6 2 OMM vs FAM 6.099048e-04 -31.6602915
## 66 Block 6 2 FAM vs CON 6.927375e-03 -0.1463324
## 67 Block 6 3 OMM vs CON 7.631619e-04 30.6019009
## 68 Block 6 3 FAM vs CON 3.721458e-04 38.5996717
## 69 Block 6 4 FAM vs CON 1.840610e-02 17.2818949
## 70 Block 6 5 OMM vs FAM 1.674826e-02 33.2912212
## 71 Block 6 5 OMM vs CON 2.347018e-05 46.6081681
## 72 Block 6 6 OMM vs CON 2.320048e-02 -26.1138486
## 73 Retention 1 OMM vs FAM 4.508165e-02 37.2758176
## 74 Retention 1 OMM vs CON 1.271498e-12 138.8451103
## 75 Retention 1 FAM vs CON 2.648589e-07 101.5692927
## 76 Retention 3 OMM vs CON 1.807403e-03 43.8259031
## 77 Retention 3 FAM vs CON 4.554302e-03 39.6534052
## 78 Retention 5 OMM vs FAM 2.197188e-03 47.6577523
## 79 Retention 5 OMM vs CON 2.506808e-04 56.8386485
## 80 Transfer 1 OMM vs FAM 4.971358e-04 128.8279246
## 81 Transfer 1 OMM vs CON 1.096730e-09 197.6400079
## 82 Transfer 1 FAM vs CON 7.235065e-03 68.8120833
## 83 Transfer 2 OMM vs CON 2.998438e-02 22.2744411
## 84 Transfer 3 OMM vs FAM 2.920290e-02 38.7690885
## 85 Transfer 4 OMM vs CON 3.120026e-02 44.4083169
## 86 Transfer 5 OMM vs FAM 1.818928e-02 59.9735688
## 87 Transfer 5 OMM vs CON 1.471898e-03 55.4082156
## 88 Transfer 6 OMM vs FAM 1.892055e-05 53.9776046
## 89 Transfer 6 OMM vs CON 1.262649e-06 61.8591475
# ==============================================================================
# 6. RT IMPROVEMENT ANALYSIS (% Change Block 1 -> 6)
# ==============================================================================
# Filter Correct Trials & Remove Outliers (2.5 SD) from raw data first
rt_clean_imp <- trial_df %>%
filter(Session %in% c(1, 6), feedback.ACC.trial == 1, feedback.RT.mean > 0) %>%
group_by(Subject, Session) %>%
mutate(mn=mean(feedback.RT.mean), sd=sd(feedback.RT.mean),
is_outlier = feedback.RT.mean > (mn+2.5*sd) | feedback.RT.mean < (mn-2.5*sd)) %>%
filter(!is_outlier) %>% ungroup()
# Calculate Means & % Difference
rt_summary <- rt_clean_imp %>%
group_by(Subject, Group, Session) %>%
summarise(Mean_RT = mean(feedback.RT.mean), .groups = "drop") %>%
pivot_wider(names_from = Session, values_from = Mean_RT, names_prefix = "Block") %>%
mutate(percentage_diff = (Block1 - Block6) / Block1 * 100) %>%
filter(!is.na(percentage_diff))
# ANOVA (Type II)
print("--- ANOVA: RT Improvement % (Block 1 to 6) ---")
## [1] "--- ANOVA: RT Improvement % (Block 1 to 6) ---"
m_improv <- lm(percentage_diff ~ Group, data = rt_summary)
print(Anova(m_improv, type="II"))
## Anova Table (Type II tests)
##
## Response: percentage_diff
## Sum Sq Df F value Pr(>F)
## Group 423.3 2 0.5197 0.5982
## Residuals 18328.0 45
# Plot (Lukas Style)
# Using model effects to get proper error bars
m_improv_eff <- allEffects(m_improv)
m_improv_df <- as.data.frame(m_improv_eff[[1]])
p_improv <- ggplot(m_improv_df, aes(x=Group, y=fit, fill=Group)) +
geom_bar(stat="identity", color="black", show.legend = FALSE, width=0.7) +
geom_errorbar(aes(ymin=lower, ymax=upper), width=.2) +
scale_fill_manual(values = custom_colors) +
labs(title="RT Improvement (Block 1 to 6)", y="% Improvement", x=NULL) +
theme_classic() + theme(legend.position="none")
print(p_improv)
# ==============================================================================
# 7. SELF-REPORT ANALYSIS (LMM + Individual Data Points)
# ==============================================================================
# 1. Prepare Data (Effort Only)
df_survey_eff <- survey_df %>%
filter(!is.na(Measure_Label)) %>%
select(Participant, Group, Measure_Label, Effort) %>%
mutate(Effort = as.numeric(Effort))
# 2. Linear Mixed Model (Group x Phase)
print("--- MODEL: SUBJECTIVE EFFORT (Group x Phase) ---")
## [1] "--- MODEL: SUBJECTIVE EFFORT (Group x Phase) ---"
m_effort <- lmer(Effort ~ Group * Measure_Label + (1|Participant),
data = df_survey_eff, REML = FALSE)
print(Anova(m_effort, type="II"))
## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Effort
## Chisq Df Pr(>Chisq)
## Group 2.8076 2 0.24566
## Measure_Label 103.1879 3 < 2e-16 ***
## Group:Measure_Label 12.9126 6 0.04445 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# 3. Exploratory Pairwise Comparisons
print("--- PAIRWISE CONTRASTS (Effort by Phase) ---")
## [1] "--- PAIRWISE CONTRASTS (Effort by Phase) ---"
emm_effort <- emmeans(m_effort, ~ Group | Measure_Label)
print(pairs(emm_effort, adjust="none"))
## Measure_Label = Baseline:
## contrast estimate SE df t.ratio p.value
## CON - FAM 0.688 1.68 180 0.410 0.6826
## CON - OMM 0.500 1.68 180 0.298 0.7661
## FAM - OMM -0.188 1.68 180 -0.112 0.9112
##
## Measure_Label = Post-Training:
## contrast estimate SE df t.ratio p.value
## CON - FAM -4.812 1.68 180 -2.867 0.0046
## CON - OMM -4.562 1.68 180 -2.718 0.0072
## FAM - OMM 0.250 1.68 180 0.149 0.8818
##
## Measure_Label = Early Learning:
## contrast estimate SE df t.ratio p.value
## CON - FAM -1.438 1.68 180 -0.856 0.3929
## CON - OMM 0.375 1.68 180 0.223 0.8235
## FAM - OMM 1.812 1.68 180 1.080 0.2817
##
## Measure_Label = Late Learning:
## contrast estimate SE df t.ratio p.value
## CON - FAM -0.562 1.68 180 -0.335 0.7379
## CON - OMM -2.250 1.68 180 -1.340 0.1818
## FAM - OMM -1.688 1.68 180 -1.005 0.3161
##
## Degrees-of-freedom method: kenward-roger
# 4. Plot (Bar Chart + Individual Points)
p_effort <- ggplot(df_survey_eff, aes(x=Measure_Label, y=Effort, fill=Group)) +
# A. The Bars
stat_summary(fun = mean, geom = "bar", position = position_dodge(0.8),
width = 0.7, color="black", alpha=0.6) +
# B. The Individual Points
geom_point(position = position_jitterdodge(jitter.width = 0.15, dodge.width = 0.8),
size = 1.5, alpha = 0.4, show.legend = FALSE) +
# C. The Error Bars (SE)
stat_summary(fun.data = mean_se, geom = "errorbar",
width = 0.2, position = position_dodge(0.8), color="black") +
# ---------------------------------------------------------
# NEW: Significance Brackets for "Post-Training" (x = 2)
# ---------------------------------------------------------
# Bracket 1: CON (1.73) vs FAM (2.0)
annotate("segment", x = 1.73, xend = 2.0, y = 22, yend = 22) +
annotate("segment", x = 1.73, xend = 1.73, y = 21, yend = 22) +
annotate("segment", x = 2.0, xend = 2.0, y = 21, yend = 22) +
annotate("text", x = 1.865, y = 22.5, label = "**", color = "red", size = 6) +
# Bracket 2: CON (1.73) vs OMM (2.27)
# (Drawn slightly higher at y=24 to stack nicely above Bracket 1)
annotate("segment", x = 1.73, xend = 2.27, y = 24, yend = 24) +
annotate("segment", x = 1.73, xend = 1.73, y = 23, yend = 24) +
annotate("segment", x = 2.27, xend = 2.27, y = 23, yend = 24) +
annotate("text", x = 2.0, y = 24.5, label = "**", color = "red", size = 6) +
# Ensure the Y-axis is tall enough to fit the brackets
scale_y_continuous(limits = c(0, 26)) +
scale_fill_manual(values = custom_colors) +
labs(title="5A) Subjective Effort across Phases",
subtitle="Bars = Mean + SE; Dots = Individual Participants",
y="NASA-TLX Effort Score", x=NULL) +
theme_classic() +
theme(legend.position="top", axis.text.x = element_text(angle=15, hjust=1))
print(p_effort)
# ==============================================================================
# 8. FINAL COMPOSITE PLOTS
# ==============================================================================
mean_ci95 <- function(x) {
m <- mean(x, na.rm = TRUE)
se <- sd(x, na.rm = TRUE) / sqrt(sum(!is.na(x)))
return(data.frame(y = m, ymin = m - 1.96 * se, ymax = m + 1.96 * se))
}
# --- 1. Prepare Participant-Level Aggregated Accuracy ---
# Calculate % accuracy per participant per block for Learning
acc_subj_learn <- trial_df %>%
filter(Session <= 6) %>%
group_by(Subject, Group, Session) %>%
summarise(Acc_Rate = mean(as.numeric(feedback.ACC.trial), na.rm = TRUE), .groups = "drop")
# Calculate % accuracy per participant per condition for Test
acc_subj_test <- trial_df %>%
filter(Session > 6) %>%
group_by(Subject, Group, Condition) %>%
summarise(Acc_Rate = mean(as.numeric(feedback.ACC.trial), na.rm = TRUE), .groups = "drop")
# --- 2. FIGURE 2A: Learning Accuracy (Boxplot + Jitter) ---
p1a <- ggplot(acc_subj_learn, aes(x = factor(Session), y = Acc_Rate, fill = Group, color = Group)) +
# Add raw participant data points
geom_point(position = position_jitterdodge(jitter.width = 0.15, dodge.width = 0.75),
alpha = 0.3, size = 1, show.legend = FALSE) +
# Add boxplots
geom_boxplot(position = position_dodge(width = 0.75),
alpha = 0.7, outlier.shape = NA, color = "black") +
# Add Significance Bracket: Block 1 CON (x=0.75) vs OMM (x=1.25)
annotate("segment", x = 0.75, xend = 1.25, y = 1.02, yend = 1.02) +
annotate("segment", x = 0.75, xend = 0.75, y = 1.00, yend = 1.02) +
annotate("segment", x = 1.25, xend = 1.25, y = 1.00, yend = 1.02) +
annotate("text", x = 1, y = 1.04, label = "**", color = "red", size = 6) +
# Set Y-axis to accommodate the bracket above 100%
scale_y_continuous(labels = scales::percent, limits = c(0, 1.1)) +
scale_fill_manual(values = custom_colors) +
scale_color_manual(values = custom_colors) +
labs(title = "2A) Accuracy (Learning)", y = "Accuracy", x = "Block") +
theme_classic() +
theme(legend.position = "top", axis.text.x = element_blank())
print(p1a)
# --- ALTERNATIVE FIGURE 2B: RT Learning (Boxplot + Jitter) ---
p1b <- ggplot(df_learn_plot, aes(x = factor(Session), y = feedback.RT.mean, fill = Group, color = Group)) +
# Add raw data points
geom_point(position = position_jitterdodge(jitter.width = 0.2, dodge.width = 0.75),
alpha = 0.3, size = 1, show.legend = FALSE) +
# Add boxplots
geom_boxplot(position = position_dodge(width = 0.75),
alpha = 0.7, outlier.shape = NA, color = "black") +
# Add Significance Bracket: Block 1 FAM (x=1) vs OMM (x=1.25)
annotate("segment", x = 1, xend = 1.25, y = 1350, yend = 1350) +
annotate("segment", x = 1, xend = 1, y = 1300, yend = 1350) +
annotate("segment", x = 1.25, xend = 1.25, y = 1300, yend = 1350) +
annotate("text", x = 1.125, y = 1400, label = "*", color = "red", size = 6) +
# Add Significance Bracket: Block 2 FAM (x=2) vs OMM (x=2.25)
annotate("segment", x = 2, xend = 2.25, y = 1150, yend = 1150) +
annotate("segment", x = 2, xend = 2, y = 1100, yend = 1150) +
annotate("segment", x = 2.25, xend = 2.25, y = 1100, yend = 1150) +
annotate("text", x = 2.125, y = 1200, label = "*", color = "red", size = 6) +
# Extend the Y-axis so the brackets don't get cut off
coord_cartesian(ylim = c(0, 1500)) +
scale_fill_manual(values = custom_colors) +
scale_color_manual(values = custom_colors) +
labs(title = "2B) Response Time (Learning)", y = "RT (ms)", x = "Block") +
theme_classic() +
theme(legend.position = "none")
print(p1b)
# --- 2A: Accuracy Test (Error Bars for Categorical) ---
# Using raw means with CI
acc_test_summ <- trial_df %>%
filter(Session > 6) %>%
group_by(Group, Condition) %>%
summarise(
M = mean(feedback.ACC.trial),
SE = sd(feedback.ACC.trial)/sqrt(n()),
CI = 1.96 * SE,
.groups = "drop"
)
p2a <- ggplot(acc_subj_test, aes(x = Condition, y = Acc_Rate, fill = Group, color = Group)) +
geom_point(position = position_jitterdodge(jitter.width = 0.15, dodge.width = 0.75),
alpha = 0.3, size = 1, show.legend = FALSE) +
geom_boxplot(position = position_dodge(width = 0.75),
alpha = 0.7, outlier.shape = NA, color = "black") +
# FIX: Use coord_cartesian to "zoom" without deleting data below 40% or above 100%
coord_cartesian(ylim = c(0.3, 1.05)) +
scale_y_continuous(labels = scales::percent) +
scale_fill_manual(values = custom_colors) +
scale_color_manual(values = custom_colors) +
labs(title = "3A) Accuracy (Test)", y = "Accuracy", x = "Condition") +
theme_classic() +
theme(legend.position = "top")
print(p2a)
# --- 2B: RT Test ---
p2b <- ggplot(df_test_plot, aes(x = Condition, y = feedback.RT.mean, fill = Group, color = Group)) +
geom_point(position = position_jitterdodge(jitter.width = 0.2, dodge.width = 0.75),
alpha = 0.3, size = 1, show.legend = FALSE) +
geom_boxplot(position = position_dodge(width = 0.75),
alpha = 0.7, outlier.shape = NA, color = "black") +
# Add Significance Bracket: Transfer FAM (x=2) vs OMM (x=2.25)
annotate("segment", x = 2, xend = 2.25, y = 1350, yend = 1350) +
annotate("segment", x = 2, xend = 2, y = 1300, yend = 1350) +
annotate("segment", x = 2.25, xend = 2.25, y = 1300, yend = 1350) +
annotate("text", x = 2.125, y = 1380, label = "*", color = "red", size = 6) +
# FIX: Extend coord_cartesian to 1500 to give the asterisk plenty of breathing room
coord_cartesian(ylim = c(0, 1500)) +
scale_fill_manual(values = custom_colors) +
scale_color_manual(values = custom_colors) +
labs(title = "3B) Response Time (Test)", y = "RT (ms)", x = "Condition") +
theme_classic() +
theme(legend.position = "none")
print(p2b)
# --- 3: Concatenation (Ribbon) ---
pd <- position_dodge(0.2)
p3a <- ggplot(df_foot_plot, aes(x=foot.position, y=feedback.RT, color=Group, fill=Group, group=Group)) +
stat_summary(fun.data = mean_ci95, geom = "ribbon", alpha = 0.2, color = NA) +
stat_summary(fun = mean, geom = "line", size = 1, position = pd) +
stat_summary(fun = mean, geom = "point", size = 1.5, position = pd) +
facet_wrap(~Panel_Label, nrow=2) +
scale_color_manual(values=custom_colors) + scale_fill_manual(values=custom_colors) +
labs(title="3. Concatenation Dynamics", subtitle="Mean RT +/- 95% CI", x="Step", y="Step RT (ms)") +
coord_cartesian(ylim=c(0, 1000)) + theme_classic() + theme(strip.background=element_blank(), legend.position="top")
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
# Alternative 3b
pd <- position_dodge(0.2)
p3b <- ggplot(df_foot_plot, aes(x = foot.position, y = feedback.RT, color = Group, fill = Group, group = Group)) +
# 1. Add a subtle shaded background explicitly highlighting Steps 4 and 5
annotate("rect", xmin = 3.8, xmax = 5.2, ymin = -Inf, ymax = Inf,
alpha = 0.1, fill = "black", color = NA) +
# 2. Drop the ribbon (since CIs are tiny) and just use clean lines and slightly larger points
stat_summary(fun = mean, geom = "line", size = 1.2, position = pd) +
stat_summary(fun = mean, geom = "point", size = 2.5, position = pd) +
facet_wrap(~Panel_Label, nrow = 2) +
scale_color_manual(values = custom_colors) +
scale_fill_manual(values = custom_colors) +
labs(title = "3. Concatenation Dynamics",
subtitle = "Shaded region highlights theoretical chunk boundaries (Steps 4 & 5)",
x = "Sequence Step", y = "Step RT (ms)") +
# Adjusted the y-limit slightly so the bump fills the space better
coord_cartesian(ylim = c(200, 1000)) +
theme_classic() +
theme(strip.background = element_blank(),
legend.position = "top",
text = element_text(size = 14))
print(p3b)
# Alternative3c.
# Create a tiny dataframe telling ggplot exactly where to put the stars
stars_df <- data.frame(
Panel_Label = c("Block 1", "Block 1", "Block 2", "Block 2", "Block 3", "Block 3", "Block 4", "Block 4", "Block 5", "Block 5", "Block 6", "Retention", "Transfer"),
foot.position = c(4, 5, 4, 5, 4, 5, 4, 5, 4, 5, 5, 5, 5),
feedback.RT = c(850, 850, 600, 600, 500, 500, 500, 500, 500, 500, 500, 500, 500), # Adjust Y heights if they overlap the lines
label = c("*", "***", "*", "***", "*", "***", "*", "***", "*", "***", "***", "***", "***")
)
# 2. Build the plot
p3c <- ggplot(df_foot_plot, aes(x = foot.position, y = feedback.RT, color = Group, group = Group)) +
stat_summary(fun = mean, geom = "line", size = 1, position = pd) +
stat_summary(fun = mean, geom = "point", size = 2, position = pd) +
# FIX: inherit.aes = FALSE prevents the 'Group not found' error
geom_text(data = stars_df, aes(x = foot.position, y = feedback.RT, label = label),
color = "red", size = 6, inherit.aes = FALSE) +
facet_wrap(~Panel_Label, nrow = 4) +
scale_color_manual(values = custom_colors) +
labs(title = "4. Concatenation Dynamics", x = "Step", y = "Step RT (ms)") +
coord_cartesian(ylim = c(0, 1000)) +
theme_classic() +
theme(strip.background = element_blank(), legend.position = "top")
print(p3c)
# Compile and Output All
grid1 <- plot_grid(p1a, p1b, ncol = 1, align = "v", rel_heights = c(0.8, 1))
grid2 <- plot_grid(p2a, p2b, ncol = 1, align = "v", rel_heights = c(0.8, 1))
print(grid1)
print(grid2)
print(p3a)
print(p3b)
print(p3c)
# ==============================================================================
# 9. CORRELATION: Does Effort Predict Learning?
# ==============================================================================
# 1. Calculate Learning Magnitude (Standard Mean Method)
# Using 2.5SD cleaning to ensure the "Change Score" isn't driven by one bad trial
rt_clean_corr <- trial_df %>%
filter(Session %in% c(1, 6), feedback.ACC.trial == 1, feedback.RT.mean > 0) %>%
group_by(Subject, Session) %>%
mutate(
mn = mean(feedback.RT.mean),
sd = sd(feedback.RT.mean),
is_outlier = feedback.RT.mean > (mn + 2.5*sd) | feedback.RT.mean < (mn - 2.5*sd)
) %>%
filter(!is_outlier) %>%
ungroup()
# Calculate % Improvement
rt_summary <- rt_clean_corr %>%
group_by(Subject, Group, Session) %>%
summarise(Mean_RT = mean(feedback.RT.mean), .groups = "drop") %>%
pivot_wider(names_from = Session, values_from = Mean_RT, names_prefix = "Block") %>%
mutate(percentage_diff = (Block1 - Block6) / Block1 * 100) %>%
filter(!is.na(percentage_diff))
# 2. Prepare Effort Data (Post-Training)
effort_m2 <- survey_df %>%
filter(Measure == 2) %>%
select(Participant, Effort) %>%
mutate(Effort = as.numeric(Effort)) %>%
rename(Subject = Participant, PostTraining_Effort = Effort)
# 3. Merge
rt_summary$Subject <- as.character(rt_summary$Subject)
effort_m2$Subject <- as.character(effort_m2$Subject)
corr_df <- inner_join(rt_summary, effort_m2, by="Subject")
# 4. Statistical Test (Spearman, One-Tailed)
# Hypothesis: Higher Effort -> Greater Improvement (Positive Correlation)
print("--- CORRELATION: RT Improvement % vs Effort ---")
## [1] "--- CORRELATION: RT Improvement % vs Effort ---"
cor_res <- cor.test(corr_df$percentage_diff, corr_df$PostTraining_Effort,
method = "spearman", alternative = "greater")
## Warning in cor.test.default(corr_df$percentage_diff,
## corr_df$PostTraining_Effort, : Cannot compute exact p-value with ties
print(cor_res)
##
## Spearman's rank correlation rho
##
## data: corr_df$percentage_diff and corr_df$PostTraining_Effort
## S = 13506, p-value = 0.03331
## alternative hypothesis: true rho is greater than 0
## sample estimates:
## rho
## 0.2669547
# 5. Visualization (Scatterplot with Regression)
p_corr <- ggplot(corr_df, aes(x=PostTraining_Effort, y=percentage_diff)) +
geom_point(aes(color=Group), size=3, alpha=0.7) +
# Add regression line to show the trend
geom_smooth(method="lm", color="black", se=TRUE, fill="grey80", alpha=0.5) +
scale_color_manual(values = custom_colors) +
labs(title="5B) Effort Predicts Motor Learning",
subtitle = paste("Spearman rho =", round(cor_res$estimate, 2), ", p =", round(cor_res$p.value, 3), "(1-tailed)"),
x="Subjective Effort (Post-Training)",
y="RT Improvement (Block 1 to 6 %)") +
theme_classic() +
theme(legend.position="top")
print(p_corr)
## `geom_smooth()` using formula = 'y ~ x'
# ==============================================================================
# 10. MECHANISM ANALYSIS: Does Group Moderate the Effort-Learning Link?
# ==============================================================================
# --- 1. Prepare Data ---
effort_m2 <- survey_df %>% filter(Measure == 2) %>%
select(Participant, Effort) %>% mutate(Effort=as.numeric(Effort)) %>%
rename(Subject = Participant, PostTraining_Effort = Effort)
rt_summary$Subject <- as.character(rt_summary$Subject)
effort_m2$Subject <- as.character(effort_m2$Subject)
mech_df <- inner_join(rt_summary, effort_m2, by="Subject")
# --- 2. The Interaction Model (Moderation) ---
# Does the relationship between Effort and Learning differ by Group?
print("--- MODEL: Learning ~ Effort * Group ---")
## [1] "--- MODEL: Learning ~ Effort * Group ---"
m_mech <- lm(percentage_diff ~ PostTraining_Effort * Group, data = mech_df)
print(Anova(m_mech, type="II"))
## Anova Table (Type II tests)
##
## Response: percentage_diff
## Sum Sq Df F value Pr(>F)
## PostTraining_Effort 568.4 1 1.4170 0.2406
## Group 209.4 2 0.2611 0.7715
## PostTraining_Effort:Group 913.4 2 1.1386 0.3300
## Residuals 16846.3 42
# --- 3. Group-Specific Correlations ---
# Let's see the 'r' value for each group individually
print("--- CORRELATIONS BY GROUP ---")
## [1] "--- CORRELATIONS BY GROUP ---"
group_cors <- mech_df %>%
group_by(Group) %>%
summarise(
r = cor(percentage_diff, PostTraining_Effort, method = "spearman"),
p = cor.test(percentage_diff, PostTraining_Effort, method = "spearman")$p.value,
n = n()
)
## Warning: There were 3 warnings in `summarise()`.
## The first warning was:
## ℹ In argument: `p = cor.test(percentage_diff, PostTraining_Effort, method =
## "spearman")$p.value`.
## ℹ In group 1: `Group = CON`.
## Caused by warning in `cor.test.default()`:
## ! Cannot compute exact p-value with ties
## ℹ Run ]8;;ide:run:dplyr::last_dplyr_warnings()dplyr::last_dplyr_warnings()]8;; to see the 2 remaining warnings.
print(group_cors)
## # A tibble: 3 × 4
## Group r p n
## <fct> <dbl> <dbl> <int>
## 1 CON 0.381 0.145 16
## 2 FAM 0.416 0.109 16
## 3 OMM -0.0488 0.858 16
# --- 4. Alternative 3-Way Interaction Model (Moderation of Learning Slope) ---
print("--- MODEL: LMM EFFORT MODERATION ON LEARNING RATE ---")
## [1] "--- MODEL: LMM EFFORT MODERATION ON LEARNING RATE ---"
# Join the effort scores to your trial-level RT learning dataframe
df_learn_effort <- df_learn %>%
inner_join(effort_m2, by="Subject")
# Run the LMM testing how Effort moderates the RT slope over Sessions
m_rt_effort <- lmer(Log_RT ~ Group * as.numeric(Session) * PostTraining_Effort + (1|Subject),
data = df_learn_effort)
print(Anova(m_rt_effort, type="II"))
## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Log_RT
## Chisq Df Pr(>Chisq)
## Group 0.6724 2 0.7145
## as.numeric(Session) 3849.7017 1 < 2.2e-16 ***
## PostTraining_Effort 0.9948 1 0.3186
## Group:as.numeric(Session) 25.2245 2 3.331e-06 ***
## Group:PostTraining_Effort 1.0389 2 0.5949
## as.numeric(Session):PostTraining_Effort 37.1745 1 1.080e-09 ***
## Group:as.numeric(Session):PostTraining_Effort 162.6103 2 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# ==============================================================================
# 11. VISUALIZING THE DISSOCIATION
# ==============================================================================
p_mech_final <- ggplot(mech_df, aes(x=PostTraining_Effort, y=percentage_diff, color=Group)) +
geom_point(size=3, alpha=0.6) +
# Separate regression lines
geom_smooth(method="lm", se=FALSE, size=1.5) +
scale_color_manual(values = custom_colors) +
labs(title="5C) Distinct Learning Regimes",
subtitle="FAM/CON rely on Effort; OMM improves via Strategy (Dissociated)",
x="Subjective Effort (Post-Training)", y="RT Improvement (%)") +
theme_classic() + theme(legend.position="top")
print(p_mech_final)
## `geom_smooth()` using formula = 'y ~ x'
# ==============================================================================
# 12. FIGURE 4: SUBJECTIVE MECHANISMS (Composite Plot)
# ==============================================================================
# --- Ensure Plots are Ready (from previous sections) ---
# p_effort (from Section 7: Bar + Raincloud)
# p_corr (from Section 8: Global Correlation)
# p_mech_final (from Section 9: Group Regimes)
# --- 1. Clean up Titles for the Composite Figure ---
# We remove the standard titles so we can use the A/B/C labels clearly
p4a <- p_effort +
#labs(title = "5A) Subjective Effort Trajectory", subtitle = NULL) +
theme(plot.margin = margin(10, 10, 10, 10))
p4b <- p_corr +
#labs(title = "5B) Global Mechanism: Effort Predicts Learning", subtitle = NULL) +
theme(plot.margin = margin(10, 10, 10, 10))
p4c <- p_mech_final +
#labs(title = "5C) Distinct Regimes: Strategy vs. Effort", subtitle = NULL) +
theme(plot.margin = margin(10, 10, 10, 10))
# --- 2. Combine Vertically using cowplot ---
# Rel_heights: Give the Bar chart (A) slightly less space than the Scatters (B/C)
figure4_grid <- plot_grid(
p4a,
p4b,
p4c,
ncol = 1,
#labels = c(),
label_size = 14,
align = "v",
rel_heights = c(0.8, 1, 1)
)
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
# --- 3. Display and Save ---
print(figure4_grid)
# Optional: Save for Manuscript (High Res)
# ggsave("Figure4_Mechanisms.png", figure4_grid, width = 6, height = 12, dpi = 300)
# ggsave("Figure4_Mechanisms.pdf", figure4_grid, width = 6, height = 12)
# Figure 2: Learning Dynamics (Tall - 2 panels)
ggsave("Figure2_Learning.png", plot = grid1, width = 6.5, height = 8, units = "in", dpi = 300)
# Figure 3: Test Phase (Tall - 2 panels)
ggsave("Figure3_Test.png", plot = grid2, width = 6.5, height = 8, units = "in", dpi = 300)
# Figure 4: Concatenation Dynamics (Wide - Faceted)
ggsave("Figure4_Concatenation.png", plot = p3c, width = 7.5, height = 9, units = "in", dpi = 300)
# Figure 5: Subjective Mechanisms (Very Tall - 3 panels)
# Note: You already built this as 'figure4_grid' in your script, but it is Figure 5 in the text
ggsave("Figure5_Mechanisms.png", plot = figure4_grid, width = 6.5, height = 9, units = "in", dpi = 300)