Supplementary File 1
Analysis Syntax
22 Jul 2024
Data & Setup
Load required packages and study data for analyses.
library(tidyverse) # Package `tidyverse` used for data manipulation.
library(haven) # Package `haven` used for data read and display.
library(gtsummary) # Package `gtsummary` used for data summary and display.
library(flextable) # Package `flextable` used for data and results display.
library(ltm) # Package `ltm` used for Cronbach's alpha analyses.
library(lavaan) # Package `lavaan` used for latent variable (factor) analyses.
library(semPlot) # Package `semPlot` used for presentation of latent variable results.
# data = read.csv("D:/zb. Publication Records/gg. Validation of ITQ in Colombian Sample/MIVIDA_T1_Complete_July2023.csv")
data <- read.csv("MIVIDA Data.csv") # Read basic datafile.
studydata <- data |> # Select study variables and view missigness patterns.
dplyr::select(pid, itq1:itq6, cptsd1:cptsd6)
naniar::vis_miss(studydata)
naniar::mcar_test(studydata)
## # A tibble: 1 × 4
## statistic df p.value missing.patterns
## <dbl> <dbl> <dbl> <int>
## 1 166. 184 0.820 19
studydata <- studydata |> # Select records containing 20% or less missingness.
dplyr::select(pid, itq1:itq6, cptsd1:cptsd6) |>
naniar::add_prop_miss() |>
dplyr::filter(prop_miss_all <= 0.20)
naniar::vis_miss(studydata)
imp_data <- missMethods::impute_EM(studydata[,-1], stochastic = FALSE) # Imputed missing values using Expectation Maximisation algorithm
studydata <- data |> dplyr::filter(pid %in% studydata$pid)
# Calculate multivariate Skewness, kurtosis, and Mahalanobis’ D-squared for primary study variables.
studydata |>
dplyr::select(itq1:itq6, cptsd1:cptsd6) |>
QuantPsyc::mult.norm()
## $mult.test
## Beta-hat kappa p-val
## Skewness 17.55496 1541.91093 0
## Kurtosis 224.07277 35.11217 0
##
## $Dsq
## [1] 7.8101663 2.3789357 25.3124837 9.0151992 4.5922404 14.6640500 14.9101819 14.2652619 11.9407174 7.2584670 3.3153556 10.5023796 3.2752811 31.2699952 14.0338503 23.4998034 23.0314206 46.9615066 21.5917045 45.1664971 5.3466227 3.2164886 17.4316553 22.8093485 8.0725499 9.1679141 4.4317530 26.1137018 19.1313479 8.4528501 2.2933530 8.6274699 15.9311870 20.6765998 23.8186598 30.9510052 23.1411681 2.2933530 24.0013370 7.2966129 11.8366161 3.2163139 8.8805025 8.9656910 19.2900048 19.3528110 17.5759675 32.9414501 4.1909718 32.8382663 2.3789357 17.8221761 19.7098975 17.3752368 13.6431963 6.9876815 14.1221589 12.0944139 5.6560329 26.3339665 10.7508520 18.3193540 15.6910734 12.4413890 14.3921608 8.0459582 8.6739987 31.7208747 10.2578199 24.1925084 31.1159191 25.7086676 20.2049943 3.7631822 15.0839029 12.9154104 13.2753601 3.9937498 7.1430047 7.3637993 11.2954127 3.4684244 18.8215692 3.3090739 10.2646221 14.3708261 1.8041156 13.3377208 2.6236537 2.6236537
## [91] 6.2198527 16.7704297 23.3247563 24.3714433 3.5515554 16.9035271 18.9036021 50.1132689 9.8872831 42.5298863 11.4767222 4.4653560 21.2303383 16.6636840 2.5870074 2.3789357 5.2353124 3.6667612 6.5655656 26.0565201 6.2742998 17.7498740 28.0667820 2.6236537 26.2804944 2.3789357 5.7772498 17.4743681 7.4633653 5.2722671 20.1713017 2.8533899 20.9837827 22.9367752 4.3892462 2.8237504 22.9367848 2.3789357 16.9607709 9.6062905 18.2429610 3.5039489 5.2752100 44.4677151 18.8286181 19.9741859 11.2059760 12.8804782 18.8674416 4.1896555 26.7835659 3.5725066 7.7276216 21.2628708 15.8671032 6.9574149 2.7659563 4.8234413 2.3789357 2.2933530 10.6910967 5.7600034 4.1701900 3.2164886 9.7560485 7.5556765 2.4572855 2.5166528 10.3148559 28.5491472 7.5660642 25.0014461 8.4258379 16.3477648 13.4969741 4.1800262 5.9115348 9.7894829 15.9851387 3.6629464 3.2195217 7.5413661 9.1010807 3.9046407 12.4627654 2.3789357 10.3864678 3.8781024 19.2150465 17.6246971
## [181] 11.6889475 4.3240824 4.8624814 14.5625862 24.7463528 6.7010791 14.2056327 4.2168524 20.1114894 17.8645319 18.0820605 14.0138054 2.6236537 7.8801731 13.2271675 18.2048520 10.6329118 12.1796007 28.6040441 2.3789357 9.7817349 11.5057322 19.0451997 16.5986183 9.8372381 3.0759502 6.7235752 26.6049453 8.8497813 6.2038206 2.3789357 17.2937285 14.1098911 11.2396282 16.1452396 6.0653050 8.5138110 13.0110813 3.5211326 8.7329076 7.4948737 9.5622261 17.0325360 4.7954968 25.2414705 2.3789357 9.9067818 7.3354404 15.3933755 5.1301824 0.6985622 6.2006744 32.4634999 31.9266938 24.8291498 2.3789357 18.9341749 13.9857712 17.8021564 4.7008169 5.9212593 17.3508736 8.5307865 9.8289054 11.3028811 16.7601241 7.8293811 2.2933530 8.7868351 18.1926046 2.2933530 23.5729739 13.5111325 25.2854123 8.9182772 11.4253604 14.5542841 4.3870867 8.0295180 2.3789357 11.4836599 20.1978683 4.2486430 22.6087057 3.2587356 10.7696991 3.8269675 4.3870867 13.8609636 16.8838116
## [271] 2.1217970 10.8173830 8.5303831 17.6882068 9.8832690 6.3760216 2.2933530 9.5499689 4.8362980 5.3380269 3.0759502 27.0914953 5.3179479 17.8991976 6.1431687 5.4014846 17.6121308 3.1915023 2.2933530 28.2306804 21.1915781 8.1648626 2.3789357 2.6236537 3.5964402 2.3789357 3.3090739 2.5450742 4.0273109 7.1276111 7.5389759 10.2312667 8.0635193 28.6901209 3.0907175 2.3789357 9.7460114 18.4089334 4.7281563 9.8225692 15.1103491 15.6937784 9.9510100 20.4618886 5.1425200 17.7552242 5.8664721 5.3069173 5.6312531 24.7220993 7.2292132 3.2432243 26.6811530 2.6236537 12.9784760 5.7313313 4.3870867 37.4255648 10.3454448 17.7351397 6.2017367 9.1845509 11.5484489 6.8997601 48.6777432 17.1403198 14.9449905 11.9582150 27.2540675 16.2053989 7.9000376 19.2584103 2.3789357 2.0786348 3.6020106 4.2486430 4.0771367 7.6741137 20.9178111 2.3789357 23.4698035 17.7340262 2.3789357 5.4522403 16.5817010 36.1374721 9.5372943 11.9683929 5.4929923 2.3789357
## [361] 8.5431963 6.2659148 4.8855789 9.7332111 5.0765551 16.3056629 3.5931629 2.3789357 14.6893845 2.3789357 0.6985622 16.5351599 15.5711784 13.1230658 3.7127995 14.3833933 4.3223171 31.1235092 4.1433902 10.0904697 5.8247663 8.0019930 2.3789357 17.1596335 22.8132331 11.7499308 28.1947374 16.2401894 3.9003819 23.9596295 12.0188065 2.3789357 25.6033285 8.0784351 11.3062205 10.2865295 6.9912080 3.2432243 12.6347592 11.8207404 38.1145566 16.6472959 12.2516974 21.9784940 3.4572649 17.9881179 18.8882255 19.8881290 3.3090739 28.1947374 15.6245793 7.8981397 7.4703466 4.5003369 6.8551533 3.7768248 17.8211694 11.4677740 18.4907351 25.5307799 2.3789357 17.4719604 11.9450228 11.2804742 37.9514375 4.6724912 6.7514993 4.9234361 3.2545845 37.1540478 9.5499689 40.5730576 2.7659563 7.3175916 13.6670152 21.1106778 5.4344646 3.5844839 17.8560069 30.1026494 10.9099712 6.1035311 6.4920741 3.9459885 22.9761497 10.3058855 4.8119095 2.3789357 2.3789357 19.6955448
## [451] 5.6312531 5.2098441 6.3245358 8.9921163 3.2164886 8.0802401 33.5797510 4.5194527 13.5264068 2.5166528 2.5166528 4.9650929 8.3115642 7.8063455 11.2452423 17.2253320 6.0815092 2.2933530 3.9542974 12.7033533 13.3284197 23.9383528 3.7222287 27.3311357 14.2448996 3.8352314 4.5289779 11.1192536 12.5449252 12.9045794 6.8964860 6.8172057 23.5762323 16.3033210 2.3789357 2.3789357 13.7340317 7.0264202 11.2386086 15.6021707 3.2432243 8.1513379 17.5720369 19.3737995 2.2933530 5.3159748 14.6559428 6.3360354 3.8863190 12.7147557 27.0523239 23.4660526 29.7074528 10.7436886 3.6667612 4.7577669 19.8476708 19.1788678 13.6165876 2.5723813 2.3789357 2.2483772 14.4989652 2.3789357 8.5973445 2.3789357 11.4391406 8.1145869 3.5725066 3.7127995 12.7393033 5.8664721 7.3076362 13.6463119 14.9556918 3.4684244 2.2906307
##
## $CriticalDsq
## [1] 28.29952Calculate sociodemographic and trauma endorsement summary statistics for Tables 1 & 2.
# Sociodemographics
table1 <-
studydata |>
dplyr::mutate(Assignment = factor(pclgroup, labels = c("Control", "Case")),
Gender = dem2.factor,
Sexuality = dem6.factor,
Martial_Status = dem4.factor,
Age = dem3,
Ethnicity = dem11.factor,
Religion = dem12.factor,
Education = dem14.factor,
Economic_Stat = dem17.factor,
Living_Area = dem8.factor,
Victim_Registry = vic1.factor,
FARC_Involved = farc.factor) |>
dplyr::mutate(Gender = recode_factor(Gender, # All demographic categories with n < 10 recorded endorsements recoded for parsimony and identity protection.
"Man" = "Male", # Recode small gender categories
"Female" = "Female",
"Transgender man" = "Other Gender Identity",
"Transgender woman" = "Other Gender Identity",
"Transgender/gender fluid" = "Other Gender Identity",
"Other" = "Other Gender Identity"),
Ethnicity = recode(Ethnicity, # Recode small ethnicity categories
"Gitano" = "Other",
"Other - Please specify" = "Other"),
Economic_Stat = recode_factor(Economic_Stat, # Collapse economic stratum.
"Stratum 1 (Very low)" = "Stratum 1-2 (Very Low - Low)",
"Stratum 2 (Low)" = "Stratum 1-2 (Very Low - Low)",
"Stratum 3 (Medium Low)" = "Stratum 3-4 (Medium Low - Medium)",
"Stratum 4 (Medium)" = "Stratum 3-4 (Medium Low - Medium)"),
Religion = recode_factor(Religion, # Recode small relgious categories
"Catholic" = "Catholic",
"Protestant" = "Other",
"Jeovah Witness" = "Other",
"Other - Please Specify" = "Other"),
Sexuality = recode_factor(Sexuality,
"Heterosexual" = "Heterosexual")) |>
dplyr::select(Gender, Sexuality, Martial_Status, Age, Ethnicity, Religion,
Education, Economic_Stat, Living_Area, Victim_Registry, FARC_Involved) |>
tbl_summary(label = list(Economic_Stat ~ "Economic Status", # Call `tbl_summary` to display demographic information and summary statistics.
Living_Area ~ "Area of Residence",
FARC_Involved ~ "Were you or a close family member (e.g. sibling or child)
forcibly recruited to be part of the revolutionary armed forces?"),
type = c(Age) ~ "continuous" ,
statistic = all_continuous() ~ "{mean} ({sd})",
digits = list(all_continuous() ~ 2, all_categorical() ~ c(0, 2)),
missing = "no")
# gt::gtsave(as_gt(table1), filename = "Table1.rtf")
# Trauma Endorsement, inc. Index
table2 <-
studydata |>
mutate(Gender = dem2.factor,
Sexuality = dem6.factor,
Age = dem3,
Ethnicity = dem11.factor,
Religion = dem12.factor,
Education = dem14.factor,
Economic_Stat = dem17.factor,
Living_Area = dem8.factor,
Victim_Registry = vic1.factor,
FARC_Involved = farc.factor) |>
dplyr::select(ntrauma, lec1.factor, lec2a.factor, lec3a.factor, lec4a.factor, # Select potentially traumatic experiences from data.
lec5a.factor, lec6a.factor, lec7a.factor, lec8a.factor,
lec9a.factor, lec10a.factor, lec11a.factor, lec12a.factor,
lec13a.factor, lec14a.factor, lec15a.factor, lec16a.factor,
lec17a.factor, lec18a.factor, lec19a.factor, lec20a.factor,
lec21a.factor, lec22a.factor, worsttrauma.factor) |>
tbl_summary( # Call summarise and relabel trauma events for display purposes.
label = list(lec1.factor = "Natural Disaster", lec2a.factor = "Fire or Explosion",
lec3a.factor = "Traffic Accident", lec4a.factor = "Serious Accident",
lec5a.factor = "Exposure to Toxic Substances", lec6a.factor = "Physical Assault",
lec7a.factor = "Armed Assault", lec8a.factor = "Sexual Assault",
lec9a.factor = "Unwanted Sexual Contact", lec10a.factor = "Combat",
lec11a.factor = "Captivity", lec12a.factor = "Life Threatening Illness or Injury",
lec13a.factor = "Severe Human Suffering", lec14a.factor = "Sudden Violent Death",
lec15a.factor = "Sudden Accidental Death", lec16a.factor = "Serious Injury, Harm, or Death Caused to Another",
lec17a.factor = "Parent or Partner Ridicule", lec18a.factor = "Physical Torture",
lec19a.factor = "Psychological Torture", lec20a.factor = "House or Property Damaged",
lec21a.factor = "Forced Displacement", lec22a.factor = "Other Stressful Event",
worsttrauma.factor = "Index Trauma"),
statistic = all_continuous() ~ "{mean} ({sd})",
digits = list(all_continuous() ~ 2, all_categorical() ~ c(0, 2)),
missing = "no")
# gt::gtsave(as_gt(table2), filename = "Table2.rtf")Calculate endorsement of CPTSD symptoms, Item response ≥ 2, and probable diagnosis according to International Trauma Questionnaire (ITQ) criteria based on these (1).
# ITQ Diagnosis Table
studydata |>
rename(itq7 = itqfi1, itq8 = itqfi2, itq9 = itqfi3,
itqdso1 = cptsd1, itqdso2 = cptsd2, itqdso3 = cptsd3,
itqdso4 = cptsd4, itqdso5 = cptsd5, itqdso6 = cptsd6,
itqdso7 = cptsdfi1, itqdso8 = cptsdfi2, itqdso9 = cptsdfi3) |>
mutate_at(vars("itq1", "itq2", "itq3", "itq4",
"itq5", "itq6", "itq7", "itq8",
"itq9", "itqdso1", "itqdso2", "itqdso3",
"itqdso4", "itqdso5", "itqdso6",
"itqdso7", "itqdso8", "itqdso9"), list(~ if_else(.x >= 2, true = 1, false = 0)))
Re <- studydata$itq1 + studydata$itq2
Av <- studydata$itq3 + studydata$itq4
Th <- studydata$itq5 + studydata$itq6
PTSDfi <- studydata$itqfi1 + studydata$itqfi2 + studydata$itqfi3
AD <- studydata$cptsd1 + studydata$cptsd2
NSC <- studydata$cptsd3 + studydata$cptsd4
DR <- studydata$cptsd5 + studydata$cptsd6
DSOfi <- studydata$cptsdfi1 + studydata$cptsdfi2 + studydata$cptsdfi3
studydata$Re_case <- ifelse(Re > 0, 1, 0)
studydata$Av_case <- ifelse(Av > 0, 1, 0)
studydata$Th_case <- ifelse(Th > 0, 1, 0)
studydata$PTSDfi_case <- ifelse(PTSDfi > 0, 1, 0)
studydata$AD_case <- ifelse(AD > 0, 1, 0)
studydata$NSC_case <- ifelse(NSC > 0, 1, 0)
studydata$DR_case <- ifelse(DR > 0, 1, 0)
studydata$DSOfi_case <- ifelse(DSOfi > 0, 1, 0)
studydata$PTSD_diag <- ifelse(studydata$Re_case + studydata$Av_case + studydata$Th_case + studydata$PTSDfi_case == 4, 1, 0)
studydata$CPTSD_diag <- ifelse(studydata$PTSD_diag + studydata$AD_case + studydata$NSC_case + studydata$DR_case + studydata$DSOfi_case == 5, 2, 0)
studydata$ITQ_diag <- ifelse(studydata$CPTSD_diag == 2, 2, ifelse(studydata$PTSD_diag == 1, 1, 0))table3 <- studydata |>
dplyr::select(Re_case, Av_case, Th_case, PTSDfi_case,
AD_case, NSC_case, DR_case, DSOfi_case,
ITQ_diag) |>
mutate(ITQ_diag = recode_factor(ITQ_diag, "0" = "No Diagnosis",
"1" = "PTSD", "2" = "CPTSD")) |>
tbl_summary(
label = list(Re_case = "Reexperiencing",
Av_case = "Avoidance",
Th_case = "Sense of Threat",
PTSDfi_case = "Functional Impairment (PTSD)",
AD_case = "Affective Dysregulation",
NSC_case = "Negative Self Concept",
DR_case = "Disturbed Relationships",
DSOfi_case = "Functional Impairment (DSO)",
ITQ_diag = "Diagnostic Category based on ICD-11 Criteria"), digits = list(all_continuous() ~ 2, all_categorical() ~ c(0, 2)),
missing = "no")
# gt::gtsave(as_gt(table3), filename = "Table3.rtf")Calculation of internal reliability (Cronbah’s α) for ITQ and subscales.
df1 <- # Select relevant study variables & relabel for clarity in primary analysis.
imp_data |>
rename(ptsd1 = itq1, ptsd2 = itq2, ptsd3 = itq3,
ptsd4 = itq4, ptsd5 = itq5, ptsd6 = itq6,
dso1 = cptsd1, dso2 = cptsd2, dso3 = cptsd3,
dso4 = cptsd4, dso5 = cptsd5, dso6 = cptsd6) |>
dplyr::select(ptsd1, ptsd2, ptsd3, ptsd4, ptsd5, ptsd6,
dso1, dso2, dso3, dso4, dso5, dso6)
ITQ <- df1 |> dplyr::select(ptsd1, ptsd2, ptsd3, ptsd4, ptsd5, ptsd6, # Select variables for scale reliability analysis.
dso1, dso2, dso3, dso4, dso5, dso6)
PTSDsubscale <- df1 |> dplyr::select(ptsd1, ptsd2, ptsd3, ptsd4, ptsd5, ptsd6) # Select variables for subscale reliability analysis.
DSOsubscale <- df1 |> dplyr::select(dso1, dso2, dso3, dso4, dso5, dso6) # Select variables for subscale reliability analysis.
df1 |> dplyr::select(ptsd1, ptsd2) |> cronbach.alpha() # Select variables for subscale reliability analysis.
##
## Cronbach's alpha for the 'dplyr::select(df1, ptsd1, ptsd2)' data-set
##
## Items: 2
## Sample units: 557
## alpha: 0.837
df1 |> dplyr::select(ptsd3, ptsd4) |> cronbach.alpha() # Select variables for subscale reliability analysis.
##
## Cronbach's alpha for the 'dplyr::select(df1, ptsd3, ptsd4)' data-set
##
## Items: 2
## Sample units: 557
## alpha: 0.845
df1 |> dplyr::select(ptsd5, ptsd6) |> cronbach.alpha() # Select variables for subscale reliability analysis.
##
## Cronbach's alpha for the 'dplyr::select(df1, ptsd5, ptsd6)' data-set
##
## Items: 2
## Sample units: 557
## alpha: 0.779
df1 |> dplyr::select(dso1, dso2) |> cronbach.alpha() # Select variables for subscale reliability analysis.
##
## Cronbach's alpha for the 'dplyr::select(df1, dso1, dso2)' data-set
##
## Items: 2
## Sample units: 557
## alpha: 0.698
df1 |> dplyr::select(dso3, dso4) |> cronbach.alpha() # Select variables for subscale reliability analysis.
##
## Cronbach's alpha for the 'dplyr::select(df1, dso3, dso4)' data-set
##
## Items: 2
## Sample units: 557
## alpha: 0.863
df1 |> dplyr::select(dso5, dso6) |> cronbach.alpha() # Select variables for subscale reliability analysis.
##
## Cronbach's alpha for the 'dplyr::select(df1, dso5, dso6)' data-set
##
## Items: 2
## Sample units: 557
## alpha: 0.855
# Internal Reliability for ITQ, PTSD, and DSO Scales
cronbach.alpha(ITQ, na.rm = T)
##
## Cronbach's alpha for the 'ITQ' data-set
##
## Items: 12
## Sample units: 557
## alpha: 0.919
cronbach.alpha(PTSDsubscale, na.rm = T)
##
## Cronbach's alpha for the 'PTSDsubscale' data-set
##
## Items: 6
## Sample units: 557
## alpha: 0.896
cronbach.alpha(DSOsubscale, na.rm = T)
##
## Cronbach's alpha for the 'DSOsubscale' data-set
##
## Items: 6
## Sample units: 557
## alpha: 0.895Confirmatory Factor Analysis
Confirmatory Factor Analysis allows for verifying the factor structure of a set of observed variables by testing the hypothesis that a relationship exists between observed/ manifest variables and specified underlying constructs (2). Confirmatory Factor Analysis approaches specifically hypothesise the structure of latent variables, i.e. specifying the number and arrangement of latent variables based on existing evidence. Fit indices, factor loadings, and inter-factor correlations of variables are estimated through specification of hypothesised models and these metrics compared to identify the model that provides the best fit to the data.
It is recommended that model fit is evaluated using a host of fit indices and the decision of optimal model fit guided by interpretation of overall goodness-of-fit across these based on established guidelines (2). The following fit indices were interpreted based on established recommendations provided by Hu and Bentler (2) indicating acceptable model fit: Comparative Fit Index (CFI) and Tucker Lewis Index (TLI) ≥ 0.95, Root Mean Square Error of Approximation (RMSEA) and Standardized Root Mean Square Residual (SRMR) ≤ 0.6. Akaike Information Criterion (AIC) and Bayesian Information Criterion (BIC) included in this Supplementary File provide additional measures of absolute model fit where lower values are considered indicative of more optimal model fit (3).
# One Factor
Model1 <- 'CPTSD =~ ptsd1 + ptsd2 + ptsd3 + ptsd4 + ptsd5 + ptsd6 +
dso1 + dso2 + dso3 + dso4 + dso5 + dso6'
cfaModel1 <- cfa(Model1, data = df1)
M1summary <- summary(cfaModel1, fit.measures=TRUE, standardized=TRUE)
# Two Factor
Model2 <- 'PTSD =~ ptsd1 + ptsd2 + ptsd3 + ptsd4 + ptsd5 + ptsd6
DSO =~ dso1 + dso2 + dso3 + dso4 + dso5 + dso6'
cfaModel2 <- cfa(Model2, data = df1)
M2summary <- summary(cfaModel2, fit.measures=TRUE, standardized=TRUE)
# Six Factor
Model3 <- 'Re =~ ptsd1 + ptsd2
Av =~ ptsd3 + ptsd4
SoT =~ ptsd5 + ptsd6
AD =~ dso1 + dso2
NSC =~ dso3 + dso4
DR =~ dso5 + dso6'
cfaModel3 <- cfa(Model3, data = df1)
M3summary <- summary(cfaModel3, fit.measures=TRUE, standardized=TRUE)
# Six Factor, One Higher Order
Model4 <- 'Re =~ ptsd1 + ptsd2
Av =~ ptsd3 + ptsd4
SoT =~ ptsd5 + ptsd6
AD =~ dso1 + dso2
NSC =~ dso3 + dso4
DR =~ dso5 + dso6
CPTSD =~ Re + Av + SoT + AD + NSC+ DR'
cfaModel4 <- cfa(Model4, data = df1)
M4summary <- summary(cfaModel4, fit.measures=TRUE, standardized=TRUE)
# Six Factor, Two Higher Order
Model5 <- 'Re =~ ptsd1 + ptsd2
Av =~ ptsd3 + ptsd4
SoT =~ ptsd5 + ptsd6
AD =~ dso1 + dso2
NSC =~ dso3 + dso4
DR =~ dso5 + dso6
PTSD =~ Re + Av + SoT
DSO =~ AD + NSC+ DR'
cfaModel5 <- cfa(Model5, data = df1)
M5summary <- summary(cfaModel5, fit.measures=TRUE, standardized=T)
# Model Comparison table
models <- list(cfaModel1, cfaModel2, cfaModel3, cfaModel4, cfaModel5)
semTable::compareLavaan(models, chidif = FALSE) |>
rownames_to_column("Model") |>
flextable()
## chisq df pvalue rmsea cfi tli srmr aic bic
## Model3 115.397 39 0 0.059 0.982 0.970 0.023 17732.14 17900.72
## Model5 169.996 47 0 0.069 0.971 0.959 0.036 17770.74 17904.74
## Model4 452.047 48 0 0.123 0.905 0.869 0.079 18050.79 18180.47
## Model2 456.546 53 0 0.117 0.905 0.882 0.052 18045.29 18153.36
## Model1 1170.549 54 0 0.193 0.737 0.678 0.100 18757.30 18861.04Model | chisq | df | pvalue | rmsea | cfi | tli | srmr | aic | bic |
|---|---|---|---|---|---|---|---|---|---|
Model3 | 115.397 | 39 | 0 | 0.059 | 0.982 | 0.970 | 0.023 | 17,732.14 | 17,900.72 |
Model5 | 169.996 | 47 | 0 | 0.069 | 0.971 | 0.959 | 0.036 | 17,770.74 | 17,904.74 |
Model4 | 452.047 | 48 | 0 | 0.123 | 0.905 | 0.869 | 0.079 | 18,050.79 | 18,180.47 |
Model2 | 456.546 | 53 | 0 | 0.117 | 0.905 | 0.882 | 0.052 | 18,045.29 | 18,153.36 |
Model1 | 1,170.549 | 54 | 0 | 0.193 | 0.737 | 0.678 | 0.100 | 18,757.30 | 18,861.04 |
Item-Factor Loadings and Inter-factor Correlations
A correlated first order model [Model 3] was found to provide best fit to these data based on the recommendations cited previously (3). Metric describing the structure and component associations in this model these are reported below. Factor loadings represent coefficients indicating the extent to which an observed/ manifest variable is associated with a latent factor. These loadings may be understood to show how well an observed variable reflects or measures the underlying latent construct (2). While cut-off values indicating acceptability are not applied to this metric increasing values are considered indicative of robust prediction of the latent variable (2,5).
semPaths(cfaModel3, what = "paths", nCharNodes = 5, residuals = FALSE)
title("Model 3", line = 1, adj = 0)
M3summary # Full model results
## lavaan 0.6.17 ended normally after 54 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 39
##
## Number of observations 557
##
## Model Test User Model:
##
## Test statistic 115.397
## Degrees of freedom 39
## P-value (Chi-square) 0.000
##
## Model Test Baseline Model:
##
## Test statistic 4307.740
## Degrees of freedom 66
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 0.982
## Tucker-Lewis Index (TLI) 0.970
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -8827.072
## Loglikelihood unrestricted model (H1) -8769.373
##
## Akaike (AIC) 17732.144
## Bayesian (BIC) 17900.724
## Sample-size adjusted Bayesian (SABIC) 17776.920
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.059
## 90 Percent confidence interval - lower 0.047
## 90 Percent confidence interval - upper 0.072
## P-value H_0: RMSEA <= 0.050 0.104
## P-value H_0: RMSEA >= 0.080 0.003
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.023
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## Re =~
## ptsd1 1.000 0.978 0.805
## ptsd2 1.165 0.053 21.846 0.000 1.139 0.895
## Av =~
## ptsd3 1.000 1.063 0.848
## ptsd4 1.065 0.047 22.556 0.000 1.132 0.864
## SoT =~
## ptsd5 1.000 1.069 0.765
## ptsd6 1.059 0.056 18.937 0.000 1.132 0.834
## AD =~
## dso1 1.000 0.650 0.602
## dso2 1.700 0.120 14.164 0.000 1.105 0.897
## NSC =~
## dso3 1.000 1.100 0.886
## dso4 0.891 0.037 24.352 0.000 0.981 0.860
## DR =~
## dso5 1.000 1.094 0.892
## dso6 0.911 0.037 24.523 0.000 0.996 0.837
##
## Covariances:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## Re ~~
## Av 0.858 0.072 11.952 0.000 0.825 0.825
## SoT 0.836 0.075 11.158 0.000 0.799 0.799
## AD 0.391 0.045 8.644 0.000 0.614 0.614
## NSC 0.601 0.063 9.590 0.000 0.559 0.559
## DR 0.576 0.062 9.361 0.000 0.539 0.539
## Av ~~
## SoT 0.921 0.080 11.521 0.000 0.810 0.810
## AD 0.401 0.047 8.494 0.000 0.580 0.580
## NSC 0.584 0.065 8.934 0.000 0.499 0.499
## DR 0.599 0.065 9.179 0.000 0.515 0.515
## SoT ~~
## AD 0.498 0.055 9.094 0.000 0.717 0.717
## NSC 0.638 0.070 9.072 0.000 0.542 0.542
## DR 0.783 0.074 10.537 0.000 0.669 0.669
## AD ~~
## NSC 0.554 0.056 9.923 0.000 0.775 0.775
## DR 0.597 0.058 10.281 0.000 0.840 0.840
## NSC ~~
## DR 1.003 0.076 13.112 0.000 0.833 0.833
##
## Variances:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .ptsd1 0.521 0.043 12.242 0.000 0.521 0.353
## .ptsd2 0.321 0.044 7.341 0.000 0.321 0.198
## .ptsd3 0.443 0.043 10.316 0.000 0.443 0.281
## .ptsd4 0.437 0.046 9.426 0.000 0.437 0.254
## .ptsd5 0.810 0.064 12.746 0.000 0.810 0.415
## .ptsd6 0.560 0.057 9.839 0.000 0.560 0.304
## .dso1 0.742 0.049 15.263 0.000 0.742 0.637
## .dso2 0.295 0.059 4.962 0.000 0.295 0.195
## .dso3 0.332 0.038 8.658 0.000 0.332 0.215
## .dso4 0.339 0.033 10.253 0.000 0.339 0.260
## .dso5 0.307 0.035 8.704 0.000 0.307 0.204
## .dso6 0.424 0.036 11.909 0.000 0.424 0.299
## Re 0.956 0.088 10.899 0.000 1.000 1.000
## Av 1.131 0.097 11.710 0.000 1.000 1.000
## SoT 1.143 0.114 10.016 0.000 1.000 1.000
## AD 0.423 0.057 7.381 0.000 1.000 1.000
## NSC 1.211 0.096 12.602 0.000 1.000 1.000
## DR 1.196 0.093 12.841 0.000 1.000 1.000
as.data.frame(inspect(cfaModel3, what="std")$lambda) |> # Factor loadings
mutate_at(vars(Re:DR), list(~recode(.,`0` = NA_real_))) |>
rownames_to_column("Variable") |>
flextable() # |> save_as_docx(path = "Table5.docx")Variable | Re | Av | SoT | AD | NSC | DR |
|---|---|---|---|---|---|---|
ptsd1 | 0.8045463 | |||||
ptsd2 | 0.8954039 | |||||
ptsd3 | 0.8477119 | |||||
ptsd4 | 0.8636058 | |||||
ptsd5 | 0.7650278 | |||||
ptsd6 | 0.8341631 | |||||
dso1 | 0.6024721 | |||||
dso2 | 0.8974392 | |||||
dso3 | 0.8857850 | |||||
dso4 | 0.8599617 | |||||
dso5 | 0.8921283 | |||||
dso6 | 0.8370257 |
Inter-factor correlations indicate the extent to which latent variables in the model correlate with others, indicative of the dimensionality in the latent structure (4). These results shown in the correlation matrix below indicate that the latent factors identified were positively correlated, and strongest correlations observed within PTSD-related factors [Experiencing, Avoidance, and Sense of Threat] and DSO-related factors [Affect Dysregulation, Negative Self Concept, and Interpersonal Disturbance].
inspect(cfaModel3, what="std")$psi # Inter-factor correlations
## Re Av SoT AD NSC DR
## Re 1.000
## Av 0.825 1.000
## SoT 0.799 0.810 1.000
## AD 0.614 0.580 0.717 1.000
## NSC 0.559 0.499 0.542 0.775 1.000
## DR 0.539 0.515 0.669 0.840 0.833 1.000
IFcorr <- as.matrix(inspect(cfaModel3, what="std")$psi)
# Correlation Matrix Plot
## Function to grab values for corr matrix
get_upper_tri <- function(cormat){
cormat[lower.tri(cormat)]<- NA
return(cormat)
}
upper_tri <- get_upper_tri(round(IFcorr, 3)) # Required inter-factor correlation figure
## Inter-Factor Correlation Matrix Figure
# jpeg("Interfactor_correlation.jpeg", width = 600, height = 300)
ggplot(data = reshape2::melt(upper_tri, na.rm = TRUE),
aes(Var2, Var1, fill = value))+
geom_tile(color = "white") +
geom_text(aes(Var2, Var1, label=sprintf("%0.2f", round(value, digits = 3))), color = "black", size = 4) +
scale_fill_gradient2(low = "#175676", high = "#ba324f", mid = "white",
midpoint = 0, limit = c(-1,1), space = "Lab",
name="Correlation") +
theme_minimal() +
theme(axis.title.x = element_blank(),
axis.title.y = element_blank(),
panel.grid.major = element_blank(),
panel.border = element_blank(),
panel.background = element_blank(),
axis.ticks = element_blank(),
legend.justification = c(1, 0),
legend.position = c(0.3, 0.7),
legend.direction = "horizontal")+
guides(fill = guide_colorbar(barwidth = 7, barheight = 1,
title.position = "top", title.hjust = 0.5))
# dev.off()Plotted figures representing models in analyses, and descriptions reference for variable labels.
# semPaths(object, what = "paths", whatLabels, style, layout = "tree",
# intercepts = TRUE, residuals = TRUE, thresholds = TRUE,
# intStyle = "multi", rotation = 1, curve, nCharNodes = 3)
# jpeg("CFA_models.jpeg", width = 1200, height = 600)
layout(matrix(c(1, 4,
1, 4,
2, 4,
2, 5,
3, 5,
3, 5), ncol=2, nrow=6, byrow=T)) # Layout matrix divided 2x6 to even spacing.
# layout.show(5) # Use `layout.show` to inspect plotting framework.
semPaths(cfaModel1, what = "paths", nCharNodes = 5, residuals = F)
title("Model 1", line = 1, adj = 0) # Title position, `adj = 0` == left aligned.
semPaths(cfaModel2, what = "paths", nCharNodes = 5, residuals = F)
title("Model 2", line = 1, adj = 0)
semPaths(cfaModel3, what = "paths", nCharNodes = 5, residuals = F)
title("Model 3", line = 3.4, adj = 0)
semPaths(cfaModel4, what = "paths", nCharNodes = 5, residuals = F)
title("Model 4", line = 1, adj = 0)
semPaths(cfaModel5, what = "paths", nCharNodes = 5, residuals = F)
title("Model 5", line = 1, adj = 0)
# dev.off()
# jpeg("Model3.jpeg", width = 1200, height = 400)
layout(1) # Reset plotting layout to single.
semPaths(cfaModel3, what = "std", edge.label.cex = 0.75,
edge.color = "#44355B", nCharNodes = 5, residuals = F)
title("Model 3", line = 1, adj = 0)
# dev.off()
## Item Labels
# PTSD1 = "Having disturbing dreams that reproduce part of the experience or are clearly related to the experience?"
# PTSD2 = "Having intense images or memories that sometimes come to your mind, in which you feel that the experience is happening again here and now?"
# PTSD3 = "Avoid thoughts, feelings, physical sensations, or other internal stimuli that remind you of the experience?"
# PTSD4 = "Avoid people, places, conversations, objects, activities, situations, or other external stimuli that remind you of the experience?"
# PTSD5 = "Be super alert, vigilant or on guard?"
# PTSD6 = "Feeling startled or easily frightened?"
# DSO1 = "When I'm upset, it takes me a long time to calm down"
# DSO2 = "I feel emotionally disconnected or emotionally off"
# DSO3 = "I feel like a failure"
# DSO4 = "I feel worthless"
# DSO5 = "I feel distant or cut-off from people"
# DSO6 = "I find it difficult to be emotionally close to people"
References
1. Cloitre M, Shevlin M, Brewin CR, Bisson JI, Roberts NP, Maercker A, et al. The International Trauma Questionnaire: development of a self-report measure of ICD-11 PTSD and complex PTSD. Acta Psychiatr Scand. 2018;138(6):536–46.
2. Brown TA, Moore MT. Confirmatory factor analysis. In: Handbook of structural equation modeling. New York, NY, US: The Guilford Press; 2012. p. 361–79.
3. Hu L t, Bentler PM. Cutoff criteria for fit indexes in covariance structure analysis: Conventional criteria versus new alternatives. Struct Equ Model Multidiscip J. 1999;6(1):1–55.
4. Schmitt TA, Sass DA, Chappelle W, Thompson W. Selecting the “Best” Factor Structure and Moving Measurement Validation Forward: An Illustration. J Pers Assess. 2018 Jul 4;100(4):345–62.
5. Kline RB. Principles and practice of structural equation modeling, 4th ed. New York, NY, US: Guilford Press; 2016. xvii, 534 p. (Principles and practice of structural equation modeling, 4th ed.).