How to analyze the accuracy and stability of a network using PRONA

This vignette gives an overview of how to perform item stability analysis, nonparametric bootstrapping, casedrop bootstrapping, and nodedrop bootstrapping on a GGM network generated with PRONA. For running these analyses, PRONA largely serves as a wrapper for the EGANet and bootnet packages. For more information on these analyses and how to interpret their results, please see Christensen & Golino, Psych, 2021 and Epskamp, Borsboom, & Fried, Behavioral Research Methods, 2018.

1. Load PRONA, download data, and construct a GGM

For this example, we will use the data from the introduction vignette to construct a GGM. Please refer to the introduction vignette for more details on where this data is from and how the network is constructed.

library(PRONA)
reduced_df <- read.csv('../../PRONA_additional_files/test_scripts/output/reduced_data.csv', stringsAsFactors = FALSE)
ptsd_network <- construct_ggm(reduced_df, normal = FALSE)

## Loading required namespace: huge

## Conducting nonparanormal (npn) transformation via skeptic....done.

plot_ggm(ptsd_network, label.size = 2.5)

2. Perform Item Stability Analysis

To perform item stability analysis, we use the run_bootEGA function. This function has 3 parameters:

df: A dataframe of symptom severity/frequency data
normal: Whether to consider all variables normally distributed. If false, conducts a non-paranormal transformation (Default: FALSE)
iter: Number of permutations to run for bootstrapping analysis (Default: 10,000)

# Replace file path with wherever you want to save the results of bootEGA
boot.ega.file <- "../../PRONA_additional_files/test_scripts/output/boot.ega.RDS"
if (!file.exists(boot.ega.file)){
  boot.ega <- run_bootEGA(reduced_df)
  saveRDS(boot.ega, boot.ega.file)}
boot.ega <- readRDS(boot.ega.file)

Let’s plot the median network structure from bootEGA.

plot_bootEGA(boot.ega, label.size = 3.5)

Next, we calculate item and dimension stability and plot the results

dimStab <- calculate_dimStab(boot.ega)

plot_item_stability(dimStab)

3. Perform nonparametric bootstrapping to analyze edge weight accuracy

To perform nonparametric bootstrapping, we use the nonparam_boot function, which has six parameters:

df: A dataframe of symptom severity/frequency data
normal: Whether to consider all variables normally distributed. If false, conducts a non-paranormal transformation (Default: FALSE)
nBoots: Number of permutations to run for bootstrapping analysis (Default: 2,500)
statistics: Vector indicating which statistics to store. See bootnet for options. (Default: c(‘edge’,‘strength’,‘closeness’,‘betweenness’)
gamma: The EBIC tuning parameter to use. Must be 0, 0.25, or 0.5. (Default: 0)
communities: If you are running bootnet on bridge centrality measures, use this parameter to set community labels (set this to ega$wc, where ega is the output from the construct_ggm function)

nonparam.boot.file <- "../../PRONA_additional_files/test_scripts/output/nonparam.boot.RDS"
if (!file.exists(nonparam.boot.file)){
  # since the construct_ggm function used a gamma = 0.5, we will do the same here
  nonparam.boot <- nonparam_boot(reduced_df, normal = FALSE, statistics = c('edge','strength','closeness','betweenness','bridgeStrength','bridgeCloseness','bridgeBetweenness'), communities = ptsd_network$wc, gamma = 0.5)
  saveRDS(nonparam.boot, nonparam.boot.file)}
nonparam.boot <- readRDS(nonparam.boot.file)

Let’s plot and summarize the results.

library(bootnet)

## Warning: package 'bootnet' was built under R version 4.2.3

## Loading required package: ggplot2

## Warning: package 'ggplot2' was built under R version 4.2.3

## This is bootnet 1.6

## For questions and issues, please see github.com/SachaEpskamp/bootnet.

plot_nonparam_boot(nonparam.boot)

summarize_nonparam_boot(nonparam.boot)

## # A tibble: 105 × 17
## # Groups:   type, node1, node2 [105]
##    type  id    node1 node2 sample    mean     sd  CIlower CIupper    q2.5  q97.5
##    <chr> <chr> <chr> <chr>  <dbl>   <dbl>  <dbl>    <dbl>   <dbl>   <dbl>  <dbl>
##  1 edge  AVOI… AVOI… BAD.… 0.0443 0.0429  0.0417 -0.0391   0.128   0      0.138 
##  2 edge  AVOI… AVOI… DIST… 0      0.0109  0.0221 -0.0441   0.0441  0      0.0778
##  3 edge  AVOI… AVOI… FEEL… 0.0514 0.0421  0.0405 -0.0296   0.132   0      0.134 
##  4 edge  AVOI… AVOI… FEEL… 0      0.0113  0.0224 -0.0449   0.0449  0      0.0797
##  5 edge  AVOI… AVOI… FEEL… 0      0.00695 0.0189 -0.0377   0.0377  0      0.0666
##  6 edge  AVOI… AVOI… HAVI… 0      0.00724 0.0204 -0.0408   0.0408  0      0.0728
##  7 edge  AVOI… AVOI… HAVI… 0      0.00690 0.0188 -0.0376   0.0376  0      0.0645
##  8 edge  AVOI… AVOI… JUMP… 0.210  0.201   0.0519  0.107    0.314   0.0969 0.302 
##  9 edge  AVOI… AVOI… LESS… 0.0971 0.0826  0.0462  0.00475  0.189   0      0.178 
## 10 edge  AVOI… AVOI… NOT.… 0      0.00147 0.0149 -0.0297   0.0297 -0.0311 0.0412
## # ℹ 95 more rows
## # ℹ 6 more variables: q2.5_non0 <dbl>, mean_non0 <dbl>, q97.5_non0 <dbl>,
## #   var_non0 <dbl>, sd_non0 <dbl>, prop0 <dbl>

4. Perform casedrop and nodedrop bootstrapping to analyze edge weight and centrality stability

To perform casedrop bootstrapping, we use the casedrop_boot function, which has the same six parameters as the nonparam_boot function.

casedrop.boot.file <- "../../PRONA_additional_files/test_scripts/output/casedrop.boot.RDS"
if (!file.exists(casedrop.boot.file)){
  casedrop.boot <- casedrop_boot(reduced_df, normal = FALSE, statistics = c('edge','strength','closeness','betweenness','bridgeStrength','bridgeCloseness','bridgeBetweenness'), communities = ptsd_network$wc, gamma = 0.5)
  saveRDS(casedrop.boot, casedrop.boot.file)}
casedrop.boot <- readRDS(casedrop.boot.file)

Plot the casedrop bootstrapping results for edge weight, centralities, and bridge centralities

plot_casedrop_boot(casedrop.boot, type = "edge")

## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
## Scale for fill is already present.
## Adding another scale for fill, which will replace the existing scale.

plot_casedrop_boot(casedrop.boot, type = "centralities")

## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
## Scale for fill is already present.
## Adding another scale for fill, which will replace the existing scale.

plot_casedrop_boot(casedrop.boot, type = "bridge centralities")

## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
## Scale for fill is already present.
## Adding another scale for fill, which will replace the existing scale.

Next, we perform correlation stability analysis on the casedrop bootstrapping results.

cor_stability_analysis(casedrop.boot)

## === Correlation Stability Analysis === 
## 
## Sampling levels tested:
##    nPerson Drop%   n
## 1       90  74.9 281
## 2      118  67.1 227
## 3      146  59.3 263
## 4      174  51.5 275
## 5      201  44.0 239
## 6      229  36.2 260
## 7      257  28.4 236
## 8      285  20.6 229
## 9      313  12.8 242
## 10     341   5.0 248
## 
## Maximum drop proportions to retain correlation of 0.7 in at least 95% of the samples:
## 
## betweenness: 0.128 
##   - For more accuracy, run bootnet(..., caseMin = 0.05, caseMax = 0.206) 
## 
## bridgeBetweenness: 0.05 (CS-coefficient is lowest level tested)
##   - For more accuracy, run bootnet(..., caseMin = 0, caseMax = 0.128) 
## 
## bridgeCloseness: 0.128 
##   - For more accuracy, run bootnet(..., caseMin = 0.05, caseMax = 0.206) 
## 
## bridgeStrength: 0.362 
##   - For more accuracy, run bootnet(..., caseMin = 0.284, caseMax = 0.44) 
## 
## closeness: 0.128 
##   - For more accuracy, run bootnet(..., caseMin = 0.05, caseMax = 0.206) 
## 
## edge: 0.44 
##   - For more accuracy, run bootnet(..., caseMin = 0.362, caseMax = 0.515) 
## 
## strength: 0.44 
##   - For more accuracy, run bootnet(..., caseMin = 0.362, caseMax = 0.515) 
## 
## Accuracy can also be increased by increasing both 'nBoots' and 'caseN'.

Finally, let’s perform nodedrop bootstrapping and plot the results.

nodedrop.boot.file <- "../../PRONA_additional_files/test_scripts/output/nodedrop.boot.RDS"
if (!file.exists(nodedrop.boot.file)){
  nodedrop.boot <- nodedrop_boot(reduced_df, normal = FALSE, statistics = c('edge','strength','closeness','betweenness'), gamma = 0.5)
  saveRDS(nodedrop.boot, nodedrop.boot.file)}
nodedrop.boot <- readRDS(nodedrop.boot.file)

plot_nodedrop_boot(nodedrop.boot, type = 'edge')

plot_nodedrop_boot(nodedrop.boot, type = 'centralities')