This is the code used in the slides for the Workshop ‘Applying Network Analysis to Psychological Data’ at the EFPSA Congress 2016. The slides can be downloaded here http://jmbh.github.io.
# load necessary pacakges
library(devtools)
#install_github('SachaEpskamp/qgraph')
#install_github('jmbh/mgm')
library(qgraph)
library(mgm)
library(httr) # downloading data from httpsConstructing Networks
# construct a network
AdjacencyMatrix <- matrix(0,4,4)
AdjacencyMatrix[1,2] <- AdjacencyMatrix[2,1] <- 1
AdjacencyMatrix[2,3] <- AdjacencyMatrix[3,2] <- 1
AdjacencyMatrix[2,4] <- AdjacencyMatrix[4,2] <- 1
AdjacencyMatrix## [,1] [,2] [,3] [,4]
## [1,] 0 1 0 0
## [2,] 1 0 1 1
## [3,] 0 1 0 0
## [4,] 0 1 0 0qgraph(AdjacencyMatrix) # visualize
Constructing Random Networks
## set up random network
p <- 20 # number of nodes
AdjMatrix <- matrix(0,p,p) #create empty matrix
set.seed(22) #set seed for reproducibility
AdjMatrix[upper.tri(AdjMatrix)] <- sample(0:1,(p*(p-1))/2,
prob=c(.9,.1),replace=TRUE)
AdjMatrix <- AdjMatrix + t(AdjMatrix) # make symmetric
AdjMatrix[1:5,1:5] # look at edges between first 5 nodes## [,1] [,2] [,3] [,4] [,5]
## [1,] 0 0 0 0 0
## [2,] 0 0 1 0 0
## [3,] 0 1 0 0 0
## [4,] 0 0 0 0 0
## [5,] 0 0 0 0 0qgraph(AdjMatrix) # visualize
Correlation Networks
url='https://jmbh.github.io/figs/efpsa_workshop/BDIdata.RDS'
GET(url, write_disk("BDIdata.RDS", overwrite=TRUE))## Response [https://jmbh.github.io/figs/efpsa_workshop/BDIdata.RDS]
## Date: 2016-05-14 09:14
## Status: 200
## Content-Type: application/octet-stream
## Size: 34.3 kB
## <ON DISK> C:\Users\jo\Dropbox\MyData\_PhD\_Talks\efpsa_congress_2016\na_workshop\BDIdata.RDSBDI_data <- readRDS('BDIdata.RDS')
# look at data
# data
BDI_data$data[1:3,1:5]## aids01 aids02 aids03 aids04 aids05
## [1,] 1 3 1 2 1
## [2,] 3 4 2 1 1
## [3,] 1 1 1 2 1#labels
BDI_data$vnames[1:5]## [1] "Falling Asleep" "Sleep During the Night"
## [3] "Waking Up Too Early" "Sleeping Too Much"
## [5] "Feeling Sad"CorMatrix <- cor(BDI_data$data)
round(CorMatrix[1:4, 1:4],2)## aids01 aids02 aids03 aids04
## aids01 1.00 0.28 0.22 0.02
## aids02 0.28 1.00 0.34 -0.06
## aids03 0.22 0.34 1.00 -0.10
## aids04 0.02 -0.06 -0.10 1.00# visualization: ring layout
qgraph(CorMatrix, nodeNames=BDI_data$vnames, legend.cex = .3, vsize=4)
# visualization: 'spring' layout (Fruchterman Reingold algorithm)
qgraph(CorMatrix, nodeNames=BDI_data$vnames,
legend.cex = .3, layout='spring', vsize=3)
Conditional Independence Networks
## BDI Data
# fit conditional independence network
fit <- mgmfit(BDI_data$data, rep('g', 28), rep(1, 28), d=2, pbar = FALSE)
round(fit$wadj[1:4, 1:4],2)## [,1] [,2] [,3] [,4]
## [1,] 0.00 0.13 0.03 0.00
## [2,] 0.13 0.00 0.22 0.00
## [3,] 0.03 0.22 0.00 0.05
## [4,] 0.00 0.00 0.05 0.00qgraph(fit$wadj, nodeNames=BDI_data$vnames,
legend.cex = .3, layout='spring', vsize=3)
## Autism data
url='http://jmbh.github.io/figs/efpsa_workshop/autism_datalist.RDS'
GET(url, write_disk("autism_datalist.RDS", overwrite=TRUE))## Response [http://jmbh.github.io/figs/efpsa_workshop/autism_datalist.RDS]
## Date: 2016-05-14 09:14
## Status: 200
## Content-Type: application/octet-stream
## Size: 99.7 kB
## <ON DISK> C:\Users\jo\Dropbox\MyData\_PhD\_Talks\efpsa_congress_2016\na_workshop\autism_datalist.RDSAutism_data <- readRDS('autism_datalist.RDS')
# fit mixed graphical model
#fit2 <- mgmfit(Autism_data$data, Autism_data$type, Autism_data$lev,
# d=2, pbar=FALSE) # commented out as it takes a while
# instead we download the fit object here:
url='http://jmbh.github.io/figs/efpsa_workshop/fitobj_mixed.RDS'
GET(url, write_disk("fitobj_mixed.RDS", overwrite=TRUE))## Response [http://jmbh.github.io/figs/efpsa_workshop/fitobj_mixed.RDS]
## Date: 2016-05-14 09:14
## Status: 200
## Content-Type: application/octet-stream
## Size: 21.2 kB
## <ON DISK> C:\Users\jo\Dropbox\MyData\_PhD\_Talks\efpsa_congress_2016\na_workshop\fitobj_mixed.RDSfit2 <- readRDS('fitobj_mixed.RDS')
# look at edges between first 10 variables
round(fit2$wadj[1:10, 1:10],2)## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
## [1,] 0 0.0 0 0 0.00 0.00 0.00 0.00 0.00 0.00
## [2,] 0 0.0 0 0 0.00 0.00 0.00 0.00 0.10 0.00
## [3,] 0 0.0 0 0 0.00 0.00 0.00 0.00 0.00 0.00
## [4,] 0 0.0 0 0 0.00 0.00 0.00 0.00 0.00 0.00
## [5,] 0 0.0 0 0 0.00 0.00 0.26 0.00 0.00 0.00
## [6,] 0 0.0 0 0 0.00 0.00 0.03 0.00 0.00 0.06
## [7,] 0 0.0 0 0 0.26 0.03 0.00 0.00 0.03 0.00
## [8,] 0 0.0 0 0 0.00 0.00 0.00 0.00 0.00 0.03
## [9,] 0 0.1 0 0 0.00 0.00 0.03 0.00 0.00 0.00
## [10,] 0 0.0 0 0 0.00 0.06 0.00 0.03 0.00 0.00# visualize
qgraph(fit2$wadj, nodeNames=Autism_data$colnames,
layout='spring', edge.color=fit2$edgecolor,
legend.cex=.3, vsize=3, legend.cex=1)
## Make a nicer graph
groups_type <- list("Demographics"=c(1,14,15,28),
"Psychological"=c(2,4,5,6,18,20,21),
"Social environment" = c(7,16,17,19,26,27),
"Medical"=c(3,8,9,10,11,12,13,22,23,24,25))
group_col <- c("#72CF53", "#53B0CF", "#FFB026", "#ED3939")# we need this graph for the layout: Q0$layout
Q0 <- qgraph(fit2$adj,
vsize=3.5, esize=2.5, layout="spring", edge.color = fit2$edgecolor,
color=group_col,
border.width=1.5,
border.color="black",
groups=groups_type,
nodeNames=Autism_data$colnames,
legend=TRUE,
legend.mode="style2",
legend.cex=.6)
# Produces the graph in the slides
qgraph(fit2$wadj,
vsize=3.5, esize=6,
layout=Q0$layout,
edge.color = fit2$edgecolor,
color=group_col,
border.width=1.5,
border.color="black",
groups=groups_type,
nodeNames=Autism_data$colnames,
legend=TRUE,
legend.mode="style2",
legend.cex=.6)