A MATLAB toolbox for Granger causal connectivity analysis

Abstract

Assessing directed functional connectivity from time series data is a key challenge in neuroscience. One approach to this problem leverages a combination of Granger causality analysis and network theory. This article describes a freely available MATLAB toolbox – ‘Granger causal connectivity analysis’ (GCCA) – which provides a core set of methods for performing this analysis on a variety of neuroscience data types including neuroelectric, neuromagnetic, functional MRI, and other neural signals. The toolbox includes core functions for Granger causality analysis of multivariate steady-state and event-related data, functions to preprocess data, assess statistical significance and validate results, and to compute and display network-level indices of causal connectivity including ‘causal density’ and ‘causal flow’. The toolbox is deliberately small, enabling its easy assimilation into the repertoire of researchers. It is however readily extensible given proficiency with the MATLAB language.

Introduction

The accelerating availability of neuroscience data is placing increasing demands on analysis methods. These data are generated at multiple levels of description, from spike trains to local field potentials to functional MRI (fMRI) BOLD signals. A key challenge when analyzing such data is to determine the functional connectivity of the underlying mechanisms. Many methods for functional connectivity analysis identify undirected connectivity; examples include synchrony (Engel and Singer, 2001) and phase coherence (Nunez et al., 2001, Doesburg et al., 2009). However, a satisfactory understanding of neural mechanisms is likely to require identification of directed functional connectivity. A powerful technique for extracting such connectivity from data is Granger causality (G-causality) (Granger, 1969, Ding et al., 2006, Seth, 2007). According to G-causality, a variable $X_1$ ‘Granger causes’ a variable $X_2$ if information in the past of $X_1$ helps predict the future of $X_2$ with better accuracy than is possible when considering only information in the past of $X_2$ itself (Granger, 1969, Seth, 2007). This paper describes a freely available software toolbox, programmed in the MATLAB (Natick, MA) environment, which allows application of a range of G-causality analyses to neuroscience data broadly construed. Taken together, the analysis methods incorporated in the toolbox form ‘Granger causal connectivity analysis’ (GCCA). A first version of the toolbox was released in 2005 (Seth, 2005); the present paper describes the first significant revision and extension of the software, as well as the theoretical infrastructure underpinning its functionality.

It should be noted that identifying directed functional connectivity is not equivalent to identifying physically instantiated causal interactions in systems. Although the two descriptions are intimately related (Seth and Edelman, 2007, Cadotte et al., 2008), physically instantiated causal structure can only be unambiguously identified by perturbing a system and observing the consequences (Pearl, 1999). Directed functional connectivity in general, and G-causality analysis in particular, is therefore best understood as a statistical relationship among observed variables that reflects but may not be identical to the underlying physical mechanism.

The importance of identifying causal structure within data, especially during exploratory analysis phases, indicates a need for easy-to-apply, transparent, and extensible software methods. Such methods are provided by the GCCA toolbox described here. The toolbox includes several different types of function. The core functions implement G-causality analysis given multivariate time series data. Other functions test whether the provided data satisfies necessary assumptions, assess the statistical significance and validity of inferred interactions, generate network-level descriptions of patterns of causal interactions, and graphically display analysis results. Functions are also included to apply various preprocessing techniques and to demonstrate the toolbox capabilities. The toolbox is intentionally small when compared to several other brain signal analysis toolboxes (e.g., Delorme and Makeig, 2004, Cui et al., 2008), in order to facilitate assimilation of the methods into the repertoire of researchers. Also, the toolbox is not targeted specifically to any particular experimental technique. G-causality analysis has been successfully applied to a wide range of data types including spike trains (Cadotte et al., 2008, Nedungadi et al., 2009), local field potentials (Kaminski et al., 2001, Brovelli et al., 2004), EEG signals (Babiloni et al., 2005, Keil et al., 2009), fMRI BOLD signals (Roebroeck et al., 2005, Sato et al., 2006), among others. Several such studies have utilized an earlier version of the present toolbox (Sridharan et al., 2008, Shanahan, 2008, Gow et al., 2008, Gow and Segawa, 2009, Gaillard et al., 2009, Chang et al., 2008, Stevens et al., 2009) or have followed the procedures it advocated (Keil et al., 2009). By November 2009, over 600 researchers from many different countries had downloaded the toolbox. The present toolbox extends the earlier version in several important ways, including (i) frequency decomposition of causal interactions, (ii) partial Granger causality, (iii) ‘Granger autonomy’, (iv) analysis of event-related data, (v) additional preprocessing and validation techniques, for example for removal of single frequency line-noise from data and for checking model consistency, and (vi) built-in functions for bootstrap and permutation statistical testing. The toolbox is also enhanced with respect to the efficiency of its core algorithms.

A screenshot of the toolbox demonstration function is provided in Fig. 1. In this figure, data are generated according to the model described in Baccalá and Sameshima (2001) (Eq. (1.1), see top-left panel).

$$\begin{array}{c}{x}_1(t)=0.95\sqrt{2}{x}_1(t-1)-0.9025x_1(t-2)+w_1(t)\hfill \\ x_2(t)=0.5x_1(t-2)+w_2(t)\hfill \\ x_3(t)=-0.4x_1(t-3)+w_3(t)\hfill \\ x_4(t)=-0.5x_1(t-2)+0.25\sqrt{2}{x}_4(t-1)+0.25\sqrt{2}{x}_5(t-1)+w_4(t)\hfill \\ x_5(t)=-0.25\sqrt{2}{x}_4(t-1)+0.25\sqrt{2}{x}_5(t-1)+w_5(t)\hfill \end{array}$$

The model generates five time-series which are shown in the bottom-right panel of the figure. $w_1-w_5$ are independent normally-distributed processes with zero mean and unit standard deviation. The remaining panels show the results of GCCA of the time-series. The top panels show G-causalities in both matrix and network form. In this case the inferred causalities are identical to the underlying physical network structure (as noted this need not be true in general). The bottom panels show network-level summary descriptions of the causal patterns, described further in Section 4. Briefly, causal flow reflects the extent to which a variable is influenced by or influences the remainder of the system; causal density expresses the overall degree of causal interactivity, either averaged across an entire network or assessed on a per-variable basis.

The remainder of this paper is organized as follows. Section 2 describes the underlying theory of G-causality and of the various extensions represented in the toolbox. Section 3 discusses the statistical assumptions involved in G-causality analysis and describes various preprocessing and validation techniques. Section 4 introduces some network-level descriptors of G-causality interaction patterns which can be useful for inferring macroscopic dynamical properties of neural systems. Section 5 discusses the important issue of filtering. Section 6 then outlines modality-specific issues involved in application of G-causality to data acquired from functional MRI (fMRI), neuroelectric and neuromagnetic signals, and spike train data. Limitations of the present instantiation of the toolbox are noted in Section 7, and Section 8 sets G-causality in the context of other frameworks for identifying directed functional connectivity. For detailed descriptions of the individual functions in the toolbox the reader is referred to the accompanying manual, available online at www.anilseth.com.