SEED <- 123; set.seed(SEED)
suppressPackageStartupMessages({
library(readxl)
library(dplyr)
library(stringr)
library(tidyr)
library(forcats)
library(lme4); library(lmerTest)
library(emmeans)
library(purrr)
library(tibble)
})
path_scales <- "C:/Users/97380/Documents/LFF_data process/Scale/scale_data_pre.xlsx"
sheet_name <- 1 # 或者 "变量最终得分"
raw <- read_excel(path_scales, sheet = sheet_name)
dat_wide <- raw %>%
mutate(
group = as.character(Group),
group = recode(group, `1`="control", `2`="training", .default = group),
group = factor(group, levels = c("control","training"))
)
# 2) 标准化 ID 与 group（1->control, 2->training）
dat_wide <- raw %>%
rename(ID = ID, Group = Group) %>%
mutate(
group = as.character(Group),
group = recode(group, `1` = "control", `2` = "training", .default = group),
group = factor(group, levels = c("control","training"))
)
# 3) 选择所有量表列：以 pre_/post_/fu_ 开头的列
scale_cols <- names(dat_wide)[str_detect(names(dat_wide), "^(pre|post|fu)_")]
# 4) 长数据：time 从列名提取，metric 为量表名，score 为数值
scales_long <- dat_wide %>%
select(ID, group, all_of(scale_cols)) %>%
pivot_longer(all_of(scale_cols), names_to = "name", values_to = "score") %>%
extract(name, into = c("time","metric"), regex = "^(pre|post|fu)_(.+)$") %>%
mutate(
time = factor(time, levels = c("pre","post","fu")),
metric = metric %>% str_replace_all("\u3000", " ") %>% str_squish(), # 清理全角空格
score = suppressWarnings(as.numeric(score))
) %>%
select(ID, group, time, metric, score) %>%
arrange(metric, ID, time)
# 选 VVIQ，准备数据
metric_name <- "VVIQ"
dat_vviq <- scales_long %>%
filter(metric == metric_name) %>%
mutate(
group = fct_relevel(group, "control", "training"),
time = fct_relevel(time, "pre", "post", "fu")
)
# 小工具：按时间窗口建模
mk_fit <- function(df, times, label){
d <- df %>% filter(time %in% times) %>% droplevels() %>%
mutate(time = fct_relevel(time, times))
fit <- lmer(score ~ group * time + (1|ID), data = d, REML = TRUE)
list(fit = fit, label = label)
}
fits <- list(
mk_fit(dat_vviq, c("pre","post"), "pre-post"),
mk_fit(dat_vviq, c("pre","fu"), "pre-fu"),
mk_fit(dat_vviq, c("post","fu"), "post-fu")
)
# ===== 固定效应：t / SE / p / d =====
fixed_tbl <- map_dfr(fits, function(x){
fit <- x$fit; lab <- x$label; sig <- sigma(fit)
co <- as.data.frame(coef(summary(fit)))
tibble(
model = lab,
term = rownames(co),
estimate = co[,"Estimate"],
std.error = co[,"Std. Error"],
df = co[,"df"],
t = co[,"t value"],
p = co[,"Pr(>|t|)"],
d = estimate / sig,
SE_d = std.error / sig,
d_low = d + qt(0.025, df) * SE_d,
d_high = d + qt(0.975, df) * SE_d
)
})
# 组内简单效应：post−pre / fu−pre
within_tbl <- map_dfr(fits[1:2], function(x){
fit <- x$fit; lab <- x$label; sig <- sigma(fit)
emm <- emmeans(fit, ~ time | group)
con <- if (lab == "pre–post") list(post_minus_pre = c(-1, 1))
else list(fu_minus_pre = c(-1, 1))
w <- summary(contrast(emm, con, by = "group", adjust = "none"))
as.data.frame(w) %>%
mutate(
model = lab,
d = estimate / sig,
SE_d = SE / sig,
d_low = d + qt(0.025, df) * SE_d,
d_high = d + qt(0.975, df) * SE_d
) %>%
select(model, group, contrast, estimate, SE, df, t.ratio, p.value, d, SE_d, d_low, d_high)
})
# 查看结果
fixed_tbl
## # A tibble: 12 × 11
## model term estimate std.error df t p d SE_d d_low
## <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 pre-post (Inter… 79.5 3.34 81.5 23.8 1.89e-38 7.35 0.309 6.74
## 2 pre-post groupt… -1.13 4.72 81.5 -0.240 8.11e- 1 -0.105 0.437 -0.974
## 3 pre-post timepo… 6.80 2.79 58.0 2.44 1.79e- 2 0.629 0.258 0.112
## 4 pre-post groupt… -17.0 3.95 58.0 -4.32 6.29e- 5 -1.58 0.365 -2.31
## 5 pre-fu (Inter… 79.5 3.37 80.4 23.6 7.65e-38 7.43 0.316 6.81
## 6 pre-fu groupt… -1.13 4.77 80.4 -0.238 8.13e- 1 -0.106 0.446 -0.994
## 7 pre-fu timefu 3.00 2.76 58.0 1.09 2.82e- 1 0.280 0.258 -0.237
## 8 pre-fu groupt… -9.74 3.90 58.0 -2.50 1.54e- 2 -0.911 0.365 -1.64
## 9 post-fu (Inter… 86.3 3.31 72.3 26.1 1.76e-38 10.1 0.387 9.33
## 10 post-fu groupt… -18.2 4.68 72.3 -3.88 2.28e- 4 -2.13 0.548 -3.22
## 11 post-fu timefu -3.80 2.21 58.0 -1.73 8.99e- 2 -0.445 0.258 -0.962
## 12 post-fu groupt… 7.29 3.12 58.0 2.34 2.28e- 2 0.854 0.365 0.123
## # ℹ 1 more variable: d_high <dbl>
within_tbl
## model group contrast estimate SE df t.ratio p.value
## 1 pre-post control fu_minus_pre 6.800000 2.791218 58 2.436212 0.0179338490
## 2 pre-post training fu_minus_pre -10.233333 2.791218 58 -3.666260 0.0005356436
## 3 pre-fu control fu_minus_pre 2.995000 2.759805 58 1.085222 0.2823136701
## 4 pre-fu training fu_minus_pre -6.744333 2.759805 58 -2.443772 0.0175975727
## d SE_d d_low d_high
## 1 0.6290273 0.2581989 0.1121860 1.1458685
## 2 -0.9466243 0.2581989 -1.4634656 -0.4297831
## 3 0.2802030 0.2581989 -0.2366382 0.7970443
## 4 -0.6309792 0.2581989 -1.1478204 -0.1141379
fixed_d_from_fit <- function(fit){
sig <- sigma(fit)
co <- as.data.frame(coef(summary(fit)))
tibble(
term = rownames(co),
estimate = co[,"Estimate"],
std.error = co[,"Std. Error"],
df = co[,"df"],
t = co[,"t value"],
p = co[,"Pr(>|t|)"],
d = estimate / sig,
SE_d = std.error / sig,
d_low = d + qt(0.025, df) * SE_d,
d_high = d + qt(0.975, df) * SE_d
)
}
within_d_from_fit <- function(fit, window){
emm <- emmeans(fit, ~ time | group)
sig <- sigma(fit)
if (identical(window, c("pre","post"))) {
w <- summary(contrast(emm, list(post_minus_pre = c(-1, 1)),
by = "group", adjust = "none"))
} else if (identical(window, c("pre","fu"))) {
w <- summary(contrast(emm, list(fu_minus_pre = c(-1, 1)),
by = "group", adjust = "none"))
} else return(NULL)
as.data.frame(w) %>%
mutate(
d = estimate / sig,
SE_d = SE / sig,
d_low = d + qt(0.025, df) * SE_d,
d_high = d + qt(0.975, df) * SE_d
) %>%
select(group, contrast, estimate, SE, df, t.ratio, p.value, d, SE_d, d_low, d_high)
}
mk_fit <- function(df, times, label){
d <- df %>% filter(time %in% times) %>% droplevels() %>%
mutate(time = fct_relevel(time, times),
group = fct_relevel(group, "control","training"))
fit <- lmer(score ~ group * time + (1|ID), data = d, REML = TRUE)
list(fit = fit, label = label, window = times)
}
###TEPS — anticipatory
metric_name <- "TEPS_anticipatory"
dat_teps_a <- scales_long |>
mutate(metric = as.character(metric)) |>
filter(metric == metric_name) |>
mutate(group = fct_relevel(group, "control","training"),
time = fct_relevel(time, "pre","post","fu"))
fits_teps_a <- list(
mk_fit(dat_teps_a, c("pre","post"), "pre–post"),
mk_fit(dat_teps_a, c("pre","fu"), "pre–fu"),
mk_fit(dat_teps_a, c("post","fu"), "post–fu")
)
TEPS_anti_fixed <- map_dfr(fits_teps_a, \(x)
fixed_d_from_fit(x$fit) |> mutate(model = x$label, .before = 1))
TEPS_anti_within <- map_dfr(fits_teps_a[1:2], \(x)
within_d_from_fit(x$fit, x$window) |> mutate(model = x$label, .before = 1))
TEPS_anti_fixed
## # A tibble: 12 × 11
## model term estimate std.error df t p d SE_d d_low
## <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 pre–post (Int… 35.6 1.22 75.3 29.2 9.25e-43 10.4 0.355 9.64
## 2 pre–post grou… -0.133 1.73 75.3 -0.0772 9.39e- 1 -0.0388 0.502 -1.04
## 3 pre–post time… -1.47 0.888 58.0 -1.65 1.04e- 1 -0.426 0.258 -0.943
## 4 pre–post grou… 4.27 1.26 58.0 3.40 1.24e- 3 1.24 0.365 0.510
## 5 pre–fu (Int… 35.6 1.19 78.3 29.8 1.64e-44 9.83 0.330 9.18
## 6 pre–fu grou… -0.133 1.69 78.3 -0.0790 9.37e- 1 -0.0368 0.466 -0.965
## 7 pre–fu time… -1.77 0.935 58.0 -1.89 6.37e- 2 -0.488 0.258 -1.00
## 8 pre–fu grou… 4.68 1.32 58.0 3.54 7.96e- 4 1.29 0.365 0.562
## 9 post–fu (Int… 34.1 1.17 74.1 29.3 1.84e-42 10.7 0.366 9.99
## 10 post–fu grou… 4.13 1.65 74.1 2.51 1.44e- 2 1.30 0.518 0.266
## 11 post–fu time… -0.300 0.822 58.0 -0.365 7.16e- 1 -0.0943 0.258 -0.611
## 12 post–fu grou… 0.413 1.16 58.0 0.355 7.24e- 1 0.130 0.365 -0.601
## # ℹ 1 more variable: d_high <dbl>
TEPS_anti_within
## model group contrast estimate SE df t.ratio p.value
## 1 pre–post control post_minus_pre -1.466667 0.8880647 58 -1.651531 0.104034848
## 2 pre–post training post_minus_pre 2.800000 0.8880647 58 3.152923 0.002559407
## 3 pre–fu control fu_minus_pre -1.766667 0.9346467 58 -1.890197 0.063731647
## 4 pre–fu training fu_minus_pre 2.912667 0.9346467 58 3.116329 0.002846887
## d SE_d d_low d_high
## 1 -0.4264235 0.2581989 -0.9432648 0.09041770
## 2 0.8140813 0.2581989 0.2972401 1.33092253
## 3 -0.4880468 0.2581989 -1.0048881 0.02879439
## 4 0.8046327 0.2581989 0.2877915 1.32147392
### TEPS — consummatory
metric_name <- "TEPS_consummatory"
dat_teps_c <- scales_long |>
mutate(metric = as.character(metric)) |>
filter(metric == metric_name) |>
mutate(group = fct_relevel(group, "control","training"),
time = fct_relevel(time, "pre","post","fu"))
fits_teps_c <- list(
mk_fit(dat_teps_c, c("pre","post"), "pre–post"),
mk_fit(dat_teps_c, c("pre","fu"), "pre–fu"),
mk_fit(dat_teps_c, c("post","fu"), "post–fu")
)
TEPS_cons_fixed <- map_dfr(fits_teps_c, \(x)
fixed_d_from_fit(x$fit) |> mutate(model = x$label, .before = 1))
TEPS_cons_within <- map_dfr(fits_teps_c[1:2], \(x)
within_d_from_fit(x$fit, x$window) |> mutate(model = x$label, .before = 1))
TEPS_cons_fixed
## # A tibble: 12 × 11
## model term estimate std.error df t p d SE_d d_low
## <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 pre–post (Inte… 45.1 1.31 78.2 34.5 4.87e-49 11.4 0.331 10.7
## 2 pre–post group… -3.43 1.85 78.2 -1.85 6.75e- 2 -0.868 0.468 -1.80
## 3 pre–post timep… -1.90 1.02 58.0 -1.86 6.80e- 2 -0.480 0.258 -0.997
## 4 pre–post group… 5.47 1.44 58.0 3.78 3.68e- 4 1.38 0.365 0.651
## 5 pre–fu (Inte… 45.1 1.32 78.7 34.1 6.55e-49 11.1 0.327 10.5
## 6 pre–fu group… -3.43 1.87 78.7 -1.83 7.05e- 2 -0.847 0.462 -1.77
## 7 pre–fu timefu -1.70 1.05 58.0 -1.62 1.10e- 1 -0.420 0.258 -0.936
## 8 pre–fu group… 5.10 1.48 58.0 3.45 1.05e- 3 1.26 0.365 0.529
## 9 post–fu (Inte… 43.2 1.25 69.5 34.7 1.39e-45 14.9 0.428 14.0
## 10 post–fu group… 2.03 1.76 69.5 1.15 2.53e- 1 0.698 0.606 -0.510
## 11 post–fu timefu 0.200 0.752 58.0 0.266 7.91e- 1 0.0687 0.258 -0.448
## 12 post–fu group… -0.363 1.06 58.0 -0.342 7.34e- 1 -0.125 0.365 -0.856
## # ℹ 1 more variable: d_high <dbl>
TEPS_cons_within
## model group contrast estimate SE df t.ratio p.value
## 1 pre–post control post_minus_pre -1.900000 1.021624 58 -1.859783 0.0679900292
## 2 pre–post training post_minus_pre 3.566667 1.021624 58 3.491172 0.0009262945
## 3 pre–fu control fu_minus_pre -1.700000 1.046183 58 -1.624954 0.1095947847
## 4 pre–fu training fu_minus_pre 3.403333 1.046183 58 3.253094 0.0019056010
## d SE_d d_low d_high
## 1 -0.4801940 0.2581989 -0.9970352 0.03664726
## 2 0.9014168 0.2581989 0.3845755 1.41825799
## 3 -0.4195613 0.2581989 -0.9364026 0.09727991
## 4 0.8399453 0.2581989 0.3231041 1.35678655
metric_name <- "BDI"
dat_bdi <- scales_long %>%
filter(metric == metric_name) %>%
mutate(group = fct_relevel(group, "control","training"),
time = fct_relevel(time, "pre","post","fu"))
# 三模型
fits_bdi <- list(
mk_fit(dat_bdi, c("pre","post"), "pre-post"),
mk_fit(dat_bdi, c("pre","fu"), "pre-fu"),
mk_fit(dat_bdi, c("post","fu"), "post-fu")
)
# 固定效应 t/SE/p/d
BDI_fixed <- map_dfr(fits_bdi, function(x){
fixed_d_from_fit(x$fit) %>%
mutate(model = x$label, .before = 1)
})
# 组内简单效应（ pre-post 与 pre-fu）
BDI_within <- map_dfr(fits_bdi[1:2], function(x){
within_d_from_fit(x$fit, x$window) %>%
mutate(model = x$label, .before = 1)
})
# 查看结果
BDI_fixed
## # A tibble: 12 × 11
## model term estimate std.error df t p d SE_d d_low
## <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 pre-post (Inte… 24.6 1.14 102. 21.6 7.15e-40 4.97 0.230 4.51
## 2 pre-post group… 0.267 1.61 102. 0.165 8.69e- 1 0.0538 0.325 -0.591
## 3 pre-post timep… -1.73 1.28 58.0 -1.35 1.81e- 1 -0.350 0.258 -0.867
## 4 pre-post group… -5.83 1.81 58.0 -3.22 2.08e- 3 -1.18 0.365 -1.91
## 5 pre-fu (Inte… 24.6 1.12 109. 21.9 9.83e-42 4.64 0.212 4.22
## 6 pre-fu group… 0.267 1.59 109. 0.168 8.67e- 1 0.0502 0.299 -0.543
## 7 pre-fu timefu -2.43 1.37 58.0 -1.77 8.13e- 2 -0.458 0.258 -0.975
## 8 pre-fu group… -7.64 1.94 58.0 -3.94 2.22e- 4 -1.44 0.365 -2.17
## 9 post-fu (Inte… 22.9 1.42 80.8 16.1 5.52e-27 5.05 0.313 4.43
## 10 post-fu group… -5.57 2.01 80.8 -2.77 6.90e- 3 -1.23 0.443 -2.11
## 11 post-fu timefu -0.700 1.17 58.0 -0.598 5.52e- 1 -0.154 0.258 -0.671
## 12 post-fu group… -1.81 1.65 58.0 -1.09 2.80e- 1 -0.399 0.365 -1.13
## # ℹ 1 more variable: d_high <dbl>
BDI_within
## model group contrast estimate SE df t.ratio
## 1 pre-post control post_minus_pre -1.733333 1.279645 58 -1.354543
## 2 pre-post training post_minus_pre -7.566667 1.279645 58 -5.913100
## 3 pre-fu control fu_minus_pre -2.433333 1.371416 58 -1.774322
## 4 pre-fu training fu_minus_pre -10.072333 1.371416 58 -7.344479
## p.value d SE_d d_low d_high
## 1 1.808167e-01 -0.3497414 0.2581989 -0.8665827 0.16709981
## 2 1.894887e-07 -1.5267558 0.2581989 -2.0435971 -1.00991461
## 3 8.125838e-02 -0.4581280 0.2581989 -0.9749693 0.05871318
## 4 7.728579e-10 -1.8963363 0.2581989 -2.4131776 -1.37949510
metric_name <- "STAIT"
dat_stait <- scales_long %>%
dplyr::filter(metric == metric_name) %>%
dplyr::mutate(
group = forcats::fct_relevel(group, "control","training"),
time = forcats::fct_relevel(time, "pre","post","fu")
)
# 三个时间窗口模型
fits_stait <- list(
mk_fit(dat_stait, c("pre","post"), "pre-post"),
mk_fit(dat_stait, c("pre","fu"), "pre-fu"),
mk_fit(dat_stait, c("post","fu"), "post-fu")
)
# 固定效应：t / SE / p / d
STAIT_fixed <- purrr::map_dfr(fits_stait, function(x){
fixed_d_from_fit(x$fit) %>%
dplyr::mutate(model = x$label, .before = 1)
})
# 组内简单效应（ pre-post 与 pre-fu）
STAIT_within <- map_dfr(fits_stait[1:2], function(x){
within_d_from_fit(x$fit, x$window) %>%
dplyr::mutate(model = x$label, .before = 1)
})
# 查看结果
STAIT_fixed
## # A tibble: 12 × 11
## model term estimate std.error df t p d SE_d d_low
## <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 pre-post (Int… 55.5 1.41 91.0 39.4 5.50e-59 10.4 0.265 9.90
## 2 pre-post grou… 0.167 1.99 91.0 0.0837 9.33e- 1 0.0313 0.374 -0.712
## 3 pre-post time… 0.633 1.37 58.0 0.461 6.47e- 1 0.119 0.258 -0.398
## 4 pre-post grou… -5.23 1.94 58.0 -2.69 9.25e- 3 -0.983 0.365 -1.71
## 5 pre-fu (Int… 55.5 1.35 92.3 41.0 4.24e-61 10.7 0.260 10.1
## 6 pre-fu grou… 0.167 1.91 92.3 0.0872 9.31e- 1 0.0320 0.367 -0.698
## 7 pre-fu time… 0.767 1.34 58.0 0.571 5.70e- 1 0.147 0.258 -0.369
## 8 pre-fu grou… -7.26 1.90 58.0 -3.82 3.25e- 4 -1.40 0.365 -2.13
## 9 post-fu (Int… 56.1 1.51 77.9 37.2 2.30e-51 12.4 0.333 11.7
## 10 post-fu grou… -5.07 2.13 77.9 -2.38 2.00e- 2 -1.12 0.471 -2.06
## 11 post-fu time… 0.133 1.17 58.0 0.114 9.10e- 1 0.0294 0.258 -0.487
## 12 post-fu grou… -2.03 1.65 58.0 -1.23 2.25e- 1 -0.447 0.365 -1.18
## # ℹ 1 more variable: d_high <dbl>
STAIT_within
## model group contrast estimate SE df t.ratio
## 1 pre-post control post_minus_pre 0.6333333 1.374264 58 0.4608527
## 2 pre-post training post_minus_pre -4.6000000 1.374264 58 -3.3472463
## 3 pre-fu control fu_minus_pre 0.7666667 1.342964 58 0.5708767
## 4 pre-fu training fu_minus_pre -6.4930000 1.342964 58 -4.8348294
## p.value d SE_d d_low d_high
## 1 0.6466277052 0.1189917 0.2581989 -0.3978496 0.6358329
## 2 0.0014374416 -0.8642553 0.2581989 -1.3810965 -0.3474140
## 3 0.5702885762 0.1473997 0.2581989 -0.3694415 0.6642410
## 4 0.0000101826 -1.2483476 0.2581989 -1.7651888 -0.7315064
metric_name <- "RRS"
dat_rrs <- scales_long %>%
dplyr::filter(metric == metric_name) %>%
dplyr::mutate(
group = forcats::fct_relevel(group, "control","training"),
time = forcats::fct_relevel(time, "pre","post","fu")
)
# 三模型
fits_rrs <- list(
mk_fit(dat_rrs, c("pre","post"), "pre-post"),
mk_fit(dat_rrs, c("pre","fu"), "pre-fu"),
mk_fit(dat_rrs, c("post","fu"), "post-fu")
)
# 固定效应 t/SE/p/d
RRS_fixed <- purrr::map_dfr(fits_rrs, function(x){
fixed_d_from_fit(x$fit) %>%
dplyr::mutate(model = x$label, .before = 1)
}) %>%
dplyr::select(model, term, t, std.error, p, d)
# 组内简单效应（pre-post 与 pre-fu）
RRS_within <- purrr::map_dfr(fits_rrs[1:2], function(x){
within_d_from_fit(x$fit, x$window) %>%
dplyr::mutate(model = x$label, .before = 1)
}) %>%
dplyr::select(model, group, contrast, t.ratio, SE, p.value, d)
RRS_fixed
## # A tibble: 12 × 6
## model term t std.error p d
## <chr> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 pre-post (Intercept) 3.37e+ 1 1.69 1.62e-52 9.13e+ 0
## 2 pre-post grouptraining 4.74e- 1 2.39 6.37e- 1 1.82e- 1
## 3 pre-post timepost -4.97e- 1 1.61 6.21e- 1 -1.28e- 1
## 4 pre-post grouptraining:timepost -1.71e+ 0 2.28 9.22e- 2 -6.25e- 1
## 5 pre-fu (Intercept) 3.48e+ 1 1.63 2.78e-57 8.34e+ 0
## 6 pre-fu grouptraining 4.90e- 1 2.31 6.25e- 1 1.66e- 1
## 7 pre-fu timefu -4.54e- 1 1.76 6.52e- 1 -1.17e- 1
## 8 pre-fu grouptraining:timefu -4.12e+ 0 2.49 1.23e- 4 -1.50e+ 0
## 9 post-fu (Intercept) 3.17e+ 1 1.77 1.40e-48 9.31e+ 0
## 10 post-fu grouptraining -1.10e+ 0 2.51 2.73e- 1 -4.59e- 1
## 11 post-fu timefu 6.25e-15 1.56 1.00e+ 0 1.61e-15
## 12 post-fu grouptraining:timefu -2.89e+ 0 2.20 5.42e- 3 -1.06e+ 0
RRS_within
## model group contrast t.ratio SE p.value d
## 1 pre-post control post_minus_pre -0.4967377 1.610508 6.212516e-01 -0.1282571
## 2 pre-post training post_minus_pre -2.9183341 1.610508 5.000843e-03 -0.7535106
## 3 pre-fu control fu_minus_pre -0.4536655 1.763414 6.517626e-01 -0.1171359
## 4 pre-fu training fu_minus_pre -6.2749497 1.763414 4.778693e-08 -1.6201850
erq_targets <- c("ERQ_cognitive reappraisal", "ERQ_expressive suppression")
ERQ_fixed <- list()
ERQ_within <- list()
for (metric_name in erq_targets) {
dat_erq <- scales_long %>%
dplyr::filter(metric == metric_name) %>%
dplyr::mutate(
group = forcats::fct_relevel(group, "control","training"),
time = forcats::fct_relevel(time, "pre","post","fu")
)
fits_erq <- list(
mk_fit(dat_erq, c("pre","post"), "pre-post"),
mk_fit(dat_erq, c("pre","fu"), "pre-fu"),
mk_fit(dat_erq, c("post","fu"), "post-fu")
)
# 固定效应：t / SE / p / d
fx <- purrr::map_dfr(fits_erq, function(x){
fixed_d_from_fit(x$fit) %>%
dplyr::mutate(metric = metric_name, model = x$label, .before = 1)
})
ERQ_fixed[[metric_name]] <- fx
# 组内简单效应（只做 pre–post / pre–fu）
wx <- purrr::map_dfr(fits_erq[1:2], function(x){
within_d_from_fit(x$fit, x$window) %>%
dplyr::mutate(metric = metric_name, model = x$label, .before = 1)
})
ERQ_within[[metric_name]] <- wx
}
# 合并&只保留关心的列
ERQ_fixed <- dplyr::bind_rows(ERQ_fixed) %>%
dplyr::select(metric, model, term, t, std.error, p, d) %>%
dplyr::rename(SE = std.error, p.value = p)
ERQ_within <- dplyr::bind_rows(ERQ_within) %>%
dplyr::select(metric, model, group, contrast, t.ratio, SE, p.value, d) %>%
dplyr::rename(t = t.ratio)
# 查看结果
ERQ_fixed
## # A tibble: 24 × 7
## metric model term t SE p.value d
## <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 ERQ_cognitive reappraisal pre-post (Intercept) 30.3 0.995 5.37e-45 10.0
## 2 ERQ_cognitive reappraisal pre-post grouptraini… -0.592 1.41 5.55e- 1 -0.277
## 3 ERQ_cognitive reappraisal pre-post timepost -0.900 0.778 3.72e- 1 -0.232
## 4 ERQ_cognitive reappraisal pre-post grouptraini… 3.06 1.10 3.35e- 3 1.12
## 5 ERQ_cognitive reappraisal pre-fu (Intercept) 32.0 0.941 5.96e-51 8.60
## 6 ERQ_cognitive reappraisal pre-fu grouptraini… -0.626 1.33 5.33e- 1 -0.238
## 7 ERQ_cognitive reappraisal pre-fu timefu 0.221 0.905 8.26e- 1 0.0571
## 8 ERQ_cognitive reappraisal pre-fu grouptraini… 1.39 1.28 1.69e- 1 0.508
## 9 ERQ_cognitive reappraisal post-fu (Intercept) 31.6 0.932 1.83e-45 11.0
## 10 ERQ_cognitive reappraisal post-fu grouptraini… 1.92 1.32 5.83e- 2 0.946
## # ℹ 14 more rows
ERQ_within
## metric model group contrast t
## 1 ERQ_cognitive reappraisal pre-post control post_minus_pre -0.89984012
## 2 ERQ_cognitive reappraisal pre-post training post_minus_pre 3.42796237
## 3 ERQ_cognitive reappraisal pre-fu control fu_minus_pre 0.22099615
## 4 ERQ_cognitive reappraisal pre-fu training fu_minus_pre 2.19007188
## 5 ERQ_expressive suppression pre-post control post_minus_pre -1.14442421
## 6 ERQ_expressive suppression pre-post training post_minus_pre 0.66758079
## 7 ERQ_expressive suppression pre-fu control fu_minus_pre -0.04595728
## 8 ERQ_expressive suppression pre-fu training fu_minus_pre -0.83963946
## SE p.value d
## 1 0.7779160 0.371926318 -0.23233772
## 2 0.7779160 0.001124862 0.88509608
## 3 0.9049931 0.825871586 0.05706096
## 4 0.9049931 0.032553728 0.56547413
## 5 0.6990415 0.257148724 -0.29548906
## 6 0.6990415 0.507048884 0.17236862
## 7 0.7253113 0.963502174 -0.01186612
## 8 0.7253113 0.404558268 -0.21679398