Exploring transient transfer entropy based on a group-wise ICA decomposition of EEG data

NeuroImage 49 (2010) 1593–1600

Contents lists available at ScienceDirect

NeuroImage j o u r n a l h o m e p a g e : w w w. e l s e v i e r. c o m / l o c a t e / y n i m g

Exploring transient transfer entropy based on a group-wise ICA decomposition of EEG data Vasily A. Vakorin a,⁎, Natasa Kovacevic a, Anthony R. McIntosh a,b a b

Rotman Research Institute of Baycrest, Canada Department of Psychology, University of Toronto, Canada

a r t i c l e

i n f o

Article history: Received 30 July 2008 Revised 17 June 2009 Accepted 11 August 2009 Available online 18 August 2009 Keywords: Information-theoretic approach Group-based independent component analysis Partial least square analysis Transfer entropy

a b s t r a c t This paper presents a data-driven pipeline for studying asymmetries in mutual interdependencies between distinct components of EEG signal. Due to volume conductance, estimating coherence between scalp electrodes may lead to spurious results. A group-based independent component analysis (ICA), which is conducted across all subjects and conditions simultaneously, is an alternative representation of the EEG measurements. Within this approach, the extracted components are independent in a global sense while short-lived or transient interdependencies may still be present between the components. In this paper, functional roles of the ICA components are specified through a partial least squares (PLS) analysis of task effects within the time course of the derived components. Functional integration is estimated within the information-theoretic approach using transfer entropy analysis based on asymmetries in mutual interdependencies of reconstructed phase dynamics. A secondary PLS analysis is performed to assess robust task-specific changes in transfer entropy estimates between functionally specific components. © 2009 Elsevier Inc. All rights reserved.

Introduction A fundamental problem in cognitive neuroscience is the question of how distributed brain areas work together to perform a specific task. Functional specialization and functional integration are two complementary principles established in the literature (Cohen and Tong, 2001). Considering the brain as a complex system, a broadly defined property characterizing a highly variable system with many parts whose behaviors strongly depend on the behavior of other parts, is one possible paradigm to study functional integration (Jirsa and McIntosh, 2007). In studies of complex systems, an analysis of the influence one system exerts over another is a part of the connectivity framework. It is important not only to reveal the existence of synchronized or temporally correlated systems but also to detect the directionality of prevailing information flows. Schreiber (2000) proposed a non-parametric measure to quantify the asymmetries in mutual, generally nonlinear interdependencies between coupled systems, designed as a deviation of transition probabilities from the Markov property. Chavez et al. (2003) applied this concept in an attempt to detect possible causal relationships in the epileptogenic networks. It can be shown that, under proper conditions, the transfer entropy is equivalent to the conditional mutual information (Palus and Vejmelka, 2007). The latter measure is the basis for the technique introduced by Palus et al. (2001) to

⁎ Corresponding author. E-mail address: [email protected] (V.A. Vakorin). 1053-8119/$ – see front matter © 2009 Elsevier Inc. All rights reserved. doi:10.1016/j.neuroimage.2009.08.027

infer directionality of coupling using information-theoretic tools. Such an approach does not need to specify a model of nonlinear interactions between signals. This attractive property of informationbased statistics has found a numerous number of applications, including those in neuroscience. Palus et al. (2001) applied their techniques in a case study of EEG recording of an epileptic patient. Hinrichs et al. (2006) estimated the direction of effective connectivity between neural activations in the brain damaged by stroke, using functional MRI. The work of Chavez et al. (2003) is an attempt to reveal active abnormal couplings in the epileptogenic networks with intracranial EEG. Unfortunately, EEG measurements do not directly represent localized brain regions in the vicinity of one electrode. Rather, due to volume conduction, the measured potentials reflect a summed signal from simultaneously active, underlying current sources (Nunez and Shrinivasan, 2005). When the signal passes the layers of cerebrospinal fluid, dura, scalp, and skull, it becomes filtered and spread out across electrodes. When calculating coherence in the space represented by electrodes, it is important to be certain that correlations between sensors are due to true physiological interactions between activated brain regions. In EEG, the potentials are measured with respect to another potential, and the amplitude of a reference electrode can cause spurious coherence estimates (Nunez et al., 1997; Fein et al., 1988). The same issue holds for estimating coherence in the phase space of the recorded signals (Guevara et al., 2005). Independent component analysis (ICA) is an alternative spatiotemporal representation of the EEG data. ICA blindly decomposes

1594

V.A. Vakorin et al. / NeuroImage 49 (2010) 1593–1600

multi-channel scalp EEG into a set of maximally independent components associated with locally coherent activity. First, singlesubject applications of ICA to EEG measurements (Makeig et al., 1996) can be traced to the times of developing early ICA algorithms. One of the optimization criteria for ICA is based on the Informax principle, which maximizes the mutual information between a set of signals extracted by the unmixing filters and the set of signal mixtures (Bell and Sejnowsky, 1995). The signal amplitude of a specific electrode at a particular time point is considered a realization of a random variable, not taking into account the order of data points within trials or the order of trials within subjects or conditions. The property of maximum information transmission leads to redundancy between signals at the output level, resulting in a set of ICA components. Extracting additive components with ICA at the level of individual subjects found numerous applications (to mention a few, see Makeig et al. (1997) and Makeig et al. (2004)). Onton et al. (2006) and Onton and Makeig (2006) provided comprehensive reviews on applying ICA to EEG data. To conduct normative research, in contrast to single-subject studies, there is a need to identify a set of similar and stable components across many participants involved in a study. The EEG literature basically identifies two approaches to address this question. The first strategy is to perform ICA, and thereafter look for similarities among independent components (IC) from different subjects in spite of differences in their projections on scalp data (Makeig et al., 2002). However, under this account, existence of reliable tools to cluster ICs across subjects remains an open issue. The second strategy to compare EEG data is to equate scalp locations, grouping channel activity across subjects. In spite of inhomogeneity in spatial contributions of source activities to a given scalp location, straightforward grouping could be a good alternative to avoid the problem related to the imprecise clustering of ICs. Calhoun et al. (2001) proposed a group ICA of functional MRI data to make group inferences. Kovacevic and McIntosh (2007) performed a group-based ICA of EEG data, exploiting similarities across individuals in cortical geometry to reveal robust spatiotemporal patterns. It should be noted that the derived components are independent in a global sense, across the entire time course for a given electrode, including trials for all subjects and conditions. However, short-lived dependencies between components may still be present. The idea of event-related coherence between EEG independent components can be traced back to early applications of ICA to EEG measurements (Makeig, 2000). For example, Delorme and Makeig (2003) reported the concomitant effects of μ-rhythm regulation on changes in phase coherence between maximally independent components accounting for posterior alpha rhythm. These types of interactions might also be referred to as transient effects. We use the term “short-lived” interdependencies to emphasize intermittent relationships which are ignored by the optimized model of global independence of the ICA components. We also use the term “transient” interdependencies to designate possible directed coherence during an interstimulus interval. This study uses the phase dynamics reconstructed from a single trial time course of the extracted group ICA components to explore asymmetries in transient transfer entropy between them. We analyzed data from a crossmodal cueing task where an auditory stimulus signaled the rule for a response to a visual stimulus (Diaconescu et al., 2008). A full description of the group ICA results from this study and their theoretical relevance are presented in the original manuscript. In the present case, we focus on the use of these data for a “proof-of-principle” demonstration of the combined use of group ICA, transfer entropy and partial least squares (PLS) to examine interdependencies. A PLS analysis (McIntosh and Lobaugh, 2004) was performed in the ICA component space to determine whether the component temporal expression varied as a function of experimental condition. In addition, a secondary PLS analysis was performed to

identify changes in the amounts of transfer entropy between the components, robustly expressed across conditions and subjects. Materials and methods Experiment EEG data were acquired using Neuroscan 4.0 with a 64-channel ElectroCap according to the standard 10/20 system. Impedances were kept below 5kΩ. EEG data were digitized at 250 Hz sampling rate and passed through a 0.01–100 Hz band-pass filter. During the recording, all electrodes were referenced to Cz but were re-referenced to an average reference for further analysis. More details on the acquisitions protocol can be found in the study of Alain et al. (2005). Each trial consisted of two stimuli (S1 and S2), each presented for 250 ms, separated by a 750 ms interstimulus interval (ISI) and followed by a response period. The time interval between the end of S2 and the beginning of S1 in the trial was 800 ms or 1200 ms (equiprobable). There were two tasks, depending on the relevance and modality of the stimuli. The visual stimulus was a square checkerboard, presented laterally, on the left or the right side of the screen. The auditory stimulus was either a low-pitch (250 Hz) or a high-pitch (4000 Hz) binaural pure tone. In the first task, S1 was an auditory cue that signaled a response rule to the upcoming S2, which was a lateralized visual target. The rule indicated a response in a manner which is either compatible or incompatible with the target presentation (button press with the hand on the same or opposite side of the target presentation). We refer to this task as the auditory– visual (AV) task. The second task was similar to the first, except for the reversed order of auditory cue and visual target. In this case, S1 was a lateralized visual target and S2 was an auditory cue. It is important to note that the cue remained relevant, and did not function as a distractor. In order to produce a correct response, participants had to remember the location of the target while waiting for the response rule signaled by the cue. We refer to this task as the visual–auditory (VA) task. Two sub-groups of subjects participated in the experiment: 11 (group A) and 12 (group B) healthy right-handed participants with normal to corrected-to-normal vision (10 males, mean age ± s.d., 21 ± 5). For the first sub-group, the high- and low-pitch tone cued subjects to produce a compatible and incompatible response, respectively. For the second sub-group, the response rule associated with the tone pitch was reversed, so that a high-pitch tone indicated incompatible response whereas a low-pitch tone indicated compatible response. Behaviorally, there were no differences between the groups in reaction to the response rule assignment (Diaconescu et al., 2008). For the purpose of estimating asymmetry in transient transfer entropy between functionally specific components in data decomposition, we combined the two sub-groups in one big group. Each task consisted of 2 blocks of 240 trials, randomized across 4 types depending on the target location (left or right) and cued response rule (compatible or incompatible). The four types of trials are denoted LC, LI, RC, and RI representing left-compatible, leftincompatible, right-compatible, and right-incompatible conditions, respectively. For details regarding preprocessing procedures see Kovacevic and McIntosh (2007). Briefly, the continuous data were band-pass filtered (0.5–45 Hz), epoched into 2.2 s long epochs timelocked to the onset of the first stimulus and 200 ms prestimulus baseline and artifact-corrected within the EEGLAB software using subject specific ICA decompositions (Delorme and Makeig, 2004). Group ICA A complete description of the group-wise ICA of EEG data can be found in Kovacevic and McIntosh (2007), which reports the group ICA

V.A. Vakorin et al. / NeuroImage 49 (2010) 1593–1600

based on the data in sub-group B. Here we will summarize the main points. Data from all subjects and conditions were combined into a single data set. The resulting data matrix X was k-by-m in size, where k is the number of channels (electrodes) and m is equal to the number of conditions times the number of subjects times the number of trials times the number of time points per trial. As a preprocessing step, spatial dimensionality of the matrix X was initially reduced with principal component analysis (PCA). ICA was subsequently performed to rotate the dimensionality reduced data by PCA into a new space spanned by temporally independent components. Regarding a criterion of independency, an Informax version of ICA was applied, as implemented in EEGLAB (Delorme and Makeig, 2004). Bayesian information criterion (BIC) was used to determine an optimal number of ICA components (Hansen et al., 2001). Bayesian probability was estimated for models with different numbers of components, and the model with maximum probability was selected for further analysis. For a given model with p components, the PCA filtering and subsequent ICA decomposition are given by the weighting matrix of size p-by-k. The matrix expressed each component as a weighted sum of all the electrodes. Single trial time series of the derived components were obtained by multiplying the electrode single trial time series by the weighting matrix.

1595

In estimating the measure D, I(Yτ, X|Y) and I(Xτ, Y|X) can be averaged over a selected range of lags τ (Palus et al., 2001), namely

DðxYyÞ =

bIðYτ ; XjYÞNτ − bIðXτ ; Y j XÞNτ bIðYτ ; XjYÞNτ + bIðXτ ; Y j XÞNτ

ð4Þ

where b…Nτ denotes averaging over τ. The averaging procedure decreases the variance of estimation error and can avoid the problem in detecting spurious high values of interdependency. In this study, we estimated transient transfer entropy between the group-wise ICA components during the interstimulus interval. The Hilbert transform was applied to extract phase dynamics from the component time courses (for more on phase estimation see Rosenblum et al. (2004)). Considering the phase dynamics as realizations of the processes X and Y, the individual and joint entropies in Eq. (2) were estimated through a generalized histogram method (Daub et al., 2004). The transfer information I(Yτ, X|Y) and I(Xτ, Y|X) was estimated for each lag τ = {1, 2…, 15} measured in data points, and then averaged. The estimated transfer information was used to calculate the statistic D according to Eq. (3). PLS analysis

Transfer entropy analysis Information-theoretic tools have been recently developed to characterize coupled systems. Within this approach, conditional mutual information is a key concept used to define asymmetries in mutual coupling (Palus et al., 2001). Using the idea of finite-order Markov processes, Schreiber (2000) proposed a measure of the coherence between systems evolving in time. The measure termed transfer entropy was based on approximate conditioning of transition probabilities. It can be shown that under proper conditions the transfer entropy is equivalent to the conditional mutual information (Palus and Vejmelka, 2007). Let X = {x1…xn} and Y = {y1…yn} be the realizations of two processes representing two interacting systems. Also, let Xτ = {x1+τ… xn+τ} and Yτ = {y1+τ…yn+τ} be the τ-future of the processes X and Y, respectively. Mutual information I(Xτ, Y) between Xτ and Y reflects the average information contained in the future of X about Y. Specifically, it can be expressed as IðXτ ; YÞ = H ðXτ Þ + H ðYÞ − HðXτ ; YÞ

ð1Þ

where H(Xτ) and H(Y) are the individual entropies of the processes, and H(Xτ, Y) is the entropy of their joint distribution. However, I(Xτ, Y) may also contain information contained in the future of X about the process X itself. The net information contained in the τfuture of X about the process Y can be refined by estimating the conditional mutual information I(Xτ, Y|X) between Xτ and Y given X as IðXτ ; Y jXÞ = H ðXτ ; XÞ + H ðY; XÞ − H ðXτ ; Y; XÞ − HðXτ Þ

ð2Þ

In a similar way, I(Yτ, X|Y) reflects the amount of information contained in the future of the signal Y about the signal X, provided that the information about Y itself is excluded. An asymmetry in the two measures, I(Xτ, Y|X) and I(Yτ, X|Y) allows one to estimate the directed coherence between two processes, X and Y. The interdependency statistic can be defined as the normalized difference between the two amounts of transfer information, namely DðxYyÞ =

IðYτ ; Xj YÞ − I ðXτ ; Y j XÞ IðYτ ; Xj YÞ + I ðXτ ; Y j XÞ

As a note, this statistic is antisymmetric: D(x → y) = − D(y → x).

ð3Þ

Partial least squares (PLS) analysis is a technique to produce a set of orthogonal factors (latent variables) to model the covariance structure between two matrices representing two sets of data (e.g., EEG and constrast, or imaging data and behavioral data). Here we give a brief description of the technique, and refer the reader to Hay et al. (2002); Lobaugh et al. (2001); Dzel et al. (2003); and McIntosh and Lobaugh (2004) for more details. A non-rotated version of PLS analysis is based on the projection of a set of a priori contrasts on the data available for analysis, for example, the original imaging data or its derivative such as connectivity estimates. The latent variables are chosen to provide a maximum correlation with dependent variables. A latent variable consists of three components. It includes: (a) a strength of differences represented by a singular value, which can be evaluated by computing the percentage of the cross-block covariance explained; (b) design saliences specifying a task contrast and indicating which tasks have different data-contrast correlations; and (c) electrode saliences specifying which time points and electrodes represent the optimal spatiotemporal relation to the identified contrast. Statistical assessment in PLS is based on resampling procedures. First, randomly permuting conditions within subjects, the significance of the singular value for each latent variable is assessed, by checking if the latent variable is statistically different from noise. The PLS results are recalculated using a relatively large number of permutations. A measure of significance is calculated as the number of times the permuted singular value is higher than the observed singular value. If the latent variable is found significant, it is subsequently tested for stability across subjects, identifying time points of robust contrast expression. The robustness is assessed through a bootstrap procedure, using resampling of subjects with replacement within conditions. A measure of stability is calculated as the ratio of the salience to the standard error of the generated bootstrap distribution. Thresholding the absolute bootstrap ratios at some level is used to identify the time points of robust contrast expression. In this paper we employed a non-rotated version of PLS analysis which is based on the covariance between a set of a priori contrasts and imaging data and its derivatives, specifically ICA component data and connectivity estimates. The significance of contrast tested and pvalue was estimated using permutation tests with 500 permutations. Spatiotemporal patterns of contrast expressions were assessed using bootstrap resampling with 500 samples. Points of robust contrast expression were identified by thresholding the absolute bootstrap

1596

V.A. Vakorin et al. / NeuroImage 49 (2010) 1593–1600

ratios to values N3.5, corresponding roughly to a 99% confidence interval.

Pipeline of the analysis A scheme shown in Fig. 1 illustrates the data-driven pipeline of the analysis performed in this paper. Two groups of subjects were involved in this study. For the complete description of experiments see the Experiment section. First, the data from both subject groups were combined and subsequently decomposed into additive components through the group-wise ICA (see Group ICA) to extract robust common spatiotemporal activity patterns across two groups of subjects. The great amount of specificity across the ICA components allowed one to suggest separate functional roles associated with each component. As an intermediary step, we focused on the interpretation of the group ICA components.

Manipulating the contrast based on the cue–target order, functional roles of the derived ICA components were specified using a non-rotated PLS analysis (see PLS analysis). The dominant effect found in two components laid the basis for the subsequent analysis aimed to compare differences in connectivity between the ICA components across conditions. Two components, which we call the evoked visual and auditory components, were identified according to their time course and its modulations in relation to the experimental manipulations. Specifically, an a priori orthogonal contrast was designed to compare two tasks: task 1 (AV) characterized by the auditory cue first–visual target second order, which includes the LC1, LI1, RC1, and RI1 experimental subconditions, and task 2 (VA) for visual target first–auditory cue second, which includes the LC2, LI2, RC2, and RI2 subconditions. The contrast was tested using a non-rotated task PLS analysis performed to collate the waveforms of the ICA components with respect to cue–target order, thus providing a pair of temporal patterns to characterize the populations. Focusing on differences between tasks AV and VA, conditions LC, LI, RC, and RI

Fig. 1. A scheme illustrating the data-driven pipeline of the analysis performed in this paper. First, a group-wise independent component analysis (ICA) was performed to decompose the group data into robust spatiotemporal activity patterns across subject groups. A partial least squares (PLS) analysis was used to specify functional roles of the derived ICA components. Transfer entropy was used as a measure of transient (short-lived) interdependencies between the components. Subsequently, a secondary PLS analysis was carried out at the level of transfer entropy estimates in search for robust differences across two types of task: auditory cue first followed by visual target, and visual target first followed by auditory cue. The primary conclusion is related to changes in asymmetry in short-lived interdependency between the components during the interstimulus interval, as a function of the sensory modality of the first stimulus.

V.A. Vakorin et al. / NeuroImage 49 (2010) 1593–1600

within the two tasks were merged for the analysis. In other words, a task-specific time course was obtained by averaging over 4 conditionspecific time courses for each subject, AV and VA. Subsequently, the transfer entropy analysis (see Transfer entropy analysis) was performed to analyze asymmetries in mutual dependencies between the waveforms of the group ICA components during interstimulus interval. Only data from the interstimulus interval were employed in the local transfer entropy analysis, specifying a common feature of four subconditions within each task, when response is not yet known. Transfer entropy was estimated for all single trials, and then averaged over trials. Let M denote the number of the ICA components chosen for the analysis. Also, let a pair of lower sub-indices {i,j}, where i,j = 1,...,M, identify the direction of dominant transfer entropy, leading from the

1597

component i to j. Regardless of the strategy applied to analyze the transfer entropy, based on averaging the lags (4) or for individual lags τ according to (3), the outcome of the entropy analysis is the M × M s;k connectivity matrices Ds;k = fα gi;j = 1; N ;M derived per subject s and per condition k. By construction, these matrices are antisymmetric for any given pair of s and k, with the diagonal elements equal to 0. The antisymmetry is related to the fact that reversing the order of how two components are put into analysis (i → j instead of j → i) should s;k lead to the sign change in the measure D, namely α s;k i;j = − α j;i for any i ≠ j. Second, the methods should not detect any difference in transfer entropy between two identical signals, that is, α s;k i;j = 0 when i = j. Finally, the mutual interdependencies between the components were tested for possible differences across conditions using a nonrotated PLS analysis (see PLS analysis), performed at the level of the

Fig. 2. PLS analysis performed to test temporal specificity of the group-based ICA components. The upper part of the figure illustrates the contrast designed to collate the waveforms of the group components in two tasks: auditory cue first (AV) and visual target first (VA). The results of the PLS analysis are arranged in rows associated with the 8 group-wise ICA components and columns related to the component's topography and waveforms. The topographic maps reflect the spatial distribution of the group ICA component signals. Plots in the second and fourth columns depict grand mean component waveforms for two groups of conditions, AV and VA. Black markers at the bottom of each plot indicate the time points of stable difference between the two tasks.

1598

V.A. Vakorin et al. / NeuroImage 49 (2010) 1593–1600

transfer entropy estimates. The transfer entropy for AV and VA tasks was analyzed in two ways: on the basis of lag-averaged data as specified in (4), and separately over a range of lag τ according to (3). The rationale for averaging across lags τ is to decrease the variance of the asymmetry statistic (Palus et al., 2001). If the data are noisy, or the transfer entropy estimation is noise sensitive, or in other words, if the distribution of flow across lags across subjects is highly variable across adjacent lags, then the average may provide a misleading result. Addressing the distribution of information flow across lags, a set of PLS analyses was performed separately for lags τ = 1,…, 15. Similar to the PLS analysis based on lag-averaging, asymmetry in transfer entropy between the components was tested with respect to the cue–target order: the conditions LC1, LI1, RC1 (task AV), and RI1 versus LC2, LI2, RC2, and RI2 (task VA). As before, our focus was task contrast, AV vs VA. Results To capture the task-related effects common to all subjects, a group ICA was performed as described in the Group ICA section. Bayesian information criterion (BIC) yielded an optimal number of ICA components equal to 8. Combined data sets from two sessions, which differ in coding the compatible and incompatible conditions with the high- and low-pitch tones, were consequently decomposed into the eight maximally independent spatial components. Fig. 2 shows the results of the non-rotated task PLS analysis performed to test functional specificity of the independent components. The top of the figure illustrates the contrast designed to compare the component waveforms: AV vs VA tasks. The contrast was statistically significant with the p-value less than 0.001. The first and third columns in Fig. 2 represent the topography of 8 group ICA components survived from BIC analysis. The topographic maps are shown as colormap images of the interpolated partial inverse of the weighting matrices that define the distribution of the component signals across all the electrodes. The interpolated scalp topographic maps were generated based on the electrode locations using spherical spline interpolation. Plots in the second and fourth columns are the corresponding grand average component waveforms. The waveforms were averaged across subjects and trials separately for two tasks differing in cue–target order (AV versus VA). The final part of the PLS analysis is a bootstrapping procedure that determines which components and at which latencies are expressing a given contrast stably across subjects. The stable points are indicated by black markers at the bottom of each plot, specifying periods of stable difference between two types of conditions. We found that the second and fifth ICA components, IC2 and IC5, differentiated strongly between visual and auditory stimulus presentation, respectively. As can be seen in Fig. 2, the topography of IC2 and IC5 resembles typical early ERP responses elicited by visual and auditory responses, respectively. Furthermore, the time course of IC2 exhibits strong evoked visual response following the onset of a visual target (intervals [1000 ms 1300 ms] in AV and [0 300 ms] in VA), observed as a strong increase in activation time-locked to the onset of presenting the visual stimulus/target. Its spatial distribution characterized by the topographic map presents an occipital positivity and is consistent with primary visual response P1. Similarly, IC5 exhibits strong evoked response time-locked to the onset of auditory stimulus: [0 ms 300 ms] in AV and [1000 ms 1300 ms] in VA. The topography of this component is consistent with primary auditory response P1. The results of the PLS analysis performed to compare localized transfer entropy between the ICA components across different conditions are illustrated in Fig. 3. The figure is based on averaging the statistic D across lags τ as defined in (4). Grand average transfer entropy matrices, D = αi,j, specifying asymmetry in localized transfer entropy between the 8 ICA components are shown separately for two tasks, AV (auditory cue first) and VA (visual target first) in the upper part (A). The PLS analysis revealed the statistically significant effect of

Fig. 3. PLS analysis performed to analyze changes in transient transfer entropy between the single trial waveforms of the group ICA components as a reaction to changes in the cue–target order. Transfer entropy is calculated for all single trials, and thereafter the connectivity is averaged across trials. The statistic for estimating asymmetry in transfer entropy was averaged across lags τ according to (4). The upper part (a) illustrates grand average causal matrices specifying local asymmetry in transfer entropy between the 8 ICA components, shown separately for two tasks, AV (auditory cue first) and VA (visual target first). The second and fifth components are associated with visual and auditory activation, respectively (see Fig. 2c). The lower part (b) depicts bootstrap ratio map. The matrices are illustrated as colormap images with the legends given on the right. The bootstrap ratio map indicates contribution of a specific pair of components across two types of conditions (AV vs VA). Positive values directly support the tested contrast, indicating the asymmetries in transfer entropy having a positive relation to the contrast due to the cue–target order: task 1 (AV) N task 2 (VA). Negative bootstrap ratio values have a negative relation to the contrast: task 1 (AV) b task 2 (VA).

the contrast ‘task 1 (AV) versus task 2(VA)’ with p-value less than 0.001. Fig. 3b illustrates the bootstrap ratio map, Q = {qi,j}. The matrix Q is depicted as a colormap image. Each element of the bootstrap ratio matrix determines the contribution of a specific connection between two components to the tested contrast. Specifically, each matrix value {qi,j}, as specified by the row i and column j, is associated with transfer entropy between the component i = 1,…, 8 and j = 1,…, 8. Positive values of the bootstrap ratio, {qi,j} N 0, directly support the contrast under investigation, task 1 (AV) versus task 2 (VA), indicating a subset of all possible combinations between components, wherein the changes in localized transfer entropy during the interstimulus interval

V.A. Vakorin et al. / NeuroImage 49 (2010) 1593–1600

due to cue–target order have a positive relation to the contrast design (task 1 N task 2). In other words, on average (across subjects), the estimates of transfer entropy interdependency αi,j, between the ith and jth components for an experiment with task 1 (AV) are higher than for an experiment with task 2 (VA). Similarly, negative bootstrap values, {qi,j} b 0, may also support the same contrast but in the reverse order. Specifically, negative bootstrap ratio values indicate a pair of components wherein differences in the transfer entropy have a negative relation to the contrast design (task 1 b task 2). Addressing the issue of the lags as specified in the Pipeline of the analysis section, we complemented the lag-averaged PLS analysis by performing a set of PLS analyses to compare asymmetry in mutual dependencies for the AV versus VA tasks separately at lags = 1,…, 15 data points. In general, the obtained results (see Fig. 4) are in accordance with previous findings (Fig. 3b), demonstrating that the tested contrast is significant, and the absolute value of the bootstrap ratio statistic for the pair of IC2 and IC5 is high. The peak in the bootstrap ratio is reached at the lag around 8–10, which corresponds approximately to 35 ms. Looking into the asymmetry in transfer entropy in AV and VA tasks, one can notice that the lag-averaged asymmetry in transfer entropy for AV tasks is negative basically due to lags around 7–12 data points, whereas the positive value of lagaveraged asymmetry in VA experiments is explained by first lags. In addition, one can see that the statistic is a smoothed function of lag τ, which supports the analysis based on averaging across lags. The key result from the PLS analysis is related to the contribution of specific pairs of components to the tested contrast within subject groups. As can be seen from Fig. 3b, the bootstrap ratio values for the asymmetry in transfer entropy between IC5, which was previously found to be associated with the evoked auditory activation, and all other components are positive: {q5,j} N 0 (red tinge for positive). This would support the idea that the auditory component has more

Fig. 4. Asymmetry in transient transfer entropy between the single trial waveforms of the group ICA components as a reaction to changes in the cue–target order: (a) subjectaveraged index D specifying local asymmetry in transfer entropy between components IC2 and IC5 for the AV task; (b) the same subject-averaged index D for the VA task; (c) bootstrap ratio values from PLS analyses performed separately at all lags τ to compare the statistic D between two tasks, AV and VA. All the statistics are shown as function of lag τ = 1,…, 15, used to construct future time series, as specified in (2).

1599

predictive power with respect to all other components in a cue first experiment than in an experiment where the cue follows the target. In other words, we are observing the relative predominant effect of the auditory component, originating from the auditory nature of the cue presented first. Further, the bootstrap ratio values for IC2, which is associated with the evoked visual activation, and the rest of the components are negative: {q2,j} b 0, as can be seen in Fig. 3b in blue. Taking into account that negative bootstrap values have a negative relation to the contrast under investigation, we can conclude that, similar to the auditory component IC5, the visual component IC2 has a relatively dominant, predictive power among the components in experiments with visual target first. In particular, it is worth noting that the most significant contribution to the contrast “task 1 (AV) versus task 2 (VA)” comes from the changes in transfer entropy between IC2 and IC5. Discussion A number of studies have recently introduced an analysis of interdependencies in ICA space with application to neuroimaging data. Londei et al. (2005) applied a transfer entropy based analysis using spatial ICA for functional MRI data. The time course of the spatially independent components was extracted and was subject to Grander causality test based on a linear regression model. Our study used temporal ICA as the basis for exploring transient changes in transfer entropy in a space spanned by the ICA components. It is crucial to realize that the derived group-based ICA components were temporally independent in a global sense over the entire time source of individual components, leaving the possibility of the existence of short-lived dependencies. In this study, we analyzed the effects of cue-driven transient transfer entropy during the interstimulus interval between functional modules represented by the waveforms of the ICA components with specific functional roles. The first PLS analysis was performed to highlight the functional representation of the derived components. To do this, we compared the averaged single trial waveforms of the derived components, analyzing task effects based on the order of visual–auditory stimuli. As specified by the experimental tasks, that is, visual stimulus first and auditory stimulus second, or vice versa, we focused our interest on two components: IC2 and IC5. These two components were found to be related to the activation of visual and auditory systems, respectively. Next, we estimated asymmetries in transfer entropy on the interstimulus interval. Subsequently, a secondary PLS was carried out to compare the patterns of transfer entropy interdependencies between the two tasks. A more complete exposition of the group ICA data is presented in Diaconescu et al. (2008), wherein the crossmodal stimulus–response compatibility paradigm was used to identify modality-independent aspects of rule processing and cued response facilitation. In the present case, we present an application of a data-driven analysis to examine interdependencies rather than to test specific predictions. In effect, this is a “proof-of-principle” paper emphasizing the method. The group ICA in essence is a non-trivial decomposition of neuronal activity. The analysis depends on a number of parameters and assumptions which may affect the final results. For example, concatenating the data from multiple subjects works on the assumption that the electrode placements and projections of sources to the surface are the same across subjects. Another issue is the model selection problem associated with estimating the numbers of components for PCA/ICA analyses. The BIC criterion is one of the techniques for estimating an asymptotic approximation of the posterior distribution of model given observed data. However, such approximation is known to have several shortcomings (Friel and Pettitt, 2008), and might not be optimal in selecting an appropriate model. Performing a non-linear analysis on top of the decomposed data, not without uncertainty with regards to their functional

1600

V.A. Vakorin et al. / NeuroImage 49 (2010) 1593–1600

specificity, makes the interpretation of the final results even less trivial. A reasonable way to show the validity of a complex data-driven pipeline is to demonstrate its performance in relatively simple, functionally unambiguous situations. That was the primary goal of this study. It is pleasing that the results based on the analysis performed in this study do follow what one may expect from changes in cue and target modality, as we discuss in the paper. We observed statistically robust changes in the transfer entropy across subjects during the interstimulus interval due to the stimulus order, specifying an increased predictive power of auditory IC5 and visual IC2 components with auditory cue first (AV) and visual target first (VA) tasks, respectively. It is interesting to look into how IC2 and IC5 are related to other components. Although statistically insignificant, the asymmetries in transfer entropy analysis revealed the increased influence of the components IC2 and IC5 in relation to all other components in the experiments with the corresponding stimulus order. Summing up, the sensory modality of the first stimulus turned out to be the contributing factor in estimating cuedriven effects of transient transfer entropy between functionally specific elements of the decomposed data. References Alain, C., Reinke, K., He, Y., Wang, C., Lobaugh, N., 2005. Hearing two things at once: neurophysiological indices of speech segregation and identification. J. Cogn. Neurosci. 17 (5), 811–818. Bell, A.J., Sejnowsky, T.J., 1995. An information-maximization approach to blind separation and blind deconvolution. Neural Comput. 7 (6), 1129–1159. Calhoun, V.D., Adali, T., Pearlson, G.D., Pekar, J.J., 2001. A method for making group inferences from functional MRI data using independent component analysis. Hum. Brain Mapp. 14 (3), 140–151. Chavez, M., Martinerie, J., Le Van Quyen, M., 2003. Statistical assessment of nonlinear causality: application to epileptic EEG signals. J. Neurosci. Methods 124 (2), 113–128. Cohen, J.D., Tong, F., 2001. Neuroscience: the face of controversy. Science 293 (5539), 2405–2407. Daub, C.O., Steuer, R., Selbig, J., Kloska, S., 2004. Estimating mutual information using b-spline functions an improved similarity measure for analysing gene data. BMC Bioinformatics 5 (18), 1–12. Delorme, A., Makeig, S., 2003. EEG changes accompanying learned regulation of 12-Hz EEG activity. IEEE Trans. Neural Syst. Rehabil. Eng. 11 (2), 133–137. Delorme, A., Makeig, S., 2004. EEGLAB: an open source toolbox for analysis of singletrial EEG dynamics including independent component analysis. J. Neurosci. Methods 134 (1), 9–21. Diaconescu, A.O., Kovacevic, N., McIntosh, A.R., 2008. Modality-independent processes in cued motor preparation revealed by cortical potentials. NeuroImage 42 (3), 1255–1265. Dzel, E., Habib, R., Schott, B., Schoenfeld, A., Lobaugh, N., McIntosh, A.R., Scholz , M., Heinze, H.J., 2003. A multivariate, spatiotemporal analysis of electromagnetic timefrequency data of recognition memory. NeuroImage 18 (2), 185–197.

Fein, G., Raz, J., Brown, F.F., Merrin, E.L., 1988. Common reference coherence data are confounded by power and phase effects. Electroencephalogr. Clin. Neurophysiol. 69 (6), 581–584. Friel, N., Pettitt, A.N., 2008. Marginal likelihood estimation via power posteriors. J. R. Stat. Soc., Ser. B 70, 589–607. Guevara, R., Velazquez, J.L., Nenadovic, V., Wennberg, R., Senjanovic, G., Dominguez, L.G., 2005. Phase synchronization measurements using electroencephalographic recordings: what can we really say about neuronal synchrony? Neuroinformatics 3 (4), 301–314. Hansen, L.K., Larsen, J., Kolenda, T., 2001. Blind detection of independent dynamic component. IEEE Int. Conf. Acoust. Speech Signal Process 5, 3197–3200. Hay, J.F., Kane, K.A., West, R., Alain, C., 2002. Event-related neural activity associated with habit and recollection. Neuropsychologia 40 (3), 260–270. Hinrichs, H., Heinze, H.J., Schoenfeld, M.A., 2006. Causal visual interactions as revealed by an information theoretic measure and fMRI. NeuroImage 31 (3), 1051–1060. Jirsa, V.K., McIntosh, A.R. (Eds.), 2007. Handbook of Brain Connectivity. Springer-Verlag, Berlin Heidelberg. Kovacevic, N., McIntosh, A.R., 2007. Groupwise independent component decomposition of EEG data and partial least square analysis. NeuroImage 35 (3), 1103–1112. Lobaugh, N.J., West, R., McIntosh, A.R., 2001. Spatiotemporal analysis of experimental differences in event-related potential data with partial least squares. Psychophysiology 38 (3), 517–530. Londei, A., D’Ausilio, A., Basso, D., Belardinelli, M.O., 2005. A new method for detecting causality in fMRI data of cognitive processing. Cogn. Process. 7 (1), 42–52. Makeig, S., Transient event-related coherence between independent EEG components. Society for Psychophysiological Research annual meeting, 2000. Makeig, S., Bell, A.J., Ghahremani, D., Sejnowski, T.J., 1996. Independent component analysis of electroencephalographic data. Adv. Neural Inf. Process Syst. 8, 145–151. Makeig, S., Jung, T.P., Bell, A.J., Ghahremani, D., Sejnowski, T.J., 1997. Blind separation of auditory event-related brain responses into independent components. Proc. Natl. Acad. Sci. U. S. A. 94 (20), 10979–10984. Makeig, S., Westerfield, M., Jung, T.P., Enghoff, S., Townsend, J., Courchesne, E., Sejnowski, T.J., 2002. Dynamic brain sources of visual evoked responses. Science 295 (5555), 690–694. Makeig, S., Debener, S., Onton, J., Delorme, A., 2004. Mining event-related brain dynamics. Trends Cogn. Sci. 8 (5), 204–210. McIntosh, A.R., Lobaugh, N.J., 2004. Partial least squares analysis of neuroimaging data: applications and advances. NeuroImage 23 (Suppl. 1), S250–S263. Nunez, P.L., Shrinivasan, R., 2005. Electric fields in the brain: the neurophysics of EEG. Oxford University Press. Nunez, P.L., Srinivasan, R., Westdorp, A.F., Wijesinghe, R.S., Tucker, D.M., Silberstein, R.B., Cadusch, P.J., 1997. EEG coherency. I: statistics, reference electrode, volume conduction, Laplacians, cortical imaging, and interpretation at multiple scales. Electroencephalogr. Clin. Neurophysiol. 103, 499–515. Onton, J., Makeig, S., 2006. Information-based modeling of event-related brain dynamics. Prog. Brain Res. 159, 99–120. Onton, J., Westerfield, M., Townsend, J., Makeig, S., 2006. Imaging human EEG dynamics using independent component analysis. Neurosci. Biobehav. Rev. 30 (6), 808–822. Palus, M., Vejmelka, M., 2007. Directionality from coupling between bivariate time series: how to avoid false causalities and missed connections. Phys. Rev., E 75, 056211. Palus, M., Komarek, V., Hrncir, Z., Sterbova, K., 2001. Synchronization as adjustment of infomation rates: detection from bivariate time series. Phys. Rev., E 63, 046211. Rosenblum, M., Pikovsky, A., Kurths, J., 2004. Synchronization approach to analysis of biological systems. Fluctuation Noise Lett. 4 (1), L53–L62. Schreiber, T., 2000. Measuring information transfer. Phys. Rev. Lett. 85 (2), 461–464.