setwd("C:/Work Files/Collaboration/Lisa Rapport")
Lisa Power Analysis
library(lavaan)
This is lavaan 0.6-19
lavaan is FREE software! Please report any bugs.
library(semTools)
###############################################################################
This is semTools 0.5-6
All users of R (or SEM) are invited to submit functions or ideas for functions.
###############################################################################
# Define the Latent Growth Model
<- '
model # Latent growth factors
i =~ 1*cog_T1 + 1*cog_T2 + 1*cog_T3 # Intercept factor (i)
s =~ 0*cog_T1 + 1*cog_T2 + 2*cog_T3 # Slope factor (s)
# Regressions for ACEs
i ~ aces # Path from ACEs to intercept
s ~ aces # Path from ACEs to slope
# Mediation through Emotion Coping
s ~ em_coping # Path from emotion coping to slope
em_coping ~ aces # Path from ACEs to emotion coping
'
# Define population parameters for the model
<- '
pop_params # Latent growth structure
i =~ 1*cog_T1 + 1*cog_T2 + 1*cog_T3
s =~ 0*cog_T1 + 1*cog_T2 + 2*cog_T3
# Paths from ACEs and Emotion Coping
i ~ 0.3*aces
s ~ 0.35*aces + 0.25*em_coping
em_coping ~ 0.3*aces
# Define means and variances for observed variables
cog_T1 ~ 0*1
cog_T2 ~ 0*1
cog_T3 ~ 0*1
aces ~ 0*1 # Mean of ACEs set to 0
em_coping ~ 0*1 # Mean of Emotion Coping set to 0
# Residual variances for observed indicators
cog_T1 ~~ 0.5*cog_T1
cog_T2 ~~ 0.5*cog_T2
cog_T3 ~~ 0.5*cog_T3
aces ~~ 1.0*aces # Variance of ACEs
em_coping ~~ 1.0*em_coping # Variance of Emotion Coping
'
# Set parameters for the Monte Carlo simulation
set.seed(123) # For reproducibility
<- 10 # Reduce to a smaller number initially for debugging
n_reps <- 150
n <- data.frame(intercept_pval = numeric(n_reps),
results slope_pval = numeric(n_reps),
indirect_pval = numeric(n_reps))
# Run the simulation
for (i in 1:n_reps) {
# Generate simulated data
<- simulateData(pop_params, sample.nobs = n)
sim_data
# Fit the model to the simulated data
<- sem(model, data = sim_data)
fit
# Extract parameter estimates and p-values
<- parameterEstimates(fit)
pvals
# Debugging output to confirm parameter names
print(pvals[, c("lhs", "op", "rhs", "pvalue")]) # Print p-values for all paths
# Extract p-values using safer indexing to avoid replacement errors
<- ifelse(any(pvals$lhs == "i" & pvals$op == "~" & pvals$rhs == "aces"),
intercept_pval $pvalue[pvals$lhs == "i" & pvals$op == "~" & pvals$rhs == "aces"], NA)
pvals<- ifelse(any(pvals$lhs == "s" & pvals$op == "~" & pvals$rhs == "aces"),
slope_pval $pvalue[pvals$lhs == "s" & pvals$op == "~" & pvals$rhs == "aces"], NA)
pvals<- ifelse(any(pvals$lhs == "em_coping" & pvals$op == "~" & pvals$rhs == "aces"),
indirect_pval $pvalue[pvals$lhs == "em_coping" & pvals$op == "~" & pvals$rhs == "aces"], NA)
pvals
$intercept_pval[i] <- intercept_pval
results$slope_pval[i] <- slope_pval
results$indirect_pval[i] <- indirect_pval
results }
lhs op rhs pvalue
1 i =~ cog_T1 NA
2 i =~ cog_T2 NA
3 i =~ cog_T3 NA
4 s =~ cog_T1 NA
5 s =~ cog_T2 NA
6 s =~ cog_T3 NA
7 i ~ aces 0.000
8 s ~ aces 0.003
9 s ~ em_coping 0.725
10 em_coping ~ aces 0.003
11 cog_T1 ~~ cog_T1 0.084
12 cog_T2 ~~ cog_T2 0.000
13 cog_T3 ~~ cog_T3 0.038
14 em_coping ~~ em_coping 0.000
15 i ~~ i 0.000
16 s ~~ s 0.000
17 i ~~ s 0.848
18 aces ~~ aces NA
lhs op rhs pvalue
1 i =~ cog_T1 NA
2 i =~ cog_T2 NA
3 i =~ cog_T3 NA
4 s =~ cog_T1 NA
5 s =~ cog_T2 NA
6 s =~ cog_T3 NA
7 i ~ aces 0.000
8 s ~ aces 0.002
9 s ~ em_coping 0.000
10 em_coping ~ aces 0.000
11 cog_T1 ~~ cog_T1 0.001
12 cog_T2 ~~ cog_T2 0.000
13 cog_T3 ~~ cog_T3 0.188
14 em_coping ~~ em_coping 0.000
15 i ~~ i 0.000
16 s ~~ s 0.000
17 i ~~ s 0.779
18 aces ~~ aces NA
lhs op rhs pvalue
1 i =~ cog_T1 NA
2 i =~ cog_T2 NA
3 i =~ cog_T3 NA
4 s =~ cog_T1 NA
5 s =~ cog_T2 NA
6 s =~ cog_T3 NA
7 i ~ aces 0.056
8 s ~ aces 0.000
9 s ~ em_coping 0.191
10 em_coping ~ aces 0.000
11 cog_T1 ~~ cog_T1 0.730
12 cog_T2 ~~ cog_T2 0.000
13 cog_T3 ~~ cog_T3 0.240
14 em_coping ~~ em_coping 0.000
15 i ~~ i 0.000
16 s ~~ s 0.000
17 i ~~ s 0.049
18 aces ~~ aces NA
lhs op rhs pvalue
1 i =~ cog_T1 NA
2 i =~ cog_T2 NA
3 i =~ cog_T3 NA
4 s =~ cog_T1 NA
5 s =~ cog_T2 NA
6 s =~ cog_T3 NA
7 i ~ aces 0.000
8 s ~ aces 0.016
9 s ~ em_coping 0.002
10 em_coping ~ aces 0.000
11 cog_T1 ~~ cog_T1 0.000
12 cog_T2 ~~ cog_T2 0.000
13 cog_T3 ~~ cog_T3 0.147
14 em_coping ~~ em_coping 0.000
15 i ~~ i 0.000
16 s ~~ s 0.000
17 i ~~ s 0.843
18 aces ~~ aces NA
lhs op rhs pvalue
1 i =~ cog_T1 NA
2 i =~ cog_T2 NA
3 i =~ cog_T3 NA
4 s =~ cog_T1 NA
5 s =~ cog_T2 NA
6 s =~ cog_T3 NA
7 i ~ aces 0.046
8 s ~ aces 0.000
9 s ~ em_coping 0.000
10 em_coping ~ aces 0.287
11 cog_T1 ~~ cog_T1 0.005
12 cog_T2 ~~ cog_T2 0.001
13 cog_T3 ~~ cog_T3 0.016
14 em_coping ~~ em_coping 0.000
15 i ~~ i 0.000
16 s ~~ s 0.000
17 i ~~ s 0.631
18 aces ~~ aces NA
lhs op rhs pvalue
1 i =~ cog_T1 NA
2 i =~ cog_T2 NA
3 i =~ cog_T3 NA
4 s =~ cog_T1 NA
5 s =~ cog_T2 NA
6 s =~ cog_T3 NA
7 i ~ aces 0.000
8 s ~ aces 0.000
9 s ~ em_coping 0.000
10 em_coping ~ aces 0.099
11 cog_T1 ~~ cog_T1 0.298
12 cog_T2 ~~ cog_T2 0.000
13 cog_T3 ~~ cog_T3 0.414
14 em_coping ~~ em_coping 0.000
15 i ~~ i 0.000
16 s ~~ s 0.000
17 i ~~ s 0.008
18 aces ~~ aces NA
lhs op rhs pvalue
1 i =~ cog_T1 NA
2 i =~ cog_T2 NA
3 i =~ cog_T3 NA
4 s =~ cog_T1 NA
5 s =~ cog_T2 NA
6 s =~ cog_T3 NA
7 i ~ aces 0.000
8 s ~ aces 0.000
9 s ~ em_coping 0.016
10 em_coping ~ aces 0.000
11 cog_T1 ~~ cog_T1 0.012
12 cog_T2 ~~ cog_T2 0.000
13 cog_T3 ~~ cog_T3 0.482
14 em_coping ~~ em_coping 0.000
15 i ~~ i 0.000
16 s ~~ s 0.000
17 i ~~ s 0.693
18 aces ~~ aces NA
lhs op rhs pvalue
1 i =~ cog_T1 NA
2 i =~ cog_T2 NA
3 i =~ cog_T3 NA
4 s =~ cog_T1 NA
5 s =~ cog_T2 NA
6 s =~ cog_T3 NA
7 i ~ aces 0.145
8 s ~ aces 0.000
9 s ~ em_coping 0.444
10 em_coping ~ aces 0.000
11 cog_T1 ~~ cog_T1 0.264
12 cog_T2 ~~ cog_T2 0.000
13 cog_T3 ~~ cog_T3 0.297
14 em_coping ~~ em_coping 0.000
15 i ~~ i 0.000
16 s ~~ s 0.000
17 i ~~ s 0.975
18 aces ~~ aces NA
lhs op rhs pvalue
1 i =~ cog_T1 NA
2 i =~ cog_T2 NA
3 i =~ cog_T3 NA
4 s =~ cog_T1 NA
5 s =~ cog_T2 NA
6 s =~ cog_T3 NA
7 i ~ aces 0.048
8 s ~ aces 0.000
9 s ~ em_coping 0.025
10 em_coping ~ aces 0.000
11 cog_T1 ~~ cog_T1 0.018
12 cog_T2 ~~ cog_T2 0.000
13 cog_T3 ~~ cog_T3 0.211
14 em_coping ~~ em_coping 0.000
15 i ~~ i 0.000
16 s ~~ s 0.000
17 i ~~ s 0.380
18 aces ~~ aces NA
lhs op rhs pvalue
1 i =~ cog_T1 NA
2 i =~ cog_T2 NA
3 i =~ cog_T3 NA
4 s =~ cog_T1 NA
5 s =~ cog_T2 NA
6 s =~ cog_T3 NA
7 i ~ aces 0.000
8 s ~ aces 0.001
9 s ~ em_coping 0.001
10 em_coping ~ aces 0.000
11 cog_T1 ~~ cog_T1 0.004
12 cog_T2 ~~ cog_T2 0.000
13 cog_T3 ~~ cog_T3 0.309
14 em_coping ~~ em_coping 0.000
15 i ~~ i 0.000
16 s ~~ s 0.000
17 i ~~ s 0.920
18 aces ~~ aces NA
# Calculate power (proportion of significant p-values below 0.05)
<- mean(results$intercept_pval < 0.05, na.rm = TRUE)
power_intercept <- mean(results$slope_pval < 0.05, na.rm = TRUE)
power_slope <- mean(results$indirect_pval < 0.05, na.rm = TRUE)
power_indirect
cat("Power for intercept path (ACES -> Intercept):", power_intercept, "\n")
Power for intercept path (ACES -> Intercept): 0.8
cat("Power for slope path (ACES -> Slope):", power_slope, "\n")
Power for slope path (ACES -> Slope): 1
cat("Power for indirect path (ACES -> Emotion Coping -> Slope):", power_indirect, "\n")
Power for indirect path (ACES -> Emotion Coping -> Slope): 0.8
To test the hypotheses that; i) the rates of improvement in cognitive functioning following a TBI vary as a function of childhood experiences of trauma (ACES) and ii) that this relationship between ACES and rates of recovery are mediated by emotion regulation, we will specify a latent growth curve model. This model assumes linear rates of cognitive improvement across 3-time points (baseline, 1 month, 3 months). We will include a direct path from ACES to the latent intercept and slope factors. We expect the path from ACES to the latent intercept to be significant and negative indicating that individuals with higher ACES scores will have more impaired cognitive functioning at discharge and that the path from ACES to the latent slope will be negative indicating that individuals with higher ACES scores will have slower and more moderate rates of improvement in cognitive functioning over time. We also expect that there will be a significant indirect effect from ACES to emotion regulation to the latent slope term indicating that the effect of ACES on rates of recovery is mediated through less adaptive emotion regulation skills.
We estimated the a priori power for this model via a Monte Carlo simulation conducted in R using the lavaan package (version 0.6-19). We specified a population model in which intercept factor loadings were set to 1 and slope factor loadings were set as equidistant loadings from 0-2, reflecting linear growth across 3 time points. The population values for the model parameters were set such that the the path from ACES to the latent intercept and latent slope factors were .3 and .35 respectively. The path from ACES to emotion regulation was set to .3 and the path from emotion regulation to the latent slope factor was set to .25. These parameter values are assumed to approximate meaningful expected values for the study. The sample size was set to n=150. Using 1,000 replications the random samples drawn from these population values resulted in an observed power of .8 or greater for all model parameters.