In epongpipat/ppi: Psychophysiological Interaction (PPI)

I believe it is helpful to think of PPI as a essentially a lot of data manipulation/wrangling followed by multi-level modeling.

This example will focus on how to perform the data manipulation portion of PPI.

This example assumes that the reader:

1. is familiar with standard task fMRI analyses and jargon
2. has a conceptual understanding of psychophysiological interaction (PPI)
3. has pre-processed the fMRI data

R Packages

First, let's load some packages.

packages <- c("tidyverse", "eepR", "furrr", "ggpubr")
xfun::pkg_attach(packages, install = T, message = F)

script_dir <- "~/Documents/GitHub/ppi/R/"
script_files <- list.files(script_dir)
script_path <- paste0(script_dir, script_files)
script_ppi <- lapply(script_path, source)

Define Parameters

Let's define some MRI parameters and an upsample factor for convolution and deconvolution.

In this example, we are going to use SPM's first gamma variate (spmg1) for the HRF. This is SPM's default HRF other current options include gam and block.

MRI

tr <- 1.5
n_volumes <- 260

Files

events_file <- "~/Box/my_mri/behavioral/nback_3tb4188_2017_Sep_13_1341_soa.csv"
out_file <- NULL

PPI Options

hrf_name <- "spmg1"
upsample_factor <- 16
detrend_factor <- 2

psy_contrast_table <- cbind(nback_vs_baseline = c(1, 1, 1, 1, -4)/5,
                                              task_vs_control = c(-3, 1, 1, 1, 0)/4,
                                              linear_task = c(0, -1, 0, 1, 0),
                                              quadratic_task = c(0, -1, 2, -1, 0)/3)
psy_contrast_table

Ideal HRF

Let's first obtain the hemodynamic response function (HRF) that we are going to use during convolution and deconvolution steps.

hrf <- get_hrf_afni(hrf_name, tr, upsample_factor)

Psychological Variables

Load Data

df_psy <- read.csv(events_file) %>%
  rename(onset = start, trial_type = block) %>%
  mutate(duration = ifelse(trial_type == "0-back", 10 * 2.5, 20 * 2.5)) %>%
  select(run, onset, duration, trial_type) %>%
  group_by(run) %>%
  nest() %>%
  mutate(data = map(data, as.data.frame))

Create

psy_var <- apply(df_psy[, "data"], 1, function(x) create_psy_var(x[[1]], psy_contrast_table, hrf, tr, n_volumes, upsample_factor))

Physiological Variable

Load Data

dir_phys <- "/Volumes/shared/KK_KR_JLBS/Wave1/MRI/FMRI/PPI/Nback_individual_covariates/PHYS_MNI_42_-42_42/3tb1780/"
files_phys <- list.files(dir_phys, "_mean.1D")
path_phys <- paste0(dir_phys, files_phys)
df_phys <- tibble(file = files_phys,
                  path = path_phys) %>%
  mutate(run = str_remove(file, "_mean.1D"),
         run = str_remove(run, "r"),
         data = future_map(path, function(x) read.csv(x, header = F, col.names = "data")))

Create

phys_var <- apply(df_phys[, "data"], 1, function(x) create_phys_var(x[[1]], detrend_factor, upsample_factor, hrf))

Psychophysiological Interaction Variable

Create

temp_psy_var <- map(psy_var, "upsample")
temp_phys_var <- map(phys_var, "deconvolve")
ppi_var <- list()
for (i in 1:3) {
  ppi_var[[i]] <- create_ppi_var(temp_psy_var[[i]], temp_phys_var[[i]], hrf, tr, n_volumes, upsample_factor)
}

Design Matrix

Initialize

temp_psy_var <- map(psy_var, "downsample")
temp_phys_var <- map(phys_var, "detrend")
temp_ppi_var <- map(ppi_var, "downsample")

df_design_mat <- tibble(run = c(1:3)) %>%
  mutate(data = future_map(run, function(x) create_design_matrix(temp_psy_var[[x]], temp_phys_var[[x]], temp_ppi_var[[x]]))) %>%
  unnest()

head(df_design_mat)

func_fig_data_heat_map <- function(data) {
  colnames(data) <- paste0(str_pad(1:ncol(data), 2, "left", "0"), "_", colnames(data))
  apply(data, 2, scale_min_max) %>%
    as.data.frame() %>%
    mutate(volume = row_number()) %>%
    gather(., variable, value, -volume) %>%
    ggplot(., aes(variable, volume, fill = value)) +
    geom_raster() + 
    scale_fill_distiller(palette = "Greys", direction = 1) + 
    theme_minimal() +
    theme(axis.text.x = element_text(hjust = 1, angle = 45)) +
    labs(x = NULL)
}

func_fig_data_ts_long <- function(data) {
  colnames(data) <- paste0(str_pad(1:ncol(data), 2, "left", "0"), "_", colnames(data))
  data %>%
  mutate(volume = row_number()) %>%
  gather(., variable, value, -volume) %>%
  ggplot(., aes(volume, value)) +
  geom_line() +
  facet_wrap(~ variable, scales = "free_y", ncol = 1) +
  theme_minimal()
}

func_fig_data_heat_map(df_design_mat)

func_fig_data_ts_long(df_design_mat)

Save

if (!is.null(out_file)) {
  write.csv(df_design_mat_final, out_file)
}

epongpipat/ppi documentation built on Jan. 31, 2024, 1:02 p.m.