--- title: "A Network Model of Putative Resilience Factors for Groups of Adolescents with and without Exposure to Childhood Adversity" author: "J. Fritz, Dr E. I. Fried, Prof I. M. Goodyer, Dr P. O. Wilkinson, & Dr A.-L. van Harmelen" date: "2018/August/22" output: word_document: default --- ## 0: Table of content ```{r} ##1: Opening data and loading packages ##2: Resilience factor (RF) and general distress comparisons for adolescents with ('CA') and adolescents without ('noCA') childhood adversity ##3: RF network models for CA adolescents (for results including the general distress variable see section 6) ##4: RF network models for noCA adolescents (for results including the general distress variable see section 7) ##5: Contrasting the CA and the noCA network structure (for results including the general distress variable see section 9) ##6: The RF network structure for CA adolescents including general distress ##7: The RF network structure for noCA adolescents including general distress ##8: Correlations and regularized partial correlations between the RFs and the general distress variable ##9: Contrasting the CA and the noCA network structure including the general distress variable ##10: Contrasting CA and noCA RF networks corrected for the general distress variable ##11: Shortest pathways between RFs and the general distress variable for CA and noCA groups ##12: correlation of the 'RF-general distress' interrelations between the two groups ##13: Methods: Sample descriptives ##14: Supplement I: Methodological Rationale ##15: Supplement II: Variable Preparation: Means and standard deviations or inter quartile ranges and frequencies for untransformed RFs and the general distress variable ##16: Supplement III: Association and concentration networks for CA and noCA adolescents ##17: Supplement IV: Adjacency matrices of the main models ##18: Supplement V: Interconnectedness of RFs ##19: Supplement VI: Accuracy and Stability of the CA network ##20: Supplement VI: Accuracy and Stability of the noCA network ##21: Supplement VII: Sensitivity analyses for the CA network ##22: Supplement VII: Sensitivity analyses for the noCA network ##23: Supplement VIII: Main models depicted with faded 'RF-RF' and 'RF-General Distress' interrelations ##24: Supplement IX: Exploring the expressive suppression variable ##25: Supplement IX: RF network of the CA group without the expressive suppression variable ##26: Supplement IX: RF network of the noCA group without the expressive suppression variable ##27: Supplement IX: Comparing the CA and the noCA network structure without the expressive suppression variable ##28: Supplement X: Node strength and expected influence for networks corrected for the general distress variable, for both CA and noCA ##29: Supplement XI: Network Pathways between the RFs and General Distress ##30: Supplement XII: Exploring the complex interplay between CA, RFs and general distress ##31: Supplement XIII: Reliability and/or Validity Information for the used Measures ##32: Supplement XIV: Individual RF Interrelation Differences between the CA and the no-CA Networks ``` ##1: Opening data and loading packages ```{r} Network_Data <- read.csv("File_directory\\Data_Name") ####A: open packages library("pastecs") library("qgraph") library("mgm") library("bootnet") library("NetworkComparisonTest") library("pastecs") library("dplyr") library("coin") library("ggplot2") library("lavaan") library("rockchalk") library("mediation") library("gvlma") ``` ##2: Resilience factor (RF) and general distress comparisons for adolescents with ('CA') and adolescents without ('noCA') childhood adversity ```{r} #CA vs no-CA table(Network_Data$CA_dic) #gender for the CA (=1) and the noCA (=0) group table(Network_Data$sex) chi_sex <- table(Network_Data$sex, Network_Data$CA_dic); chi_sex #age stat.desc(Network_Data$age_srqt1) #a: friendship support stat.desc(Network_Data$friendsupport_tr, norm = TRUE) t.test(Network_Data$friendsupport_tr ~ Network_Data$CA_dic, var.equal = FALSE) #b: family support stat.desc(Network_Data$familysupport_tr, norm = TRUE) t.test(Network_Data$familysupport_tr ~ Network_Data$CA_dic, var.equal = FALSE) #c: family cohesion stat.desc(Network_Data$familycohesion_tr, norm = TRUE) t.test(Network_Data$familycohesion_tr ~ Network_Data$CA_dic, var.equal = FALSE) #d: negative self-esteem stat.desc(Network_Data$negativeselfesteem_tr, norm = TRUE) t.test(Network_Data$negativeselfesteem_tr ~ Network_Data$CA_dic, var.equal = FALSE) #e: positive self-esteem stat.desc(Network_Data$positiveselfesteem_tr, norm = TRUE) t.test(Network_Data$positiveselfesteem_tr ~ Network_Data$CA_dic, var.equal = FALSE) #f: brooding stat.desc(Network_Data$brooding_tr, norm = TRUE) t.test(Network_Data$brooding_tr ~ Network_Data$CA_dic, var.equal = FALSE) #g: reflection stat.desc(Network_Data$reflection_tr, norm = TRUE) t.test(Network_Data$reflection_tr ~ Network_Data$CA_dic, var.equal = FALSE) #h: distress tolerance stat.desc(Network_Data$distresstolerance_tr, norm = TRUE) t.test(Network_Data$distresstolerance_tr ~ Network_Data$CA_dic, var.equal = FALSE) #i: aggression table(Network_Data$aggressionbin) chi_agg <- table(Network_Data$aggressionbin, Network_Data$CA_dic); chi_agg chisq.test(chi_agg) #j: expressive suppression table(Network_Data$expressivesuppressionbin) chi_exs <- table(Network_Data$expressivesuppressionbin, Network_Data$CA_dic); chi_exs chisq.test(chi_exs) #k: general distress stat.desc(Network_Data$generaldistress_tr, norm = TRUE) t.test(Network_Data$generaldistress_tr ~ Network_Data$CA_dic, var.equal = FALSE) ``` ##3: RF network models for CA adolescents (for results including the general distress variable see section 6) ```{r} ####A: Subsetting the data in CA and noCA CA_Network_Data <- subset(Network_Data, Network_Data$CA_dic==1) head(CA_Network_Data); nrow(CA_Network_Data) ####B: Making the CA network data frame CA_network_c <- data.frame(CA_Network_Data$friendsupport_tr, CA_Network_Data$familysupport_tr, CA_Network_Data$familycohesion_tr, CA_Network_Data$positiveselfesteem_tr, CA_Network_Data$negativeselfesteem_tr, CA_Network_Data$brooding_tr, CA_Network_Data$reflection_tr, CA_Network_Data$distresstolerance_tr, CA_Network_Data$aggressionbin, CA_Network_Data$expressivesuppressionbin) ####C: Making the correlation matrix CA_network_cor <- cor_auto(CA_network_c); CA_network_cor #Variables detected as ordinal: CA_Network_Data.aggressionbin; CA_Network_Data.expressivesuppressionbin max(CA_network_cor[upper.tri(CA_network_cor,diag=FALSE)]) # 0.686106 ####D: Main model with average layout (as presented in paper) #a: average layout adopted from models including general distress, to enhance the visual comparability Layout <- structure(c(-0.126424547383095, 0.799405165379372, 0.41431789943595, 1, 0.826992489872837, 0.0167398316516412, -0.526596792250546, -1, -0.365765573385925, -0.532477942540536, -0.432974713716992, -0.790641633253371, -0.943141335210814, 0.0462132496070069, 0.580670902412863, 1, 0.945560460165707, -0.178658845633278, 0.183925750895849, -1), .Dim = c(10L, 2L)) #b: regularized partial correlation network for CA #pdf(file = "CA_RF_glasso.pdf") CA_RF_glasso <- qgraph(CA_network_cor, layout = Layout, graph = "glasso", details = TRUE, nodeNames = colnames(CA_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(CA_network_c)) #dev.off() ####E: Centrality analysis #a: node strength and node expected influence CA_centrality <- centrality_auto(CA_RF_glasso); CA_centrality #b: node predictability CA_network_mgm <- na.omit(CA_network_c); nrow(CA_network_mgm) # 508 RFtype<- c(rep("g", 8),rep("c", 2)); RFtype RFlev<- c(rep(1,8),rep(2,2)); RFlev set.seed(1) fit_obj <- mgm(data = CA_network_mgm, type = RFtype, lev = RFlev, rule.reg = "AND", k = 2, binarySign = TRUE, scale = TRUE) fit_obj$pairwise$wadj fit_obj$pairwise$signs pred_obj <- predict(fit_obj, data = CA_network_mgm, errorCon = c("RMSE", "R2"), errorCat = c("CC", "nCC")) pred_obj$errors errors <- c(0.14, 0.48, 0.54, 0.42, 0.57, 0.61, 0.44, 0.11, 0.00, 0.07); errors pieColor <- c(rep("#90B4D4",8),rep("#EB9446",2)) #pdf(file = "CA_RF_pred.pdf") CA_RF_pred <- qgraph(fit_obj$pairwise$wadj, edge.color = fit_obj$pairwise$edgecolor, pie = errors, pieColor = pieColor, layout=Layout, nodeNames = colnames(CA_network_mgm), label.cex = .9, theme = "classic", curveAll = FALSE, fade = FALSE, legend = TRUE, cut = 0) #dev.off() CA_am_mgm <- fit_obj$pairwise$wadj; CA_am_mgm CA_am_mgm_signs <- fit_obj$pairwise$signs; CA_am_mgm_signs ``` ##4: RF network models for noCA adolescents (for results including the general distress variable see section 7) ```{r} ####A: Subsetting the data in CA and noCA noCA_Network_Data <- subset(Network_Data, Network_Data$CA_dic==0) head(noCA_Network_Data); nrow(noCA_Network_Data) ####B: Making the noCA network data frame noCA_network_c <- data.frame(noCA_Network_Data$friendsupport_tr, noCA_Network_Data$familysupport_tr, noCA_Network_Data$familycohesion_tr, noCA_Network_Data$positiveselfesteem_tr, noCA_Network_Data$negativeselfesteem_tr, noCA_Network_Data$brooding_tr, noCA_Network_Data$reflection_tr, noCA_Network_Data$distresstolerance_tr, noCA_Network_Data$aggressionbin, noCA_Network_Data$expressivesuppressionbin) ####C: Making the correlation matrix noCA_network_cor <- cor_auto(noCA_network_c); noCA_network_cor #Variables detected as ordinal: noCA_Network_Data.aggressionbin; noCA_Network_Data.expressivesuppressionbin max(noCA_network_cor[upper.tri(noCA_network_cor,diag=FALSE)]) # 0.6609587 ####D: Main model with average layout (as presented in paper) #a: regularized partial correlation network for noCA #pdf(file = "noCA_RF_glasso.pdf") noCA_RF_glasso <- qgraph(noCA_network_cor, layout = Layout, graph = "glasso", details = TRUE, nodeNames = colnames(noCA_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(noCA_network_c)) #dev.off() ####E: Centrality analysis #a: node strength and node expected influence noCA_centrality <- centrality_auto(noCA_RF_glasso); noCA_centrality #b: node predictability noCA_network_mgm <- na.omit(noCA_network_c); nrow(noCA_network_mgm) # 443 RFtype<- c(rep("g", 8),rep("c", 2)); RFtype RFlev<- c(rep(1,8),rep(2,2)); RFlev set.seed(1) fit_obj <- mgm(data = noCA_network_mgm, type = RFtype, lev = RFlev, rule.reg = "AND", k = 2, binarySign = TRUE, scale = TRUE) fit_obj$pairwise$wadj fit_obj$pairwise$signs pred_obj <- predict(fit_obj, data = noCA_network_mgm, errorCon = c("RMSE", "R2"), errorCat = c("CC", "nCC")) pred_obj$errors errors <- c(0.18, 0.47, 0.48, 0.42, 0.52, 0.53, 0.34, 0.08, 0.02, 0.02); errors pieColor <- c(rep("#90B4D4",8),rep("#EB9446",2)) #pdf(file = "noCA_RF_pred.pdf") noCA_RF_pred <- qgraph(fit_obj$pairwise$wadj, edge.color = fit_obj$pairwise$edgecolor, pie = errors, pieColor = pieColor, layout=Layout, nodeNames = colnames(noCA_network_mgm), label.cex = .9, theme = "classic", curveAll = FALSE, fade = FALSE, legend = TRUE, cut = 0) #dev.off() noCA_am_mgm <- fit_obj$pairwise$wadj; noCA_am_mgm noCA_am_mgm_signs <- fit_obj$pairwise$signs; noCA_am_mgm_signs ``` ##5: Contrasting the CA and the noCA network structure (for results including the general distress variable see section 9) ```{r} ####A: correlate CA with noCA RF interrelation matrix adjmat_CA <- getWmat(CA_RF_glasso); adjmat_CA cor_CA <- adjmat_CA[upper.tri(adjmat_CA,diag=FALSE)]; cor_CA adjmat_noCA <- getWmat(noCA_RF_glasso); adjmat_noCA cor_noCA <- adjmat_noCA[upper.tri(adjmat_noCA,diag=FALSE)]; cor_noCA cor_pears <- cor(cor_CA,cor_noCA, method = "pearson"); cor_pears # 0.9362756 cor_spear <- cor(cor_CA,cor_noCA, method = "spearman"); cor_spear # 0.6290396 ####B: Network comparison test (NCT; permutation test) approach for the comparison of the CA vs noCA networks #Reference: van Borkulo, C. D. et al. Comparing network structures on three aspects: A permutation test. (Manuscript submitted for publication, 2017). #RUN THE BELONGING NCT SCRIPT WITH HOLM-BONFERRONI CORRECTION, WHICH CAN BE DERIVED FROM: https://github.com/cvborkulo/NetworkComparisonTest/blob/master/R/NCT_cor_auto.R #a: rename variables so that they are equal for the two groups CA_network_NCT <- CA_network_c CA_network_NCT <- rename(CA_network_NCT, friendsupport = CA_Network_Data.friendsupport_tr) CA_network_NCT <- rename(CA_network_NCT, familysupport = CA_Network_Data.familysupport_tr) CA_network_NCT <- rename(CA_network_NCT, familycohesion = CA_Network_Data.familycohesion_tr) CA_network_NCT <- rename(CA_network_NCT, positiveselfesteem = CA_Network_Data.positiveselfesteem_tr) CA_network_NCT <- rename(CA_network_NCT, negativeselfesteem = CA_Network_Data.negativeselfesteem_tr) CA_network_NCT <- rename(CA_network_NCT, expressivesuppression = CA_Network_Data.expressivesuppressionbin) CA_network_NCT <- rename(CA_network_NCT, distresstolerance = CA_Network_Data.distresstolerance_tr) CA_network_NCT <- rename(CA_network_NCT, aggression = CA_Network_Data.aggressionbin) CA_network_NCT <- rename(CA_network_NCT, brooding = CA_Network_Data.brooding_tr) CA_network_NCT <- rename(CA_network_NCT, reflection = CA_Network_Data.reflection_tr) head(CA_network_NCT);nrow(CA_network_NCT) noCA_network_NCT <- noCA_network_c noCA_network_NCT <- rename(noCA_network_NCT, friendsupport = noCA_Network_Data.friendsupport_tr) noCA_network_NCT <- rename(noCA_network_NCT, familysupport = noCA_Network_Data.familysupport_tr) noCA_network_NCT <- rename(noCA_network_NCT, familycohesion = noCA_Network_Data.familycohesion_tr) noCA_network_NCT <- rename(noCA_network_NCT, positiveselfesteem = noCA_Network_Data.positiveselfesteem_tr) noCA_network_NCT <- rename(noCA_network_NCT, negativeselfesteem = noCA_Network_Data.negativeselfesteem_tr) noCA_network_NCT <- rename(noCA_network_NCT, expressivesuppression = noCA_Network_Data.expressivesuppressionbin) noCA_network_NCT <- rename(noCA_network_NCT, distresstolerance = noCA_Network_Data.distresstolerance_tr) noCA_network_NCT <- rename(noCA_network_NCT, aggression = noCA_Network_Data.aggressionbin) noCA_network_NCT <- rename(noCA_network_NCT, brooding = noCA_Network_Data.brooding_tr) noCA_network_NCT <- rename(noCA_network_NCT, reflection = noCA_Network_Data.reflection_tr) head(noCA_network_NCT);nrow(noCA_network_NCT) #b:perform NCT test set.seed(12) CA_vs_noCA <- NCT_cor_auto(CA_network_NCT, noCA_network_NCT, it=5000, binary.data = FALSE, paired = FALSE, weighted=TRUE, test.edges = TRUE, edges = list(c(1,2),c(1,3),c(1,4),c(1,5),c(1,6),c(1,7),c(1,8),c(1,9),c(1,10),c(2,3),c(2,4),c(2,5),c(2,6),c(2,7),c(2,8), c(2,9),c(2,10),c(3,4),c(3,5),c(3,6),c(3,7),c(3,8),c(3,9),c(3,10),c(4,5),c(4,6),c(4,7),c(4,8),c(4,9),c(4,10), c(5,6),c(5,7),c(5,8),c(5,9),c(5,10),c(6,7),c(6,8),c(6,9),c(6,10),c(7,8),c(7,9),c(7,10),c(8,9),c(8,10),c(9,10)), progressbar=TRUE) #c: descriptives summary(CA_vs_noCA) plot(CA_vs_noCA, what="network") plot(CA_vs_noCA, what="strength") plot(CA_vs_noCA, what="edge") ####C:To compare the Global Network Expected Influence, please run the function in the following scipt: "2_b_2_Analysis_Script_NCT_EI-auto_JF-EF-IG-PW-ALvH_2018.07July.03_v1". The script is adapted from: van Borkulo, C. D. et al. Comparing network structures on three aspects: A permutation test (Manuscript submitted for publication, 2017). The adapted script can be requested from the first author jf585@cam.ac.uk. #Please note, what reads now as 'Global strength per group' contains the values of the global network expected influence. #a:perform NCT test set.seed(12) EI_CA_vs_noCA <- NCT_cor_auto_EI(CA_network_NCT, noCA_network_NCT, it=5000, binary.data = FALSE, paired = FALSE, weighted=TRUE, test.edges = TRUE, edges = list(c(1,2),c(1,3),c(1,4),c(1,5),c(1,6),c(1,7),c(1,8),c(1,9),c(1,10),c(2,3),c(2,4),c(2,5),c(2,6),c(2,7),c(2,8), c(2,9),c(2,10),c(3,4),c(3,5),c(3,6),c(3,7),c(3,8),c(3,9),c(3,10),c(4,5),c(4,6),c(4,7),c(4,8),c(4,9),c(4,10), c(5,6),c(5,7),c(5,8),c(5,9),c(5,10),c(6,7),c(6,8),c(6,9),c(6,10),c(7,8),c(7,9),c(7,10),c(8,9),c(8,10),c(9,10)), progressbar=TRUE) summary(EI_CA_vs_noCA) ####D: To reveal the results with the uncorrected p-values, the "method = 'holm'" commands in the NCT script need to be replaced by "method = 'none'". ``` ##6: The RF network structure for CA adolescents including general distress ```{r} ####A: Making the network data frame CA_PP_network_c <- data.frame(CA_Network_Data$friendsupport_tr, CA_Network_Data$familysupport_tr, CA_Network_Data$familycohesion_tr, CA_Network_Data$positiveselfesteem_tr, CA_Network_Data$negativeselfesteem_tr, CA_Network_Data$brooding_tr, CA_Network_Data$reflection_tr, CA_Network_Data$distresstolerance_tr, CA_Network_Data$aggressionbin, CA_Network_Data$expressivesuppressionbin, CA_Network_Data$generaldistress_tr) ####B: Making the correlation matrix CA_PP_network_cor <- cor_auto(CA_PP_network_c); CA_PP_network_cor #Variables detected as ordinal: CA_Network_Data.aggressionbin; CA_Network_Data.expressivesuppressionbin max(CA_PP_network_cor[upper.tri(CA_PP_network_cor,diag=FALSE)]) # 0.686106 ####C: Main model with average layout (as presented in paper) #a: see syntax for retrieving layout in section VII PP_Layout <- structure(c(-0.126424547383095, 0.799405165379372, 0.41431789943595, 1, 0.826992489872837, 0.0167398316516412, -0.526596792250546, -1, -0.365765573385925, -0.532477942540536, 0.298338406009852, -0.432974713716992, -0.790641633253371, -0.943141335210814, 0.0462132496070069, 0.580670902412863, 1, 0.945560460165707, -0.178658845633278, 0.183925750895849, -1, 0.483519659486572), .Dim = c(11L, 2L)) #b: regularized partial correlation network for CA with the general distress variable #pdf(file = "CA_PP_RF_glasso.pdf") CA_PP_RF_glasso <- qgraph(CA_PP_network_cor, layout = PP_Layout, graph = "glasso", details = TRUE, nodeNames = colnames(CA_PP_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(CA_PP_network_c)) #dev.off() ####D: Checking Robustness and Stability for centrality analyses #a: interrelation accuracy set.seed(3) CA_PP_accuracy <- bootnet(CA_PP_network_c, default = "EBICglasso", type = "nonparametric", nCores = 8, nBoots = 2000, statistics = c("edge", "strength","expectedInfluence")) plot(CA_PP_accuracy$sample, layout = PP_Layout, fade = FALSE) plot(CA_PP_accuracy, labels = TRUE, order = "sample") plot(CA_PP_accuracy, labels = FALSE, order = "sample") #b: network centrality coefficient stability set.seed(4) CA_PP_stability <- bootnet(CA_PP_network_c, default = "EBICglasso", type = "case", nCores = 8, caseN = 25, nBoots =2000, statistics = c("strength","expectedInfluence")) plot(CA_PP_stability$sample, layout = PP_Layout, fade = FALSE) plot(CA_PP_stability, statistics = c("strength", "expectedInfluence")) corStability(CA_PP_stability, statistics = c("strength", "expectedInfluence"), cor = 0.7) ####E: Centrality analysis #a: node strength and node expected influence CA_PP_centrality <- centrality_auto(CA_PP_RF_glasso); CA_PP_centrality #b: node predictability CA_PP_network_mgm <- na.omit(CA_PP_network_c); nrow(CA_PP_network_mgm) # 508 RFtype<- c(rep("g", 8),rep("c", 2),rep("g", 1)); RFtype RFlev<- c(rep(1,8),rep(2,2),rep(1,1)); RFlev set.seed(1) fit_obj <- mgm(data = CA_PP_network_mgm, type = RFtype, lev = RFlev, rule.reg = "AND", k = 2, binarySign = TRUE, scale = TRUE) fit_obj$pairwise$wadj fit_obj$pairwise$signs pred_obj <- predict(fit_obj, data = CA_PP_network_mgm, errorCon = c("RMSE", "R2"), errorCat = c("CC", "nCC")) pred_obj$errors errors <- c(0.17, 0.47, 0.54, 0.43, 0.64, 0.65, 0.44, 0.12, 0.00, 0.02, 0.69); errors pieColor <- c(rep("#90B4D4",8),rep("#EB9446",2),rep("#90B4D4",1)) #pdf(file = "CA_PP_RF_pred.pdf") CA_PP_RF_pred <- qgraph(fit_obj$pairwise$wadj, edge.color = fit_obj$pairwise$edgecolor, pie = errors, pieColor = pieColor, layout=PP_Layout, nodeNames = colnames(CA_PP_network_mgm), label.cex = .9, theme = "classic", curveAll = FALSE, fade = FALSE, legend = TRUE, cut = 0) #dev.off() CA_PP_am_mgm <- fit_obj$pairwise$wadj; CA_PP_am_mgm CA_PP_am_mgm_signs <- fit_obj$pairwise$signs; CA_PP_am_mgm_signs ``` ##7: The RF network structure for noCA adolescents including general distress ```{r} ####A: Making the network data frame noCA_PP_network_c <- data.frame(noCA_Network_Data$friendsupport_tr, noCA_Network_Data$familysupport_tr, noCA_Network_Data$familycohesion_tr, noCA_Network_Data$positiveselfesteem_tr, noCA_Network_Data$negativeselfesteem_tr, noCA_Network_Data$brooding_tr, noCA_Network_Data$reflection_tr, noCA_Network_Data$distresstolerance_tr, noCA_Network_Data$aggressionbin, noCA_Network_Data$expressivesuppressionbin, noCA_Network_Data$generaldistress_tr) ####B: Making the correlation matrix noCA_PP_network_cor <- cor_auto(noCA_PP_network_c); noCA_PP_network_cor #Variables detected as ordinal: noCA_Network_Data.aggressionbin; noCA_Network_Data.expressivesuppressionbin max(noCA_PP_network_cor[upper.tri(noCA_PP_network_cor,diag=FALSE)]) # 0.6609587 ####C: Main model with average layout (as presented in paper) ##a1: regularized partial correlation network CA without average layout #set.seed(99) #CA_PP_RF_glasso_lay <- qgraph(CA_PP_network_cor, layout = "spring", graph = "glasso", details = TRUE, # nodeNames = colnames(CA_PP_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, # theme = "colorblind", fade = FALSE, sampleSize = nrow(CA_PP_network_c)) ##a2: regularized partial correlation network noCA without average layout #set.seed(99) #noCA_PP_RF_glasso_lay <- qgraph(noCA_PP_network_cor, layout = "spring", graph = "glasso", details = TRUE, # nodeNames = colnames(noCA_PP_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, # theme = "colorblind", fade = FALSE, sampleSize = nrow(noCA_PP_network_c)) #set.seed(99) #PP_Layout <- averageLayout(CA_PP_RF_glasso_lay,noCA_PP_RF_glasso_lay) #dput(PP_Layout) #PP_Layout <- structure(c(-0.126424547383095, 0.799405165379372, 0.41431789943595, #1, 0.826992489872837, 0.0167398316516412, -0.526596792250546, #-1, -0.365765573385925, -0.532477942540536, 0.298338406009852, #-0.432974713716992, -0.790641633253371, -0.943141335210814, 0.0462132496070069, #0.580670902412863, 1, 0.945560460165707, -0.178658845633278, #0.183925750895849, -1, 0.483519659486572), .Dim = c(11L, 2L)) #b: regularized partial correlation network for noCA with the general distress variable #pdf(file = "noCA_PP_RF_glasso.pdf") noCA_PP_RF_glasso <- qgraph(noCA_PP_network_cor, layout = PP_Layout, graph = "glasso", details = TRUE, nodeNames = colnames(noCA_PP_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(noCA_PP_network_c)) #dev.off() ####D: Checking Robustness and Stability for centrality analyses #a: interrelation accuracy set.seed(3) noCA_PP_accuracy <- bootnet(noCA_PP_network_c, default = "EBICglasso", type = "nonparametric", nCores = 8, nBoots = 2000, statistics = c("edge", "strength","expectedInfluence")) plot(noCA_PP_accuracy$sample, layout = PP_Layout, fade = FALSE) plot(noCA_PP_accuracy, labels = TRUE, order = "sample") plot(noCA_PP_accuracy, labels = FALSE, order = "sample") #b: network centrality coefficient stability set.seed(4) noCA_PP_stability <- bootnet(noCA_PP_network_c, default = "EBICglasso", type = "case", nCores = 8, caseN = 25, nBoots =2000, statistics = c("strength","expectedInfluence")) plot(noCA_PP_stability$sample, layout = PP_Layout, fade = FALSE) plot(noCA_PP_stability, statistics = c("strength", "expectedInfluence")) corStability(noCA_PP_stability, statistics = c("strength", "expectedInfluence"), cor = 0.7) ####E: Centrality analysis #a: node strength and node expected influence noCA_PP_centrality <- centrality_auto(noCA_PP_RF_glasso); noCA_PP_centrality #b: node predictability noCA_PP_network_mgm <- na.omit(noCA_PP_network_c); nrow(noCA_PP_network_mgm) # 442 RFtype<- c(rep("g", 8),rep("c", 2),rep("g", 1)); RFtype RFlev<- c(rep(1,8),rep(2,2),rep(1,1)); RFlev set.seed(1) fit_obj <- mgm(data = noCA_PP_network_mgm, type = RFtype, lev = RFlev, rule.reg = "AND", k = 2, binarySign = TRUE, scale = TRUE) fit_obj$pairwise$wadj fit_obj$pairwise$signs pred_obj <- predict(fit_obj, data = noCA_PP_network_mgm, errorCon = c("RMSE", "R2"), errorCat = c("CC", "nCC")) pred_obj$errors errors <- c(0.19, 0.46, 0.49, 0.43, 0.60, 0.60, 0.35, 0.09, 0.00, 0.02, 0.68); errors pieColor <- c(rep("#90B4D4",8),rep("#EB9446",2),rep("#90B4D4",1)) #pdf(file = "noCA_PP_RF_pred.pdf") noCA_PP_RF_pred <- qgraph(fit_obj$pairwise$wadj, edge.color = fit_obj$pairwise$edgecolor, pie = errors, pieColor = pieColor, layout=PP_Layout, nodeNames = colnames(noCA_PP_network_mgm), label.cex = .9, theme = "classic", curveAll = FALSE, fade = FALSE, legend = TRUE, cut = 0) #dev.off() noCA_PP_am_mgm <- fit_obj$pairwise$wadj; noCA_PP_am_mgm noCA_PP_am_mgm_signs <- fit_obj$pairwise$signs; noCA_PP_am_mgm_signs ``` ##8: Correlations and regularized partial correlations between the RFs and the general distress variable ```{r} ####A: for CA #a: correlation CA_PP_RF_cor <- qgraph(CA_PP_network_cor, layout = PP_Layout, graph = "cor", details = TRUE, nodeNames = colnames(CA_PP_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(CA_PP_network_c)) x <- getWmat(CA_PP_RF_cor); x[,11] #b: regularized partial correlation y <- getWmat(CA_PP_RF_glasso); y[,11] ####B: for noCA #a: correlation noCA_PP_RF_cor <- qgraph(noCA_PP_network_cor, layout = PP_Layout, graph = "cor", details = TRUE, nodeNames = colnames(noCA_PP_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(noCA_PP_network_c)) z <- getWmat(noCA_PP_RF_cor); z[,11] #b: regularized partial correlation xx <- getWmat(noCA_PP_RF_glasso); xx[,11] ``` ##9: Contrasting the CA and the noCA network structure including the general distress variable ```{r} ####A: correlate CA with noCA interrelation matrix including the general distress variable adjmat_CA_PP <- getWmat(CA_PP_RF_glasso); adjmat_CA_PP cor_CA_PP <- adjmat_CA_PP[upper.tri(adjmat_CA_PP,diag=FALSE)]; cor_CA_PP adjmat_noCA_PP <- getWmat(noCA_PP_RF_glasso); adjmat_noCA_PP cor_noCA_PP <- adjmat_noCA_PP[upper.tri(adjmat_noCA_PP,diag=FALSE)]; cor_noCA_PP cor_pears_PP <- cor(cor_CA_PP,cor_noCA_PP, method = "pearson"); cor_pears_PP # 0.9094228 cor_spear_PP <- cor(cor_CA_PP,cor_noCA_PP, method = "spearman"); cor_spear_PP # 0.6775919 ####B: NCT approach for the comparison of the CA vs noCA networks with the general distress variable #a: rename variables so that they are equal for the two groups CA_PP_network_NCT <- CA_PP_network_c CA_PP_network_NCT <- rename(CA_PP_network_NCT, friendsupport = CA_Network_Data.friendsupport_tr) CA_PP_network_NCT <- rename(CA_PP_network_NCT, familysupport = CA_Network_Data.familysupport_tr) CA_PP_network_NCT <- rename(CA_PP_network_NCT, familycohesion = CA_Network_Data.familycohesion_tr) CA_PP_network_NCT <- rename(CA_PP_network_NCT, positiveselfesteem = CA_Network_Data.positiveselfesteem_tr) CA_PP_network_NCT <- rename(CA_PP_network_NCT, negativeselfesteem = CA_Network_Data.negativeselfesteem_tr) CA_PP_network_NCT <- rename(CA_PP_network_NCT, expressivesuppression = CA_Network_Data.expressivesuppressionbin) CA_PP_network_NCT <- rename(CA_PP_network_NCT, distresstolerance = CA_Network_Data.distresstolerance_tr) CA_PP_network_NCT <- rename(CA_PP_network_NCT, aggression = CA_Network_Data.aggressionbin) CA_PP_network_NCT <- rename(CA_PP_network_NCT, brooding = CA_Network_Data.brooding_tr) CA_PP_network_NCT <- rename(CA_PP_network_NCT, reflection = CA_Network_Data.reflection_tr) CA_PP_network_NCT <- rename(CA_PP_network_NCT, generaldistress = CA_Network_Data.generaldistress_tr) head(CA_PP_network_NCT);nrow(CA_PP_network_NCT) noCA_PP_network_NCT <- noCA_PP_network_c noCA_PP_network_NCT <- rename(noCA_PP_network_NCT, friendsupport = noCA_Network_Data.friendsupport_tr) noCA_PP_network_NCT <- rename(noCA_PP_network_NCT, familysupport = noCA_Network_Data.familysupport_tr) noCA_PP_network_NCT <- rename(noCA_PP_network_NCT, familycohesion = noCA_Network_Data.familycohesion_tr) noCA_PP_network_NCT <- rename(noCA_PP_network_NCT, positiveselfesteem = noCA_Network_Data.positiveselfesteem_tr) noCA_PP_network_NCT <- rename(noCA_PP_network_NCT, negativeselfesteem = noCA_Network_Data.negativeselfesteem_tr) noCA_PP_network_NCT <- rename(noCA_PP_network_NCT, expressivesuppression = noCA_Network_Data.expressivesuppressionbin) noCA_PP_network_NCT <- rename(noCA_PP_network_NCT, distresstolerance = noCA_Network_Data.distresstolerance_tr) noCA_PP_network_NCT <- rename(noCA_PP_network_NCT, aggression = noCA_Network_Data.aggressionbin) noCA_PP_network_NCT <- rename(noCA_PP_network_NCT, brooding = noCA_Network_Data.brooding_tr) noCA_PP_network_NCT <- rename(noCA_PP_network_NCT, reflection = noCA_Network_Data.reflection_tr) noCA_PP_network_NCT <- rename(noCA_PP_network_NCT, generaldistress = noCA_Network_Data.generaldistress_tr) head(noCA_PP_network_NCT);nrow(noCA_PP_network_NCT) #b:perform NCT test set.seed(4) CA_vs_noCA_PP <- NCT_cor_auto(CA_PP_network_NCT, noCA_PP_network_NCT, it=5000, binary.data = FALSE, paired = FALSE, weighted=TRUE, test.edges = TRUE, edges = list(c(1,2),c(1,3),c(1,4),c(1,5),c(1,6),c(1,7),c(1,8),c(1,9),c(1,10),c(1,11),c(2,3),c(2,4),c(2,5),c(2,6),c(2,7),c(2,8), c(2,9),c(2,10),c(2,11),c(3,4),c(3,5),c(3,6),c(3,7),c(3,8),c(3,9),c(3,10),c(3,11),c(4,5),c(4,6),c(4,7),c(4,8),c(4,9),c(4,10),c(4,11), c(5,6),c(5,7),c(5,8),c(5,9),c(5,10),c(5,11),c(6,7),c(6,8),c(6,9),c(6,10),c(6,11),c(7,8),c(7,9),c(7,10),c(7,11),c(8,9),c(8,10),c(8,11),c(9,10),c(9,11),c(10,11)), progressbar=TRUE) #f: descriptives summary(CA_vs_noCA_PP) plot(CA_vs_noCA_PP, what="network") plot(CA_vs_noCA_PP, what="strength") plot(CA_vs_noCA_PP, what="edge") ####C:To compare the Global Network Expected Influence, please run the function in the following scipt: "2_b_2_Analysis_Script_NCT_EI-auto_JF-EF-IG-PW-ALvH_2018.07July.03_v1". The script is adapted from: van Borkulo, C. D. et al. Comparing network structures on three aspects: A permutation test (Manuscript submitted for publication, 2017). The adapted script can be requested from the first author jf585@cam.ac.uk. #Please note, what reads now as 'Global strength per group' contains the values of the global network expected influence. #a:perform NCT test set.seed(4) EI_CA_vs_noCA_PP <- NCT_cor_auto_EI(CA_PP_network_NCT, noCA_PP_network_NCT, it=5000, binary.data = FALSE, paired = FALSE, weighted=TRUE, test.edges = TRUE, edges = list(c(1,2),c(1,3),c(1,4),c(1,5),c(1,6),c(1,7),c(1,8),c(1,9),c(1,10),c(1,11),c(2,3),c(2,4),c(2,5),c(2,6),c(2,7),c(2,8), c(2,9),c(2,10),c(2,11),c(3,4),c(3,5),c(3,6),c(3,7),c(3,8),c(3,9),c(3,10),c(3,11),c(4,5),c(4,6),c(4,7),c(4,8),c(4,9),c(4,10),c(4,11), c(5,6),c(5,7),c(5,8),c(5,9),c(5,10),c(5,11),c(6,7),c(6,8),c(6,9),c(6,10),c(6,11),c(7,8),c(7,9),c(7,10),c(7,11),c(8,9),c(8,10),c(8,11),c(9,10),c(9,11),c(10,11)), progressbar=TRUE) summary(EI_CA_vs_noCA_PP) ####D: To reveal the results with the uncorrected p-values, the "method = 'holm'" commands in the NCT script need to be replaced by "method = 'none'". ``` ##10: Contrasting CA and noCA RF networks corrected for the general distress variable ```{r} ####A: comparing RF-RF interrelations of the CA and noCA networks with correction for the general distress variable # CA with correction CA_PP_RF_glasso_mat <- as.data.frame(getWmat(CA_PP_RF_glasso)); CA_PP_RF_glasso_mat CA_PP_RF_glasso_mat2 <- CA_PP_RF_glasso_mat[-c(11),-c(11)]; CA_PP_RF_glasso_mat2 cor_CA_PP_RF_glasso_mat <- CA_PP_RF_glasso_mat2[upper.tri(CA_PP_RF_glasso_mat2, diag=FALSE)]; cor_CA_PP_RF_glasso_mat # noCA with correction noCA_PP_RF_glasso_mat <- as.data.frame(getWmat(noCA_PP_RF_glasso)); noCA_PP_RF_glasso_mat noCA_PP_RF_glasso_mat2 <- noCA_PP_RF_glasso_mat[-c(11),-c(11)]; noCA_PP_RF_glasso_mat2 cor_noCA_PP_RF_glasso_mat <- noCA_PP_RF_glasso_mat2[upper.tri(noCA_PP_RF_glasso_mat2, diag=FALSE)]; cor_noCA_PP_RF_glasso_mat #Comparing RF-RF interrelations of the networks without and with correction for the general distress variable cor(cor_CA_PP_RF_glasso_mat, cor_noCA_PP_RF_glasso_mat, method = "pearson") # 0.8846993 cor(cor_CA_PP_RF_glasso_mat, cor_noCA_PP_RF_glasso_mat, method = "spearman") # 0.5387118 ####B: NCT approach for the comparison of the CA vs noCA networks corrected for the general distress variable #b:perform NCT test #Please run the function in the following scipt: "2_c_1_Analysis_Script_NCT-indep_JF-EF-IG-PW-ALvH_2018.06June.20_v1". The script is adapted from: van Borkulo, C. D. et al. Comparing network structures on three aspects: A permutation test (Manuscript submitted for publication, 2017). The adapted script can be requested from the first author jf585@cam.ac.uk. set.seed(10) CA_vs_noCA_corPP <- NCT_cor_auto_indep(CA_PP_network_NCT, noCA_PP_network_NCT, it=5000, binary.data = FALSE, paired = FALSE, weighted=TRUE, test.edges = TRUE, edges = list(c(1,2),c(1,3),c(1,4),c(1,5),c(1,6),c(1,7),c(1,8),c(1,9),c(1,10),c(2,3),c(2,4),c(2,5),c(2,6),c(2,7),c(2,8), c(2,9),c(2,10),c(3,4),c(3,5),c(3,6),c(3,7),c(3,8),c(3,9),c(3,10),c(4,5),c(4,6),c(4,7),c(4,8),c(4,9),c(4,10), c(5,6),c(5,7),c(5,8),c(5,9),c(5,10),c(6,7),c(6,8),c(6,9),c(6,10),c(7,8),c(7,9),c(7,10),c(8,9),c(8,10),c(9,10)), progressbar=TRUE) summary(CA_vs_noCA_corPP) plot(CA_vs_noCA_corPP, what="network") plot(CA_vs_noCA_corPP, what="strength") plot(CA_vs_noCA_corPP, what="edge") ####C:To compare the Global Network Expected Influence, please run the function in the following scipt: "2_c_2_Analysis_Script_NCT_EI-indep_JF-EF-IG-PW-ALvH_2018.07July.03_v1". The script is adapted from: van Borkulo, C. D. et al. Comparing network structures on three aspects: A permutation test (Manuscript submitted for publication, 2017). The adapted script can be requested from the first author jf585@cam.ac.uk. #Please note, what reads now as 'Global strength per group' contains the values of the global network expected influence. #a:perform NCT test set.seed(10) EI_CA_vs_noCA_corPP <- NCT_cor_auto_indep_EI(CA_PP_network_NCT, noCA_PP_network_NCT, it=5000, binary.data = FALSE, paired = FALSE, weighted=TRUE, test.edges = TRUE, edges = list(c(1,2),c(1,3),c(1,4),c(1,5),c(1,6),c(1,7),c(1,8),c(1,9),c(1,10),c(2,3),c(2,4),c(2,5),c(2,6),c(2,7),c(2,8), c(2,9),c(2,10),c(3,4),c(3,5),c(3,6),c(3,7),c(3,8),c(3,9),c(3,10),c(4,5),c(4,6),c(4,7),c(4,8),c(4,9),c(4,10), c(5,6),c(5,7),c(5,8),c(5,9),c(5,10),c(6,7),c(6,8),c(6,9),c(6,10),c(7,8),c(7,9),c(7,10),c(8,9),c(8,10),c(9,10)), progressbar=TRUE) summary(EI_CA_vs_noCA_corPP) ####D: To reveal the results with the uncorrected p-values, the "method = 'holm'" commands in the NCT script need to be replaced by "method = 'none'". ``` ##11: Shortest pathways between RFs and the general distress variable for CA and noCA groups ```{r} ####A: for CA #a: renaming the variables in the data frame CA_PP_network_c <- data.frame(CA_Network_Data$friendsupport_tr, CA_Network_Data$familysupport_tr, CA_Network_Data$familycohesion_tr, CA_Network_Data$positiveselfesteem_tr, CA_Network_Data$negativeselfesteem_tr, CA_Network_Data$brooding_tr, CA_Network_Data$reflection_tr, CA_Network_Data$distresstolerance_tr, CA_Network_Data$aggressionbin, CA_Network_Data$expressivesuppressionbin, CA_Network_Data$generaldistress_tr) CA_PP_network_c <- rename(CA_PP_network_c, frn = CA_Network_Data.friendsupport_tr) CA_PP_network_c <- rename(CA_PP_network_c, fms = CA_Network_Data.familysupport_tr) CA_PP_network_c <- rename(CA_PP_network_c, fmc = CA_Network_Data.familycohesion_tr) CA_PP_network_c <- rename(CA_PP_network_c, pst = CA_Network_Data.positiveselfesteem_tr) CA_PP_network_c <- rename(CA_PP_network_c, ngt = CA_Network_Data.negativeselfesteem_tr) CA_PP_network_c <- rename(CA_PP_network_c, exp = CA_Network_Data.expressivesuppressionbin) CA_PP_network_c <- rename(CA_PP_network_c, dst = CA_Network_Data.distresstolerance_tr) CA_PP_network_c <- rename(CA_PP_network_c, agg = CA_Network_Data.aggressionbin) CA_PP_network_c <- rename(CA_PP_network_c, brd = CA_Network_Data.brooding_tr) CA_PP_network_c <- rename(CA_PP_network_c, rfl = CA_Network_Data.reflection_tr) CA_PP_network_c <- rename(CA_PP_network_c, GD = CA_Network_Data.generaldistress_tr) head(CA_PP_network_c);nrow(CA_PP_network_c) #b: making the correlation matrix CA_PP_network_cor <- cor_auto(CA_PP_network_c); CA_PP_network_cor max(CA_PP_network_cor[upper.tri(CA_PP_network_cor,diag=FALSE)]) # 0.686106 #c: making the network CA_PP_RF_glasso <- qgraph(CA_PP_network_cor, layout = PP_Layout, graph = "glasso", details = FALSE, nodeNames = colnames(CA_PP_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(CA_PP_network_c)) ####B: for no-CA #a: renaming the variables in the data frame noCA_PP_network_c <- data.frame(noCA_Network_Data$friendsupport_tr, noCA_Network_Data$familysupport_tr, noCA_Network_Data$familycohesion_tr, noCA_Network_Data$positiveselfesteem_tr, noCA_Network_Data$negativeselfesteem_tr, noCA_Network_Data$brooding_tr, noCA_Network_Data$reflection_tr, noCA_Network_Data$distresstolerance_tr, noCA_Network_Data$aggressionbin, noCA_Network_Data$expressivesuppressionbin, noCA_Network_Data$generaldistress_tr) head(noCA_PP_network_c) noCA_PP_network_c <- rename(noCA_PP_network_c, frn = noCA_Network_Data.friendsupport_tr) noCA_PP_network_c <- rename(noCA_PP_network_c, fms = noCA_Network_Data.familysupport_tr) noCA_PP_network_c <- rename(noCA_PP_network_c, fmc = noCA_Network_Data.familycohesion_tr) noCA_PP_network_c <- rename(noCA_PP_network_c, pst = noCA_Network_Data.positiveselfesteem_tr) noCA_PP_network_c <- rename(noCA_PP_network_c, ngt = noCA_Network_Data.negativeselfesteem_tr) noCA_PP_network_c <- rename(noCA_PP_network_c, exp = noCA_Network_Data.expressivesuppressionbin) noCA_PP_network_c <- rename(noCA_PP_network_c, dst = noCA_Network_Data.distresstolerance_tr) noCA_PP_network_c <- rename(noCA_PP_network_c, agg = noCA_Network_Data.aggressionbin) noCA_PP_network_c <- rename(noCA_PP_network_c, brd = noCA_Network_Data.brooding_tr) noCA_PP_network_c <- rename(noCA_PP_network_c, rfl = noCA_Network_Data.reflection_tr) noCA_PP_network_c <- rename(noCA_PP_network_c, GD = noCA_Network_Data.generaldistress_tr) head(noCA_PP_network_c);nrow(noCA_PP_network_c) #b: making the correlation matrix noCA_PP_network_cor <- cor_auto(noCA_PP_network_c); noCA_PP_network_cor max(noCA_PP_network_cor[upper.tri(noCA_PP_network_cor,diag=FALSE)]) # 0.6609587 #c:making the network noCA_PP_RF_glasso <- qgraph(noCA_PP_network_cor, layout = PP_Layout, graph = "glasso", details = FALSE, nodeNames = colnames(noCA_PP_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(noCA_PP_network_c)) ####C: Differing shortest pathways for the CA and the no-CA group #pdf(file = "Differing_SPL.pdf", width = 8.27, height = 11.69) #layout(matrix(c(1,2,3,4,5,6,7,8,9,10),5,2, byrow = TRUE), # widths=c(lcm(5.5),lcm(5.5)), # heights=c(lcm(5.5),lcm(5.5),lcm(5.5),lcm(5.5),lcm(5.5))) #friendship support #pdf(file = "CA_FRN_SPL.pdf") pathways(CA_PP_RF_glasso,from=1,to=11) #dev.off() #pdf(file = "noCA_FRN_SPL.pdf") pathways(noCA_PP_RF_glasso,from=1,to=11) #dev.off() #family support #pdf(file = "CA_FMS_SPL.pdf") pathways(CA_PP_RF_glasso,from=2,to=11) #dev.off() #pdf(file = "noCA_FMS_SPL.pdf") pathways(noCA_PP_RF_glasso,from=2,to=11) #dev.off() #family cohesion #pdf(file = "CA_FMC_SPL.pdf") pathways(CA_PP_RF_glasso,from=3,to=11) #dev.off() #pdf(file = "noCA_FMC_SPL.pdf") pathways(noCA_PP_RF_glasso,from=3,to=11) #dev.off() #distress tolerance #pdf(file = "CA_DST_SPL.pdf") pathways(CA_PP_RF_glasso,from=8,to=11) #dev.off() #pdf(file = "noCA_DST_SPL.pdf") pathways(noCA_PP_RF_glasso,from=8,to=11) #dev.off() #expressive suppression #pdf(file = "CA_EXP_SPL.pdf") pathways(CA_PP_RF_glasso,from=10,to=11) #dev.off() #pdf(file = "noCA_EXP_SPL.pdf") pathways(noCA_PP_RF_glasso,from=10,to=11) #dev.off() ####D: Equivalent shortest pathways for the CA and the no-CA group #positive self-esteem pathways(CA_PP_RF_glasso,from=4,to=11) pathways(noCA_PP_RF_glasso,from=4,to=11) #negative self-esteem pathways(CA_PP_RF_glasso,from=5,to=11) pathways(noCA_PP_RF_glasso,from=5,to=11) #brooding pathways(CA_PP_RF_glasso,from=6,to=11) pathways(noCA_PP_RF_glasso,from=6,to=11) #reflection pathways(CA_PP_RF_glasso,from=7,to=11) pathways(noCA_PP_RF_glasso,from=7,to=11) #aggression pathways(CA_PP_RF_glasso,from=9,to=11) pathways(noCA_PP_RF_glasso,from=9,to=11) #pdf(file = "Equal_SPL.pdf", width = 3, height = 1.5) #layout(matrix(c(1,2),1,2, byrow = TRUE), # widths=c(lcm(3.5),lcm(3.5)), # heights=c(lcm(3.5))) #jpeg(file = "Equal_SPL.jpg", width = 7680, height = 3840) layout(matrix(c(1,2),1,2, byrow = TRUE)) pathways(CA_PP_RF_glasso,from=c(4,5,6,7,9),to=c(11)) pathways(noCA_PP_RF_glasso,from=c(4,5,6,7,9),to=c(11)) #dev.off() ``` ##12: correlation of the 'RF-general distress' interrelations between the two groups ```{r} #a: correlation for correlational interrelations cor_pears_RFGD <- cor(x[,11],z[,11], method = "pearson"); cor_pears_RFGD # 0.96063 cor_spear_RFGD <- cor(x[,11],z[,11], method = "spearman"); cor_spear_RFGD # 0.9272727 #b: correlation for regularized partial correlational interrelation cor_pears_RFGD <- cor(y[,11],xx[,11], method = "pearson"); cor_pears_RFGD # 0.9241573 cor_spear_RFGD <- cor(y[,11],xx[,11], method = "spearman"); cor_spear_RFGD # 0.880771 ``` ##13: Methods: Sample descriptives ```{r} #a: CA vs no-CA table(Network_Data$CA_dic) #b: gender table(Network_Data$sex) chi_sex <- table(Network_Data$sex, Network_Data$CA_dic); chi_sex chisq.test(chi_sex) #c: age age <- na.omit(data.frame(Network_Data$age_srqt1, Network_Data$CA_dic)) stat.desc(age$Network_Data.age_srqt1, norm = TRUE) stat.desc(CA_Network_Data$age_srqt1, norm = TRUE) stat.desc(noCA_Network_Data$age_srqt1, norm = TRUE) t.test(Network_Data$age_srqt1 ~ Network_Data$CA_dic, var.equal = FALSE) #d: Socio-economit status (SES) chi_SES <- table(Network_Data$acorn_ses, Network_Data$CA_dic); chi_SES #As SES was split in 5 ordered categories, we applied the Asymptotic Cochran-Armitage test chi_SES2 <- matrix( c(30, 77, 11, 36, 105, 170, 41, 37, 314, 318), byrow = TRUE, ncol = 2, dimnames = list( "SES" = c("hardpressed", "moderatemeans", "comfoff", "urbanprosperity", "wealthyachievers"), "CA" = c("0", "1") ) ) chi_SES2 <- as.table(chi_SES2);chi_SES2 chi_SES2 <- chisq_test(chi_SES2, scores = list("SES" = c(1,2,3,4,5)));chi_SES2 #e: prior psychiatric history table(Network_Data$banydiag) Network_Data$banydiag_dic <- if_else(Network_Data$banydiag != "none", "any", "none") chi_diag <- table(Network_Data$banydiag_dic, Network_Data$CA_dic); chi_diag chisq.test(chi_diag) #f: mood stat.desc(CA_Network_Data$t1_mfqtotal_4cat, norm = TRUE) stat.desc(noCA_Network_Data$t1_mfqtotal_4cat, norm = TRUE) t.test(Network_Data$t1_mfqtotal_4cat ~ Network_Data$CA_dic, var.equal = FALSE) #g: anxiety stat.desc(CA_Network_Data$t1_rcmastotal_4cat, norm = TRUE) stat.desc(noCA_Network_Data$t1_rcmastotal_4cat, norm = TRUE) t.test(Network_Data$t1_rcmastotal_4cat ~ Network_Data$CA_dic, var.equal = FALSE) ``` ##14: Supplement I: Methodological Rationale ```{r} #This supplement does not contain analyses ``` ##15: Supplement II: Variable Preparation: Means and standard deviations or inter quartile ranges and frequencies for untransformed RFs and the general distress variable ```{r} ####A: Numbers for Table 1 #a: friend support stat.desc(CA_Network_Data$friendsupport, norm = TRUE) stat.desc(noCA_Network_Data$friendsupport, norm = TRUE) #b: family support stat.desc(CA_Network_Data$familysupport, norm = TRUE) stat.desc(noCA_Network_Data$familysupport, norm = TRUE) #c: family cohesion stat.desc(CA_Network_Data$familycohesion, norm = TRUE) stat.desc(noCA_Network_Data$familycohesion, norm = TRUE) #d: negative self-esteem stat.desc(CA_Network_Data$negativeselfesteem, norm = TRUE) stat.desc(noCA_Network_Data$negativeselfesteem, norm = TRUE) #e: positive self-esteem stat.desc(CA_Network_Data$positiveselfesteem, norm = TRUE) stat.desc(noCA_Network_Data$positiveselfesteem, norm = TRUE) #f: brooding stat.desc(CA_Network_Data$brooding, norm = TRUE) stat.desc(noCA_Network_Data$brooding, norm = TRUE) #g: reflection stat.desc(CA_Network_Data$reflection, norm = TRUE) stat.desc(noCA_Network_Data$reflection, norm = TRUE) #h: distress tolerance stat.desc(CA_Network_Data$distresstolerance, norm = TRUE) stat.desc(noCA_Network_Data$distresstolerance, norm = TRUE) #i: aggression stat.desc(CA_Network_Data$aggression, norm = TRUE) stat.desc(noCA_Network_Data$aggression, norm = TRUE) #j: expressive suppression stat.desc(CA_Network_Data$expressivesuppression, norm = TRUE) IQR(na.omit(CA_Network_Data$expressivesuppression)) stat.desc(noCA_Network_Data$expressivesuppression, norm = TRUE) IQR(na.omit(noCA_Network_Data$expressivesuppression)) #k: general distress stat.desc(CA_Network_Data$b_f_gen, norm = TRUE) stat.desc(noCA_Network_Data$b_f_gen, norm = TRUE) ####B: Supplement Figure 1 #a: make the data frame names=c(rep("frn", 1238) , rep("fms", 1238) , rep("fmc", 1238), rep("ngt", 1238) , rep("pst", 1238) , rep("rfl", 1238), rep("brd", 1238), rep("dst", 1238), rep("agg", 1238), rep("GD", 1238)) value=c(Network_Data$friendsupport, Network_Data$familysupport, Network_Data$familycohesion, Network_Data$negativeselfesteem, Network_Data$positiveselfesteem, Network_Data$reflection, Network_Data$brooding, Network_Data$distresstolerance, Network_Data$aggression, Network_Data$b_f_gen) group=as.factor(c(Network_Data$CA_dic, Network_Data$CA_dic, Network_Data$CA_dic, Network_Data$CA_dic, Network_Data$CA_dic, Network_Data$CA_dic, Network_Data$CA_dic, Network_Data$CA_dic, Network_Data$CA_dic, Network_Data$CA_dic)) data=na.omit(data.frame(names,value,group)) data$group <- as.factor(data$group) #b: make the graph #pdf(file = "Distribution.pdf") plot <- ggplot(data, aes(x = group, y = value, group = group)) + geom_point(aes(colour = group), alpha=0.3, position = "jitter", size = 0.5) + scale_colour_manual(values=c("#009900","magenta4"), breaks = c("0", "1"), labels = c("No-CA", "CA"))+ geom_boxplot(alpha=0, outlier.colour="black", outlier.shape=16, outlier.size=0.5, notch=FALSE) + facet_wrap( ~ names, ncol = 5) + theme_minimal() + theme(axis.title.x = element_blank(), axis.title.y = element_blank(), legend.position = "right", panel.grid = element_blank(), legend.title = element_blank(), axis.text.x = element_blank())+ guides(colour = guide_legend(override.aes = list(size=2.5))) plot #dev.off() ``` ##16: Supplement III: Association and concentration networks for CA and noCA adolescents ```{r} ####A: Models for CA #a: association #pdf(file = "CA_RF_cor.pdf") CA_RF_cor <- qgraph(CA_network_cor, layout = Layout, graph = "cor", details = TRUE, nodeNames = colnames(CA_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(CA_network_c)) #dev.off() #pdf(file = "CA_RF_cor_fade.pdf") CA_RF_cor_fade <- qgraph(CA_network_cor, layout = Layout, graph = "cor", details = TRUE, nodeNames = colnames(CA_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = TRUE, sampleSize = nrow(CA_network_c)) #dev.off() #b: concentration #pdf(file = "CA_RF_pcor.pdf") CA_RF_pcor <- qgraph(CA_network_cor, layout = Layout, graph = "concentration", details = TRUE, nodeNames = colnames(CA_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(CA_network_c)) #dev.off() #pdf(file = "CA_RF_pcor_fade.pdf") CA_RF_pcor_fade <- qgraph(CA_network_cor, layout = Layout, graph = "concentration", details = TRUE, nodeNames = colnames(CA_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = TRUE, sampleSize = nrow(CA_network_c)) #dev.off() ####B: Models for noCA #a: association #pdf(file = "noCA_RF_cor.pdf") noCA_RF_cor <- qgraph(noCA_network_cor, layout = Layout, graph = "cor", details = TRUE, nodeNames = colnames(noCA_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(noCA_network_c)) #dev.off() #pdf(file = "noCA_RF_cor_fade.pdf") noCA_RF_cor_fade <- qgraph(noCA_network_cor, layout = Layout, graph = "cor", details = TRUE, nodeNames = colnames(noCA_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = TRUE, sampleSize = nrow(noCA_network_c)) #dev.off() #b: concentration #pdf(file = "noCA_RF_pcor.pdf") noCA_RF_pcor <- qgraph(noCA_network_cor, layout = Layout, graph = "concentration", details = TRUE, nodeNames = colnames(noCA_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(noCA_network_c)) #dev.off() #pdf(file = "noCA_RF_pcor_fade.pdf") noCA_RF_pcor_fade <- qgraph(noCA_network_cor, layout = Layout, graph = "concentration", details = TRUE, nodeNames = colnames(noCA_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade =TRUE, sampleSize = nrow(noCA_network_c)) #dev.off() ``` ##17: Supplement IV: Adjacency matrices of the main models ```{r} ####A: Models for CA #a: RF-RF getWmat(CA_RF_glasso) #b: RF-RF + general dsitress getWmat(CA_PP_RF_glasso) ####B: Models for noCA #a: RF-RF getWmat(noCA_RF_glasso) #b: RF-RF + general dsitress getWmat(noCA_PP_RF_glasso) ``` ##18: Supplement V: Interconnectedness of RFs ```{r} #a: Models for CA: node strength and node expected influence CA_centrality <- centrality_auto(CA_RF_glasso); CA_centrality$node.centrality #b: Models for noCA: node strength and node expected influence noCA_centrality <- centrality_auto(noCA_RF_glasso); noCA_centrality$node.centrality ``` ##19: Supplement VI: Accuracy and Stability of the CA network ```{r} #a: estimate main model network with bootnet CA_stability_RF_glasso <- estimateNetwork(CA_network_c, default = "EBICglasso") #b: check wehther model has the same characteristics as main model network with qgraph plot(CA_stability_RF_glasso, layout = Layout, labels = TRUE, fade = FALSE) CA_cent <- centrality_auto(CA_stability_RF_glasso); CA_cent$node.centrality #c: RF interrelation (edge weight) accuracy set.seed(3) CA_accuracy <- bootnet(CA_network_c, default = "EBICglasso", type = "nonparametric", nCores = 8, nBoots = 2000, statistics = c("edge", "strength","expectedInfluence")) plot(CA_accuracy$sample, layout = Layout, fade = FALSE) #pdf(file = "bootnet_accuracy_CA.pdf") plot(CA_accuracy, labels = TRUE, order = "sample") #dev.off() #pdf(file = "bootnet_accuracy_CA_withoutLabels.pdf") plot(CA_accuracy, labels = FALSE, order = "sample") #dev.off() #d: network characteristic coefficient stability set.seed(4) CA_stability <- bootnet(CA_network_c, default = "EBICglasso", type = "case", nCores = 8, caseN = 25, nBoots = 2000, statistics = c("strength","expectedInfluence")) plot(CA_stability$sample, layout = Layout, fade = FALSE) CA_stability$sample$graph #pdf(file = "bootnet_stability_CA.pdf") plot(CA_stability, statistics = c("strength", "expectedInfluence")) #dev.off() #e: stability (CS) coefficient corStability(CA_stability, statistics = c("strength", "expectedInfluence"), cor = 0.7) #f: difference tests of RF interrelations (edge weights) for all pairs of RFs #pdf(file = "bootnet_edge_difference_CA.pdf") CA_plotacc <- plot(CA_accuracy, "edge", plot = "difference", onlyNonZero = TRUE, order = "sample"); CA_plotacc #dev.off() #pdf(file = "bootnet_edge_difference_CA_withoutLabels.pdf") CA_plotacc <- plot(CA_accuracy, "edge", plot = "difference", onlyNonZero = TRUE, order = "sample", labels = FALSE); CA_plotacc #dev.off() #g: difference test for node strength #pdf(file = "bootnet_strength_difference_CA.pdf") CA_plotstrength <- plot(CA_accuracy, "strength"); CA_plotstrength #dev.off() #pdf(file = "bootnet_expectedInfluence_difference_CA.pdf") CA_plotexpectedInfluence <- plot(CA_accuracy, "expectedInfluence"); CA_plotexpectedInfluence #dev.off() ``` ##20: Supplement VI: Accuracy and Stability of the noCA network ```{r} #a: estimate main model network with bootnet noCA_stability_RF_glasso <- estimateNetwork(noCA_network_c, default = "EBICglasso") #b: check wehther model has the same characteristics as main model network with qgraph plot(noCA_stability_RF_glasso, layout = Layout, labels = TRUE, fade = FALSE) noCA_cent <- centrality_auto(noCA_stability_RF_glasso); noCA_cent$node.centrality #c: RF interrelation (edge weight) accuracy set.seed(3) noCA_accuracy <- bootnet(noCA_network_c, default = "EBICglasso", type = "nonparametric", nCores = 8, nBoots = 2000, statistics = c("edge", "strength","expectedInfluence")) noCA_accuracy$sample$graph plot(noCA_accuracy$sample, layout = Layout, fade = FALSE) #pdf(file = "bootnet_accuracy_noCA.pdf") plot(noCA_accuracy, labels = TRUE, order = "sample") #dev.off() #pdf(file = "bootnet_accuracy_noCA_withoutLabels.pdf") plot(noCA_accuracy, labels = FALSE, order = "sample") #dev.off() #d: network characteristic coefficient stability set.seed(4) noCA_stability <- bootnet(noCA_network_c, default = "EBICglasso", type = "case", nCores = 8, caseN = 25, nBoots =2000, statistics = c("strength","expectedInfluence")) plot(noCA_stability$sample, layout = Layout, fade = FALSE) #pdf(file = "bootnet_stability_noCA.pdf") plot(noCA_stability, statistics = c("strength", "expectedInfluence")) #dev.off() #e: stability (CS) coefficient corStability(noCA_stability, statistics = c("strength", "expectedInfluence"), cor = 0.7) #f: difference tests of RF interrelations (edge weights) for all pairs of RFs #pdf(file = "bootnet_edge_difference_noCA.pdf") noCA_plotacc <- plot(noCA_accuracy, "edge", plot = "difference", onlyNonZero = TRUE, order = "sample"); noCA_plotacc #dev.off() #pdf(file = "bootnet_edge_difference_noCA_withoutLabels.pdf") noCA_plotacc <- plot(noCA_accuracy, "edge", plot = "difference", onlyNonZero = TRUE, order = "sample", labels = FALSE); noCA_plotacc #dev.off() #g: difference test for node strength #pdf(file = "bootnet_strength_difference_noCA.pdf") noCA_plotstrength <- plot(noCA_accuracy, "strength"); noCA_plotstrength #dev.off() #pdf(file = "bootnet_expectedInfluence_difference_noCA.pdf") noCA_plotexpectedInfluence <- plot(noCA_accuracy, "expectedInfluence"); noCA_plotexpectedInfluence #dev.off() ``` ##21: Supplement VII: Sensitivity analyses for the CA network ```{r} ####A: Making the network data frames for sensitivity analyses for the CA group #b: main model without transformation CA_network_withouttr <- data.frame(CA_Network_Data$friendsupport, CA_Network_Data$familysupport, CA_Network_Data$familycohesion, CA_Network_Data$positiveselfesteem, CA_Network_Data$negativeselfesteem, CA_Network_Data$brooding, CA_Network_Data$reflection, CA_Network_Data$distresstolerance, CA_Network_Data$aggressionbin, CA_Network_Data$expressivesuppressionbin) #c: mean score model CA_network_mean <- data.frame(CA_Network_Data$friendsupport_M_tr, CA_Network_Data$familysupport_M_tr, CA_Network_Data$familycohesion_M_tr, CA_Network_Data$positiveselfesteem_M_tr, CA_Network_Data$negativeselfesteem_M_tr, CA_Network_Data$brooding_M_tr, CA_Network_Data$reflection_M_tr, CA_Network_Data$distresstolerance_M_tr, CA_Network_Data$aggressionbin, CA_Network_Data$expressivesuppressionbin) #d: mean score model without transformation CA_network_mean_withouttr <- data.frame(CA_Network_Data$friendsupport_M, CA_Network_Data$familysupport_M, CA_Network_Data$familycohesion_M, CA_Network_Data$positiveselfesteem_M, CA_Network_Data$negativeselfesteem_M, CA_Network_Data$brooding_M, CA_Network_Data$reflection_M, CA_Network_Data$distresstolerance_M, CA_Network_Data$aggressionbin, CA_Network_Data$expressivesuppressionbin) ####B: Making the correlation matrices #b: cor matrix for main model without transformation CA_network_withouttr_cor <- cor_auto(CA_network_withouttr); CA_network_withouttr_cor #Variables detected as ordinal: CA_Network_Data.CA_aggressionbin; CA_Network_Data.CA_expressivesuppressionbin max(CA_network_withouttr_cor[upper.tri(CA_network_withouttr_cor,diag=FALSE)]) # 0.6874924 #c: cor matrix for mean score model CA_network_mean_cor <- cor_auto(CA_network_mean); CA_network_mean_cor #Variables detected as ordinal: CA_Network_Data.CA_aggressionbin_M; CA_Network_Data.CA_expressivesuppressionbin_M max(CA_network_mean_cor[upper.tri(CA_network_mean_cor,diag=FALSE)]) # 0.7079723 #d: cor matrix for mean score model without transformation CA_network_mean_withouttr_cor <- cor_auto(CA_network_mean_withouttr); CA_network_mean_withouttr_cor #Variables detected as ordinal: CA_Network_Data.CA_aggressionbin_M; CA_Network_Data.CA_expressivesuppressionbin_M max(CA_network_mean_withouttr_cor[upper.tri(CA_network_mean_withouttr_cor,diag=FALSE)]) # 0.7062286 ####C: Establish adjacency matrices for sensitivity analyses for the CA group #a: main model am_ggmfull <- getWmat(CA_RF_glasso); am_ggmfull cor1 <- am_ggmfull[upper.tri(am_ggmfull,diag=FALSE)]; cor1 #b: main model without transformation CA_RF_glasso_withouttr <- qgraph(CA_network_withouttr_cor, layout = Layout, graph = "glasso", details = TRUE, nodeNames = colnames(CA_network_withouttr_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(CA_network_withouttr)) am_ggmfull_withouttr <- getWmat(CA_RF_glasso_withouttr); am_ggmfull_withouttr cor2 <- am_ggmfull_withouttr[upper.tri(am_ggmfull_withouttr,diag=FALSE)]; cor2 #c: mean score model CA_RF_glasso_mean <- qgraph(CA_network_mean_cor, layout = Layout, graph = "glasso", details = TRUE, nodeNames = colnames(CA_network_mean_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(CA_network_mean)) am_ggmmeanfull <- getWmat(CA_RF_glasso_mean); am_ggmmeanfull cor3 <- am_ggmmeanfull[upper.tri(am_ggmmeanfull,diag=FALSE)]; cor3 #d: mean score model without transformation CA_RF_glasso_mean_withouttr <- qgraph(CA_network_mean_withouttr_cor, layout = Layout, graph = "glasso", details = TRUE, nodeNames = colnames(CA_network_mean_withouttr_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(CA_network_mean_withouttr)) am_ggmmeanfull_withouttr <- getWmat(CA_RF_glasso_mean_withouttr); am_ggmmeanfull_withouttr cor4 <- am_ggmmeanfull_withouttr[upper.tri(am_ggmmeanfull_withouttr,diag=FALSE)]; cor4 #e: MGM model am_mgmdat_CA <- as.data.frame(CA_am_mgm) CA_am_mgm_signs write.csv(am_mgmdat_CA, file = "am_mgmdat_CA_16.08.2018.csv")#then enter signs into excel file am_mgm_CA <- read.csv("File_directory\\am_mgmdat_CA_16.08.2018.csv"); am_mgm_CA am_mgm_CA <- as.matrix(am_mgm_CA); am_mgm_CA cor5 <- am_mgm_CA[upper.tri(am_mgm_CA, diag=FALSE)]; cor5 #f: main model for complete sample CA_network_cor_mgm <- cor_auto(CA_network_mgm); CA_network_cor_mgm #Variables detected as ordinal: CA_Network_Data.CA_aggressionbin; CA_Network_Data.CA_expressivesuppressionbin max(CA_network_cor_mgm[upper.tri(CA_network_cor_mgm,diag=FALSE)]) # 0.6861858 CA_RF_glasso_mgm <- qgraph(CA_network_cor_mgm, layout = Layout, graph = "glasso", details = TRUE, nodeNames = colnames(CA_network_cor_mgm), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(CA_network_mgm)) am_ggm <- getWmat(CA_RF_glasso_mgm); am_ggm cor6 <- am_ggm[upper.tri(am_ggm,diag=FALSE)]; cor6 ####D: Sensitivity analyses for the CA group #a: main model (qgraph full sample) vs main model without transformation (qgraph full sample) cor(cor1, cor2, method = "pearson") # 0.9992543 cor(cor1, cor2, method = "spearman") # 0.9893359 #b: main model (qgraph full sample) vs mean score model (qgraph full sample) cor(cor1, cor3, method = "pearson") # 0.9930884 cor(cor1, cor3, method = "spearman") # 0.9692338 #c: main model (qgraph full sample) vs mean score model without transformation (qgraph full sample) cor(cor1, cor4, method = "pearson") # 0.9918727 cor(cor1, cor4, method = "spearman") # 0.9680878 #d: main model (qgraph full sample) vs mgm model (mgm subsample) cor(cor1, cor5, method = "pearson") # 0.9350924 cor(cor1, cor5, method = "spearman") # 0.8686555 #e: main model (qgraph full sample) vs complete sample model (qgraph subsample) cor(cor1, cor6, method = "pearson") # 0.9930838 cor(cor1, cor6, method = "spearman") # 0.9496684 #f: mgm model (mgm subsample) vs complete sample model (qgraph subsample) cor(cor5, cor6, method = "pearson") # 0.9160366 cor(cor5, cor6, method = "spearman") # 0.8419105 ``` ##22: Supplement VII: Sensitivity analyses for the noCA network ```{r} ####A: Making the network data frames for sensitivity analyses for the noCA group #b: main model without transformation noCA_network_withouttr <- data.frame(noCA_Network_Data$friendsupport, noCA_Network_Data$familysupport, noCA_Network_Data$familycohesion, noCA_Network_Data$positiveselfesteem, noCA_Network_Data$negativeselfesteem, noCA_Network_Data$brooding, noCA_Network_Data$reflection, noCA_Network_Data$distresstolerance, noCA_Network_Data$aggressionbin, noCA_Network_Data$expressivesuppressionbin) #c: mean score model noCA_network_mean <- data.frame(noCA_Network_Data$friendsupport_M_tr, noCA_Network_Data$familysupport_M_tr, noCA_Network_Data$familycohesion_M_tr, noCA_Network_Data$positiveselfesteem_M_tr, noCA_Network_Data$negativeselfesteem_M_tr, noCA_Network_Data$brooding_M_tr, noCA_Network_Data$reflection_M_tr, noCA_Network_Data$distresstolerance_M_tr, noCA_Network_Data$aggressionbin, noCA_Network_Data$expressivesuppressionbin) #d: mean score model without transformation noCA_network_mean_withouttr <- data.frame(noCA_Network_Data$friendsupport_M, noCA_Network_Data$familysupport_M, noCA_Network_Data$familycohesion_M, noCA_Network_Data$positiveselfesteem_M, noCA_Network_Data$negativeselfesteem_M, noCA_Network_Data$brooding_M, noCA_Network_Data$reflection_M, noCA_Network_Data$distresstolerance_M, noCA_Network_Data$aggressionbin, noCA_Network_Data$expressivesuppressionbin) ####B: Making the correlation matrices #b: cor matrix for main model without transformation noCA_network_withouttr_cor <- cor_auto(noCA_network_withouttr); noCA_network_withouttr_cor #Variables detected as ordinal: noCA_Network_Data.noCA_aggressionbin; noCA_Network_Data.noCA_expressivesuppressionbin max(noCA_network_withouttr_cor[upper.tri(noCA_network_withouttr_cor,diag=FALSE)]) # 0.6631588 #c: cor matrix for mean score model noCA_network_mean_cor <- cor_auto(noCA_network_mean); noCA_network_mean_cor #Variables detected as ordinal: noCA_Network_Data.noCA_aggressionbin_M; noCA_Network_Data.noCA_expressivesuppressionbin_M max(noCA_network_mean_cor[upper.tri(noCA_network_mean_cor,diag=FALSE)]) # 0.6875651 #d: cor matrix for mean score model without transformation noCA_network_mean_withouttr_cor <- cor_auto(noCA_network_mean_withouttr); noCA_network_mean_withouttr_cor #Variables detected as ordinal: noCA_Network_Data.noCA_aggressionbin_M; noCA_Network_Data.noCA_expressivesuppressionbin_M max(noCA_network_mean_withouttr_cor[upper.tri(noCA_network_mean_withouttr_cor,diag=FALSE)]) # 0.6831912 ####C: Establish adjacency matrices for sensitivity analyses for the noCA group #a: main model am_ggmfull <- getWmat(noCA_RF_glasso); am_ggmfull cor1 <- am_ggmfull[upper.tri(am_ggmfull,diag=FALSE)]; cor1 #b: main model without transformation noCA_RF_glasso_withouttr <- qgraph(noCA_network_withouttr_cor, layout = Layout, graph = "glasso", details = TRUE, nodeNames = colnames(noCA_network_withouttr_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(noCA_network_withouttr)) am_ggmfull_withouttr <- getWmat(noCA_RF_glasso_withouttr); am_ggmfull_withouttr cor2 <- am_ggmfull_withouttr[upper.tri(am_ggmfull_withouttr,diag=FALSE)]; cor2 #c: mean score model noCA_RF_glasso_mean <- qgraph(noCA_network_mean_cor, layout = Layout, graph = "glasso", details = TRUE, nodeNames = colnames(noCA_network_mean_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(noCA_network_mean)) am_ggmmeanfull <- getWmat(noCA_RF_glasso_mean); am_ggmmeanfull cor3 <- am_ggmmeanfull[upper.tri(am_ggmmeanfull,diag=FALSE)]; cor3 #d: mean score model without transformation noCA_RF_glasso_mean_withouttr <- qgraph(noCA_network_mean_withouttr_cor, layout = Layout, graph = "glasso", details = TRUE, nodeNames = colnames(noCA_network_mean_withouttr_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(noCA_network_mean_withouttr)) am_ggmmeanfull_withouttr <- getWmat(noCA_RF_glasso_mean_withouttr); am_ggmmeanfull_withouttr cor4 <- am_ggmmeanfull_withouttr[upper.tri(am_ggmmeanfull_withouttr,diag=FALSE)]; cor4 #e: MGM model am_mgmdat_noCA <- as.data.frame(noCA_am_mgm) noCA_am_mgm_signs write.csv(am_mgmdat_noCA, file = "am_mgmdat_noCA_16.08.2018.csv")# then enter signs into excel file am_mgm_noCA <- read.csv("File_directory\\am_mgmdat_noCA_16.08.2018.csv"); am_mgm_noCA am_mgm_noCA <- as.matrix(am_mgm_noCA); am_mgm_noCA cor5 <- am_mgm_noCA[upper.tri(am_mgm_noCA, diag=FALSE)]; cor5 #f: main model for complete sample noCA_network_cor_mgm <- cor_auto(noCA_network_mgm); noCA_network_cor_mgm #Variables detected as ordinal: noCA_Network_Data.noCA_aggressionbin; noCA_Network_Data.noCA_expressivesuppressionbin max(noCA_network_cor_mgm[upper.tri(noCA_network_cor_mgm,diag=FALSE)]) # 0.6651351 noCA_RF_glasso_mgm <- qgraph(noCA_network_cor_mgm, layout = Layout, graph = "glasso", details = TRUE, nodeNames = colnames(noCA_network_cor_mgm), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(noCA_network_mgm)) am_ggm <- getWmat(noCA_RF_glasso_mgm); am_ggm cor6 <- am_ggm[upper.tri(am_ggm,diag=FALSE)]; cor6 ####D: Sensitivity analyses for the noCA group #a: main model (qgraph full sample) vs main model without transformation (qgraph full sample) cor(cor1, cor2, method = "pearson") # 0.9991662 cor(cor1, cor2, method = "spearman") # 0.9953488 #b: main model (qgraph full sample) vs mean score model (qgraph full sample) cor(cor1, cor3, method = "pearson") # 0.9875692 cor(cor1, cor3, method = "spearman") # 0.949658 #c: main model (qgraph full sample) vs mean score model without transformation (qgraph full sample) cor(cor1, cor4, method = "pearson") # 0.9844638 cor(cor1, cor4, method = "spearman") # 0.9233595 #d: main model (qgraph full sample) vs mgm model (mgm subsample) cor(cor1, cor5, method = "pearson") # 0.9737368 cor(cor1, cor5, method = "spearman") # 0.897561 #e: main model (qgraph full sample) vs complete sample model (qgraph subsample) cor(cor1, cor6, method = "pearson") # 0.9966776 cor(cor1, cor6, method = "spearman") # 0.9845693 #f: mgm model (mgm subsample) vs complete sample model (qgraph subsample) cor(cor5, cor6, method = "pearson") # 0.9817638 cor(cor5, cor6, method = "spearman") # 0.9291292 ``` ##23: Supplement VIII: Main models depicted with faded 'RF-RF' and 'RF-General Distress' interrelations ```{r} ####A: Models for CA #a: RF-RF #pdf(file = "CA_RF_glasso_fade.pdf") CA_RF_glasso_fade <- qgraph(CA_network_cor, layout = Layout, graph = "glasso", details = TRUE, nodeNames = colnames(CA_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = TRUE, sampleSize = nrow(CA_network_c)) #dev.off() #b: RF-RF + general distress #pdf(file = "CA_PP_RF_glasso_fade.pdf") CA_PP_RF_glasso_fade <- qgraph(CA_PP_network_cor, layout = PP_Layout, graph = "glasso", details = TRUE, nodeNames = colnames(CA_PP_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = TRUE, sampleSize = nrow(CA_PP_network_c)) #dev.off() ####B: Models for noCA #a: RF-RF #pdf(file = "noCA_RF_glasso_fade.pdf") noCA_RF_glasso_fade <- qgraph(noCA_network_cor, layout = Layout, graph = "glasso", details = TRUE, nodeNames = colnames(noCA_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = TRUE, sampleSize = nrow(noCA_network_c)) #dev.off() #b: RF-RF + general distress #pdf(file = "noCA_PP_RF_glasso_fade.pdf") noCA_PP_RF_glasso_fade <- qgraph(noCA_PP_network_cor, layout = PP_Layout, graph = "glasso", details = TRUE, nodeNames = colnames(noCA_PP_network_cor), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = TRUE, sampleSize = nrow(noCA_PP_network_c)) #dev.off() ``` ##24: Supplement IX: Exploring the expressive suppression variable ```{r} ####A: moderation model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$expressivesuppressionbin) mod <- na.omit(model) maineffects <- lm(mod$Network_Data.generaldistress_tr ~ mod$Network_Data.CA_dic + mod$Network_Data.expressivesuppressionbin); summary(maineffects) moderation <- lm(mod$Network_Data.generaldistress_tr ~ mod$Network_Data.CA_dic*mod$Network_Data.expressivesuppressionbin, data=mod); summary(moderation) ####B: mediation mediation <-'Network_Data.generaldistress_tr ~ c*Network_Data.CA_dic Network_Data.expressivesuppressionbin ~ a*Network_Data.CA_dic Network_Data.generaldistress_tr ~ b*Network_Data.expressivesuppressionbin ab := a*b total := c + (a*b)' fit_mediation <- sem(mediation, data = mod, estimator = "MLR") summary(fit_mediation) parameterEstimates(fit_mediation, zstat = TRUE, pvalue = TRUE, ci = TRUE, level = 0.95, standardized = TRUE, rsquare = TRUE) ``` ##25: Supplement IX: RF network of the CA group without the expressive suppression variable ```{r} ####A: Making the network data frame without the expressive suppression variable CA_network_c_exp <- data.frame(CA_Network_Data$friendsupport_tr, CA_Network_Data$familysupport_tr, CA_Network_Data$familycohesion_tr, CA_Network_Data$positiveselfesteem_tr, CA_Network_Data$negativeselfesteem_tr, CA_Network_Data$brooding_tr, CA_Network_Data$reflection_tr, CA_Network_Data$distresstolerance_tr, CA_Network_Data$aggressionbin) ####B: Making the correlation matrix CA_network_cor_exp <- cor_auto(CA_network_c_exp); CA_network_cor_exp #Variables detected as ordinal: CA_Network_Data.CA_aggressionbin max(CA_network_cor_exp[upper.tri(CA_network_cor_exp,diag=FALSE)]) # 0.686106 ####C: Main model with average layout (as presented in paper) #a: see syntax for retrieving layout in section 7 exp_Layout <- structure(c(-0.126424547383095, 0.799405165379372, 0.41431789943595, 1, 0.826992489872837, 0.0167398316516412, -0.526596792250546, -1, -0.365765573385925, -0.432974713716992, -0.790641633253371, -0.943141335210814, 0.0462132496070069, 0.580670902412863, 1, 0.945560460165707, -0.178658845633278, 0.183925750895849), .Dim = c(9L, 2L)) #b: regularized partial correlation network for CA without the expressive suppression variable #pdf(file = "CA_RF_glasso_exp.pdf") CA_RF_glasso_exp <- qgraph(CA_network_cor_exp, layout = exp_Layout, graph = "glasso", details = TRUE, nodeNames = colnames(CA_network_cor_exp), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(CA_network_c_exp)) #dev.off() #pdf(file = "CA_RF_glasso_exp_fade.pdf") CA_RF_glasso_exp_fade <- qgraph(CA_network_cor_exp, layout = exp_Layout, graph = "glasso", details = TRUE, nodeNames = colnames(CA_network_cor_exp), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = TRUE, sampleSize = nrow(CA_network_c_exp)) #dev.off() ####D: Checking Robustness and Stability for centrality analyses #a: interrelation accuracy set.seed(3) CA_accuracy_exp <- bootnet(CA_network_c_exp, default = "EBICglasso", type = "nonparametric", nCores = 8, nBoots = 2000, statistics = c("edge", "strength","expectedInfluence")) plot(CA_accuracy_exp$sample, layout = exp_Layout, fade = FALSE) plot(CA_accuracy_exp, labels = TRUE, order = "sample") plot(CA_accuracy_exp, labels = FALSE, order = "sample") #b: network centrality coefficients stability set.seed(4) CA_stability_exp <- bootnet(CA_network_c_exp, default = "EBICglasso", type = "case", nCores = 8, caseN = 25, nBoots =2000, statistics = c("strength","expectedInfluence")) plot(CA_stability_exp$sample, layout = exp_Layout, fade = FALSE) plot(CA_stability_exp, statistics = c("strength", "expectedInfluence")) corStability(CA_stability_exp, statistics = c("strength", "expectedInfluence"), cor = 0.7) ####E: Centrality analysis #a: node strength and node expected influence CA_centrality_exp <- centrality_auto(CA_RF_glasso_exp); CA_centrality_exp #b: node predictability CA_network_mgm_exp <- na.omit(CA_network_c_exp); nrow(CA_network_mgm_exp) # 512 RFtype<- c(rep("g", 8),rep("c", 1)); RFtype RFlev<- c(rep(1,8),rep(2,1)); RFlev set.seed(1) fit_obj <- mgm(data = CA_network_mgm_exp, type = RFtype, lev = RFlev, rule.reg = "AND", k = 2, binarySign = TRUE, scale = TRUE) fit_obj$pairwise$wadj fit_obj$pairwise$signs pred_obj <- predict(fit_obj, data = CA_network_mgm_exp, errorCon = c("RMSE", "R2"), errorCat = c("CC", "nCC")) pred_obj$errors ``` ##26: Supplement IX: RF network of the noCA group without the expressive suppression variable ```{r} ####A: Making the network data frame noCA_network_c_exp <- data.frame(noCA_Network_Data$friendsupport_tr, noCA_Network_Data$familysupport_tr, noCA_Network_Data$familycohesion_tr, noCA_Network_Data$positiveselfesteem_tr, noCA_Network_Data$negativeselfesteem_tr, noCA_Network_Data$brooding_tr, noCA_Network_Data$reflection_tr, noCA_Network_Data$distresstolerance_tr, noCA_Network_Data$aggressionbin) ####B: Making the correlation matrix noCA_network_cor_exp <- cor_auto(noCA_network_c_exp); noCA_network_cor_exp #Variables detected as ordinal: noCA_Network_Data.noCA_aggressionbin max(noCA_network_cor_exp[upper.tri(noCA_network_cor_exp,diag=FALSE)]) # 0.6609587 ####C: Main model with average layout (as presented in paper) #b: regularized partial correlation network for noCA without the expressive suppression variable #pdf(file = "noCA_RF_glasso_exp.pdf") noCA_RF_glasso_exp <- qgraph(noCA_network_cor_exp, layout = exp_Layout, graph = "glasso", details = TRUE, nodeNames = colnames(noCA_network_cor_exp), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = FALSE, sampleSize = nrow(noCA_network_c_exp)) #dev.off() #pdf(file = "noCA_RF_glasso_exp_fade.pdf") noCA_RF_glasso_exp_fade <- qgraph(noCA_network_cor_exp, layout = exp_Layout, graph = "glasso", details = TRUE, nodeNames = colnames(noCA_network_cor_exp), maximum = 0.71, minimum = 0, cut = 0, legend = FALSE, theme = "colorblind", fade = TRUE, sampleSize = nrow(noCA_network_c_exp)) #dev.off() ####D: Checking Robustness and Stability for centrality analyses #a: interrelation accuracy set.seed(3) noCA_accuracy_exp <- bootnet(noCA_network_c_exp, default = "EBICglasso", type = "nonparametric", nCores = 8, nBoots = 2000, statistics = c("edge", "strength","expectedInfluence")) plot(noCA_accuracy_exp$sample, layout = exp_Layout, fade = FALSE) plot(noCA_accuracy_exp, labels = TRUE, order = "sample") plot(noCA_accuracy_exp, labels = FALSE, order = "sample") #b: network centrality coefficients stability set.seed(4) noCA_stability_exp <- bootnet(noCA_network_c_exp, default = "EBICglasso", type = "case", nCores = 8, caseN = 25, nBoots =2000, statistics = c("strength","expectedInfluence")) plot(noCA_stability_exp$sample, layout = exp_Layout, fade = FALSE) plot(noCA_stability_exp, statistics = c("strength", "expectedInfluence")) corStability(noCA_stability_exp, statistics = c("strength", "expectedInfluence"), cor = 0.7) ####E: Centrality analyses #a: node strength and node expected influence noCA_centrality_exp <- centrality_auto(noCA_RF_glasso_exp); noCA_centrality_exp #b: node predictability noCA_network_mgm_exp <- na.omit(noCA_network_c_exp); nrow(noCA_network_mgm_exp) # 444 RFtype<- c(rep("g", 8),rep("c", 1)); RFtype RFlev<- c(rep(1,8),rep(2,1)); RFlev set.seed(1) fit_obj <- mgm(data = noCA_network_mgm_exp, type = RFtype, lev = RFlev, rule.reg = "AND", k = 2, binarySign = TRUE, scale = TRUE) fit_obj$pairwise$wadj fit_obj$pairwise$signs pred_obj <- predict(fit_obj, data = noCA_network_mgm_exp, errorCon = c("RMSE", "R2"), errorCat = c("CC", "nCC")) pred_obj$errors ``` ##27: Supplement IX: Comparing the CA and the noCA network structure without the expressive suppression variable ```{r} ####A: correlate CA with noCA RF interrelation matrix without the expressive suppression variable adjmat_CA_exp <- getWmat(CA_RF_glasso_exp); adjmat_CA_exp cor_CA_exp <- adjmat_CA_exp[upper.tri(adjmat_CA_exp,diag=FALSE)]; cor_CA_exp adjmat_noCA_exp <- getWmat(noCA_RF_glasso_exp); adjmat_noCA_exp cor_noCA_exp <- adjmat_noCA_exp[upper.tri(adjmat_noCA_exp,diag=FALSE)]; cor_noCA_exp cor_pears_exp <- cor(cor_CA_exp,cor_noCA_exp, method = "pearson"); cor_pears_exp # 0.9524978 cor_spear_exp <- cor(cor_CA_exp,cor_noCA_exp, method = "spearman"); cor_spear_exp # 0.7011594 ####B: NCT approach for the comparison of the CA vs noCA networks without the expressive suppression variable #a: rename variables so that they are equal for the two groups CA_network_NCT_exp <- CA_network_c_exp CA_network_NCT_exp <- rename(CA_network_NCT_exp, friendsupport = CA_Network_Data.friendsupport_tr) CA_network_NCT_exp <- rename(CA_network_NCT_exp, familysupport = CA_Network_Data.familysupport_tr) CA_network_NCT_exp <- rename(CA_network_NCT_exp, familycohesion = CA_Network_Data.familycohesion_tr) CA_network_NCT_exp <- rename(CA_network_NCT_exp, positiveselfesteem = CA_Network_Data.positiveselfesteem_tr) CA_network_NCT_exp <- rename(CA_network_NCT_exp, negativeselfesteem = CA_Network_Data.negativeselfesteem_tr) CA_network_NCT_exp <- rename(CA_network_NCT_exp, distresstolerance = CA_Network_Data.distresstolerance_tr) CA_network_NCT_exp <- rename(CA_network_NCT_exp, aggression = CA_Network_Data.aggressionbin) CA_network_NCT_exp <- rename(CA_network_NCT_exp, brooding = CA_Network_Data.brooding_tr) CA_network_NCT_exp <- rename(CA_network_NCT_exp, reflection = CA_Network_Data.reflection_tr) head(CA_network_NCT_exp);nrow(CA_network_NCT_exp) noCA_network_NCT_exp <- noCA_network_c_exp noCA_network_NCT_exp <- rename(noCA_network_NCT_exp, friendsupport = noCA_Network_Data.friendsupport_tr) noCA_network_NCT_exp <- rename(noCA_network_NCT_exp, familysupport = noCA_Network_Data.familysupport_tr) noCA_network_NCT_exp <- rename(noCA_network_NCT_exp, familycohesion = noCA_Network_Data.familycohesion_tr) noCA_network_NCT_exp <- rename(noCA_network_NCT_exp, positiveselfesteem = noCA_Network_Data.positiveselfesteem_tr) noCA_network_NCT_exp <- rename(noCA_network_NCT_exp, negativeselfesteem = noCA_Network_Data.negativeselfesteem_tr) noCA_network_NCT_exp <- rename(noCA_network_NCT_exp, distresstolerance = noCA_Network_Data.distresstolerance_tr) noCA_network_NCT_exp <- rename(noCA_network_NCT_exp, aggression = noCA_Network_Data.aggressionbin) noCA_network_NCT_exp <- rename(noCA_network_NCT_exp, brooding = noCA_Network_Data.brooding_tr) noCA_network_NCT_exp <- rename(noCA_network_NCT_exp, reflection = noCA_Network_Data.reflection_tr) head(noCA_network_NCT_exp);nrow(noCA_network_NCT_exp) #b:perform NCT test set.seed(7) CA_vs_noCA_exp <- NCT_cor_auto(CA_network_NCT_exp, noCA_network_NCT_exp, it=5000, binary.data = FALSE, paired = FALSE, weighted=TRUE, test.edges = TRUE, edges = list(c(1,2),c(1,3),c(1,4),c(1,5),c(1,6),c(1,7),c(1,8),c(1,9),c(2,3),c(2,4),c(2,5),c(2,6),c(2,7),c(2,8), c(2,9),c(3,4),c(3,5),c(3,6),c(3,7),c(3,8),c(3,9),c(4,5),c(4,6),c(4,7),c(4,8),c(4,9), c(5,6),c(5,7),c(5,8),c(5,9),c(6,7),c(6,8),c(6,9),c(7,8),c(7,9),c(8,9)), progressbar=TRUE) #c: descriptives summary(CA_vs_noCA_exp) plot(CA_vs_noCA_exp, what="network") plot(CA_vs_noCA_exp, what="strength") plot(CA_vs_noCA_exp, what="edge") ####C:To compare the Global Network Expected Influence, please run the function in the following scipt: "2_b_2_Analysis_Script_NCT_EI-auto_JF-EF-IG-PW-ALvH_2018.07July.03_v1". The script is adapted from: van Borkulo, C. D. et al. Comparing network structures on three aspects: A permutation test (Manuscript submitted for publication, 2017). The adapted script can be requested from the first author jf585@cam.ac.uk. #Please note, what reads now as 'Global strength per group' contains the values of the global network expected influence. #a:perform NCT test set.seed(7) EI_CA_vs_noCA_exp <- NCT_cor_auto_EI(CA_network_NCT_exp, noCA_network_NCT_exp, it=5000, binary.data = FALSE, paired = FALSE, weighted=TRUE, test.edges = TRUE, edges = list(c(1,2),c(1,3),c(1,4),c(1,5),c(1,6),c(1,7),c(1,8),c(1,9),c(2,3),c(2,4),c(2,5),c(2,6),c(2,7),c(2,8), c(2,9),c(3,4),c(3,5),c(3,6),c(3,7),c(3,8),c(3,9),c(4,5),c(4,6),c(4,7),c(4,8),c(4,9), c(5,6),c(5,7),c(5,8),c(5,9),c(6,7),c(6,8),c(6,9),c(7,8),c(7,9),c(8,9)), progressbar=TRUE) summary(EI_CA_vs_noCA_exp) ``` ##28: Supplement X: Node strength and expected influence for networks corrected for the general distress variable, for both CA and noCA ```{r} ####A: CA group: Node strength and expected influence #a: CA network without the general distress variable CA_strength <- data.frame(colSums(abs(getWmat(CA_RF_glasso)))); CA_strength CA_expInf <- data.frame(colSums(getWmat(CA_RF_glasso))); CA_expInf CA_NCC <- cbind(CA_strength, CA_expInf); CA_NCC #b: CA network corrected for the general distress variable x <- getWmat(CA_PP_RF_glasso); x CA_correctPP_RF_glasso <- x[-c(11),-c(11)]; CA_correctPP_RF_glasso CA_PP_strength <- data.frame(colSums(abs(CA_correctPP_RF_glasso))); CA_PP_strength CA_PP_expInf <- data.frame(colSums(CA_correctPP_RF_glasso)); CA_PP_expInf CA_PP_NCC <- cbind(CA_PP_strength, CA_PP_expInf); CA_PP_NCC #Comparing #node strength cor(CA_NCC$colSums.abs.getWmat.CA_RF_glasso..., CA_PP_NCC$colSums.abs.CA_correctPP_RF_glasso.., method = "pearson") # 0.748537 cor(CA_NCC$colSums.abs.getWmat.CA_RF_glasso..., CA_PP_NCC$colSums.abs.CA_correctPP_RF_glasso.., method = "spearman") # 0.7333333 #expected influence cor(CA_NCC$colSums.getWmat.CA_RF_glasso.., CA_PP_NCC$colSums.CA_correctPP_RF_glasso., method = "pearson") # 0.9267678 cor(CA_NCC$colSums.getWmat.CA_RF_glasso.., CA_PP_NCC$colSums.CA_correctPP_RF_glasso., method = "spearman") # 0.8060606 ####B: noCA group: Node strength and expected influence #a: noCA network without the general distress variable noCA_strength <- data.frame(colSums(abs(getWmat(noCA_RF_glasso)))); noCA_strength noCA_expInf <- data.frame(colSums(getWmat(noCA_RF_glasso))); noCA_expInf noCA_NCC <- cbind(noCA_strength, noCA_expInf); noCA_NCC #b: noCA network corrected for the general distress variable y <- getWmat(noCA_PP_RF_glasso); y noCA_correctPP_RF_glasso <- y[-c(11),-c(11)]; noCA_correctPP_RF_glasso noCA_PP_strength <- data.frame(colSums(abs(noCA_correctPP_RF_glasso))); noCA_PP_strength noCA_PP_expInf <- data.frame(colSums(noCA_correctPP_RF_glasso)); noCA_PP_expInf noCA_PP_NCC <- cbind(noCA_PP_strength, noCA_PP_expInf); noCA_PP_NCC #Comparing #node strength cor(noCA_NCC$colSums.abs.getWmat.noCA_RF_glasso..., noCA_PP_NCC$colSums.abs.noCA_correctPP_RF_glasso.., method = "pearson") # 0.7862676 cor(noCA_NCC$colSums.abs.getWmat.noCA_RF_glasso..., noCA_PP_NCC$colSums.abs.noCA_correctPP_RF_glasso.., method = "spearman") # 0.830303 #expected influence cor(noCA_NCC$colSums.getWmat.noCA_RF_glasso.., noCA_PP_NCC$colSums.noCA_correctPP_RF_glasso., method = "pearson") # 0.8400943 cor(noCA_NCC$colSums.getWmat.noCA_RF_glasso.., noCA_PP_NCC$colSums.noCA_correctPP_RF_glasso., method = "spearman") # 0.830303 ``` ##29: Supplement XI: Network Pathways between the RFs and General Distress ```{r} #All analyses that are reported in Supplement XI were already performed in section 11 ("Shortest pathways between RFs and the general distress variable for CA and noCA groups"). ``` ##30: Supplement XII: Exploring the complex interplay between CA, RFs and general distress ```{r} ####A: Moderation and mediation example figures #a: moderation figure: brooding model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$brooding_tr) mod <- na.omit(model) maineffects <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.brooding_tr, data=mod) moderation <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic*Network_Data.brooding_tr, data=mod) #pdf(file = "moderation_figure_brooding.pdf") attach(mtcars) par(mfrow=c(1,3)) x <- lm(generaldistress_tr ~ brooding_tr, data=CA_Network_Data); summary(x) psx <- plotSlopes(x, plotx="brooding_tr", xlab = "ruminative brooding CA group", ylab = "general distress", plotLegend = FALSE, plotxRange = c(-3, 2), plotyRange = c(-3, 3.5)) y <- lm(generaldistress_tr ~ brooding_tr, data=noCA_Network_Data); summary(y) psy <- plotSlopes(y, plotx="brooding_tr", xlab = "ruminative brooding no-CA group", ylab = "general distress", plotLegend = FALSE, plotxRange = c(-3, 2), plotyRange = c(-3, 3.5)) ps <- plotSlopes(moderation, plotx="Network_Data.brooding_tr", modx="Network_Data.CA_dic", xlab = "ruminative brooding", ylab = "general distress", col = c("gray70","black"), plotxRange = c(-3, 2), plotyRange = c(-3, 3.5)) #dev.off() #moderation figure: aggression model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$aggressionbin); nrow(model) mod <- na.omit(model); nrow(mod) maineffects <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.aggressionbin, data=mod) moderation <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic*Network_Data.aggressionbin, data=mod) #pdf(file = "moderation_figure_aggression.pdf") ps <- plotSlopes(moderation, plotx="Network_Data.aggressionbin", modx="Network_Data.CA_dic", xlab = "aggression", ylab = "general distress", col = c("gray70","black"), plotxRange = c(-3, 2), plotyRange = c(-3, 3.5)) #dev.off() #b: mediation figure: brooding model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$brooding_tr) mod <- na.omit(model) x<-mod$Network_Data.CA_dic m<-mod$Network_Data.brooding_tr y<-mod$Network_Data.generaldistress_tr #To compute the mediation plot, please run the function in the following scipt: "2_d_Analysis_Script_Mediation-Plot_JF-EF-IG-PW-ALvH_2018.06June.22_v1" The script is adapted from Fritz, M. S. & MacKinnon, D. P. A graphical representation of the mediated effect. Behav. Res. Methods 40, 55-60 (2008). https://doi.org/10.3758/BRM.40.1.55; #the original scripts can be found in the belonging supplement material: Fritz, M. S. & MacKinnon, D. P. med_plots.R. Retrieved June 21, 2018 from Psychonomic Society Web Archive: www.psychonomic.org/archive (2008). or at: https://link.springer.com/article/10.3758/BRM.40.1.55#SupplementaryMaterial. The adapted script can be requested from the first author jf585@cam.ac.uk. #jpeg(file = "mediation_plot.jpg", width = 2000, height = 2000, res = 300) #pdf(file = "mediation_plot.pdf") a<-med_plots(cont=0,int=0);a #dev.off() ####B: Moderation and mediation analyses for all RFs #friendship ####: moderation model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$friendsupport_tr); nrow(model) mod <- na.omit(model); nrow(mod) maineffects <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.friendsupport_tr, data=mod); summary(maineffects) moderation <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic*Network_Data.friendsupport_tr, data=mod); summary(moderation) ####: mediation mediation <-'Network_Data.generaldistress_tr ~ c*Network_Data.CA_dic Network_Data.friendsupport_tr ~ a*Network_Data.CA_dic Network_Data.generaldistress_tr ~ b*Network_Data.friendsupport_tr ab := a*b total := c + (a*b)' fit_mediation <- sem(mediation, data = mod, estimator = "MLR") summary(fit_mediation) parameterEstimates(fit_mediation, zstat = TRUE, pvalue = TRUE, ci = TRUE, level = 0.95, standardized = TRUE, rsquare = TRUE) #familysupport ####: moderation model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$familysupport_tr); nrow(model) mod <- na.omit(model); nrow(mod) maineffects <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.familysupport_tr, data=mod); summary(maineffects) moderation <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic*Network_Data.familysupport_tr, data=mod); summary(moderation) ####: mediation mediation <-'Network_Data.generaldistress_tr ~ c*Network_Data.CA_dic Network_Data.familysupport_tr ~ a*Network_Data.CA_dic Network_Data.generaldistress_tr ~ b*Network_Data.familysupport_tr ab := a*b total := c + (a*b)' fit_mediation <- sem(mediation, data = mod, estimator = "MLR") summary(fit_mediation) parameterEstimates(fit_mediation, zstat = TRUE, pvalue = TRUE, ci = TRUE, level = 0.95, standardized = TRUE, rsquare = TRUE) #familycohesion ####: moderation model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$familycohesion_tr); nrow(model) mod <- na.omit(model); nrow(mod) maineffects <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.familycohesion_tr, data=mod); summary(maineffects) moderation <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic*Network_Data.familycohesion_tr, data=mod); summary(moderation) ####: mediation mediation <-'Network_Data.generaldistress_tr ~ c*Network_Data.CA_dic Network_Data.familycohesion_tr ~ a*Network_Data.CA_dic Network_Data.generaldistress_tr ~ b*Network_Data.familycohesion_tr ab := a*b total := c + (a*b)' fit_mediation <- sem(mediation, data = mod, estimator = "MLR") summary(fit_mediation) parameterEstimates(fit_mediation, zstat = TRUE, pvalue = TRUE, ci = TRUE, level = 0.95, standardized = TRUE, rsquare = TRUE) #positiveselfesteem ####: moderation model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$positiveselfesteem_tr); nrow(model) mod <- na.omit(model); nrow(mod) maineffects <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.positiveselfesteem_tr, data=mod); summary(maineffects) moderation <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic*Network_Data.positiveselfesteem_tr, data=mod); summary(moderation) ####: mediation mediation <-'Network_Data.generaldistress_tr ~ c*Network_Data.CA_dic Network_Data.positiveselfesteem_tr ~ a*Network_Data.CA_dic Network_Data.generaldistress_tr ~ b*Network_Data.positiveselfesteem_tr ab := a*b total := c + (a*b)' fit_mediation <- sem(mediation, data = mod, estimator = "MLR") summary(fit_mediation) parameterEstimates(fit_mediation, zstat = TRUE, pvalue = TRUE, ci = TRUE, level = 0.95, standardized = TRUE, rsquare = TRUE) #negativeselfesteem ####: moderation model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$negativeselfesteem_tr); nrow(model) mod <- na.omit(model); nrow(mod) maineffects <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.negativeselfesteem_tr, data=mod); summary(maineffects) moderation <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic*Network_Data.negativeselfesteem_tr, data=mod); summary(moderation) ####: mediation mediation <-'Network_Data.generaldistress_tr ~ c*Network_Data.CA_dic Network_Data.negativeselfesteem_tr ~ a*Network_Data.CA_dic Network_Data.generaldistress_tr ~ b*Network_Data.negativeselfesteem_tr ab := a*b total := c + (a*b)' fit_mediation <- sem(mediation, data = mod, estimator = "MLR") summary(fit_mediation) parameterEstimates(fit_mediation, zstat = TRUE, pvalue = TRUE, ci = TRUE, level = 0.95, standardized = TRUE, rsquare = TRUE) #reflection ####: moderation model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$reflection_tr); nrow(model) mod <- na.omit(model); nrow(mod) maineffects <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.reflection_tr, data=mod); summary(maineffects) moderation <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic*Network_Data.reflection_tr, data=mod); summary(moderation) ####: mediation mediation <-'Network_Data.generaldistress_tr ~ c*Network_Data.CA_dic Network_Data.reflection_tr ~ a*Network_Data.CA_dic Network_Data.generaldistress_tr ~ b*Network_Data.reflection_tr ab := a*b total := c + (a*b)' fit_mediation <- sem(mediation, data = mod, estimator = "MLR") summary(fit_mediation) parameterEstimates(fit_mediation, zstat = TRUE, pvalue = TRUE, ci = TRUE, level = 0.95, standardized = TRUE, rsquare = TRUE) #brooding ####: moderation model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$brooding_tr); nrow(model) mod <- na.omit(model); nrow(mod) maineffects <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.brooding_tr, data=mod); summary(maineffects) moderation <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic*Network_Data.brooding_tr, data=mod); summary(moderation) ####: mediation mediation <-'Network_Data.generaldistress_tr ~ c*Network_Data.CA_dic Network_Data.brooding_tr ~ a*Network_Data.CA_dic Network_Data.generaldistress_tr ~ b*Network_Data.brooding_tr ab := a*b total := c + (a*b)' fit_mediation <- sem(mediation, data = mod, estimator = "MLR") summary(fit_mediation) parameterEstimates(fit_mediation, zstat = TRUE, pvalue = TRUE, ci = TRUE, level = 0.95, standardized = TRUE, rsquare = TRUE) #distress tolerance ####: moderation model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$distresstolerance_tr); nrow(model) mod <- na.omit(model); nrow(mod) maineffects <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.distresstolerance_tr, data=mod); summary(maineffects) moderation <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic*Network_Data.distresstolerance_tr, data=mod); summary(moderation) ####: mediation mediation <-'Network_Data.generaldistress_tr ~ c*Network_Data.CA_dic Network_Data.distresstolerance_tr ~ a*Network_Data.CA_dic Network_Data.generaldistress_tr ~ b*Network_Data.distresstolerance_tr ab := a*b total := c + (a*b)' fit_mediation <- sem(mediation, data = mod, estimator = "MLR") summary(fit_mediation) parameterEstimates(fit_mediation, zstat = TRUE, pvalue = TRUE, ci = TRUE, level = 0.95, standardized = TRUE, rsquare = TRUE) #additional approach fitM <- lm(Network_Data.distresstolerance_tr ~ Network_Data.CA_dic, data=mod) fitY <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.distresstolerance_tr, data=mod) fitMed <- mediate(fitM, fitY, treat="Network_Data.CA_dic", mediator="Network_Data.distresstolerance_tr") summary(fitMed) fitMedBoot <- mediate(fitM, fitY, boot=TRUE, sims=1000, treat="Network_Data.CA_dic", mediator="Network_Data.distresstolerance_tr") summary(fitMedBoot) #aggression ####: moderation model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$aggressionbin); nrow(model) mod <- na.omit(model); nrow(mod) maineffects <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.aggressionbin, data=mod); summary(maineffects) moderation <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic*Network_Data.aggressionbin, data=mod); summary(moderation) ####: mediation mediation <-'Network_Data.generaldistress_tr ~ c*Network_Data.CA_dic Network_Data.aggressionbin ~ a*Network_Data.CA_dic Network_Data.generaldistress_tr ~ b*Network_Data.aggressionbin ab := a*b total := c + (a*b)' fit_mediation <- sem(mediation, data = mod, estimator = "MLR") summary(fit_mediation) parameterEstimates(fit_mediation, zstat = TRUE, pvalue = TRUE, ci = TRUE, level = 0.95, standardized = TRUE, rsquare = TRUE) #additional approach fitM <- lm(Network_Data.aggressionbin ~ Network_Data.CA_dic, data=mod) fitY <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.aggressionbin, data=mod) fitMed <- mediate(fitM, fitY, treat="Network_Data.CA_dic", mediator="Network_Data.aggressionbin") summary(fitMed) fitMedBoot <- mediate(fitM, fitY, boot=TRUE, sims=1000, treat="Network_Data.CA_dic", mediator="Network_Data.aggressionbin") summary(fitMedBoot) #expressive suppression ####: moderation model <- data.frame(Network_Data$generaldistress_tr, Network_Data$CA_dic, Network_Data$expressivesuppressionbin); nrow(model) mod <- na.omit(model); nrow(mod) maineffects <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic + Network_Data.expressivesuppressionbin, data=mod); summary(maineffects) moderation <- lm(Network_Data.generaldistress_tr ~ Network_Data.CA_dic*Network_Data.expressivesuppressionbin, data=mod); summary(moderation) ####: mediation mediation <-'Network_Data.generaldistress_tr ~ c*Network_Data.CA_dic Network_Data.expressivesuppressionbin ~ a*Network_Data.CA_dic Network_Data.generaldistress_tr ~ b*Network_Data.expressivesuppressionbin ab := a*b total := c + (a*b)' fit_mediation <- sem(mediation, data = mod, estimator = "MLR") summary(fit_mediation) parameterEstimates(fit_mediation, zstat = TRUE, pvalue = TRUE, ci = TRUE, level = 0.95, standardized = TRUE, rsquare = TRUE) ``` ##31: Supplement XIII: Reliability and/or Validity Information for the used Measures ```{r} #This supplement does not contain analyses ``` ##32: Supplement XIV: Individual RF Interrelation Differences between the CA and the no-CA Networks ```{r} #All analyses that are reported in Supplement XIV were already performed in section 5 ("Contrasting the CA and the noCA network structure"), section 9 ("Contrasting the CA and the noCA network structure including the general distress variable"), and section 10 ("Contrasting CA and noCA RF networks corrected for the general distress variable"). ```