# =============================================================================
# EMOTIONAL AWARENESS AND EXPRESSION THERAPY (EAET) vs CBT
# Emotional Processing Mediation Analysis
# Does EAET reduce pain by increasing emotional processing?
# =============================================================
# =============================================================================
# STEP 1: DATA IMPORT AND SETUP
# =============================================================================
# Load required packages
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.3 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ ggplot2 3.4.4 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.0
## ✔ purrr 1.0.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(readxl)
library(psych)
##
## Attaching package: 'psych'
##
## The following objects are masked from 'package:ggplot2':
##
## %+%, alpha
library(lavaan) # For SEM mediation
## This is lavaan 0.6-16
## lavaan is FREE software! Please report any bugs.
##
## Attaching package: 'lavaan'
##
## The following object is masked from 'package:psych':
##
## cor2cov
library(effectsize) # For Cohen's d
##
## Attaching package: 'effectsize'
##
## The following object is masked from 'package:psych':
##
## phi
options(scipen=999)
# Set working directory and import data
setwd("/Users/melissalagunas/Desktop/Yarns Research Project")
df_raw <- read_excel("EAET_vs_CBT_Study 2023-11-07.xlsx")
## New names:
## • `currentpain_b` -> `currentpain_b...37`
## • `currentpain_b` -> `currentpain_b...42`
# CRITICAL: Filter to only RANDOMIZED participants
df_raw <- df_raw %>%
filter(!is.na(eaet))
cat("=== DATA IMPORT COMPLETE ===\n")
## === DATA IMPORT COMPLETE ===
cat("Total randomized participants:", nrow(df_raw), "\n\n")
## Total randomized participants: 126
# =============================================================================
# STEP 2: CREATE CORE VARIABLES
# =============================================================================
cat("=== CREATING ANALYSIS VARIABLES ===\n\n")
## === CREATING ANALYSIS VARIABLES ===
df <- df_raw %>%
mutate(
# Treatment variable
treatment = factor(eaet,
levels = c(0, 1),
labels = c("CBT", "EAET")),
# Race/Ethnicity (White=0 [reference], BIPOC=1)
race_white_vs_bipoc = case_when(
race_ethnicity == 2 ~ 0, # White = 0
race_ethnicity %in% c(1, 3, 4, 5, 6) ~ 1, # BIPOC = 1
TRUE ~ NA_real_
),
race_group = factor(race_white_vs_bipoc,
levels = c(0, 1),
labels = c("White", "BIPOC")),
# Completion status
completed_6m = if_else(dropout == 4, 1, 0)
)
# Create 6-month analysis dataset
df_6m <- df %>%
filter(!is.na(meanpain_6m)) # Has 6-month outcome data
cat("Sample sizes:\n")
## Sample sizes:
cat(" Total randomized: N =", nrow(df), "\n")
## Total randomized: N = 126
cat(" With 6-month data: N =", nrow(df_6m), "\n")
## With 6-month data: N = 104
cat(" EAET: N =", sum(df_6m$eaet == 1), "\n")
## EAET: N = 54
cat(" CBT: N =", sum(df_6m$eaet == 0), "\n\n")
## CBT: N = 50
# =============================================================================
# STEP 3: CREATE EMOTIONAL PROCESSING VARIABLES
# =============================================================================
cat("=== CREATING EMOTIONAL PROCESSING CHANGE SCORES ===\n\n")
## === CREATING EMOTIONAL PROCESSING CHANGE SCORES ===
df_emotional <- df_6m %>%
mutate(
# -------------------------------------------------------------------------
# 1. EMOTIONAL APPROACH COPING (EAC) - Change from baseline to post
# -------------------------------------------------------------------------
# Expression subscale (4 items, range 4-16)
# Higher = more emotional expression (GOOD)
EAC_Expression_change = EAC_pEX - EAC_bEX,
# Experience subscale (4 items, range 4-16)
# Higher = more emotional experiencing (GOOD)
EAC_Experience_change = EAC_pEP - EAC_bEP,
# Total EAC (8 items, range 8-32)
# Higher = better emotional approach coping (GOOD)
EAC_Total_change = EAC_pT - EAC_bT,
# -------------------------------------------------------------------------
# 2. AMBIVALENCE OVER EMOTIONS (AEQ) - Change from baseline to post
# -------------------------------------------------------------------------
# Positive subscale (8 items, range 8-40)
# LOWER scores = less ambivalence about positive emotions (GOOD)
AEQ_Positive_change = AEQ_pPos - AEQ_bPos,
# Entitlement subscale (6 items, range 6-30)
# LOWER scores = less feeling unentitled to emotions (GOOD)
AEQ_Entitlement_change = AEQ_pEnt - AEQ_bEnt,
# Total AEQ (14 items, range 14-70)
# LOWER total = less overall emotional ambivalence (GOOD)
AEQ_Total_change = AEQ_pT - AEQ_bT,
# -------------------------------------------------------------------------
# 3. SURVEY OF PAIN ATTITUDES (SOPA) - Emotion subscales
# -------------------------------------------------------------------------
# Emotion subscale (2 items, range 0-8)
# Higher = stronger belief that emotions affect pain
SOPA_Emotion_change = SOPA_pEmot - SOPA_bEmot,
# Second emotion subscale (5 items, range 0-20)
SOPA_Emotion2_change = SOPA_pEmot2 - SOPA_bEmot2,
# Third emotion subscale (3 items, range 0-12)
# Higher = stronger awareness of emotion-pain connection
SOPA_Emotion3_change = SOPA_pEmot3 - SOPA_bEmot3
)
# Check data availability
cat("Data availability:\n")
## Data availability:
cat(" Total with 6m outcomes:", nrow(df_emotional), "\n")
## Total with 6m outcomes: 104
cat(" With EAC data:", sum(!is.na(df_emotional$EAC_Total_change)), "\n")
## With EAC data: 104
cat(" With AEQ data:", sum(!is.na(df_emotional$AEQ_Total_change)), "\n")
## With AEQ data: 104
cat(" With SOPA emotion data:", sum(!is.na(df_emotional$SOPA_Emotion3_change)), "\n\n")
## With SOPA emotion data: 104
# =============================================================================
# STEP 4: DESCRIPTIVE STATISTICS - EMOTIONAL PROCESSING BY TREATMENT
# =============================================================================
cat("=== EMOTIONAL PROCESSING CHANGE BY TREATMENT GROUP ===\n\n")
## === EMOTIONAL PROCESSING CHANGE BY TREATMENT GROUP ===
# Descriptive statistics
emotional_descriptives <- df_emotional %>%
group_by(treatment) %>%
summarise(
N = n(),
# EAC - Higher is better
EAC_Expression_M = mean(EAC_Expression_change, na.rm = TRUE),
EAC_Expression_SD = sd(EAC_Expression_change, na.rm = TRUE),
EAC_Experience_M = mean(EAC_Experience_change, na.rm = TRUE),
EAC_Experience_SD = sd(EAC_Experience_change, na.rm = TRUE),
EAC_Total_M = mean(EAC_Total_change, na.rm = TRUE),
EAC_Total_SD = sd(EAC_Total_change, na.rm = TRUE),
# AEQ - Lower (negative change) is better
AEQ_Total_M = mean(AEQ_Total_change, na.rm = TRUE),
AEQ_Total_SD = sd(AEQ_Total_change, na.rm = TRUE),
AEQ_Positive_M = mean(AEQ_Positive_change, na.rm = TRUE),
AEQ_Positive_SD = sd(AEQ_Positive_change, na.rm = TRUE),
# SOPA - Higher is better (more awareness)
SOPA_Emotion3_M = mean(SOPA_Emotion3_change, na.rm = TRUE),
SOPA_Emotion3_SD = sd(SOPA_Emotion3_change, na.rm = TRUE),
.groups = "drop"
)
print(emotional_descriptives)
## # A tibble: 2 × 14
## treatment N EAC_Expression_M EAC_Expression_SD EAC_Experience_M
## <fct> <int> <dbl> <dbl> <dbl>
## 1 CBT 50 -0.34 2.52 -0.0600
## 2 EAET 54 -0.0741 2.73 0.204
## # ℹ 9 more variables: EAC_Experience_SD <dbl>, EAC_Total_M <dbl>,
## # EAC_Total_SD <dbl>, AEQ_Total_M <dbl>, AEQ_Total_SD <dbl>,
## # AEQ_Positive_M <dbl>, AEQ_Positive_SD <dbl>, SOPA_Emotion3_M <dbl>,
## # SOPA_Emotion3_SD <dbl>
# =============================================================================
# STEP 5: TEST IF EAET INCREASES EMOTIONAL PROCESSING MORE THAN CBT
# =============================================================================
cat("\n\n=== DOES EAET INCREASE EMOTIONAL PROCESSING MORE THAN CBT? ===\n")
##
##
## === DOES EAET INCREASE EMOTIONAL PROCESSING MORE THAN CBT? ===
cat("(Testing the 'a' path of mediation)\n\n")
## (Testing the 'a' path of mediation)
# Test each emotional processing measure
cat("1. EAC EXPRESSION (higher change = more expression)\n")
## 1. EAC EXPRESSION (higher change = more expression)
cat(" Expected: EAET > CBT\n")
## Expected: EAET > CBT
t_eac_exp <- t.test(EAC_Expression_change ~ treatment, data = df_emotional)
d_eac_exp <- cohens_d(EAC_Expression_change ~ treatment, data = df_emotional)$Cohens_d
cat(sprintf(" t(%.1f) = %.2f, p = %.3f, d = %.2f %s\n\n",
t_eac_exp$parameter, t_eac_exp$statistic, t_eac_exp$p.value, d_eac_exp,
ifelse(t_eac_exp$p.value < 0.05, "***", ifelse(t_eac_exp$p.value < 0.10, "*", ""))))
## t(102.0) = -0.52, p = 0.607, d = -0.10
cat("2. EAC EXPERIENCE (higher change = more experiencing)\n")
## 2. EAC EXPERIENCE (higher change = more experiencing)
cat(" Expected: EAET > CBT\n")
## Expected: EAET > CBT
t_eac_exp2 <- t.test(EAC_Experience_change ~ treatment, data = df_emotional)
d_eac_exp2 <- cohens_d(EAC_Experience_change ~ treatment, data = df_emotional)$Cohens_d
cat(sprintf(" t(%.1f) = %.2f, p = %.3f, d = %.2f %s\n\n",
t_eac_exp2$parameter, t_eac_exp2$statistic, t_eac_exp2$p.value, d_eac_exp2,
ifelse(t_eac_exp2$p.value < 0.05, "***", ifelse(t_eac_exp2$p.value < 0.10, "*", ""))))
## t(102.0) = -0.48, p = 0.633, d = -0.09
cat("3. EAC TOTAL (higher change = better emotional approach coping)\n")
## 3. EAC TOTAL (higher change = better emotional approach coping)
cat(" Expected: EAET > CBT\n")
## Expected: EAET > CBT
t_eac_total <- t.test(EAC_Total_change ~ treatment, data = df_emotional)
d_eac_total <- cohens_d(EAC_Total_change ~ treatment, data = df_emotional)$Cohens_d
cat(sprintf(" t(%.1f) = %.2f, p = %.3f, d = %.2f %s\n\n",
t_eac_total$parameter, t_eac_total$statistic, t_eac_total$p.value, d_eac_total,
ifelse(t_eac_total$p.value < 0.05, "***", ifelse(t_eac_total$p.value < 0.10, "*", ""))))
## t(101.8) = -1.05, p = 0.295, d = -0.21
cat("4. AEQ TOTAL (NEGATIVE change = reduced ambivalence = better)\n")
## 4. AEQ TOTAL (NEGATIVE change = reduced ambivalence = better)
cat(" Expected: EAET shows MORE NEGATIVE change than CBT\n")
## Expected: EAET shows MORE NEGATIVE change than CBT
t_aeq <- t.test(AEQ_Total_change ~ treatment, data = df_emotional)
d_aeq <- cohens_d(AEQ_Total_change ~ treatment, data = df_emotional)$Cohens_d
cat(sprintf(" t(%.1f) = %.2f, p = %.3f, d = %.2f %s\n\n",
t_aeq$parameter, t_aeq$statistic, t_aeq$p.value, d_aeq,
ifelse(t_aeq$p.value < 0.05, "***", ifelse(t_aeq$p.value < 0.10, "*", ""))))
## t(102.0) = -2.13, p = 0.036, d = -0.42 ***
cat("5. SOPA EMOTION3 (higher change = more awareness of emotion-pain link)\n")
## 5. SOPA EMOTION3 (higher change = more awareness of emotion-pain link)
cat(" Expected: EAET > CBT\n")
## Expected: EAET > CBT
t_sopa <- t.test(SOPA_Emotion3_change ~ treatment, data = df_emotional)
d_sopa <- cohens_d(SOPA_Emotion3_change ~ treatment, data = df_emotional)$Cohens_d
cat(sprintf(" t(%.1f) = %.2f, p = %.3f, d = %.2f %s\n\n",
t_sopa$parameter, t_sopa$statistic, t_sopa$p.value, d_sopa,
ifelse(t_sopa$p.value < 0.05, "***", ifelse(t_sopa$p.value < 0.10, "*", ""))))
## t(101.1) = -1.21, p = 0.229, d = -0.24
cat("*** p < .05, * p < .10\n\n")
## *** p < .05, * p < .10
# =============================================================================
# STEP 6: MULTIVARIATE MEDIATION MODEL
# Treatment → Emotional Processing → Outcomes
# =============================================================================
cat("\n")
cat("=============================================================================\n")
## =============================================================================
cat("MULTIVARIATE EMOTIONAL PROCESSING MEDIATION MODEL\n")
## MULTIVARIATE EMOTIONAL PROCESSING MEDIATION MODEL
cat("Testing if EAET works THROUGH increasing emotional processing\n")
## Testing if EAET works THROUGH increasing emotional processing
cat("=============================================================================\n\n")
## =============================================================================
# Define the full mediation model
emotional_mediation_model <- '
# -------------------------------------------------------------------------
# MEDIATORS: Emotional Processing Measures (POST-treatment change)
# -------------------------------------------------------------------------
# EAC subscales (higher = better emotional processing)
EAC_Expression_change ~ a1*eaet + meanpain_b + age + female + race_white_vs_bipoc
EAC_Experience_change ~ a2*eaet + meanpain_b + age + female + race_white_vs_bipoc
# AEQ Total (LOWER = less ambivalence = better)
AEQ_Total_change ~ a3*eaet + meanpain_b + age + female + race_white_vs_bipoc
# SOPA Emotion (higher = stronger emotion-pain connection awareness)
SOPA_Emotion3_change ~ a4*eaet + meanpain_b + age + female + race_white_vs_bipoc
# -------------------------------------------------------------------------
# OUTCOMES (6-month follow-up)
# -------------------------------------------------------------------------
# Pain Severity at 6 months
meanpain_6m ~ c1_pain*eaet +
b1_pain*EAC_Expression_change +
b2_pain*EAC_Experience_change +
b3_pain*AEQ_Total_change +
b4_pain*SOPA_Emotion3_change +
meanpain_b + age + female + race_white_vs_bipoc
# Pain Interference at 6 months
painint_6m ~ c1_painint*eaet +
b1_painint*EAC_Expression_change +
b2_painint*EAC_Experience_change +
b3_painint*AEQ_Total_change +
b4_painint*SOPA_Emotion3_change +
painint_b + age + female + race_white_vs_bipoc
# Depression at 6 months
depress_6m ~ c1_depress*eaet +
b1_depress*EAC_Expression_change +
b2_depress*EAC_Experience_change +
b3_depress*AEQ_Total_change +
b4_depress*SOPA_Emotion3_change +
depress_b + age + female + race_white_vs_bipoc
# -------------------------------------------------------------------------
# INDIRECT EFFECTS - Pain Severity
# -------------------------------------------------------------------------
ind_EAC_Expression_pain := a1 * b1_pain
ind_EAC_Experience_pain := a2 * b2_pain
ind_AEQ_pain := a3 * b3_pain
ind_SOPA_pain := a4 * b4_pain
ind_total_pain := ind_EAC_Expression_pain + ind_EAC_Experience_pain + ind_AEQ_pain + ind_SOPA_pain
total_pain := c1_pain + ind_total_pain
direct_pain := c1_pain
# -------------------------------------------------------------------------
# INDIRECT EFFECTS - Pain Interference
# -------------------------------------------------------------------------
ind_EAC_Expression_painint := a1 * b1_painint
ind_EAC_Experience_painint := a2 * b2_painint
ind_AEQ_painint := a3 * b3_painint
ind_SOPA_painint := a4 * b4_painint
ind_total_painint := ind_EAC_Expression_painint + ind_EAC_Experience_painint + ind_AEQ_painint + ind_SOPA_painint
total_painint := c1_painint + ind_total_painint
direct_painint := c1_painint
# -------------------------------------------------------------------------
# INDIRECT EFFECTS - Depression
# -------------------------------------------------------------------------
ind_EAC_Expression_depress := a1 * b1_depress
ind_EAC_Experience_depress := a2 * b2_depress
ind_AEQ_depress := a3 * b3_depress
ind_SOPA_depress := a4 * b4_depress
ind_total_depress := ind_EAC_Expression_depress + ind_EAC_Experience_depress + ind_AEQ_depress + ind_SOPA_depress
total_depress := c1_depress + ind_total_depress
direct_depress := c1_depress
# -------------------------------------------------------------------------
# COVARIANCES
# -------------------------------------------------------------------------
# Mediators likely correlate
EAC_Expression_change ~~ EAC_Experience_change
EAC_Expression_change ~~ AEQ_Total_change
EAC_Expression_change ~~ SOPA_Emotion3_change
EAC_Experience_change ~~ AEQ_Total_change
EAC_Experience_change ~~ SOPA_Emotion3_change
AEQ_Total_change ~~ SOPA_Emotion3_change
# Outcomes likely correlate
meanpain_6m ~~ painint_6m
meanpain_6m ~~ depress_6m
painint_6m ~~ depress_6m
'
cat("Fitting multivariate emotional processing mediation model...\n")
## Fitting multivariate emotional processing mediation model...
cat("Using bootstrap with 5000 iterations for confidence intervals\n")
## Using bootstrap with 5000 iterations for confidence intervals
cat("(This will take several minutes - please be patient)\n\n")
## (This will take several minutes - please be patient)
# Fit the model
emotional_fit <- sem(emotional_mediation_model,
data = df_emotional,
se = "bootstrap",
bootstrap = 5000,
missing = "fiml")
# =============================================================================
# STEP 7: MODEL FIT AND RESULTS
# =============================================================================
cat("\n=== MODEL FIT INDICES ===\n")
##
## === MODEL FIT INDICES ===
fit_emotional <- fitMeasures(emotional_fit, c("chisq", "df", "pvalue",
"cfi", "tli", "rmsea",
"rmsea.ci.lower", "rmsea.ci.upper", "srmr"))
print(fit_emotional)
## chisq df pvalue cfi tli
## 18.323 14.000 0.192 0.982 0.909
## rmsea rmsea.ci.lower rmsea.ci.upper srmr
## 0.054 0.000 0.116 0.025
cat("\nInterpretation:\n")
##
## Interpretation:
cat(" CFI/TLI > 0.95 = excellent, > 0.90 = acceptable\n")
## CFI/TLI > 0.95 = excellent, > 0.90 = acceptable
cat(" RMSEA < 0.05 = excellent, < 0.08 = acceptable\n")
## RMSEA < 0.05 = excellent, < 0.08 = acceptable
cat(" SRMR < 0.08 = good fit\n\n")
## SRMR < 0.08 = good fit
cat("=== FULL MODEL SUMMARY ===\n")
## === FULL MODEL SUMMARY ===
summary(emotional_fit, fit.measures = TRUE, standardized = TRUE, rsquare = TRUE)
## lavaan 0.6.16 ended normally after 233 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 70
##
## Number of observations 104
## Number of missing patterns 1
##
## Model Test User Model:
##
## Test statistic 18.323
## Degrees of freedom 14
## P-value (Chi-square) 0.192
##
## Model Test Baseline Model:
##
## Test statistic 307.946
## Degrees of freedom 70
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 0.982
## Tucker-Lewis Index (TLI) 0.909
##
## Robust Comparative Fit Index (CFI) 0.982
## Robust Tucker-Lewis Index (TLI) 0.909
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -1956.660
## Loglikelihood unrestricted model (H1) -1947.499
##
## Akaike (AIC) 4053.321
## Bayesian (BIC) 4238.428
## Sample-size adjusted Bayesian (SABIC) 4017.298
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.054
## 90 Percent confidence interval - lower 0.000
## 90 Percent confidence interval - upper 0.116
## P-value H_0: RMSEA <= 0.050 0.412
## P-value H_0: RMSEA >= 0.080 0.290
##
## Robust RMSEA 0.054
## 90 Percent confidence interval - lower 0.000
## 90 Percent confidence interval - upper 0.116
## P-value H_0: Robust RMSEA <= 0.050 0.412
## P-value H_0: Robust RMSEA >= 0.080 0.290
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.025
##
## Parameter Estimates:
##
## Standard errors Bootstrap
## Number of requested bootstrap draws 5000
## Number of successful bootstrap draws 5000
##
## Regressions:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## EAC_Expression_change ~
## eaet (a1) 0.257 0.535 0.479 0.632 0.257 0.049
## mnpn_ 0.013 0.134 0.097 0.923 0.013 0.009
## age -0.063 0.047 -1.357 0.175 -0.063 -0.137
## femal 1.219 0.869 1.403 0.161 1.219 0.131
## rc___ -0.487 0.629 -0.775 0.438 -0.487 -0.085
## EAC_Experience_change ~
## eaet (a2) 0.300 0.556 0.540 0.589 0.300 0.054
## mnpn_ 0.070 0.146 0.481 0.630 0.070 0.045
## age -0.033 0.048 -0.683 0.495 -0.033 -0.066
## femal 0.475 1.156 0.411 0.681 0.475 0.048
## rc___ -0.022 0.764 -0.029 0.977 -0.022 -0.004
## AEQ_Total_change ~
## eaet (a3) 4.697 2.078 2.260 0.024 4.697 0.228
## mnpn_ 0.479 0.592 0.809 0.419 0.479 0.082
## age -0.032 0.180 -0.178 0.859 -0.032 -0.018
## femal -4.855 3.350 -1.449 0.147 -4.855 -0.132
## rc___ -0.463 2.442 -0.190 0.850 -0.463 -0.021
## SOPA_Emotion3_change ~
## eaet (a4) 0.626 0.650 0.962 0.336 0.626 0.097
## mnpn_ -0.324 0.189 -1.711 0.087 -0.324 -0.179
## age -0.008 0.061 -0.132 0.895 -0.008 -0.014
## femal 0.342 1.269 0.269 0.788 0.342 0.030
## rc___ 0.103 0.620 0.166 0.868 0.103 0.015
## meanpain_6m ~
## eaet (c1_pan) -1.148 0.371 -3.097 0.002 -1.148 -0.284
## EAC_E (b1_pan) 0.121 0.066 1.833 0.067 0.121 0.157
## EAC_E (b2_pan) -0.125 0.067 -1.849 0.064 -0.125 -0.172
## AEQ_T (b3_pan) 0.014 0.020 0.716 0.474 0.014 0.073
## SOPA_ (b4_pan) 0.010 0.055 0.182 0.856 0.010 0.016
## mnpn_ 0.535 0.107 4.994 0.000 0.535 0.469
## age -0.020 0.031 -0.663 0.507 -0.020 -0.057
## femal -0.485 0.806 -0.601 0.548 -0.485 -0.067
## rc___ 0.002 0.444 0.005 0.996 0.002 0.000
## painint_6m ~
## eaet (c1_pnn) -0.699 1.404 -0.498 0.619 -0.699 -0.045
## EAC_E (b1_pnn) 0.419 0.263 1.597 0.110 0.419 0.139
## EAC_E (b2_pnn) -0.230 0.267 -0.864 0.388 -0.230 -0.082
## AEQ_T (b3_pnn) 0.069 0.069 0.995 0.320 0.069 0.090
## SOPA_ (b4_pnn) 0.048 0.198 0.241 0.809 0.048 0.020
## pnnt_ 0.630 0.096 6.531 0.000 0.630 0.619
## age 0.046 0.120 0.381 0.703 0.046 0.033
## femal -0.989 2.357 -0.420 0.675 -0.989 -0.035
## rc___ 1.440 1.697 0.848 0.396 1.440 0.084
## depress_6m ~
## eaet (c1_d) -2.255 1.100 -2.050 0.040 -2.255 -0.156
## EAC_E (b1_d) 0.283 0.230 1.233 0.217 0.283 0.103
## EAC_E (b2_d) -0.387 0.209 -1.850 0.064 -0.387 -0.150
## AEQ_T (b3_d) 0.016 0.049 0.329 0.742 0.016 0.023
## SOPA_ (b4_d) 0.216 0.172 1.254 0.210 0.216 0.096
## dprs_ 0.634 0.057 11.042 0.000 0.634 0.689
## age -0.072 0.115 -0.624 0.533 -0.072 -0.056
## femal -1.089 2.236 -0.487 0.626 -1.089 -0.042
## rc___ 0.141 1.188 0.118 0.906 0.141 0.009
##
## Covariances:
## Estimate Std.Err z-value P(>|z|) Std.lv
## .EAC_Expression_change ~~
## .EAC_Exprnc_chn 2.664 0.752 3.542 0.000 2.664
## .AEQ_Total_chng -4.139 2.778 -1.490 0.136 -4.139
## .SOPA_Emtn3_chn -1.439 0.781 -1.843 0.065 -1.439
## .EAC_Experience_change ~~
## .AEQ_Total_chng -0.697 2.540 -0.275 0.784 -0.697
## .SOPA_Emtn3_chn -0.283 0.970 -0.291 0.771 -0.283
## .AEQ_Total_change ~~
## .SOPA_Emtn3_chn 1.378 2.616 0.527 0.598 1.378
## .meanpain_6m ~~
## .painint_6m 4.958 1.455 3.408 0.001 4.958
## .depress_6m 2.938 0.862 3.407 0.001 2.938
## .painint_6m ~~
## .depress_6m 13.499 3.251 4.152 0.000 13.499
## Std.all
##
## 0.376
## -0.162
## -0.179
##
## -0.025
## -0.032
##
## 0.044
##
## 0.514
## 0.367
##
## 0.452
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .EAC_Exprssn_ch 4.363 3.578 1.220 0.223 4.363 1.671
## .EAC_Exprnc_chn 1.820 3.930 0.463 0.643 1.820 0.652
## .AEQ_Total_chng -3.327 13.989 -0.238 0.812 -3.327 -0.323
## .SOPA_Emtn3_chn 2.963 4.700 0.630 0.528 2.963 0.922
## .meanpain_6m 4.179 2.372 1.762 0.078 4.179 2.067
## .painint_6m 3.813 9.721 0.392 0.695 3.813 0.486
## .depress_6m 11.372 8.508 1.337 0.181 11.372 1.578
##
## Variances:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .EAC_Exprssn_ch 6.521 0.747 8.725 0.000 6.521 0.957
## .EAC_Exprnc_chn 7.683 1.025 7.496 0.000 7.683 0.987
## .AEQ_Total_chng 99.537 17.129 5.811 0.000 99.537 0.937
## .SOPA_Emtn3_chn 9.868 1.858 5.311 0.000 9.868 0.956
## .meanpain_6m 2.585 0.417 6.201 0.000 2.585 0.633
## .painint_6m 36.028 5.412 6.657 0.000 36.028 0.585
## .depress_6m 24.716 3.539 6.984 0.000 24.716 0.476
##
## R-Square:
## Estimate
## EAC_Exprssn_ch 0.043
## EAC_Exprnc_chn 0.013
## AEQ_Total_chng 0.063
## SOPA_Emtn3_chn 0.044
## meanpain_6m 0.367
## painint_6m 0.415
## depress_6m 0.524
##
## Defined Parameters:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## ind_EAC_Exprs_ 0.031 0.078 0.398 0.691 0.031 0.008
## ind_EAC_Exprn_ -0.037 0.078 -0.480 0.631 -0.037 -0.009
## ind_AEQ_pain 0.068 0.102 0.661 0.509 0.068 0.017
## ind_SOPA_pain 0.006 0.050 0.126 0.900 0.006 0.002
## ind_total_pain 0.068 0.141 0.478 0.633 0.068 0.017
## total_pain -1.080 0.337 -3.207 0.001 -1.080 -0.267
## direct_pain -1.148 0.371 -3.097 0.002 -1.148 -0.284
## ind_EAC_Exprs_ 0.108 0.297 0.362 0.717 0.108 0.007
## ind_EAC_Exprn_ -0.069 0.210 -0.329 0.742 -0.069 -0.004
## ind_AEQ_painnt 0.322 0.358 0.902 0.367 0.322 0.021
## ind_SOPA_pannt 0.030 0.186 0.161 0.872 0.030 0.002
## ind_total_pnnt 0.391 0.491 0.796 0.426 0.391 0.025
## total_painint -0.308 1.250 -0.247 0.805 -0.308 -0.020
## direct_painint -0.699 1.405 -0.498 0.619 -0.699 -0.045
## ind_EAC_Exprs_ 0.073 0.208 0.350 0.726 0.073 0.005
## ind_EAC_Exprn_ -0.116 0.251 -0.463 0.643 -0.116 -0.008
## ind_AEQ_deprss 0.075 0.249 0.303 0.762 0.075 0.005
## ind_SOPA_dprss 0.135 0.203 0.666 0.505 0.135 0.009
## ind_totl_dprss 0.167 0.390 0.429 0.668 0.167 0.012
## total_depress -2.088 1.037 -2.013 0.044 -2.088 -0.145
## direct_depress -2.255 1.100 -2.050 0.040 -2.255 -0.156
# Extract all parameters
emotional_params <- parameterEstimates(emotional_fit, boot.ci.type = "perc")
# =============================================================================
# STEP 8: TREATMENT → EMOTIONAL PROCESSING (A PATHS)
# =============================================================================
cat("\n\n")
cat("=============================================================================\n")
## =============================================================================
cat("PATH A: DOES EAET INCREASE EMOTIONAL PROCESSING?\n")
## PATH A: DOES EAET INCREASE EMOTIONAL PROCESSING?
cat("=============================================================================\n\n")
## =============================================================================
a_paths <- emotional_params[emotional_params$label %in% c("a1", "a2", "a3", "a4"),
c("label", "est", "se", "ci.lower", "ci.upper", "pvalue")]
a_paths_table <- data.frame(
Mediator = c("EAC Expression (↑ better)",
"EAC Experience (↑ better)",
"AEQ Total (↓ better)",
"SOPA Emotion3 (↑ better)"),
Estimate = sprintf("%.3f", a_paths$est),
SE = sprintf("%.3f", a_paths$se),
CI_95 = sprintf("[%.3f, %.3f]", a_paths$ci.lower, a_paths$ci.upper),
p_value = sprintf("%.3f", a_paths$pvalue),
Sig = ifelse(a_paths$pvalue < 0.05, "***",
ifelse(a_paths$pvalue < 0.10, "*", ""))
)
print(a_paths_table)
## Mediator Estimate SE CI_95 p_value Sig
## 1 EAC Expression (↑ better) 0.257 0.535 [-0.759, 1.350] 0.632
## 2 EAC Experience (↑ better) 0.300 0.556 [-0.783, 1.379] 0.589
## 3 AEQ Total (↓ better) 4.697 2.078 [0.596, 8.738] 0.024 ***
## 4 SOPA Emotion3 (↑ better) 0.626 0.650 [-0.634, 1.919] 0.336
cat("\n*** p < .05, * p < .10\n")
##
## *** p < .05, * p < .10
cat("\nInterpretation: Positive estimates for EAC and SOPA mean EAET increases these.\n")
##
## Interpretation: Positive estimates for EAC and SOPA mean EAET increases these.
cat(" Negative estimates for AEQ mean EAET reduces ambivalence (good!).\n\n")
## Negative estimates for AEQ mean EAET reduces ambivalence (good!).
# =============================================================================
# STEP 9: EMOTIONAL PROCESSING → OUTCOMES (B PATHS)
# =============================================================================
cat("\n")
cat("=============================================================================\n")
## =============================================================================
cat("PATH B: DOES EMOTIONAL PROCESSING PREDICT OUTCOMES?\n")
## PATH B: DOES EMOTIONAL PROCESSING PREDICT OUTCOMES?
cat("=============================================================================\n\n")
## =============================================================================
# Pain Severity
cat("--- PAIN SEVERITY (6 months) ---\n")
## --- PAIN SEVERITY (6 months) ---
b_pain <- emotional_params[emotional_params$label %in% c("b1_pain", "b2_pain", "b3_pain", "b4_pain"),
c("label", "est", "se", "ci.lower", "ci.upper", "pvalue")]
b_pain_table <- data.frame(
Mediator = c("EAC Expression", "EAC Experience", "AEQ Total", "SOPA Emotion3"),
Estimate = sprintf("%.3f", b_pain$est),
SE = sprintf("%.3f", b_pain$se),
CI_95 = sprintf("[%.3f, %.3f]", b_pain$ci.lower, b_pain$ci.upper),
p_value = sprintf("%.3f", b_pain$pvalue),
Sig = ifelse(b_pain$pvalue < 0.05, "***",
ifelse(b_pain$pvalue < 0.10, "*", ""))
)
print(b_pain_table)
## Mediator Estimate SE CI_95 p_value Sig
## 1 EAC Expression 0.121 0.066 [-0.000, 0.259] 0.067 *
## 2 EAC Experience -0.125 0.067 [-0.253, 0.012] 0.064 *
## 3 AEQ Total 0.014 0.020 [-0.026, 0.053] 0.474
## 4 SOPA Emotion3 0.010 0.055 [-0.102, 0.116] 0.856
cat("\nInterpretation: Negative estimates = increases in emotional processing → less pain\n\n")
##
## Interpretation: Negative estimates = increases in emotional processing → less pain
# Pain Interference
cat("--- PAIN INTERFERENCE (6 months) ---\n")
## --- PAIN INTERFERENCE (6 months) ---
b_painint <- emotional_params[emotional_params$label %in% c("b1_painint", "b2_painint", "b3_painint", "b4_painint"),
c("label", "est", "se", "ci.lower", "ci.upper", "pvalue")]
b_painint_table <- data.frame(
Mediator = c("EAC Expression", "EAC Experience", "AEQ Total", "SOPA Emotion3"),
Estimate = sprintf("%.3f", b_painint$est),
SE = sprintf("%.3f", b_painint$se),
CI_95 = sprintf("[%.3f, %.3f]", b_painint$ci.lower, b_painint$ci.upper),
p_value = sprintf("%.3f", b_painint$pvalue),
Sig = ifelse(b_painint$pvalue < 0.05, "***",
ifelse(b_painint$pvalue < 0.10, "*", ""))
)
print(b_painint_table)
## Mediator Estimate SE CI_95 p_value Sig
## 1 EAC Expression 0.419 0.263 [-0.064, 0.966] 0.110
## 2 EAC Experience -0.230 0.267 [-0.789, 0.255] 0.388
## 3 AEQ Total 0.069 0.069 [-0.059, 0.215] 0.320
## 4 SOPA Emotion3 0.048 0.198 [-0.360, 0.443] 0.809
# Depression
cat("\n--- DEPRESSION (6 months) ---\n")
##
## --- DEPRESSION (6 months) ---
b_depress <- emotional_params[emotional_params$label %in% c("b1_depress", "b2_depress", "b3_depress", "b4_depress"),
c("label", "est", "se", "ci.lower", "ci.upper", "pvalue")]
b_depress_table <- data.frame(
Mediator = c("EAC Expression", "EAC Experience", "AEQ Total", "SOPA Emotion3"),
Estimate = sprintf("%.3f", b_depress$est),
SE = sprintf("%.3f", b_depress$se),
CI_95 = sprintf("[%.3f, %.3f]", b_depress$ci.lower, b_depress$ci.upper),
p_value = sprintf("%.3f", b_depress$pvalue),
Sig = ifelse(b_depress$pvalue < 0.05, "***",
ifelse(b_depress$pvalue < 0.10, "*", ""))
)
print(b_depress_table)
## Mediator Estimate SE CI_95 p_value Sig
## 1 EAC Expression 0.283 0.230 [-0.181, 0.704] 0.217
## 2 EAC Experience -0.387 0.209 [-0.778, 0.046] 0.064 *
## 3 AEQ Total 0.016 0.049 [-0.084, 0.108] 0.742
## 4 SOPA Emotion3 0.216 0.172 [-0.114, 0.566] 0.210
# =============================================================================
# STEP 10: INDIRECT EFFECTS (THE MEDIATION TEST)
# =============================================================================
cat("\n\n")
cat("=============================================================================\n")
## =============================================================================
cat("INDIRECT EFFECTS: Treatment → Emotional Processing → Outcomes\n")
## INDIRECT EFFECTS: Treatment → Emotional Processing → Outcomes
cat("(This is the KEY mediation test)\n")
## (This is the KEY mediation test)
cat("=============================================================================\n\n")
## =============================================================================
# Extract all indirect effects
indirect_effects <- emotional_params[grepl("^ind_", emotional_params$label) &
!grepl("total|direct", emotional_params$label),
c("label", "est", "se", "ci.lower", "ci.upper", "pvalue")]
# Pain Severity
cat("╔════════════════════════════════════════════════════════════════╗\n")
## ╔════════════════════════════════════════════════════════════════╗
cat("║ PAIN SEVERITY (6 months) ║\n")
## ║ PAIN SEVERITY (6 months) ║
cat("╚════════════════════════════════════════════════════════════════╝\n\n")
## ╚════════════════════════════════════════════════════════════════╝
ind_pain <- indirect_effects[grepl("pain$", indirect_effects$label) & !grepl("painint", indirect_effects$label), ]
ind_pain_table <- data.frame(
Pathway = c("Via EAC Expression", "Via EAC Experience", "Via AEQ Total", "Via SOPA Emotion3"),
Indirect_Effect = sprintf("%.3f", ind_pain$est),
CI_95 = sprintf("[%.3f, %.3f]", ind_pain$ci.lower, ind_pain$ci.upper),
p_value = sprintf("%.3f", ind_pain$pvalue),
Sig = ifelse(ind_pain$pvalue < 0.05, "***",
ifelse(ind_pain$pvalue < 0.10, "*", ""))
)
print(ind_pain_table)
## Pathway Indirect_Effect CI_95 p_value Sig
## 1 Via EAC Expression 0.031 [-0.110, 0.222] 0.691
## 2 Via EAC Experience -0.037 [-0.202, 0.121] 0.631
## 3 Via AEQ Total 0.068 [-0.147, 0.272] 0.509
## 4 Via SOPA Emotion3 0.006 [-0.082, 0.132] 0.900
# Total indirect and direct effects
total_ind_pain <- emotional_params[emotional_params$label == "ind_total_pain",
c("est", "ci.lower", "ci.upper", "pvalue")]
direct_pain <- emotional_params[emotional_params$label == "direct_pain",
c("est", "ci.lower", "ci.upper", "pvalue")]
total_pain <- emotional_params[emotional_params$label == "total_pain",
c("est", "ci.lower", "ci.upper", "pvalue")]
cat(sprintf("\n%-25s: %6.3f [%6.3f, %6.3f], p = %.3f %s\n",
"Total Indirect Effect",
total_ind_pain$est, total_ind_pain$ci.lower,
total_ind_pain$ci.upper, total_ind_pain$pvalue,
ifelse(total_ind_pain$pvalue < 0.05, "***",
ifelse(total_ind_pain$pvalue < 0.10, "*", ""))))
##
## Total Indirect Effect : 0.068 [-0.190, 0.370], p = 0.633
cat(sprintf("%-25s: %6.3f [%6.3f, %6.3f], p = %.3f %s\n",
"Direct Effect (c')",
direct_pain$est, direct_pain$ci.lower,
direct_pain$ci.upper, direct_pain$pvalue,
ifelse(direct_pain$pvalue < 0.05, "***", "")))
## Direct Effect (c') : -1.148 [-1.898, -0.458], p = 0.002 ***
cat(sprintf("%-25s: %6.3f [%6.3f, %6.3f], p = %.3f %s\n",
"Total Effect (c)",
total_pain$est, total_pain$ci.lower,
total_pain$ci.upper, total_pain$pvalue,
ifelse(total_pain$pvalue < 0.05, "***", "")))
## Total Effect (c) : -1.080 [-1.772, -0.445], p = 0.001 ***
# Proportion mediated (only if both indirect and direct are significant)
if(total_ind_pain$pvalue < 0.10 && direct_pain$pvalue < 0.10) {
prop_med <- abs(total_ind_pain$est) / abs(total_pain$est)
cat(sprintf("\n%-25s: %.1f%% of treatment effect explained by emotional processing\n",
"Proportion Mediated", prop_med * 100))
}
# Pain Interference
cat("\n\n╔════════════════════════════════════════════════════════════════╗\n")
##
##
## ╔════════════════════════════════════════════════════════════════╗
cat("║ PAIN INTERFERENCE (6 months) ║\n")
## ║ PAIN INTERFERENCE (6 months) ║
cat("╚════════════════════════════════════════════════════════════════╝\n\n")
## ╚════════════════════════════════════════════════════════════════╝
ind_painint <- indirect_effects[grepl("painint", indirect_effects$label), ]
ind_painint_table <- data.frame(
Pathway = c("Via EAC Expression", "Via EAC Experience", "Via AEQ Total", "Via SOPA Emotion3"),
Indirect_Effect = sprintf("%.3f", ind_painint$est),
CI_95 = sprintf("[%.3f, %.3f]", ind_painint$ci.lower, ind_painint$ci.upper),
p_value = sprintf("%.3f", ind_painint$pvalue),
Sig = ifelse(ind_painint$pvalue < 0.05, "***",
ifelse(ind_painint$pvalue < 0.10, "*", ""))
)
print(ind_painint_table)
## Pathway Indirect_Effect CI_95 p_value Sig
## 1 Via EAC Expression 0.108 [-0.366, 0.863] 0.717
## 2 Via EAC Experience -0.069 [-0.615, 0.273] 0.742
## 3 Via AEQ Total 0.322 [-0.283, 1.134] 0.367
## 4 Via SOPA Emotion3 0.030 [-0.300, 0.500] 0.872
total_ind_painint <- emotional_params[emotional_params$label == "ind_total_painint",
c("est", "ci.lower", "ci.upper", "pvalue")]
direct_painint <- emotional_params[emotional_params$label == "direct_painint",
c("est", "ci.lower", "ci.upper", "pvalue")]
cat(sprintf("\n%-25s: %6.3f [%6.3f, %6.3f], p = %.3f %s\n",
"Total Indirect Effect",
total_ind_painint$est, total_ind_painint$ci.lower,
total_ind_painint$ci.upper, total_ind_painint$pvalue,
ifelse(total_ind_painint$pvalue < 0.05, "***",
ifelse(total_ind_painint$pvalue < 0.10, "*", ""))))
##
## Total Indirect Effect : 0.391 [-0.428, 1.519], p = 0.426
cat(sprintf("%-25s: %6.3f [%6.3f, %6.3f], p = %.3f %s\n",
"Direct Effect",
direct_painint$est, direct_painint$ci.lower,
direct_painint$ci.upper, direct_painint$pvalue,
ifelse(direct_painint$pvalue < 0.05, "***", "")))
## Direct Effect : -0.699 [-3.650, 1.928], p = 0.619
# Depression
cat("\n\n╔════════════════════════════════════════════════════════════════╗\n")
##
##
## ╔════════════════════════════════════════════════════════════════╗
cat("║ DEPRESSION (6 months) ║\n")
## ║ DEPRESSION (6 months) ║
cat("╚════════════════════════════════════════════════════════════════╝\n\n")
## ╚════════════════════════════════════════════════════════════════╝
ind_depress <- indirect_effects[grepl("depress", indirect_effects$label), ]
ind_depress_table <- data.frame(
Pathway = c("Via EAC Expression", "Via EAC Experience", "Via AEQ Total", "Via SOPA Emotion3"),
Indirect_Effect = sprintf("%.3f", ind_depress$est),
CI_95 = sprintf("[%.3f, %.3f]", ind_depress$ci.lower, ind_depress$ci.upper),
p_value = sprintf("%.3f", ind_depress$pvalue),
Sig = ifelse(ind_depress$pvalue < 0.05, "***",
ifelse(ind_depress$pvalue < 0.10, "*", ""))
)
print(ind_depress_table)
## Pathway Indirect_Effect CI_95 p_value Sig
## 1 Via EAC Expression 0.073 [-0.282, 0.622] 0.726
## 2 Via EAC Experience -0.116 [-0.745, 0.305] 0.643
## 3 Via AEQ Total 0.075 [-0.483, 0.535] 0.762
## 4 Via SOPA Emotion3 0.135 [-0.187, 0.632] 0.505
total_ind_depress <- emotional_params[emotional_params$label == "ind_total_depress",
c("est", "ci.lower", "ci.upper", "pvalue")]
direct_depress <- emotional_params[emotional_params$label == "direct_depress",
c("est", "ci.lower", "ci.upper", "pvalue")]
total_depress <- emotional_params[emotional_params$label == "total_depress",
c("est", "ci.lower", "ci.upper", "pvalue")]
cat(sprintf("\n%-25s: %6.3f [%6.3f, %6.3f], p = %.3f %s\n",
"Total Indirect Effect",
total_ind_depress$est, total_ind_depress$ci.lower,
total_ind_depress$ci.upper, total_ind_depress$pvalue,
ifelse(total_ind_depress$pvalue < 0.05, "***",
ifelse(total_ind_depress$pvalue < 0.10, "*", ""))))
##
## Total Indirect Effect : 0.167 [-0.655, 0.930], p = 0.668
cat(sprintf("%-25s: %6.3f [%6.3f, %6.3f], p = %.3f %s\n",
"Direct Effect",
direct_depress$est, direct_depress$ci.lower,
direct_depress$ci.upper, direct_depress$pvalue,
ifelse(direct_depress$pvalue < 0.05, "***", "")))
## Direct Effect : -2.255 [-4.328, 0.005], p = 0.040 ***
cat(sprintf("%-25s: %6.3f [%6.3f, %6.3f], p = %.3f %s\n",
"Total Effect",
total_depress$est, total_depress$ci.lower,
total_depress$ci.upper, total_depress$pvalue,
ifelse(total_depress$pvalue < 0.05, "***", "")))
## Total Effect : -2.088 [-4.072, -0.035], p = 0.044 ***
if(total_ind_depress$pvalue < 0.10 && direct_depress$pvalue < 0.10) {
prop_med_dep <- abs(total_ind_depress$est) / abs(total_depress$est)
cat(sprintf("\n%-25s: %.1f%% of treatment effect explained by emotional processing\n",
"Proportion Mediated", prop_med_dep * 100))
}
# =============================================================================
# STEP 11: VARIANCE EXPLAINED (R²)
# =============================================================================
cat("\n\n")
cat("=============================================================================\n")
## =============================================================================
cat("VARIANCE EXPLAINED (R²)\n")
## VARIANCE EXPLAINED (R²)
cat("=============================================================================\n\n")
## =============================================================================
rsq <- inspect(emotional_fit, "r2")
cat("Emotional Processing Mediators:\n")
## Emotional Processing Mediators:
cat(sprintf(" EAC Expression: R² = %.3f (%.1f%% variance explained)\n",
rsq["EAC_Expression_change"], rsq["EAC_Expression_change"]*100))
## EAC Expression: R² = 0.043 (4.3% variance explained)
cat(sprintf(" EAC Experience: R² = %.3f (%.1f%% variance explained)\n",
rsq["EAC_Experience_change"], rsq["EAC_Experience_change"]*100))
## EAC Experience: R² = 0.013 (1.3% variance explained)
cat(sprintf(" AEQ Total: R² = %.3f (%.1f%% variance explained)\n",
rsq["AEQ_Total_change"], rsq["AEQ_Total_change"]*100))
## AEQ Total: R² = 0.063 (6.3% variance explained)
cat(sprintf(" SOPA Emotion3: R² = %.3f (%.1f%% variance explained)\n",
rsq["SOPA_Emotion3_change"], rsq["SOPA_Emotion3_change"]*100))
## SOPA Emotion3: R² = 0.044 (4.4% variance explained)
cat("\nOutcome Models:\n")
##
## Outcome Models:
cat(sprintf(" Pain Severity: R² = %.3f (%.1f%% variance explained)\n",
rsq["meanpain_6m"], rsq["meanpain_6m"]*100))
## Pain Severity: R² = 0.367 (36.7% variance explained)
cat(sprintf(" Pain Interference: R² = %.3f (%.1f%% variance explained)\n",
rsq["painint_6m"], rsq["painint_6m"]*100))
## Pain Interference: R² = 0.415 (41.5% variance explained)
cat(sprintf(" Depression: R² = %.3f (%.1f%% variance explained)\n",
rsq["depress_6m"], rsq["depress_6m"]*100))
## Depression: R² = 0.524 (52.4% variance explained)
# =============================================================================
# SUMMARY AND INTERPRETATION
# =============================================================================
cat("\n\n")
cat("╔════════════════════════════════════════════════════════════════════════╗\n")
## ╔════════════════════════════════════════════════════════════════════════╗
cat("║ ANALYSIS COMPLETE - SUMMARY ║\n")
## ║ ANALYSIS COMPLETE - SUMMARY ║
cat("╚════════════════════════════════════════════════════════════════════════╝\n\n")
## ╚════════════════════════════════════════════════════════════════════════╝
cat("KEY FINDINGS:\n\n")
## KEY FINDINGS:
cat("1. TREATMENT → EMOTIONAL PROCESSING (a paths):\n")
## 1. TREATMENT → EMOTIONAL PROCESSING (a paths):
for(i in 1:nrow(a_paths_table)) {
cat(sprintf(" %s: %s %s\n",
a_paths_table$Mediator[i],
a_paths_table$Estimate[i],
ifelse(a_paths_table$Sig[i] != "", paste0("(", a_paths_table$Sig[i], ")"), "")))
}
## EAC Expression (↑ better): 0.257
## EAC Experience (↑ better): 0.300
## AEQ Total (↓ better): 4.697 (***)
## SOPA Emotion3 (↑ better): 0.626
cat("\n2. EMOTIONAL PROCESSING → PAIN (b paths):\n")
##
## 2. EMOTIONAL PROCESSING → PAIN (b paths):
for(i in 1:nrow(b_pain_table)) {
cat(sprintf(" %s: %s %s\n",
b_pain_table$Mediator[i],
b_pain_table$Estimate[i],
ifelse(b_pain_table$Sig[i] != "", paste0("(", b_pain_table$Sig[i], ")"), "")))
}
## EAC Expression: 0.121 (*)
## EAC Experience: -0.125 (*)
## AEQ Total: 0.014
## SOPA Emotion3: 0.010
cat("\n3. MEDIATION (indirect effects for PAIN):\n")
##
## 3. MEDIATION (indirect effects for PAIN):
for(i in 1:nrow(ind_pain_table)) {
cat(sprintf(" %s: %s %s\n",
ind_pain_table$Pathway[i],
ind_pain_table$Indirect_Effect[i],
ifelse(ind_pain_table$Sig[i] != "", paste0("(", ind_pain_table$Sig[i], ")"), "")))
}
## Via EAC Expression: 0.031
## Via EAC Experience: -0.037
## Via AEQ Total: 0.068
## Via SOPA Emotion3: 0.006
cat("\n\nINTERPRETATION GUIDE:\n")
##
##
## INTERPRETATION GUIDE:
cat(" *** = Significant mediation (p < .05)\n")
## *** = Significant mediation (p < .05)
cat(" * = Marginal significance (p < .10)\n")
## * = Marginal significance (p < .10)
cat("\n For mediation to be significant:\n")
##
## For mediation to be significant:
cat(" - Treatment must affect the mediator (a path)\n")
## - Treatment must affect the mediator (a path)
cat(" - Mediator must predict the outcome (b path)\n")
## - Mediator must predict the outcome (b path)
cat(" - Their product (indirect effect) must be significant\n\n")
## - Their product (indirect effect) must be significant
cat("✓ Emotional processing mediation analysis complete!\n")
## ✓ Emotional processing mediation analysis complete!
cat("=============================================================================\n")
## =============================================================================