# =============================================================================
# PAIN CATASTROPHIZING AS MEDIATOR OF EAET EFFECTS
# Research Question: Does EAET reduce pain by reducing catastrophizing?
# =============================================================================
library(tidyverse)
library(mediation)
library(psych)
library(effectsize)
##
## Attaching package: 'effectsize'
## The following object is masked from 'package:mvtnorm':
##
## standardize
## The following object is masked from 'package:psych':
##
## phi
library(ggplot2)
# Load your prepared data
df_6m <- readRDS("data_post_prepared.rds") %>%
filter(!is.na(meanpain_6m)) # Keep only those with 6-month data
# =============================================================================
# STEP 1: CREATE CATASTROPHIZING CHANGE VARIABLE
# =============================================================================
df_6m <- df_6m %>%
mutate(
# Catastrophizing change from baseline to post-treatment
PCS_change_post = PCS_pT - PCS_bT,
# Also calculate baseline to 6-month change for sensitivity analysis
PCS_change_6m = PCS_6mT - PCS_bT
)
cat("\n=== PAIN CATASTROPHIZING MEDIATION ANALYSIS ===\n\n")
##
## === PAIN CATASTROPHIZING MEDIATION ANALYSIS ===
cat("Sample size:", nrow(df_6m), "participants with 6-month data\n")
## Sample size: 104 participants with 6-month data
cat("Complete PCS data:", sum(!is.na(df_6m$PCS_change_post)), "participants\n\n")
## Complete PCS data: 104 participants
# =============================================================================
# STEP 2: DESCRIPTIVE STATISTICS - CATASTROPHIZING BY TREATMENT
# =============================================================================
cat("=== DESCRIPTIVE STATISTICS ===\n\n")
## === DESCRIPTIVE STATISTICS ===
# Baseline catastrophizing by treatment
cat("BASELINE Pain Catastrophizing:\n")
## BASELINE Pain Catastrophizing:
df_6m %>%
group_by(treatment) %>%
summarise(
N = sum(!is.na(PCS_bT)),
Mean = mean(PCS_bT, na.rm = TRUE),
SD = sd(PCS_bT, na.rm = TRUE),
.groups = "drop"
) %>%
print()
## # A tibble: 2 × 4
## treatment N Mean SD
## <fct> <int> <dbl> <dbl>
## 1 CBT 50 23.2 12.6
## 2 EAET 54 22.1 13.1
# Post-treatment catastrophizing by treatment
cat("\nPOST-TREATMENT Pain Catastrophizing:\n")
##
## POST-TREATMENT Pain Catastrophizing:
df_6m %>%
group_by(treatment) %>%
summarise(
N = sum(!is.na(PCS_pT)),
Mean = mean(PCS_pT, na.rm = TRUE),
SD = sd(PCS_pT, na.rm = TRUE),
.groups = "drop"
) %>%
print()
## # A tibble: 2 × 4
## treatment N Mean SD
## <fct> <int> <dbl> <dbl>
## 1 CBT 50 21.0 13.3
## 2 EAET 54 18.4 11.6
# Change in catastrophizing by treatment
cat("\nCHANGE in Pain Catastrophizing (Baseline → Post):\n")
##
## CHANGE in Pain Catastrophizing (Baseline → Post):
df_6m %>%
group_by(treatment) %>%
summarise(
N = sum(!is.na(PCS_change_post)),
Mean = mean(PCS_change_post, na.rm = TRUE),
SD = sd(PCS_change_post, na.rm = TRUE),
.groups = "drop"
) %>%
print()
## # A tibble: 2 × 4
## treatment N Mean SD
## <fct> <int> <dbl> <dbl>
## 1 CBT 50 -2.2 7.50
## 2 EAET 54 -3.70 11.3
# =============================================================================
# STEP 3: TEST PATH A - Does EAET Reduce Catastrophizing?
# =============================================================================
cat("\n\n=== PATH A: TREATMENT → CATASTROPHIZING CHANGE ===\n\n")
##
##
## === PATH A: TREATMENT → CATASTROPHIZING CHANGE ===
# Independent samples t-test
t_pcs <- t.test(PCS_change_post ~ treatment, data = df_6m)
d_pcs <- cohens_d(PCS_change_post ~ treatment, data = df_6m)$Cohens_d
cat(sprintf("Does EAET reduce catastrophizing more than CBT?\n"))
## Does EAET reduce catastrophizing more than CBT?
cat(sprintf("t(%.1f) = %.2f, p = %.3f, Cohen's d = %.2f\n\n",
t_pcs$parameter, t_pcs$statistic, t_pcs$p.value, d_pcs))
## t(92.6) = 0.80, p = 0.424, Cohen's d = 0.16
if(t_pcs$p.value < 0.05) {
cat("✓ SIGNIFICANT: EAET reduces catastrophizing more than CBT!\n")
cat(" This is the 'a path' - necessary for mediation.\n\n")
} else if(t_pcs$p.value < 0.10) {
cat("~ MARGINAL: EAET shows trend toward reducing catastrophizing.\n")
cat(" Continue to test full mediation model.\n\n")
} else {
cat("✗ NON-SIGNIFICANT: No difference in catastrophizing reduction.\n")
cat(" Mediation is unlikely, but test to be thorough.\n\n")
}
## ✗ NON-SIGNIFICANT: No difference in catastrophizing reduction.
## Mediation is unlikely, but test to be thorough.
# =============================================================================
# STEP 4: MEDIATION ANALYSIS - OUTCOME 1: PAIN SEVERITY AT 6 MONTHS
# =============================================================================
cat("\n")
cat("╔════════════════════════════════════════════════════════════════╗\n")
## ╔════════════════════════════════════════════════════════════════╗
cat("║ MEDIATION ANALYSIS: PAIN SEVERITY (6 MONTHS) ║\n")
## ║ MEDIATION ANALYSIS: PAIN SEVERITY (6 MONTHS) ║
cat("╚════════════════════════════════════════════════════════════════╝\n\n")
## ╚════════════════════════════════════════════════════════════════╝
cat("Pathway: Treatment → Catastrophizing Change → 6-month Pain\n\n")
## Pathway: Treatment → Catastrophizing Change → 6-month Pain
# Step 1: Treatment → Catastrophizing (Path a)
cat("Step 1: Treatment → Catastrophizing Change (Path a)\n")
## Step 1: Treatment → Catastrophizing Change (Path a)
model_tx_pcs_pain <- lm(PCS_change_post ~ eaet + PCS_bT + meanpain_b +
age + female + race_white_vs_bipoc,
data = df_6m)
summary(model_tx_pcs_pain)
##
## Call:
## lm(formula = PCS_change_post ~ eaet + PCS_bT + meanpain_b + age +
## female + race_white_vs_bipoc, data = df_6m)
##
## Residuals:
## Min 1Q Median 3Q Max
## -33.281 -5.552 1.406 5.815 18.967
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 9.95014 12.98251 0.766 0.445284
## eaet -1.36945 1.75848 -0.779 0.438012
## PCS_bT -0.33360 0.08243 -4.047 0.000104 ***
## meanpain_b 0.31404 0.60711 0.517 0.606149
## age -0.08498 0.16540 -0.514 0.608592
## female -6.13583 3.14109 -1.953 0.053652 .
## race_white_vs_bipoc 0.10290 2.04858 0.050 0.960041
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.821 on 97 degrees of freedom
## Multiple R-squared: 0.2155, Adjusted R-squared: 0.167
## F-statistic: 4.441 on 6 and 97 DF, p-value: 0.0005189
# Step 2: Treatment + Catastrophizing → Pain Outcome (Paths b and c')
cat("\n\nStep 2: Treatment + Catastrophizing → 6-month Pain (Paths b and c')\n")
##
##
## Step 2: Treatment + Catastrophizing → 6-month Pain (Paths b and c')
model_full_pain <- lm(meanpain_6m ~ eaet + PCS_change_post + PCS_bT + meanpain_b +
age + female + race_white_vs_bipoc,
data = df_6m)
summary(model_full_pain)
##
## Call:
## lm(formula = meanpain_6m ~ eaet + PCS_change_post + PCS_bT +
## meanpain_b + age + female + race_white_vs_bipoc, data = df_6m)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.2902 -0.9112 0.0320 1.0184 3.8838
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.96785 2.44175 1.625 0.107441
## eaet -1.03711 0.33077 -3.135 0.002277 **
## PCS_change_post 0.03764 0.01904 1.977 0.050899 .
## PCS_bT 0.04247 0.01671 2.541 0.012641 *
## meanpain_b 0.40902 0.11400 3.588 0.000527 ***
## age -0.01950 0.03106 -0.628 0.531647
## female -0.27931 0.60047 -0.465 0.642879
## race_white_vs_bipoc -0.07366 0.38414 -0.192 0.848344
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.654 on 96 degrees of freedom
## Multiple R-squared: 0.3861, Adjusted R-squared: 0.3413
## F-statistic: 8.625 on 7 and 96 DF, p-value: 0.0000000356
# Bootstrap mediation test
cat("\n\nBootstrap Mediation Analysis (5000 iterations):\n")
##
##
## Bootstrap Mediation Analysis (5000 iterations):
cat("This will take 1-2 minutes...\n\n")
## This will take 1-2 minutes...
set.seed(123) # For reproducibility
med_pcs_pain <- mediate(model_tx_pcs_pain, model_full_pain,
treat = "eaet", mediator = "PCS_change_post",
boot = TRUE, sims = 5000, boot.ci.type = "perc")
## Running nonparametric bootstrap
summary(med_pcs_pain)
##
## Causal Mediation Analysis
##
## Nonparametric Bootstrap Confidence Intervals with the Percentile Method
##
## Estimate 95% CI Lower 95% CI Upper p-value
## ACME -0.051549 -0.217971 0.085647 0.4956
## ADE -1.037112 -1.693850 -0.406765 0.0036 **
## Total Effect -1.088661 -1.756496 -0.424975 0.0028 **
## Prop. Mediated 0.047351 -0.097561 0.212622 0.4952
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Sample Size Used: 104
##
##
## Simulations: 5000
# Save results
pcs_pain_results <- data.frame(
Outcome = "Pain Severity (6m)",
Mediator = "Pain Catastrophizing Change",
Path_a = coef(model_tx_pcs_pain)["eaet"],
Path_a_p = summary(model_tx_pcs_pain)$coefficients["eaet", "Pr(>|t|)"],
Path_b = coef(model_full_pain)["PCS_change_post"],
Path_b_p = summary(model_full_pain)$coefficients["PCS_change_post", "Pr(>|t|)"],
ACME = med_pcs_pain$z0,
ACME_CI_lower = med_pcs_pain$z0.ci[1],
ACME_CI_upper = med_pcs_pain$z0.ci[2],
ACME_p = med_pcs_pain$z0.p,
Direct_Effect = med_pcs_pain$d0,
Direct_p = med_pcs_pain$d0.p,
Total_Effect = med_pcs_pain$tau.coef,
Total_p = med_pcs_pain$tau.p,
Prop_Mediated = med_pcs_pain$n0
)
# =============================================================================
# STEP 5: MEDIATION ANALYSIS - OUTCOME 2: PAIN INTERFERENCE AT 6 MONTHS
# =============================================================================
cat("\n\n")
cat("╔════════════════════════════════════════════════════════════════╗\n")
## ╔════════════════════════════════════════════════════════════════╗
cat("║ MEDIATION ANALYSIS: PAIN INTERFERENCE (6 MONTHS) ║\n")
## ║ MEDIATION ANALYSIS: PAIN INTERFERENCE (6 MONTHS) ║
cat("╚════════════════════════════════════════════════════════════════╝\n\n")
## ╚════════════════════════════════════════════════════════════════╝
# Models
model_tx_pcs_painint <- lm(PCS_change_post ~ eaet + PCS_bT + painint_b +
age + female + race_white_vs_bipoc,
data = df_6m)
model_full_painint <- lm(painint_6m ~ eaet + PCS_change_post + PCS_bT + painint_b +
age + female + race_white_vs_bipoc,
data = df_6m)
cat("Bootstrap Mediation Analysis:\n")
## Bootstrap Mediation Analysis:
set.seed(124)
med_pcs_painint <- mediate(model_tx_pcs_painint, model_full_painint,
treat = "eaet", mediator = "PCS_change_post",
boot = TRUE, sims = 5000, boot.ci.type = "perc")
## Running nonparametric bootstrap
summary(med_pcs_painint)
##
## Causal Mediation Analysis
##
## Nonparametric Bootstrap Confidence Intervals with the Percentile Method
##
## Estimate 95% CI Lower 95% CI Upper p-value
## ACME -0.20534 -0.72705 0.40006 0.5128
## ADE -0.04097 -2.59480 2.17408 0.8944
## Total Effect -0.24632 -2.82982 2.05898 0.7912
## Prop. Mediated 0.83367 -2.95851 3.11043 0.7856
##
## Sample Size Used: 104
##
##
## Simulations: 5000
pcs_painint_results <- data.frame(
Outcome = "Pain Interference (6m)",
Mediator = "Pain Catastrophizing Change",
Path_a = coef(model_tx_pcs_painint)["eaet"],
Path_a_p = summary(model_tx_pcs_painint)$coefficients["eaet", "Pr(>|t|)"],
Path_b = coef(model_full_painint)["PCS_change_post"],
Path_b_p = summary(model_full_painint)$coefficients["PCS_change_post", "Pr(>|t|)"],
ACME = med_pcs_painint$z0,
ACME_CI_lower = med_pcs_painint$z0.ci[1],
ACME_CI_upper = med_pcs_painint$z0.ci[2],
ACME_p = med_pcs_painint$z0.p,
Direct_Effect = med_pcs_painint$d0,
Direct_p = med_pcs_painint$d0.p,
Total_Effect = med_pcs_painint$tau.coef,
Total_p = med_pcs_painint$tau.p,
Prop_Mediated = med_pcs_painint$n0
)
# =============================================================================
# STEP 6: MEDIATION ANALYSIS - OUTCOME 3: DEPRESSION AT 6 MONTHS
# =============================================================================
cat("\n\n")
cat("╔════════════════════════════════════════════════════════════════╗\n")
## ╔════════════════════════════════════════════════════════════════╗
cat("║ MEDIATION ANALYSIS: DEPRESSION (6 MONTHS) ║\n")
## ║ MEDIATION ANALYSIS: DEPRESSION (6 MONTHS) ║
cat("╚════════════════════════════════════════════════════════════════╝\n\n")
## ╚════════════════════════════════════════════════════════════════╝
# Models
model_tx_pcs_depress <- lm(PCS_change_post ~ eaet + PCS_bT + depress_b +
age + female + race_white_vs_bipoc,
data = df_6m)
model_full_depress <- lm(depress_6m ~ eaet + PCS_change_post + PCS_bT + depress_b +
age + female + race_white_vs_bipoc,
data = df_6m)
cat("Bootstrap Mediation Analysis:\n")
## Bootstrap Mediation Analysis:
set.seed(125)
med_pcs_depress <- mediate(model_tx_pcs_depress, model_full_depress,
treat = "eaet", mediator = "PCS_change_post",
boot = TRUE, sims = 5000, boot.ci.type = "perc")
## Running nonparametric bootstrap
## Warning in predict.lm(new.fit.M, type = "response", newdata = pred.data.t):
## prediction from rank-deficient fit; attr(*, "non-estim") has doubtful cases
## Warning in predict.lm(new.fit.M, type = "response", newdata = pred.data.c):
## prediction from rank-deficient fit; attr(*, "non-estim") has doubtful cases
## Warning in predict.lm(new.fit.Y, type = "response", newdata = pred.data.t):
## prediction from rank-deficient fit; attr(*, "non-estim") has doubtful cases
## Warning in predict.lm(new.fit.Y, type = "response", newdata = pred.data.c):
## prediction from rank-deficient fit; attr(*, "non-estim") has doubtful cases
## Warning in predict.lm(new.fit.Y, type = "response", newdata = pred.data.t):
## prediction from rank-deficient fit; attr(*, "non-estim") has doubtful cases
## Warning in predict.lm(new.fit.Y, type = "response", newdata = pred.data.c):
## prediction from rank-deficient fit; attr(*, "non-estim") has doubtful cases
## Warning in predict.lm(new.fit.Y, type = "response", newdata = pred.data.t):
## prediction from rank-deficient fit; attr(*, "non-estim") has doubtful cases
## Warning in predict.lm(new.fit.Y, type = "response", newdata = pred.data.c):
## prediction from rank-deficient fit; attr(*, "non-estim") has doubtful cases
## Warning in predict.lm(new.fit.Y, type = "response", newdata = pred.data.t):
## prediction from rank-deficient fit; attr(*, "non-estim") has doubtful cases
## Warning in predict.lm(new.fit.Y, type = "response", newdata = pred.data.c):
## prediction from rank-deficient fit; attr(*, "non-estim") has doubtful cases
summary(med_pcs_depress)
##
## Causal Mediation Analysis
##
## Nonparametric Bootstrap Confidence Intervals with the Percentile Method
##
## Estimate 95% CI Lower 95% CI Upper p-value
## ACME -0.249232 -1.050528 0.294656 0.4608
## ADE -1.775088 -3.782688 0.187286 0.0800 .
## Total Effect -2.024321 -4.070380 0.010346 0.0524 .
## Prop. Mediated 0.123119 -0.326884 0.862479 0.4644
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Sample Size Used: 104
##
##
## Simulations: 5000
pcs_depress_results <- data.frame(
Outcome = "Depression (6m)",
Mediator = "Pain Catastrophizing Change",
Path_a = coef(model_tx_pcs_depress)["eaet"],
Path_a_p = summary(model_tx_pcs_depress)$coefficients["eaet", "Pr(>|t|)"],
Path_b = coef(model_full_depress)["PCS_change_post"],
Path_b_p = summary(model_full_depress)$coefficients["PCS_change_post", "Pr(>|t|)"],
ACME = med_pcs_depress$z0,
ACME_CI_lower = med_pcs_depress$z0.ci[1],
ACME_CI_upper = med_pcs_depress$z0.ci[2],
ACME_p = med_pcs_depress$z0.p,
Direct_Effect = med_pcs_depress$d0,
Direct_p = med_pcs_depress$d0.p,
Total_Effect = med_pcs_depress$tau.coef,
Total_p = med_pcs_depress$tau.p,
Prop_Mediated = med_pcs_depress$n0
)
# =============================================================================
# STEP 8: VISUALIZATIONS
# =============================================================================
cat("\n\n=== CREATING VISUALIZATIONS ===\n\n")
##
##
## === CREATING VISUALIZATIONS ===
# Plot 1: Catastrophizing change by treatment
p1 <- ggplot(df_6m, aes(x = treatment, y = PCS_change_post, fill = treatment)) +
geom_boxplot(alpha = 0.7, outlier.shape = NA) +
geom_jitter(width = 0.2, alpha = 0.3, size = 2) +
stat_summary(fun = mean, geom = "point", shape = 23, size = 4,
fill = "red", color = "darkred") +
geom_hline(yintercept = 0, linetype = "dashed", color = "gray50") +
scale_fill_manual(values = c("CBT" = "#4393C3", "EAET" = "#D6604D")) +
labs(title = "Change in Pain Catastrophizing by Treatment",
subtitle = sprintf("EAET: M=%.2f, CBT: M=%.2f (d=%.2f, p=%.3f)",
mean(df_6m$PCS_change_post[df_6m$treatment == "EAET"], na.rm=T),
mean(df_6m$PCS_change_post[df_6m$treatment == "CBT"], na.rm=T),
d_pcs, t_pcs$p.value),
x = "Treatment",
y = "Change in PCS (Baseline → Post)\nNegative = Improvement",
caption = "Diamond = mean, box = median and IQR") +
theme_minimal(base_size = 14) +
theme(legend.position = "none",
plot.title = element_text(face = "bold", size = 16),
plot.subtitle = element_text(size = 12))
print(p1)
ggsave("catastrophizing_change_by_treatment.png", p1, width = 8, height = 6, dpi = 300)
# Plot 2: Mediation diagram for pain severity (if significant)
if(pcs_pain_results$ACME_p < 0.10) {
library(DiagrammeR)
# Create mediation diagram
med_diagram <- grViz("
digraph mediation {
graph [layout = dot, rankdir = LR]
node [shape = box, style = filled, fillcolor = lightblue]
Treatment [label = 'Treatment\n(EAET vs CBT)']
PCS [label = 'Catastrophizing\nChange\n(Baseline→Post)']
Pain [label = 'Pain Severity\n(6 months)']
Treatment -> PCS [label = <<b>a</b> = @PATH_A@<br/>p = @PA_P@>]
PCS -> Pain [label = <<b>b</b> = @PATH_B@<br/>p = @PB_P@>]
Treatment -> Pain [label = <<b>c'</b> = @DIRECT@<br/>p = @DIR_P@>]
{rank = same; Treatment; Pain}
}
" %>%
gsub("@PATH_A@", sprintf("%.2f", pcs_pain_results$Path_a), .) %>%
gsub("@PA_P@", sprintf("%.3f", pcs_pain_results$Path_a_p), .) %>%
gsub("@PATH_B@", sprintf("%.2f", pcs_pain_results$Path_b), .) %>%
gsub("@PB_P@", sprintf("%.3f", pcs_pain_results$Path_b_p), .) %>%
gsub("@DIRECT@", sprintf("%.2f", pcs_pain_results$Direct_Effect), .) %>%
gsub("@DIR_P@", sprintf("%.3f", pcs_pain_results$Direct_p), .)
)
cat("\nMediation diagram created (view in RStudio Viewer)\n")
}
##
## Mediation diagram created (view in RStudio Viewer)
# Plot 3: Scatter plot - Catastrophizing change vs Pain change
df_6m_plot <- df_6m %>%
mutate(pain_change_6m = meanpain_6m - meanpain_b)
p3 <- ggplot(df_6m_plot, aes(x = PCS_change_post, y = pain_change_6m, color = treatment)) +
geom_point(alpha = 0.6, size = 3) +
geom_smooth(method = "lm", se = TRUE, linewidth = 1.2) +
geom_hline(yintercept = 0, linetype = "dashed", color = "gray50") +
geom_vline(xintercept = 0, linetype = "dashed", color = "gray50") +
scale_color_manual(values = c("CBT" = "#4393C3", "EAET" = "#D6604D")) +
labs(title = "Association: Catastrophizing Reduction → Pain Reduction",
subtitle = "Both negative = improvement (lower-left quadrant)",
x = "Change in Catastrophizing (Baseline → Post)\nNegative = Less Catastrophizing",
y = "Change in Pain (Baseline → 6m)\nNegative = Less Pain",
color = "Treatment") +
theme_minimal(base_size = 14) +
theme(plot.title = element_text(face = "bold", size = 16),
legend.position = "top")
print(p3)
## `geom_smooth()` using formula = 'y ~ x'
ggsave("catastrophizing_pain_association.png", p3, width = 10, height = 6, dpi = 300)
## `geom_smooth()` using formula = 'y ~ x'
# =============================================================================
# STEP 9: SENSITIVITY ANALYSIS - Test 6-month catastrophizing
# =============================================================================
cat("\n\n=== SENSITIVITY ANALYSIS ===\n")
##
##
## === SENSITIVITY ANALYSIS ===
cat("Testing if 6-month catastrophizing (instead of post-treatment) mediates effects\n\n")
## Testing if 6-month catastrophizing (instead of post-treatment) mediates effects
# Only test for pain severity outcome
model_tx_pcs6m <- lm(PCS_change_6m ~ eaet + PCS_bT + meanpain_b +
age + female + race_white_vs_bipoc,
data = df_6m)
model_full_pain6m <- lm(meanpain_6m ~ eaet + PCS_change_6m + PCS_bT + meanpain_b +
age + female + race_white_vs_bipoc,
data = df_6m)
set.seed(126)
med_pcs6m_pain <- mediate(model_tx_pcs6m, model_full_pain6m,
treat = "eaet", mediator = "PCS_change_6m",
boot = TRUE, sims = 5000, boot.ci.type = "perc")
## Running nonparametric bootstrap
cat("Using 6-month catastrophizing change:\n")
## Using 6-month catastrophizing change:
summary(med_pcs6m_pain)
##
## Causal Mediation Analysis
##
## Nonparametric Bootstrap Confidence Intervals with the Percentile Method
##
## Estimate 95% CI Lower 95% CI Upper p-value
## ACME -0.087118 -0.328309 0.106712 0.4064
## ADE -1.001543 -1.601808 -0.373601 0.0016 **
## Total Effect -1.088661 -1.746532 -0.422695 0.0008 ***
## Prop. Mediated 0.080023 -0.140076 0.306031 0.4072
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Sample Size Used: 104
##
##
## Simulations: 5000
# =============================================================================
# STEP 10: FINAL INTERPRETATION & NEXT STEPS
# =============================================================================
cat("\n\n")
cat("╔════════════════════════════════════════════════════════════════════════╗\n")
## ╔════════════════════════════════════════════════════════════════════════╗
cat("║ FINAL SUMMARY & NEXT STEPS ║\n")
## ║ FINAL SUMMARY & NEXT STEPS ║
cat("╚════════════════════════════════════════════════════════════════════════╝\n\n")
## ╚════════════════════════════════════════════════════════════════════════╝
# Count significant findings
n_sig <- sum(c(pcs_pain_results$ACME_p,
pcs_painint_results$ACME_p,
pcs_depress_results$ACME_p) < 0.05)
n_marginal <- sum(c(pcs_pain_results$ACME_p,
pcs_painint_results$ACME_p,
pcs_depress_results$ACME_p) >= 0.05 &
c(pcs_pain_results$ACME_p,
pcs_painint_results$ACME_p,
pcs_depress_results$ACME_p) < 0.10)
cat("RESULTS SUMMARY:\n")
## RESULTS SUMMARY:
cat("================\n")
## ================
cat(sprintf("Significant mediation effects (p < .05): %d of 3\n", n_sig))
## Significant mediation effects (p < .05): 1 of 3
cat(sprintf("Marginal mediation effects (p < .10): %d of 3\n\n", n_marginal))
## Marginal mediation effects (p < .10): 1 of 3