```{r, include=FALSE} knitr::opts_chunk$set(echo = FALSE) ``` ## Necessary Software - R: https://www.r-project.org - RStudio: https://www.rstudio.com Download Code (optional): http://jmbh.github.io/data/CS_WS_Utrecht_March15.Rmd ## Install & Load Packages ```{r, echo = TRUE, message=FALSE} # For Estimation # install.packages('mgm') library(mgm) # For Visualization # install.packages('qgraph') library(qgraph) ``` ## The Data Subset of 16 continuous variables of larger ESM dataset: - One individual diagnosed with Major Depression - Up to 10 measurements a day - Over 238 consecutive days - In total 1476 Measurement points Openly available: https://osf.io/j4fg8/ We select 16 mood related variables: Feeling relaxed, feeling down, feeling irritated, feeling lonely, feeling anxious, feeling satisfied, etc. ## The Data
## Load Data ```{r, echo = TRUE} data <- read.csv('http://jmbh.github.io/data/mood_data.csv', header = TRUE) data[1:3, 1:8] data_mood <- data[, -c(1:3)] # delete id variables dim(data_mood) p <- ncol(data_mood) ``` ## Estimate Mixed Graphical Model ```{r, echo = TRUE, message = FALSE} mgm_obj <- mgmfit(data = data_mood, type = rep('g', p), lev = rep(1, p), d = 1, lambda.sel = 'CV', pbar = FALSE) # Subset of weighted adjacency matrix round(mgm_obj$wadj, 2)[1:4, 1:4] ``` ## Get Signs of Pairwise Interactions ```{r, echo = TRUE} # Add signs to weighted adjacency matrix wadj <- mgm_obj$wadj * mgm_obj$signs wadj[is.na(wadj)] <- 0 round(wadj, 2)[1:4, 1:4] ``` ## Visualize as Network ```{r, echo = TRUE} Q <- qgraph(wadj, layout = 'spring', labels = colnames(data_mood)) ``` ## Compute Prediction Error ```{r, echo = TRUE} pred_obj <- predict(mgm_obj, data_mood, error.continuous = 'VarExpl') pred_obj$error[1:10,] ``` ## Visualize Network with Prediction Error ```{r, echo = TRUE} qgraph(wadj, layout = 'spring', labels = colnames(data_mood), pie = pred_obj$error$Error, pieColor = 'darkgrey') ``` ## Estimate lagged (VAR) model ```{r, echo = TRUE, messages=FALSE} var_obj <- var.mgm(data = data_mood, type = rep('g', p), lev = rep(1, p), lags = 1, lambda.sel = 'CV', pbar = FALSE) round(var_obj$wadj, 2)[1:4, 1:4] ``` ## Get Signs of Lag-1 effects ```{r, echo = TRUE} # Add signs to lagged effects evar <- var_obj$wadj * var_obj$signs evar[is.na(evar)] <- 0 round(evar, 2)[1:4, 1:4] ``` ## Visualize VAR Network ```{r, echo = TRUE} qgraph(evar, layout = 'spring', labels = colnames(data_mood)) ``` ## Assumption of Stationarity
## Assumption of Stationarity?
## Specification of Time Vector ```{r, echo = TRUE, fig.width=4, fig.height=3} # Calculate Time Vector dates <- as.character(data$date) times <- as.character(data$beeptime) x <- paste(dates, times) date_time <- strptime(x, "%d/%m/%y %H:%M") timeFromZero <- date_time - date_time[1] timeFromZero <- as.numeric(timeFromZero) plot(timeFromZero, type = 'l') # how irregular? ``` ## Estimating time-varying MGM ```{r, echo = TRUE} tv_mgm <- tv.mgmfit(data = data_mood, type = rep('g', p), lev = rep(1, p), timepoints = timeFromZero, tsteps = 15, bandwidth = .1, d = 1, pbar = FALSE) # Add signs tvwadj <- tv_mgm$wadj * tv_mgm$sign tvwadj[is.na(tvwadj)] <- 0 ``` ## Visualize tv MGM at t = 1, 8, 15 ```{r, echo = TRUE, fig.width=8, fig.height=3} cln <- colnames(data_mood) par(mfrow=c(1,3)) qgraph(tvwadj[, , 1], layout = Q$layout, labels = cln, title = 'Time = 1') qgraph(tvwadj[, , 8], layout = Q$layout, labels = cln, title = 'Time = 8') qgraph(tvwadj[, , 15], layout = Q$layout, labels = cln, title = 'Time = 15') ``` ## Visualize tv MGM: Specific edges ```{r, echo = TRUE, fig.width=7, fig.height=4} plot(tvwadj[2, 5, ], type = 'l', ylim = c(-.2,.8), xlab = 'Time', ylab = 'Parameter value') lines(tvwadj[5, 15, ], col = 'red') legend('topleft', c('Down - Lonely', 'Lonely - Worrying'), lwd=c(1, 1), col = c('black', 'red')) ``` ## Other features - Mixed Variables - Higher-order interactions - Time-lagged (VAR) models with lag k > 1 More Info: - Paper & Tutorials: http://jmbh.github.io - Package: https://cran.r-project.org/package=mgm