…and me leaving academia
I am no longer working in academia. It was time for a change. 😉 To tackle most of the questions I got, I wanted to update the NCT package and add example code of the new functionalities in the help file.
You can now use NCT in combination with the estimateNetwork() function. This enables using various network estimation methods (e.g., mgm). Also, you can find example code for assessing differences in various centrality measures (e.g., bridge centrality). You can find the most recent version on GitHub. This version will soon be uploaded to CRAN. I copied the example code below for your convenience.
For now, I will still be available for questions in case you get errors or something. If your question needs a more involved answer, I might refer you to another expert. 😉
library("IsingSampler") library("IsingFit") library("bootnet") ### Simulate binary datasets under null hypothesis: ### underlying network structures are similar # Input: N <- 6 # Number of nodes nSample <- 500 # Number of samples # Ising parameters: set.seed(123) Graph <- matrix(sample(0:1,N^2,TRUE,prob = c(0.8, 0.2)),N,N) * runif(N^2,0.5,2) Graph <- pmax(Graph,t(Graph)) Graph[4,1] <- Graph[4,1]*-1 Graph[1,4] <- Graph[1,4]*-1 Graph[5,1] <- Graph[5,1]*-1 Graph[1,5] <- Graph[1,5]*-1 Graph[6,1] <- Graph[6,1]*-1 Graph[1,6] <- Graph[1,6]*-1 diag(Graph) <- 0 Thresh <- -rowSums(Graph) / 2 # Simulate: data1 <- IsingSampler(nSample, Graph, Thresh) data2 <- IsingSampler(nSample, Graph, Thresh) colnames(data1) <- colnames(data2) <- c('V1', 'V2', 'V3', 'V4', 'V5', 'V6') ### Compare networks of data sets using NCT ### ## Networks can be compared by either (1) feeding the data directly into NCT (whereby ## you need to specify arguments such as "gamma" and "binary.data") or (2) by using ## estimateNetwork() (bootnet package) and feeding that output into NCT. For the latter ## option, we refer to the help file of estimateNetwork() for its usage. Below, both ## options are illustrated. We recommend using estimateNetwork(), since this function ## has implemented many network estimation methods. ## gamma = 0 (in estimateNetwork this hyperparameter is called "tuning"; to illustrate # how to specify a different value than the default) ## iterations (it) set to 10 to save time ## Note: Low number of iterations can give unreliable results; should be 1000 at least ## Testing whether there are differences in the three aspects that are validated # (network invariance, global strength, edge weight) ## 2 edges are tested here: between variable 1 and 2, and between 3 and 6 (can be # "list(c(2,1),c(6,3))" as well) ## (1) Feeding data directly into NCT set.seed(123) NCT_a <- NCT(data1, data2, gamma=0, it=10, binary.data = TRUE, test.edges=TRUE, edges=list(c(1,2),c(3,6))) summary(NCT_a) ## Plot results of global strength invariance test (not reliable with only 10 # permutations!) plot(NCT_a, what="strength") ## (2) Feeding the estimateNetwork() output into NCT est_1 <- estimateNetwork(data1, default = "IsingFit", tuning = 0) est_2 <- estimateNetwork(data2, default = "IsingFit", tuning = 0) ## When using estimateNetwork() output, there is no need to specify gamma and binary.data ## This yields similar output as NCT_a set.seed(123) NCT_b <- NCT(est_1, est_2, it=10, test.edges=TRUE, edges=list(c(1,2),c(3,6))) summary(NCT_b) ## Next, an example of testing whether there are differences in node strength # when data is paired (e.g., a group which is measured pre- and post-treatement). # Also, here you can see how to specify that you want to take the sign of node strength # into account (by default, the absolute value is taken and, therefore, the sign is # ignored). ## abs = FALSE set.seed(123) NCT_c = NCT(est_1, est_2, paired = TRUE, abs = FALSE, test.edges = TRUE, edges = list(c(1,2),c(3,6)), test.centrality = TRUE, centrality = c("strength"), nodes = "all", it=10) summary(NCT_c) ## Finally, an example how to test for differences in centrality (e.g., expectedInfluence) set.seed(123) NCT_d = NCT(est_1, est_2, paired = TRUE, abs = FALSE, test.edges = TRUE, edges = list(c(1,2),c(3,6)), test.centrality = TRUE, centrality = c("expectedInfluence"), nodes = "all", it=10) summary(NCT_d)