## Combining independent component analysis and Granger causality to investigate brain network dynamics with fNIRS measurements |

Biomedical Optics Express, Vol. 4, Issue 11, pp. 2629-2643 (2013)

http://dx.doi.org/10.1364/BOE.4.002629

Acrobat PDF (3963 KB)

### Abstract

In this study a new strategy that combines Granger causality mapping (GCM) and independent component analysis (ICA) is proposed to reveal complex neural network dynamics underlying cognitive processes using functional near infrared spectroscopy (fNIRS) measurements. The GCM-ICA algorithm implements the following two procedures: (i) extraction of the region of interests (ROIs) of cortical activations by ICA, and (ii) estimation of the direct causal influences in local brain networks using Granger causality among voxels of ROIs. Our results show that the use of GCM in conjunction with ICA is able to effectively identify the directional brain network dynamics in time-frequency domain based on fNIRS recordings.

© 2013 Optical Society of America

## 1. Introduction

1. F. F. Jöbsis, “Noninvasive, infrared monitoring of cerebral and myocardial oxygen sufficiency and circulatory parameters,” Science **198**(4323), 1264–1267 (1977). [CrossRef] [PubMed]

9. C. B. Akgül, A. Akin, and B. Sankur, “Extraction of cognitive activity-related waveforms from functional near-infrared spectroscopy signals,” Med. Biol. Eng. Comput. **44**(11), 945–958 (2006). [CrossRef] [PubMed]

*HbO*

_{2}) and deoxyhemoglobin (

*HbR*), which plays an important role in the

*in vivo*study of cognitive processing in the human brain. Now advances in fNIRS are undergoing a transition from mapping sites of cortical activations towards identifying the brain networks that connect these sites together into dynamic systems in time-frequency domain. However, the algorithms for quantifying brain networks with fNIRS measurements are mainly confined to the temporal cross-correlation and the coherence spectrum analysis [7

7. R. C. Mesquita, M. A. Franceschini, and D. A. Boas, “Resting state functional connectivity of the whole head with near-infrared spectroscopy,” Biomed. Opt. Express **1**(1), 324–336 (2010). [CrossRef] [PubMed]

10. J. Cui, L. Xu, S. L. Bressler, M. Z. Ding, and H. L. Liang, “BSMART: A MATLAB/C toolbox for analysis of multichannel neural time series,” Neural Netw. **21**(8), 1094–1104 (2008). [CrossRef] [PubMed]

11. M. Ding, S. L. Bressler, W. Yang, and H. Liang, “Short-window spectral analysis of cortical event-related potentials by adaptive multivariate autoregressive modeling: data preprocessing, model validation, and variability assessment,” Biol. Cybern. **83**(1), 35–45 (2000). [CrossRef] [PubMed]

## 2. Methods

### 2.1 The Theory Model for fNIRS

1. F. F. Jöbsis, “Noninvasive, infrared monitoring of cerebral and myocardial oxygen sufficiency and circulatory parameters,” Science **198**(4323), 1264–1267 (1977). [CrossRef] [PubMed]

*HbO*and

_{2}*HbR*at time

*t*with wavelength

*λ,*in which

*OD*is the optical density as determined from the negative log ratio of the detected intensity of light with respect to the incident intensity of light using CW measurements,

*r*,

*DPF*(

*r*) is the unitless differential path length factor,

*l*(r) (mm) is the distance between the source and the detector,

*i*th chromophore at wavelength

*λ*of laser source,

*HbO*and

_{2}*HbR*is written asin which

**(**

*Q**r*,

*t*) vectors are the product of the inversion matrix of the extinction coefficients and the optical density change vectors. Similar operational procedures could be extended to

*n*th wavelength case based on regularization methods:where the matrix

**is the extinction coefficient matrix and**

*E***is defined as the a priori estimate of the covariance of the measurement error. The change of total hemoglobin concentration**

*R*### 2.2 ICA

13. A. J. Bell and T. J. Sejnowski, “An information-maximization approach to blind separation and blind deconvolution,” Neural Comput. **7**(6), 1129–1159 (1995). [CrossRef] [PubMed]

13. A. J. Bell and T. J. Sejnowski, “An information-maximization approach to blind separation and blind deconvolution,” Neural Comput. **7**(6), 1129–1159 (1995). [CrossRef] [PubMed]

*K*channels with

*M*dimensional random time vector for each channel, the measurement matrix for

**is then decomposed by ICA, estimating the optimal inverse of the mixing matrix**

Δ H b O 2 **, and a set of source time courses**

*A***. In terms of ICA estimation,**

*S*12. O. Demirci, M. C. Stevens, N. C. Andreasen, A. Michael, J. Liu, T. White, G. D. Pearlson, V. P. Clark, and V. D. Calhoun, “Investigation of relationships between fMRI brain networks in the spectral domain using ICA and Granger causality reveals distinct differences between schizophrenia patients and healthy controls,” Neuroimage **46**(2), 419–431 (2009). [CrossRef] [PubMed]

**is the**

*A**K-by-N*mixing matrix,

*N*is the number of unmixed sources and

**are the**

*S**N*-by-

*M*component time courses. Typically we utilize

*K*≥

*N*so that

**is usually of full rank. The goal of ICA is to estimate an unmixing matrix**

*A*^{}

**is a good approximation to the ‘true’ neural sources**

*X***and the inversion of weighted matrix**

*S***will extract the brain activity maps of the unmixed spatial sources. The infomax ICA algorithm uses the gradient ascent iteration algorithm to compute**

*W***by maximizing the entropy of the output of a single layer neural network. The resulting updated equation for the algorithm to calculate**

*W***is,**

*W**t*represents a given approximation step, is a general function that specifies the sizes of the steps for the unmixing matrix updates(usually an exponential function or a constant),

**is the identity matrix,**

*I**T*is the transposition operator, and is the nonlinearity in the neural network.

### 2.3 GCM Method

*X*reduces the prediction error of

_{2}*X*in a linear regression model of

_{1}*X*and

_{1}*X*, as compared to a model which includes only previous observations of

_{2}*X*,

_{1}*X*will cause

_{2}*X*. To illustrate G-C, suppose that the temporal dynamics of two time series

_{1}*X*(

_{1}*t*) and

*X*(

_{2}*t*) can be described by a bivariate autoregressive model [16

16. A. K. Seth, “A MATLAB toolbox for Granger causal connectivity analysis,” J. Neurosci. Methods **186**(2), 262–273 (2010). [CrossRef] [PubMed]

*p*is the maximum number of lagged observations included in the model,

**defines the coefficients of the model, and**

*A**E*

_{1}and

*E*

_{2}are the variances for each time series. If the variance of

*E*

_{1}(or

*E*

_{2}) is reduced by the inclusion of the

*X*(or

_{2}*X*) terms in Eq. (8) (or Eq. (9)), we define

_{1}*X*(or

_{2}*X*) causes

_{1}*X*(or

_{1}*X*). Assuming that

_{2}*X*

_{1}and

*X*

_{2}are covariance stationarity, the magnitude of this interaction can be measured by the log ratio of the prediction error variances for the restricted (

*R*) and unrestricted (

*U*) models [16

16. A. K. Seth, “A MATLAB toolbox for Granger causal connectivity analysis,” J. Neurosci. Methods **186**(2), 262–273 (2010). [CrossRef] [PubMed]

*E*

_{1}

_{R}_{(12)}is derived from the model omitting the

*A*

_{12,}

*(for all*

_{j}*j*) coefficients in Eq. (8) and

*E*

_{1}

*is derived from the full model in Eq. (8).*

_{U}*X*

_{2}on

*X*

_{1}is described in the context of multiple additional variables

*X*

_{3}. . .

*X*

_{n}. In this case,

*X*

_{2}causes

*X*

_{1}if knowing

*X*

_{2}reduces the variance in

*X*

_{1}’s prediction error when all other variables

*X*

_{3}. . .

*X*

_{n}are also included in the regression model.

*E*

_{i}

*are estimated from the autoregressive model including all variables. For*

_{U}*n*variables there are

*n*restricted models, with each restricted model omitting a different predictor variable. For example, the noise covariance matrix of the restricted model omitting variable 2, iswhere all

*E*

_{i}

*are estimated from the autoregressive model omitting variable 2. The G-C from variable 2 to variable 1, conditioned on variable 3, is given by*

_{R}**A**

*matrices and utilizes the Yule-Walker equations of the model to calculate the covariance matrix*

_{m}**of the noise vector**

*V***(t) [10**

*E*10. J. Cui, L. Xu, S. L. Bressler, M. Z. Ding, and H. L. Liang, “BSMART: A MATLAB/C toolbox for analysis of multichannel neural time series,” Neural Netw. **21**(8), 1094–1104 (2008). [CrossRef] [PubMed]

*j*to channel

*i*is given by for a specific frequency

*f*[10

10. J. Cui, L. Xu, S. L. Bressler, M. Z. Ding, and H. L. Liang, “BSMART: A MATLAB/C toolbox for analysis of multichannel neural time series,” Neural Netw. **21**(8), 1094–1104 (2008). [CrossRef] [PubMed]

### 2.4 Behavior Tasks and fNIRS Recordings

*DPF*(

*r*)

*=*4. During the task periods, subject was instructed to perform a finger flexion and extension action repeatedly.

*HbO*distributions, runs were imaged in (bottom-to-top) order of their occurrence during the experiment. For example, Fig. 3 shows their distributions for an arbitrary selected channel 36 that is close to the motor cortex [5].

_{2}/HbR/HbT## 3. Results

*HbO*measurements (without filtering or pre-processing). To display the contribution of independent components (ICs) to the

_{2}*HbO*recordings, we plot in Fig. 4(a) the components that contribute the most to

_{2}*HbO*measurements in terms of the most

_{2}*HbO*variance of all the 48 ICs. We also capture the components that contribute the most to the raw

_{2}*HbR and HbT*measurements, which are provided in Fig. 4(b) and Fig. 4(c), respectively. It should be noted in Figs. 4(a)-4(c) that the black thick lines indicate the data envelope (i.e. minimum and maximum of all channels at every time point) and the colored ones show the

*HbO*components with the largest contributions to the measurements.

_{2}/HbR/HbT*HbO*, six nodes for

_{2}*HbR*, and five nodes for

*HbT*, which are displayed on the fourth column of Figs. 5(a)-5(c), respectively. The nodes within the ROIs are considered as the primary cortex areas directly activated by the stimuli and will be employed for network analysis with GCM method. The VAR model order utilized in this study is 10, the window length are 80 time points and the frequency interval for fast Fourier transformation is between 1 and 100Hz. To investigate the brain networks underlying right finger tapping, the

*HbO*,

_{2}*HbR*and

*HbT*signals from all the nodes located within the ROIs are processed to calculate the time-frequency causal influences.

*HbO*is provided in Fig. 6, in which Fig. 6(a) in the left column plots the time-frequency Granger causality influences among different nodes and Fig. 6 (b) in the right column shows the single spectra at 25 seconds (stimuli onset time). Likewise the computed time-frequency Granger causality influences for the

_{2}*HbR*and

*HbT*are displayed in Fig. 7 and Fig. 8, respectively.

*P <*0.05. Based on this, the Granger casual connectivity networks are reconstructed for arbitrarily selected time windows center at 25s (stimuli onset time) and 35s (stimulus processing time) for the results in Figs. 6-8, which are plotted in Fig. 9 for the three chromophores. The frequency selected for constructing the brain networks is 3 Hz, in which the brain activation patterns show significant correlations with the right finger tapping tasks, as displayed in Figs. 6-8. The edges of the brain networks are labeled with blue lines when they are estimated from

*HbO*signals, where the arrow represents the direction of Granger causality influences between two channels while the width of the blue lines represents the strength of Granger causality influences.

_{2}/ HbR/ HbT## 4. Discussion

*HbR*, in which strong causal influences are identified during the onset time around 25s while weak causality influences are revealed during the stimulus processing, as displayed in Fig. 7.

*HbO*networks in the first column of Fig. 9 that strong causal influences mainly occur in the left motor cortex (LMC) including Left primary motor cortex and left primary somatosensory cortex, which validate that stimuli in the right finger tapping tasks will yield the activations in LMC. In addition, we found in the first column of Fig. 9(a) that Granger causality influences are also observed from parietal cortex (PC) lobe region to LMC and from LMC to the supplementary motor area (SMA) during the onset time. These findings are consistent with the organization of motor cortex in its early stages. Interestingly during the stimulus processing, besides the strong Granger causality influences within LMC, increased causality influences among LMC, SMA and PC are observed from the second network of

_{2}*HbO*at 35s, which demonstrates motor cortex information flows further to the motor association cortex such as PC. Further, we found during the stimulus processing SMA also exerts Granger causality influences on LMC. The generated effective networks for

_{2}*HbO*correlate well with the structural networks and anatomical pathways in Fig. 10, in which the activation areas of motor stimuli are mainly located in the SMA and LMC for the right finger tapping tasks. In addition, brain activity is also observed in the motor association cortex such as PC during the stimulus processing, as displayed in Figs. 9 and 10. Interestingly we didn’t identify any neural activity from premotor cortex which generally involves the imaging of movement, but not the real motor stimuli such as finger tapping. The structural networks for finger tapping tasks could be found from the website (http://neuroscience.uth.tmc.edu/s3/chapter03.html).

_{2}*HbR*networks in the second column Fig. 9(a) that during the onset time strong causal influences mainly occur in LMC, but causal influences from LMC to SMA and PC are also identified. Moreover, we can see from the second column of Fig. 9(b) during the stimulus processing the influence within LMC is significantly reduced while influences from LMC to SMA are increased for

*HbR*. In particular, the right primary motor cortex (RMC) also exerts influences on the LMC. These features are quite different from the findings in

*HbO*networks.

_{2}*HbT*networks in the third column of Fig. 9 with those from

*HbO*and

_{2}*HbR*, we found during both onset and stimulus processing, LMC exerts strong influences on SMC. And during the stimulus processing, strong Granger causality influence from LMC to PC is confirmed as well, as shown in the third column of Fig. 9(b).

12. O. Demirci, M. C. Stevens, N. C. Andreasen, A. Michael, J. Liu, T. White, G. D. Pearlson, V. P. Clark, and V. D. Calhoun, “Investigation of relationships between fMRI brain networks in the spectral domain using ICA and Granger causality reveals distinct differences between schizophrenia patients and healthy controls,” Neuroimage **46**(2), 419–431 (2009). [CrossRef] [PubMed]

17. W. Liu, L. Forrester, and J. Whitall, “A note on time-frequency analysis of finger tapping,” J. Mot. Behav. **38**(1), 18–28 (2006). [CrossRef] [PubMed]

*HbO*and

_{2}, HbR*HbT*measurements though they show different activation patterns and network dynamics in temporal domain. It should be noted here that the ICs that reveal body movement and physiological noise such as ICs 16 and 44 in Fig. 11 are not involving the GCM analysis since only identified ROIs are employed for network analysis.

**21**(8), 1094–1104 (2008). [CrossRef] [PubMed]

## Acknowledgments

## References and links

1. | F. F. Jöbsis, “Noninvasive, infrared monitoring of cerebral and myocardial oxygen sufficiency and circulatory parameters,” Science |

2. | B. Yodh and B. Chance, “Spectroscopy and imaging with diffusing light,” Phys. Today |

3. | J. C. Ye, S. Tak, K. E. Jang, J. Jung, and J. Jang, “NIRS-SPM: statistical parametric mapping for near-infrared spectroscopy,” Neuroimage |

4. | Z. Yuan, Q. Zhang, E. S. Sobel, and H. Jiang, “Image-guided optical spectroscopy in diagnosis of osteoarthritis: A clinical study,” Biomed. Opt. Express |

5. | Z. Yuan, “Spatiotemporal and time-frequency analysis of functional near-infrared spectroscopy brain signals using independent component analysis,” J. Biomed. Opt. |

6. | I. Schelkanova and V. Toronov, “Independent component analysis of broadband near-infrared spectroscopy data acquired on adult human head,” Biomed. Opt. Express |

7. | R. C. Mesquita, M. A. Franceschini, and D. A. Boas, “Resting state functional connectivity of the whole head with near-infrared spectroscopy,” Biomed. Opt. Express |

8. | T. J. Huppert, S. G. Diamond, M. A. Franceschini, and D. A. Boas, “Homer: a review of time-series analysis methods for near-infrared spectroscopy of the brain,” Appl. Opt. |

9. | C. B. Akgül, A. Akin, and B. Sankur, “Extraction of cognitive activity-related waveforms from functional near-infrared spectroscopy signals,” Med. Biol. Eng. Comput. |

10. | J. Cui, L. Xu, S. L. Bressler, M. Z. Ding, and H. L. Liang, “BSMART: A MATLAB/C toolbox for analysis of multichannel neural time series,” Neural Netw. |

11. | M. Ding, S. L. Bressler, W. Yang, and H. Liang, “Short-window spectral analysis of cortical event-related potentials by adaptive multivariate autoregressive modeling: data preprocessing, model validation, and variability assessment,” Biol. Cybern. |

12. | O. Demirci, M. C. Stevens, N. C. Andreasen, A. Michael, J. Liu, T. White, G. D. Pearlson, V. P. Clark, and V. D. Calhoun, “Investigation of relationships between fMRI brain networks in the spectral domain using ICA and Granger causality reveals distinct differences between schizophrenia patients and healthy controls,” Neuroimage |

13. | A. J. Bell and T. J. Sejnowski, “An information-maximization approach to blind separation and blind deconvolution,” Neural Comput. |

14. | A. Delorme and S. Makeig, “EEGLAB: an open source toolbox for analysis of single-trial EEG dynamics including independent component analysis,” J. Neurosci. Methods |

15. | Z. Yuan, J. Ye, “Fusion of fNIRS and fMRI data: Identifying when and where hemodynamic signal are changing in human brains,” Front. Hum. Neurosci., doi: 10.3389/ fnhum.2013.00676(2013). |

16. | A. K. Seth, “A MATLAB toolbox for Granger causal connectivity analysis,” J. Neurosci. Methods |

17. | W. Liu, L. Forrester, and J. Whitall, “A note on time-frequency analysis of finger tapping,” J. Mot. Behav. |

**OCIS Codes**

(170.0170) Medical optics and biotechnology : Medical optics and biotechnology

(300.0300) Spectroscopy : Spectroscopy

(100.4996) Image processing : Pattern recognition, neural networks

**ToC Category:**

Image Processing

**History**

Original Manuscript: September 16, 2013

Revised Manuscript: October 20, 2013

Manuscript Accepted: October 21, 2013

Published: October 25, 2013

**Citation**

Zhen Yuan, "Combining independent component analysis and Granger causality to investigate brain network dynamics with fNIRS measurements," Biomed. Opt. Express **4**, 2629-2643 (2013)

http://www.opticsinfobase.org/boe/abstract.cfm?URI=boe-4-11-2629

Sort: Year | Journal | Reset

### References

- F. F. Jöbsis, “Noninvasive, infrared monitoring of cerebral and myocardial oxygen sufficiency and circulatory parameters,” Science198(4323), 1264–1267 (1977). [CrossRef] [PubMed]
- B. Yodh and B. Chance, “Spectroscopy and imaging with diffusing light,” Phys. Today48(3), 34–40 (1995). [CrossRef]
- J. C. Ye, S. Tak, K. E. Jang, J. Jung, and J. Jang, “NIRS-SPM: statistical parametric mapping for near-infrared spectroscopy,” Neuroimage44(2), 428–447 (2009). [CrossRef] [PubMed]
- Z. Yuan, Q. Zhang, E. S. Sobel, and H. Jiang, “Image-guided optical spectroscopy in diagnosis of osteoarthritis: A clinical study,” Biomed. Opt. Express1(1), 74–86 (2010). [CrossRef] [PubMed]
- Z. Yuan, “Spatiotemporal and time-frequency analysis of functional near-infrared spectroscopy brain signals using independent component analysis,” J. Biomed. Opt.18(10), 16011 (2013).
- I. Schelkanova and V. Toronov, “Independent component analysis of broadband near-infrared spectroscopy data acquired on adult human head,” Biomed. Opt. Express3(1), 64–74 (2012). [CrossRef] [PubMed]
- R. C. Mesquita, M. A. Franceschini, and D. A. Boas, “Resting state functional connectivity of the whole head with near-infrared spectroscopy,” Biomed. Opt. Express1(1), 324–336 (2010). [CrossRef] [PubMed]
- T. J. Huppert, S. G. Diamond, M. A. Franceschini, and D. A. Boas, “Homer: a review of time-series analysis methods for near-infrared spectroscopy of the brain,” Appl. Opt.48(10), D280–D298 (2009). [CrossRef] [PubMed]
- C. B. Akgül, A. Akin, and B. Sankur, “Extraction of cognitive activity-related waveforms from functional near-infrared spectroscopy signals,” Med. Biol. Eng. Comput.44(11), 945–958 (2006). [CrossRef] [PubMed]
- J. Cui, L. Xu, S. L. Bressler, M. Z. Ding, and H. L. Liang, “BSMART: A MATLAB/C toolbox for analysis of multichannel neural time series,” Neural Netw.21(8), 1094–1104 (2008). [CrossRef] [PubMed]
- M. Ding, S. L. Bressler, W. Yang, and H. Liang, “Short-window spectral analysis of cortical event-related potentials by adaptive multivariate autoregressive modeling: data preprocessing, model validation, and variability assessment,” Biol. Cybern.83(1), 35–45 (2000). [CrossRef] [PubMed]
- O. Demirci, M. C. Stevens, N. C. Andreasen, A. Michael, J. Liu, T. White, G. D. Pearlson, V. P. Clark, and V. D. Calhoun, “Investigation of relationships between fMRI brain networks in the spectral domain using ICA and Granger causality reveals distinct differences between schizophrenia patients and healthy controls,” Neuroimage46(2), 419–431 (2009). [CrossRef] [PubMed]
- A. J. Bell and T. J. Sejnowski, “An information-maximization approach to blind separation and blind deconvolution,” Neural Comput.7(6), 1129–1159 (1995). [CrossRef] [PubMed]
- A. Delorme and S. Makeig, “EEGLAB: an open source toolbox for analysis of single-trial EEG dynamics including independent component analysis,” J. Neurosci. Methods134(1), 9–21 (2004). [CrossRef] [PubMed]
- Z. Yuan, J. Ye, “Fusion of fNIRS and fMRI data: Identifying when and where hemodynamic signal are changing in human brains,” Front. Hum. Neurosci., doi: 10.3389/ fnhum.2013.00676(2013).
- A. K. Seth, “A MATLAB toolbox for Granger causal connectivity analysis,” J. Neurosci. Methods186(2), 262–273 (2010). [CrossRef] [PubMed]
- W. Liu, L. Forrester, and J. Whitall, “A note on time-frequency analysis of finger tapping,” J. Mot. Behav.38(1), 18–28 (2006). [CrossRef] [PubMed]

## Cited By |
Alert me when this paper is cited |

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.

« Previous Article | Next Article »

OSA is a member of CrossRef.