An Optimized Method for Bayesian Connectivity Change Point Model

Abstract

The brain undergoes functional dynamic changes at all times. Investigating functional dynamics has been recently verified to be helpful for detecting psychological conditions and powerful for analyzing disease-related abnormalities of the brain. This article aims to detect functional dynamics. Specifically, we focus on how to effectively distinguish corresponding functional connectivity and change points from functional magnetic resonance imaging (fMRI) data. By combining Bayesian connectivity change point model (BCCPM), a modified genetic algorithm (GA) is presented to optimize the evolutionary procedure toward the most probable distributions of real change points in fMRI. We randomly initialize different binary indicator vectors to represent different distributions of change points. Each indicator vector represents an individual in GA, and together they form an initial population. Then we calculate Bayesian posterior probability and use it as the fitness of each individual. Finally, we evolve individuals of current generation toward the next higher fitness generation by a series of modified genetic operators. After several evolutionary procedures, individuals in the final generation may have outstanding fitness and the one with highest fitness can represent the most likely change point distribution in the corresponding fMRI data. Furthermore, the most probable change point distribution could be resolved. We test the optimized method for BCCPM on several synthesized data sets, and the experimental results verify that the proposed model produces higher accuracy results with lower time consumption. Also, we apply the new model to real block-designed task-based fMRI data set and excellent results are obtained.

1. Introduction

Understanding functional localizations and exploring functional interactions within the brain are ongoing challenges in the area of neuroscience (Passingham et al., 2002; Lv et al., 2015). Among several methods, neuroimaging is an efficient way to achieve this task. Functional magnetic resonance imaging (fMRI) is a functional neuroimaging method (Wikipedia, 2017). By quantifying blood flow using MRI technology, fMRI data can be used to measure human brain activities (Friston et al., 1998; Logothetis et al., 2001; Wikipedia, 2017).

In recent years, multiple neuroscience researches on neuronal network-level activities using fMRI data set have invoked increasing number of attentions (Zhang et al., 2014; Ou et al., 2015b). Modeling functional connectivity and abrupt boundaries among regions of interest (ROIs) in fMRI data has been generally considered as a powerful way to investigate brain functional interactions (Malonek and Grinvald, 1996; Ou et al., 2013, 2014; Zhang et al., 2014).

Inspired by the successes of signal analysis technologies, sliding time window, short-time Fourier transform, and wavelet transform-based schemes are proposed to detect brain functional dynamics in recent literature (Friston et al., 1994; Bassett et al., 2006; Allen et al., 2014; Zhang et al., 2014). Unfortunately, these frameworks are limited due to some difficulties such as how to decide the length of the sliding window. Sparse representation is widely used in applications such as image processing, audio processing, and document analysis (Passingham et al., 2002; Mairal et al., 2008; Wright et al., 2009; Plumbley et al., 2010; Zhao et al., 2010; Lyu et al., 2016). Recently, there are several researches that utilize sparse representation for fMRI signal analysis, and some excellent results are achieved. However, sparse representation is based on the typical assumptions that the neural integration of those components is linear and the components of each voxel's fMRI signal are sparse (Passingham et al., 2002). In general, these basic assumptions might not be easy to reach for fMRI data analysis. From a technical perspective, several other very popular tools for signal analysis can also be developed for exploring brain functional interactions, such as Markov random field models, principal component analysis, independent component analysis, and autoregressive spatial models (Passingham et al., 2002; Ou et al., 2013, 2015a; Lyu et al., 2016).

In some recent studies, several Bayesian inference-based methods have been proposed for exploring brain functional dynamics. These methods are less sensitive to fMRI noises and more reliable in estimating functional interactions (Stevenson et al., 2009). Lian et al. (2014a) proposed a Bayesian connectivity change point model (BCCPM) to detect change points (defined as abrupt boundaries of functional interactions in brain networks) in fMRI. It uses Bayesian inference method to calculate the probability of multivariate time series while some certain time points are assumed as change points among fMRI data. It obtains significant performance. However, it is based on Markov Chain and Monte Carlo (MCMC) strategy and not easy to be speeded up with multiprocessor computers or graphics processing unit (GPU).

In this article, to develop a higher accuracy algorithm that can be parallelized in multiprocessor computers or GPU in the future realization, we investigate most Bayesian inference-based methods for exploring brain functional interactions (Stevenson et al., 2009; Lian et al., 2014a,–d; Zhang et al., 2014; Guo et al., 2016; Xiao et al., 2017). By combining BCCPM, we propose a modified genetic algorithm (GA) to optimize the evolutionary procedure to improve the detection accuracy and decrease the time consumption. We test GA-based BCCPM on several synthesized data sets, and the experimental results verify that our new model produces results with higher accuracy in less time. Also, we apply the new model to real block-designed task-based fMRI data set, and excellent results are obtained. In this section, we mainly describe the modified GA and BCCPM. The modified GA will be utilized as an optimization strategy and Bayesian inference theory will be used to calculate fitness of each individual in GA.

2. Modified GA and Bayesian Inference Theory

In this section, we mainly describe the modified GA and BCCPM. The modified GA will be utilized as an optimization strategy and Bayesian inference theory will be used to calculate fitness of each individual in GA.

2.1. Modified GA

As aforementioned, we need an optimization method to achieve better performance such as higher accuracy and lower time consumption. Moreover, the potential parallelized realization is an important factor we may concern.

GA is a kind of evolutionary method inspired by the process of natural selection (Liu et al., 2016; Pei et al., 2016). Commonly, we represent the possible solutions of a problem with indicator vectors called individuals. Then, we evolve the initial individuals to the next generation using selection, crossover, and mutation operators.

Actually, GA may utilize very different selection, crossover, and mutation strategies for different problems. To effectively explore brain functional interactions, we design a series of special selection, crossover, and mutation operators as follows.

Selection and crossover:

Step 1. Randomly produce an integer n ∈ 1∼n0.

Step 2. Randomly produce two different integers a1 and a2 ∈ 1∼n, and select the a1-th and a2-th individuals in the sorted i-th generation.

Step 3. Randomly produce two different integers b1 and b2 ∈ 1∼N as the selected positions.

Step 4. Crossover the selected individuals at the selected positions.

Mutation:

Step 5. Randomly produce a floating point number u ∈ [0,1].

Step 6. if u>u0, go to step 8.

Step 7. Randomly produce integers c1,c2,…,cs ∈ [0,N].

Step 8. Change the c1-th, c2-th,…, and cs-th positions from 1 to 0 or from 0 to 1.

Step 9. If all the individuals have been generated, stop; otherwise, go to Step 1.

Figure 1 illustrates the structure and flowchart of the modified GA. The main procedures of the algorithm in Figure 1 can be summarized as follows. First, we randomly initialize some binary indicator vectors to represent different distributions of change points. The dimensions of these indicator vectors are the same as the number of time points in the fMRI data. From the view of GA, these indicator vectors can be regarded as individuals of the initial population. We then calculate Bayesian posterior probability and use it as the fitness of everyone. Finally, we evolve individuals of the current generation toward the next higher fitness generation by a series of modified genetic operators, that is, special selection, special crossover, and special mutation strategies. After several evolutionary procedures, individuals in the final generation may have outstanding fitness, and the highest fitness individual can represent the most likely change point distribution in the corresponding fMRI data. Thus, the most probable change point distribution could be resolved.

FIG. 1.

Flowchart of the proposed GA. BCCPM, Bayesian connectivity change point model; fMRI, functional magnetic resonance imaging; GA, genetic algorithm.

Moreover, if we carefully inspect every step illustrated in Figure 1, we can see that the modified GA is different from a traditional algorithm. Besides typical selection, crossover, and mutation operators, it has a sorting operation and a duplication operation. These two operations exert very important influence on all other operators, and may directly improve the performance.

2.2. Bayesian connectivity change point model

In this subsection, we mainly introduce BCCPM, which was first proposed by Lian et al. (2014a) in the literature. We will use Bayesian inference theory to calculate fitness of each individual in GA.

Given an \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{R}} \times { \rm{T}}$$ \end{document} data set \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$X = ( {x_1} , {x_2} , \cdots , {x_T} )$$ \end{document} , in which T is the number of observations and R is the number of ROIs, we are interested in whether there are some differences in the joint probabilities within these ROIs between different time periods.

We define a block indicator vector as \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} \begin{align*} \vec I = ( {I_1} , {I_2} , \cdots , {I_T} ) , \tag{1} \end{align*} \end{document}

Now, suppose a set of vectors \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${x_1} , {x_2} , \cdots , {x_T}$$ \end{document} is i.i.d. (independent and identically distributed) from R-dimensional multivariate normal distribution, that is, \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{x}} \sim N ( \mu , \sum )$$ \end{document} t = 1, 2, …, T, where T denotes the number of vectors, R denotes the dimension of vectors \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${x_1} , {x_2} , \cdots , {x_T}$$ \end{document} , μ denotes the R-dimensional mean vector, and ∑ denotes the R × R covariance matrix. The conjugate prior distribution of (μ, ∑) is N-Inv-Wishart \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$( { \mu _{ \rm{0}}} , { \Lambda _0} / { \kappa _0} , { \nu _0} , { \Lambda _0} )$$ \end{document} (Gelman et al., 2014), and the posterior distribution of (μ, ∑) based on the data \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$X = ( {x_1} , {x_2} , \cdots , {x_T} )$$ \end{document} is also N-Inv-Wishart \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$( { \mu _{ \rm{T}}} , { \Lambda _{ \rm{T}}} / { \kappa _{ \rm{T}}} , { \nu _{ \rm{T}}} , { \Lambda _{\rm T}} )$$ \end{document} (Gelman et al., 2014). Therefore, we can calculate the probability of \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${x_1} , {x_2} , \cdots , {x_T}$$ \end{document} as follows: \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} \begin{align*} p (x_1,x_2,\ldots, x_{{\rm T}}) = \frac{p(x_1, x_2, \ldots, x_S; \mu, \sum)}{p (\mu, \sum \vert x_1, x_2, \ldots, x_{{\rm T}})} = (\frac{1}{2 \pi})^{{\rm R} T/2} (\frac{K_0}{K_T})^{R/2} \frac{\Gamma_{{\rm R} (\frac{V_T}{2})(\det({\Lambda_0})) {V_0}/2}}{\Gamma_{{\rm R}} (\frac{{{V_0}}}{2})(\det ({\Lambda_T})){V_T}/2}{2^{{{\rm R}}T/2}}, \tag { 2 } \end{align*} \end{document}

Consider the block indicator vector in Equation (1), the likelihood of the data matrix \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $$X = ( {x_1} , {x_2} , \cdots , {x_T} )$$ \end{document} is \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} \begin{align*} p ( X \vert { \vec I} ) = \Pi_{b = 1}^{\sum I_{{\rm k}}} p ( {X_b}), \tag{4} \end{align*} \end{document}

where \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${X_{ \rm{b}}}$$ \end{document} is the temporal observations that belong to b-th block and \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} $${ \rm{p}} ( {X_{ \rm{b}}} )$$ \end{document} can be calculated according to Equation (2). The temporal blocks are independent from each other; therefore, the posterior distribution of the configuration is \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\usepackage{upgreek}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document} \begin{align*} p ( \vec I \left\vert X \right. ) \propto p ( \vec I ) p ( X \left\vert { \vec I} \right. ) , \ { \rm where} \ p ( \vec I ) = \mathop \Pi \limits_{t = 1}^T p ( { \vec I_t} ) \ { \rm and} \ p ( { \vec I_t} ) \ { \rm is \ Bern} \ ( 0.5 ). \end{align*} \end{document}

It is worth noting that Equation (4) will be regarded as the fitness function of the proposed GA to calculate the fitness of every new individual generated by the evolutionary operators.

3. Simulation Study and Real Data Application

In this section, intensive simulation study is carried out to evaluate and validate MCMC-based BCCPM (Lian et al., 2014a) and GA-based BCCPM on several simulation data sets. Furthermore, the application results on a real emotion processing task-based fMRI data set are compared and summarized.

3.1. Simulation study

In the simulation study, eight data sets are generated based on eight different structures of dynamic networks (Fig. 2a–h) to testify that GA-based BCCPM can effectively detect the change points. There are various types of interaction patterns in the eight different structures of dynamic networks to make sure the applicability of the model. Multiple experiments are completed to validate the proposed model. The designed change points of eight different dynamic networks are shown in Figure 2 by using a vertical solid line with a position number below it.

FIG. 2.

Eight different structures of dynamic networks (a–h) and their change point distributions.

For each of the eight structures of dynamic networks, we repeat the experiment 50 times and all the results are saved to calculate their average performance. The number of repetitions of GA-based BCCPM is set as 100 and that of MCMC-based BCCPM is set as 20,000 to achieve good convergence and change point detection results for both methods. Note that all parameters of BCCPM are set at the same values for the sake of fairness.

The simulation results are summarized in Figures 3 and 4 to compare the convergence curves and change point detection results. In Figure 3, the convergence curves of the two methods are displayed on the left and right columns, respectively, from network (a) to network (h). We can observe that there are some vibrations within the highest peaks of the convergence curves in MCMC-based BCCPM. However, thanks to the advantage of the proposed modified GA, which copies the subjects with highest fitness directly to the next generation, GA-based BCCPM performs better as it always climbs to its peak monotonically. In Figure 4, the change point detection results of the two methods are shown on the left and right columns, respectively, from networks (a) to (h). We can observe that the GA-based BCCPM has better performance than the MCMC-based BCCPM as GA-based BCCPM detects all the designed change points precisely, whereas MCMC-based BCCPM makes some minor mistakes in networks (b), (c), (d), (f), and (h).

FIG. 3.

Convergence curves for networks (a–h). Left column shows results from MCMC-based BCCPM. Right column shows results from GA-based BCCPM. MCMC, Markov Chain and Monte Carlo.

FIG. 4.

Change point detection results for networks (a–h). Left column shows results from MCMC-based BCCPM. Right column shows results from GA-based BCCPM.

Moreover, the computational cost of each method is also taken into consideration. The running times (in milliseconds) of MCMC-based BCCPM and the proposed method are recorded and summarized in Table 1. We can easily see that the former has lower computational cost in all eight networks. Note that our computing environment is Windows 10 Pro; processor: Intel^® Core™ i7-4770K CPU @ 3.50 GHz; installed RAM: 16.0 GB; system type: 64-bit operating system, x64-based processor.

Table 1.

Comparison of Computational Cost (Time in Milliseconds)

	MCMC-based BCCPM	GA-based BCCPM
Network (a)	3572.2	3168.6
Network (b)	5284.6	3618.8
Network (c)	8687.4	4468.6
Network (d)	5000.2	2790.4
Network (e)	6009.4	3206.4
Network (f)	8784.6	4525.0
Network (g)	7890.2	5706.0
Network (h)	8218.8	5728.2

BCCPM, Bayesian connectivity change point model; GA, genetic algorithm; MCMC, Markov Chain and Monte Carlo.

3.2. Real data application

In this study, we applied the GA-based BCCPM on a real emotion processing task-based fMRI data and compared the results with MCMC-based BCCPM. These nine fMRI data sets are acquired from an emotion processing task on a Siemens Skyra 3T scanner (Zhu et al., 2013). The total scan length was 176 frames (i.e., time points). The data sets were preprocessed by Lian et al. (2014a) using the open source DICCCOL tools in Human Connectome Project (2013). From each of the scanned subject's brain, fMRI time series are extracted from 358 DICCCOL ROIs. Details of the data processing can be found in Human Connectome Project (2013) and Lian et al. (2014a). The resulting size of each data set is 358*176, where 358 is the total number of ROIs and 176 is the total number of frames (i.e., time points). The application results are summarized in Table 2. Similar to Lian et al. (2014a), we define change points that are detected within ±5 time points as aligned with boundaries, and change points that are detected within ±10 time points as partially aligned. From the results, we can see that our proposed GA-based BCCPM detects change points more precisely than MCMC-based BCCPM among most of the nine subjects. On average, out of the seven expected change points, MCMC-based BCCPM detects 75% (aligned and partially aligned) and GA-based BCCPM detects 81% (aligned and partially aligned).

Table 2.

Results on Real Data Application

Subject	Model	No. of change points detected	No. of change points aligned with boundaries	No. of change points partially aligned with boundaries	No. of change points not aligned with any boundaries
1	MCMC-based BCCPM	7	1	2	4
1	GA-based BCCPM	7	5	1	1
2	MCMC-based BCCPM	7	3	2	2
2	GA-based BCCPM	8	5	2	1
3	MCMC-based BCCPM	8	5	1	2
3	GA-based BCCPM	7	2	2	3
4	MCMC-based BCCPM	7	4	1	2
4	GA-based BCCPM	7	2	3	2
6	MCMC-based BCCPM	8	4	2	2
6	GA-based BCCPM	7	6	1	0
7	MCMC-based BCCPM	8	2	3	3
7	GA-based BCCPM	8	4	3	1
8	MCMC-based BCCPM	8	4	2	2
8	GA-based BCCPM	7	2	2	3
9	MCMC-based BCCPM	8	4	2	2
9	GA-based BCCPM	8	2	4	2
10	MCMC-based BCCPM	8	3	2	3
10	GA-based BCCPM	7	2	3	2

In addition, change points that are not aligned with any boundaries are also shown in Table 2. The average number of not aligned change points detected by MCMC-based BCCPM is 2.44, whereas that detected by GA-based BCCPM is 1.67. It also indicates that the GA-based BCCPM has better performance than MCMC-based BCCPM. Furthermore, these not aligned change points could be used to infer the dynamic changes happened within the designed blocks.

To better understand the functional interactions between the ROIs within the temporal blocks found by our method, the PC (named after its authors with their first names, Peter and Clark) algorithm (Spirtes and Glymour, 1991) is applied to search the patterns that can represent the underlying causal structure from the data. Since subject 6 has the most aligned change points detected by our method, we used subject 6 as an example to show the functional interaction patterns searched by the PC algorithm. In Figure 5, the top graph shows that the change points detected by GA-based BCCPM align with designed task blocks and only the last one is not aligned per our definition. From these temporal blocks, three of them are labeled as I, II, and III, which represent a face block, a shape block, and a fixation block in the experiment, respectively. The bottom row in Figure 5 shows the functional interaction patterns searched by the PC algorithm for I, face block; II, shape block; and III, fixation block. It is obvious that the interaction patterns are different from each other in the two task blocks and the fixation block.

FIG. 5.

Top: Detected change points shown in blue ridges with designed blocks split by red dotted lines of subject 6. Bottom: Functional interaction patterns searched by PC algorithm from three temporal blocks: I, face; II, shape; III, fixation.

4. Conclusion and Future Work

Finally, we conclude our method and discuss further work in this section. In this article, an optimized method for BCCPM is presented. We apply it to detect change points in synthesized and real fMRI data, and excellent results are obtained.

In the future, we will combine our method with other Bayesian inference models, such as BMCPM, DBVPM, and BCBCPM. Furthermore, for the potential parallel realization of the GA, GA-based BCCPM could be easily implemented in a parallel mode and run efficiently within GPU or multiprocessor computers, thus, we will also investigate its GPU or multiprocessor version.

Footnotes

Acknowledgments

This research is supported by the Brains-Behavior Seed grant and Molecular Basis of Disease (MBD) from Georgia State University.

Author Disclosure Statement

No competing financial interests exist.

References

Allen

E.A.

, Damaraju

, Plis

S.M.

, et al. 2014. Tracking whole-brain connectivity dynamics in the resting state. Cereb. Cortex. 24, 663–676.

Bassett

D.S.

, Meyer-Lindenberg

, Achard

, et al. 2006. Adaptive reconfiguration of fractal small-world human brain functional networks. Proc. Natl. Acad. Sci. 103, 19518–19523.

Friston

K.J.

, Fletcher

, Josephs

, et al. 1998. Event-related fMRI: Characterizing differential responses. Neuroimage. 7, 30–40.

Friston

K.J.

, Jezzard

, and Turner

1994. Analysis of functional MRI time-series. Hum. Brain Mapp. 1, 153–171.

Wikipedia. Functional magnetic resonance imaging. Available at: https://en.wikipedia.org/wiki/Functional_magnetic_resonance_imaging. Accessed June 1, 2017.

Gelman

, Carlin

J.B.

, Stern

H.S.

, et al. 2014. Bayesian Data Analysis, 2nd ed. Boca Raton, FL: Chapman and Hall/CRC.

Guo

, Liu

, Chen

, et al. 2016. Bayesian inference for functional dynamics exploring in fMRI data. Comput. Math. Methods Med. 2016, 9.

Human Connectome Project. 2013. HCP Protocols. Available at: http://protocols.humanconnectome.org/HCP/3T/task-fMRI-protocol-details.html. Accessed July 26, 2017.

Lian

, Li

, Xing

, et al. 2014a. Exploring functional brain dynamics via a Bayesian connectivity change point model. 2014 IEEE 11th International Symposium on Biomedical Imaging (ISBI), Beijing, China, 600–603.

10.

Lian

, Lv

, Xing

, et al. 2014b. Generalized fMRI activation detection via Bayesian magnitude change point model. 2014 IEEE 11th International Symposium on Biomedical Imaging (ISBI), Beijing, China, 21–24.

11.

Lian

, Li

, Young

, et al. 2014c. Dynamic network partition via Bayesian connectivity bi-partition change point model. 2014 IEEE 11th International Symposium on Biomedical Imaging (ISBI), Beijing, China, 545–548.

12.

Lian

, Li

, Zhang

, et al. 2014d. Detecting cell assembly interaction patterns via Bayesian based change-point detection and graph inference model. 2014 IEEE 11th International Symposium on Biomedical Imaging (ISBI), Beijing, China, 17–20.

13.

Liu

, Kong

, and Su

2016. An improved genetic algorithm based on the shortest path problem. 2016 IEEE International Conference on Information and Automation (ICIA), 328–332.

14.

Logothetis

N.K.

, Pauls

, Augath

, et al. 2001. Neurophysiological investigation of the basis of the fMRI signal. Nature. 412, 150.

15.

, Jiang

, Li

, et al. 2015. Assessing effects of prenatal alcohol exposure using group-wise sparse representation of fMRI data. Psychiatry Res. Neuroimaging. 233, 254–268.

16.

Lyu

, Li

, Lv

, et al. 2016. Identifying group-wise consistent sub-networks via spatial sparse representation of natural stimulus FMRI data. 2016 IEEE 13th International Symposium on Biomedical Imaging (ISBI), 62–65.

17.

Mairal

, Elad

, and Sapiro

2008. Sparse representation for color image restoration. IEEE Trans. Image Process. 17, 53–69.

18.

Malonek

, and Grinvald

1996. Interactions between electrical activity and cortical microcirculation revealed by imaging spectroscopy: Implications for functional brain mapping. Science. 272, 551.

19.

, Lian

, Xie

, et al. 2014. Atomic dynamic functional interaction patterns for characterization of ADHD. Hum. Brain Mapp. 35, 5262–5278.

20.

, Xie

, Jin

, et al. 2015a. Characterizing and differentiating brain state dynamics via hidden Markov models. Brain Topogr. 28, 666–679.

21.

, Xie

, Li

, et al. 2015b. Atomic connectomics signatures for characterization and differentiation of mild cognitive impairment. Brain Imaging Behav. 9, 663–677.

22.

, Xie

, Wang

, et al. 2013. Modeling brain functional dynamics via hidden Markov models. I 2013 6th International IEEE/EMBS Conference on Neural Engineering (NER), 569–572.

23.

Passingham

R.E.

, Stephan

K.E.

, and Kötter

2002. The anatomical basis of functional localization in the cortex. Nat. Rev. Neurosci. 3, 606–616.

24.

Pei

, Wang

, Cui

, et al. 2016. Codelet scheduling by genetic algorithm. Trustcom/BigDataSE/ISPA, 2016 IEEE, 1492–1499.

25.

Plumbley

M.D.

, Blumensath

, Daudet

, et al. 2010. Sparse representations in audio and music: From coding to source separation. Proc. IEEE. 98, 995–1005.

26.

Spirtes

, and Glymour

1991. An algorithm for fast recovery of sparse causal graphs. Soc. Sci. Comput. Rev. 9, 62–72.

27.

Stevenson

I.H.

, Rebesco

J.M.

, Hatsopoulos

N.G.

, et al. 2009. Bayesian inference of functional connectivity and network structure from spikes. IEEE Trans. Neural. Syst. Rehabil. Eng. 17, 203–213.

28.

Wright

, Yang

A.Y.

, Ganesh

, et al. 2009. Robust face recognition via sparse representation. IEEE Trans. Pattern Anal. Mach. Intell. 31, 210–227.

29.

Xiao

, Liu

, Zhang

, et al. 2017. Detecting change points in fMRI data via Bayesian inference and genetic algorithm model, 314–324. In Cai

, Daescu

, Li

, eds. Bioinformatics Research and Applications. ISBRA 2017. Lecture Notes in Computer Science, Vol. 10330. Springer, Cham.

30.

Zhang

, Li

, et al. 2014. Inferring functional interaction and transition patterns via dynamic Bayesian variable partition models. Hum. Brain Mapp. 35, 3314–3331.

31.

Zhao

, Li

, and Kwok

2010. Text detection in images using sparse representation with discriminative dictionaries. Image Vis Comput. 28, 1590–1599.

32.

Zhu

, Li

, Guo

, et al. 2013. DICCCOL: Dense individualized and common connectivity-based cortical landmarks. Cereb. Cortex. 23, 786–800.