A stimulus-locked vector autoregressive model for slow event-related fMRI designs

Abstract

Neuroscientists have become increasingly interested in exploring dynamic relationships among brain regions. Such a relationship, when directed from one region toward another, is denoted by “effective connectivity.” An fMRI experimental paradigm which is well-suited for examination of effective connectivity is the slow event-related design. This design presents stimuli at sufficient temporal spacing for determining within-trial trajectories of BOLD activation, allowing for the analysis of stimulus-locked temporal covariation of brain responses in multiple regions. This may be especially important for emotional stimuli processing, which can evolve over the course of several seconds, if not longer. However, while several methods have been devised for determining fMRI effective connectivity, few are adapted to event-related designs, which include nonstationary BOLD responses and multiple levels of nesting. We propose a model tailored for exploring effective connectivity of multiple brain regions in event-related fMRI designs — a semi-parametric adaptation of vector autoregressive (VAR) models, termed “stimulus-locked VAR” (SloVAR). Connectivity coefficients vary as a function of time relative to stimulus onset, are regularized via basis expansions, and vary randomly across subjects. SloVAR obtains flexible, data-driven estimates of effective connectivity and hence is useful for building connectivity models when prior information on dynamic regional relationships is sparse. Indices derived from the coefficient estimates can also be used to relate effective connectivity estimates to behavioral or clinical measures. We demonstrate the SloVAR model on a sample of clinically depressed and normal controls, showing that early but not late cortico–amygdala connectivity appears crucial to emotional control and early but not late cortico–cortico connectivity predicts depression severity in the depressed group, relationships that would have been missed in a more traditional VAR analysis.

Introduction

Interest among neuroscientists has increasingly centered on the determination of dynamic relationships among specialized brain regions in processing information. The determination of dynamic relationships among brain regions has been termed “connectivity” analysis (Friston, 1994a, Friston, 1994b, Horwitz, 2003). In particular, an “effective connectivity” analysis attempts to determine directed influences from one brain region to another. However, researchers often have little prior knowledge about the nature of effective connectivity among regions, including the strength or direction (excitatory or inhibitory) of such relationships. Moreover, even less may be known about how these effective connectivity relationships may vary across time after the application of an experimental stimulus. We propose an exploratory model for effective connectivity analysis which is specifically tailored to elucidate directed dynamic relationships among multiple brain regions in slow event-related fMRI experimental designs, while making minimal assumptions about the nature of these relationships. We were particularly motivated by observations of unpredicted lagged and non-linear relationships between regions (Wagner et al., 2001, Siegle et al., 2007) as well as frequently observed responses that do not match a canonical hemodynamic model (e.g., Siegle et al., 2002, Siegle et al., 2007) to develop a more exploratory approach. This exploratory approach for effective connectivity is in contradistinction to structural equation models (SEM; McIntosh and Gonzales-Lima, 1994) and dynamic causal models (DCM; Friston et al., 2003), which generally require a priori the specification of effective connectivity relationships. Our exploratory model's estimates of effective connectivity can be used to assist in building confirmatory parametric models (including SEMs or DCMs). Additionally, because we explicitly allow for variation in effective connectivity across groups and subjects, indices derived from the connectivity coefficients can be used to predict individual variation in other variables, such as behavioral or clinical measures. We provide an example of this in the analysis of the application dataset below.

Vector autoregressive (VAR) models have proven popular for effective connectivity analysis over the last several years (Harrison et al., 2003, Goebel et al., 2003, Ringo Ho et al., 2005; Bhattacharya et al., 2006, Rykhlevskaia and Gratton, 2006). One reason for this popularity is that in VAR modeling, the direction and valence of effective connectivity relationships do not need to be pre-specified. Specifically, the standard VAR model relates the time-dependent structure of fMRI activation of multiple regions through a multivariate regression of prior levels of activation on current levels. Directed connectivity coefficients between regions are derived from the off-diagonal elements of the autoregressive matrices. The VAR approach has been utilized for the construction of graphical models for effective connectivity through the use of Granger causality analysis (Goebel et al., 2003, Eichler, 2005, Roebroeck et al., 2005, Valdes-Sosa et al., 2005). Granger causality (Granger, 1969) is based on the common-sense idea that causes precede effects. With this in mind, activation in one region is said to “Granger cause” activation in another if prior levels of activation in the first region improve prediction of current values in the second; for stationary Gaussian VAR processes this is equivalent to testing whether the off-diagonal elements of the autoregressive matrices corresponding to these two regions are nonzero (Eichler, 2005). However, as noted by a referee, one must be careful not to over-interpret Granger causal relationships as truly causal; among other difficulties, VAR modeling is open to confounding from unmeasured relationships, as is any statistical modeling procedure (including SEM and DCM). Nonetheless, VAR models have proven a useful methodology for the effective connectivity analysis of neuroimaging data.

To date, however, VAR models have not been adapted to exploit specific characteristics of event-related fMRI experimental designs, which include nonstationary BOLD response curves and potentially dynamic effective connectivity relationships among brain regions. Adaptation of VAR models to slow event-related designs is particularly desirable because these designs are well-suited for exploring temporal dynamics of activity in regions involved in complex cognitive and emotional tasks, which may take several seconds to process. In the slow event-related experimental paradigm, brief stimuli are separated by an intertrial interval usually on the order of 10 to 20 s. This duration is generally sufficient to allow a hemodynamic response associated with the stimulus to occur and fully decay (Glover, 1999), allowing for the examination of post-stimulus temporal dynamics among brain regions. Experimental stimuli are usually presented in a randomized order and mean stimulus-locked hemodynamic response curves can be computed by averaging over repetitions of trials of the same type. Locking activation curves to stimulus onset enables exploration of BOLD response trajectories and dynamic interrelationships among regions which can be contrasted across stimulus types, individual subjects, and/or experimental subgroups.

An example of this approach is presented in Fig. 1; data in this figure are group mean stimulus-locked response trajectories for three regions of interest (ROIs) for 24 mentally-healthy controls and 32 acutely depressed subjects. As described in more detail below, these subjects were presented with emotionally-negative words (200 ms) with an inter-stimulus interval of 10.8 s. The mean fMRI response depends both on group and, as is to be expected, time since stimulus. These data thus exhibit first-order nonstationarity.

While curve averaging within subjects or trial types provides a simple method for studying stimulus-locked first-order nonstationarities in BOLD responses, it is less readily apparent how to handle potential nonstationarities in effective connectivity relationships. Such nonstationarities arise if effective connectivity varies in a systematic fashion across the time course of a given trial type. For example, two regions may have a directed relationship near the beginning of a trial but not toward the end or an initially excitatory relationship may become less so or even inhibitory. Several studies have shown that connectivity is time and context dependent (Aertson and Preissl, 1991, Friston, 1994a, Friston, 1994b, Sato et al., 2006, Bhattacharya et al., 2006). Thus, current VAR methods which assume stationarity and otherwise do not explicitly incorporate the structure of event-related designs do a poor job at modeling important aspects of regional dynamics following application of a stimulus.

Empirical evidence also indicates that both activation and effective connectivity varies strongly from person to person (Mechelli et al., 2002) and from group to group (e.g., Anand et al., 2005). Indeed, testing for individual or group differences in connectivity among regions is often of considerable interest in and of itself (Mayberg et al., 1999, Siegle et al., 2007). Hence, another desirable property for a effective connectivity analysis is the ability to explicitly account for between-subject and between-group variability in effective connectivity relationships.

We therefore propose an adaptation of VAR models tailored for the exploration of effective connectivity in multi-subject slow event-related designs. This model, termed the “stimulus-locked VAR” (SloVAR) model, allows for fixed-effect differences in effective connectivity across stimulus types or experimental subgroups and random effect differences in connectivity across subjects within groups. Stochastic variation in BOLD response trajectories across scans and trials is modeled via multivariate Gaussian inputs, while effective connectivity is modeled via a generalized VAR framework. This generalized framework allows for a nonstationary VAR covariance structure on the regional BOLD responses by allowing the autoregression parameters to vary as a function of lag and time since last stimulus onset. This highly parameterized effective connectivity model is regularized via a low-rank basis expansion using a bivariate tensor product of cubic B-splines. We implement this using a mixed model formulation (Wand, 2003) with roughness penalty parameters treated as components of variance and estimated from the data. The mixed model is implemented within a Bayesian framework by specifying priors on the parameters, including roughness penalty parameters and using a Markov Chain Monte Carlo (MCMC) Gibbs sampling algorithm.

The output of SloVAR modeling is thus a set of surfaces representing the directed relationships of each structure in the model, at each lag, to each other structure. We envision two primary uses for this information. First, visual inspection of observed relationships facilitates exploration of non-trivial time-dependent relationships among structures. These relationships could otherwise be missed or not accounted for in analyses. Following an inspection, parameters capturing relevant aspects of variation can be extracted for testing in parametric models (e.g., the strength of association of early activity in one structure with later activity in another structure). This usage follows common practice in using nonparametric models as aides in building parametric ones (Ruppert et al., 2003). Formal inference from such empirically-determined models can be sought in the subsequent analysis of data obtained from additional subjects or by data splitting (Piccard and Berk, 1990). Second, posterior estimates of connectivity parameters can be related to known non-fMRI cognitive or clinical information to better understand the relationships of fMRI-connectivity to these variables. If quantitative relationships are observed, features of the connectivity surfaces that generated the relationships can be explored post-hoc, again, necessitating extraction of relevant parameters to allow hypothesis testing.

Following a demonstration of the utility of the SloVAR model with simulated data, each of these applications is illustrated in the Results section below, using a sample psychiatric neuroimaging dataset in which we had previously examined functional connectivity using more conventional measures. As shown in the Results, use of SloVAR increases both intuitive understanding of the data as well as its clinical applicability. These findings suggest that early (within the first 2 s) but not late within-trial cortical to amygdala connectivity appears crucial to emotional control. This result matches theoretical explanations for mutually inhibitory influences of the prefrontal cortex on the amygdala and vice-versa but for the first time adds temporal specificity to the theoretical picture. Moreover, temporal specificity is important for relating connectivity to depression severity within the depressed group: early but not late cortico–cortico connectivity predicts depression severity on a clinical measure. This predictive relationship would have been missed in more traditional VAR analyses which assume stationarity and hence do not account for the stimulus-locked temporal variation in effective connectivity observed in these data across the course of a trial.

A MATLAB toolbox implementing the SloVAR model is available from the first author upon request.