# Dynamic Functional Connectivity Change-Point Detection With Random Matrix Theory Inference

^{1}Department of Statistics, Duksung Women's University, Seoul, South Korea^{2}College of Sungsim General Education, Youngsan University, Gyeongnam, South Korea^{3}Department of Neurosurgery, Seoul National University Hospital, Seoul, South Korea^{4}Department of Brain and Cognitive Sciences, Seoul National University, Seoul, South Korea

To study the dynamic nature of brain activity, functional magnetic resonance imaging (fMRI) data is useful including some temporal dependencies between the corresponding neural activity estimates. Recent studies have shown that the functional connectivity (FC) varies according to time and location which should be incorporated into the model. Modeling this dynamic FC (DFC) requires time-varying measures of spatial region of interest (ROI) sets. To know about the DFC, change-point detection in FC is of particular interest. In this paper, we propose a method of detecting a change-point based on the maximum of eigenvalues via random matrix theory (RMT). From covariance matrices for FC of all ROI's, the temporal change-point of FC is decided by an RMT approach. Simulation results show that our proposed method can detect meaningful FC change-points. We also illustrate the effectiveness of our FC detection approach by applying our method to epilepsy data where change-points detected are explained by the changes in memory capacity. Our study shows the possibility of RMT based approach in DFC change-point problem and in studying the complex dynamic pattern of functional brain interactions.

## 1. Introduction

Our brain is a complex network which consists of spatially distributed, but functionally linked regions that continuously share information with each other. Recent advanced functional neuroimaging studies have provided new tools to measure and to examine functional interactions between brain regions. These make the examination of DFC in the human brain and lead a development of theory associated with the complex FC structures. The fMRI data contains information about brain activity with complex spatial-temporal correlations. FC in the brain may be represented by these spatial dependence patterns and may reveal important characteristics of brain function and individual variations in cognition and behavior. FC is fundamentally a statistical concept, and is generally defined and assessed using statistical measures such as correlation (Biswall et al., 1995), and cross-coherence (Sun et al., 2004) in the temporal evolution of neural activity between distinct brain locations. Hutchison et al. (2013a,b) researched promise, issues, and interpretations as related with DFC.

Resting-state fMRI measures spontaneous fluctuations in the blood-oxygenation-level-dependent (BOLD) signal in gray matter regions. The fMRI connectivity approaches estimate statistical dependencies between gray matter activity arising from different regions (Biswall et al., 1995; Sato et al., 2006; Friston, 2011). For example, FC alterations are linked to depression, schizophrenia, and other illnesses (Mayberg et al., 1999; Jafria et al., 2008). Koike et al. (2011) provided that the default mode network (DMN) and subsystems might serve to integrate brain regions with performing function specific to each level of arousal.

Functional connectivity modelings are addressed with independent component analysis (ICA) (McKeown et al., 1998; Calhoun et al., 2001), pairwise correlation analysis (Supekar et al., 2008), and partial correlation analysis (Salvador et al., 2005). A key assumption of the aforementioned approaches is that functional connectivity remains stationary throughout the scanning session. However, many recent researches provide evidence of the dynamic nature of the brain's functional organization (Hellyer et al., 2014). Functional networks appear to fluctuate on a time scale and to reorganize under tasks during the scanning session (Chang and Glover, 2010).

Time-varying and dynamic nature of connectivity is considered in many recent neuroimaging studies. A well-known method for exploring time-varying connectivity is applying a sliding window of pre-specified length. This method calculates the correlations between distinct locations throughout the duration of the window at each step. Limitations with this sliding window approach, that the window length strongly influences the temporal smoothness or variability observed in the resulting correlations are discussed (Lindquist et al., 2014).

Identifying change-points is an important subject in neuroimaging study as it shows properties of brain networks as they relate to experimental stimuli and disease processes. Many change-point methods for fMRI signals have been proposed (Lindquist et al., 2007; Robinsona et al., 2010; Cribben et al., 2013; Xu and Lindquist, 2015). Xu and Lindquist (2015) assume that the timing of a subject's activation onset and duration are random variables, and their distributions are used to approximate the probability that a voxel/region is activated as a function of time with multi-subject change-point modeling. Bayesian approaches (Cheon and Kim, 2010; Kim and Cheon, 2010) and a non-parametric Fourier functional approach (Kim and Hart, 2011) are proposed to estimate change-points with distributional parameter changes which can be extended to the covariance change-point problems. Shao and Zhang (2010) study tests for change-point in time series, which is applicable to each individual fMRI series. Recently change-point detection for the multivariate covariance matrix are proposed for brain networks (Kim, 2019; Chenouri et al., 2020).

There have been many studies on testing for distributional structure changes for a sequence of independent and identically distributed (IID) random variables, including Csörgő and Horváth (1997). However test statistics developed for change-point detection in the IID context may not work with the time series observed in neuroimaging as they do not take into account the temporal dependence in the data. Aue et al. (2009) proposed break detection in covariance based on the *vech* operator to summarize information in the covariance matrix. Here *vech* is the vectorization operator which converts the matrix into a column vector. Aue and Horváth (2013) reviewed break points estimation in covariance structure in time series. Inclán and Tiao (1994) proposed variance change detection with cumulative sums for independent series. Petersen and Muller (2016) studied functional data covariance and applied to Alzheimer fMRI data. For detecting correlation change, Wied and Krömer (2012) considered the cumulated sums of empirical correlations to make a change test.

In the present study, we aim to deal with spatial time series data to detect FC change-points based on random matrix theory. Our approach does not make any parametric assumptions. It is especially designed for detecting one change-point but is extendable to multiple change-points in FC. In addition, we applied our developed method to fMRI data of post-surgical epilepsy patients. Previous studies have shown that activity and FC of the medial temporal lobe (MTL) during the post-encoding offline period are closely related to the memory consolidation (Tambini et al., 2010; Staresina et al., 2013; Tambini and Davachi, 2013). However, how the brain supports the memory consolidation process in the absence of one MTL structure has not been investigated. Here, we applied our developed method to the post-encoding rest period to evaluate its usefulness in characterization of memory consolidation networks in MTL resected brains.

This paper is organized as follows. In section 2, we begin by looking at the covariance change-point problem with the proposed statistics. Simulation studies are shown in section 3. Application of the proposed method to the experimental epilepsy data is shown in section 4. The paper concludes with a discussion.

## 2. Methods

### 2.1. Statistical Problems for Dynamic Functional Connectivity

A natural unit of fMRI analysis is multi-voxel regions of the brain that change their level of activation. The fMRI data can be widely used to detect and delineate regions of the brain that change in response to specific stimuli and tasks. We represent a single scan from a subject as a 3-D rectangular lattice consisting of volume elements (voxels). From one scanning session, each voxel contains serial measures of localized brain activity called blood-oxygen-level dependent (BOLD) fMRI responses. Through the hemodynamic response (HR) process, blood releases oxygen to neurons at a greater rate than to inactive neurons. This causes a change of the relative levels of oxyhemoglobin and deoxyhemoglobin that is detectable due to magnetic susceptibility. Since the number of intracranial voxels is too large to estimate a global spatial correlation matrix including all voxel pairs, we partition the voxels into mutually exclusive neuroanatomical regions and select ROIs in order to make the ROI-level inference.

Let *y*_{ivt} be an observation measured at voxel or region * v* from the serial fMRI BOLD responses for an

*th subject. Let*

**i***i*= 1, …,

*N*represent subjects,

*v*= 1, …,

*voxels (or ROIs), and*

**d***t*= 1, …,

*T*time points. For the

*i*th subject,

*d*×

*T*observed data matrix is obtained. It makes a multivariate time series for each subject and its covariance change is of interest. Tests that assess the structural stability of volatilities or cross-volatilities for multivariate non-linear time series such as fMRI data and change-point estimation is of great importance for understanding the underlying dynamics.

Let's fix the *i*th subject for the notational simplicity. Consider the observed **y**_{t} is a * d* × 1 random vector at the time point

*t*with $E\left[\Vert {y}_{t}{\Vert}^{2}\right]<\infty $. Here ‖ · ‖ denotes the Euclidean norm in ℝ

^{d}. Based on the observations of the random vector

**y**

_{1}, …,

**y**

_{T}, consider a test to identify covariance matrix change for each subject.

The null hypothesis is

which indicates the constancy of the covariances during the observational time period. Against this null hypothesis the alternative hypothesis with at most one change in the covariance is

where *τ* is the change-point. Let **Σ**_{t} = *Cov*(**y**_{t}).

### 2.2. Dynamic Functional Connectivity Test With the Largest Eigenvalues

We summarize the covariance structure using maximum eigenvalues, defined as the largest eigenvalues of the covariance matrix as a functional representative measure of the matrix.

Let λ_{1}(*t*) = max{λ_{1}(*t*), …, λ_{d}(*t*)} denote the largest eigenvalues of the matrix **Σ**_{t}. The null hypothesis in (1) can be considered as that

which indicates the constancy of the maximum eigenvalues of covariances during the observation time period. The alternative hypothesis can be written with the maximum eigenvalue measure as

where *τ* is denoted as the change-point.

Suppose that the data are from *d*−dimensional Gaussian distribution with mean zero and positive definite covariance matrix **Σ**. The sample covariance matrix (SCM) formed with these *T* samples as

has the (central) Wishart distribution (Goodman, 1963; Raj Rao et al., 2008).

Estimating the eigenvalues of a population covariance matrix from a sample covariance matrix is a problem of fundamental importance in multivariate statistics since the eigenvalues of covariance matrices play a key role in many widely used techniques with a suitable theory (Goodman, 1963; Okamoto, 1973; Gupta and Nagar, 1999; Mohsen, 2013). The amount of variance explained is measured by the eigenvalues of the population covariance, **Σ**_{d×d}. Random matrix theory predicts that the eigenvalues of a sample covariance matrix **S**_{T} are not good estimators of the eigenvalues of the population covariance when the sample size *T* and the number of variables *d* are large. El Karoui (2008) proposed a better estimator for the eigenvalues of large dimensional covariance matrices. Under the weak assumption that the Marčenko and Pastur (1967) equation holds, the El Karoui (2008) estimator can be thought as a shrinkage in a non-linear fashion.

Let λ_{1} ≥ λ_{2} ≥ ⋯ ≥ λ_{d} be the eigenvalues of the population covariance matrix. Accordingly let *l*_{1} ≥ *l*_{2} ≥ ⋯ ≥*l*_{d} be the eigenvalues of the sample covariance matrix. Anderson (1963) showed that the case where **Σ** = **I**_{d×d} with all the population eigenvalues one, under some moment conditions, the maximum eigenvalue is not a consistent estimator of λ_{1} and can be extremely biased when *T* and *d* are both large. To de-bias extreme sample eigenvalues, we use population spectral distribution, a probability measure that characterizes the population eigenvalues (El Karoui, 2007). The limiting distribution of the sample eigenvalues are dependent on data dependency and their covariance matrix. Johnstone (2001) derived the limiting distribution of the largest sample eigenvalue with independent Gaussian case and Bai and Silverstein (2010) dealt with this problem the sub-Gaussian case. In terms of the spectral decomposition of a covariance matrix, an estimator is rotation-equivariant if and only if the eigenvectors are the same as those of the sample covariance matrix (Mohsen, 2013). Thus, the differences between two such estimators appear only in their eigenvalues.

### 2.3. Dynamic Functional Connectivity Test With Tracy-Widom Distributions

Assume that **y**_{u}'s are independent from * d*−dimensional distribution,

*u*= 1,…,

*T*. Let the sample covariance matrix up to

*t*and after

*t*be

where

are the sample mean vectors separated by *t*. If **y**_{u}'s are independent from *d*−dimensional normal distribution, **A**_{t} and **B**_{t} are *d* × *d* matrices following Wishart distributions denoted by **A**_{t} ~ *W*_{d}(**Σ**, *t*) and **B**_{t} ~ *W*_{d}(**Σ**, *T* − *t*) respectively under H_{0}.

Consider that *t*, *T* − *t* ≥ *d*. The scale matrix **Σ** has no effect on the distribution of the eigenvalues, and so consider the largest eigenvalue λ_{1}(*t*) of ${({\text{A}}_{t}+{\text{B}}_{t})}^{-1}{\text{A}}_{t}$ as the greatest root statistic. If there exists change and change-point occurs at time point t, there is change in *A*_{t} and *B*_{t}. This change is reflected in the test statistic *G*_{t}.

Under *H*_{0} in (3), it holds that λ_{1}(*t**) = ⋯ = λ_{1}(*T* − *t**). For example, we can set *t** = 30 to have enough degrees of freedom. We consider the following statistic to provide information about change and change-points in the covariance matrix.

Under *H*_{0}, *E*[*G*_{t}] = 0. Since **A**_{t} and **B**_{t} are positive definite, 0 < *ϕ* < 1. Equivalently *ϕ* is the largest root of the determinantal equation

Johnstone (2009) provided the appropriate centering and scaling, and the logit transform *W*_{t} = *logit*(*ϕ*) = log(*ϕ*/(1 − *ϕ*)) is approximately Tracy-Widom distributed:

The distribution *F*_{1} was found by Tracy and Widom (1996) as the limiting law of the largest eigenvalue of *d*×*d* Gaussian symmetric matrix. The centering and scaling parameters are given by

where the angle parameters *φ* = *φ*_{t}, *γ* = *γ*_{t} are defined by

We apply Tracy and Widom scaling and centering for each maximum eigenvalues of ${({\text{A}}_{t}+{\text{B}}_{t})}^{-1}{\text{A}}_{t}$ and then $maxeigen\left[{({\text{A}}_{t}+{\text{B}}_{t})}^{-1}{\text{B}}_{t}\right]$ respectively.

Therefore our test is modified by Tracy and Widom correction to de-bias of the eigenvalues as follows

Here max *eigen*^{TW} denotes that the maximum eigenvalue is modified by Tracy and Widom correction. We propose the test statistic for change as

which rejects H_{0} if Λ_{T} is large. Accordingly the proposed change-point estimator is

Permutation tests have the flexible class of tests for both non-parametric and parametric (model-based). The proposed statistics can be used regardless of the data distributions using this permutation approach. The threshold calculation based on block permutation can be applied for testing for change-points (Strasser and Weber, 1999; Husková, 2004). Exchangeability of the errors might be too strong of an assumption for time series applications due to the dependence structure of the observations. Block permutation is an applicable alternative. Block permutation principles are suitable for testing autoregressive series and they are refined (Kirch and Steinebach, 2006; Kirch, 2007; Zeileis and Hothorn, 2013). Approximations of critical values of the change test statistics can be obtained through block permutation methods (Husková, 2004; Luger, 2006; Kirch, 2007).

To do our testing procedure, the approximate critical values are required. For each subject with dependent time series data, each critical value can be obtained according to the block permutation tests. This block permutation tests provide asymptotically valid results in the presence of dependency.

## 3. Simulations

### 3.1. Simulation Study

In this section we set up the simulation in order to assess and compare the performance of the proposed estimator for change-point estimation. The simulation study is focused on change-point estimation with multi-subject data from the multivariate vector autoregressive (MVAR) models. The objective of each simulation is to detect the functional connectivity change-point.

Before change-point estimation, testing for change is done. Since it is not possible to obtain the exact distribution of the test statistic Λ_{T} analytically, it is determined by simulation. The data are synthetically generated from multivariate normal distribution *N*(**0, Σ**_{d}) for the null distribution with the given dimension *d* as the number of ROIs, and the number of time points *T*. Several types of **Σ**_{d} are designed. For the comparison we also consider the Aue et al. (2009) method.

Aue et al. (2009) proposed the test for covariance change based on

Their change test statistic is

where ${\widehat{\Sigma}}_{T}$ is a consistent covariance estimator. Aue et al. (2009) change-point estimator is

We use ${\widehat{\Sigma}}_{T}$ as the sampled version of the moment estimator. The asymptotic critical values of the test statistics are computed from the block permutation based empirical distribution.

We simulate data according to MVAR(1) and MVAR(2) models. with the regional dimension *d* = 5 and *T* = 200 time points in 1, 000 repetitions. Since the brain network is high-dimensional, we also consider *d* = 20. The empirical critical values under the null distribution are obtained from no change-point model using block permutations for each subject. The fraction of one change-point is considered as *θ* = *τ*/*T* = 0.5 and 0.7 in the simulated data. We calculate the performance measure such as mean, median, square root of MSE (mean squared error), and the proportion around the true change-point as $p5=\widehat{P}(|\widehat{\tau}-\tau |<5)$ with 95 % confidence level.

Data are generated from the following model, **y**_{u} ~ MVAR(p) as

where **B**_{k} is *k*th parameter coefficient matrix of *MVAR*(*p*), *k* = 1, …, *p*. And **Σ**_{1} and **Σ**_{2} are covariance matrix of **ϵ**_{u}'s before and after the change-point respectively. We consider *MVAR*(1) models in cases (i) and (ii), and *MVAR*(2) models in (iii) and (iv) in the followings.

I. MVAR(1) with coefficient matrix **B**_{1}

where **I** is the identity matrix and **J** is a matrix whose components are all 1's.

(i) MVAR(1) with |**Σ**_{1}| < |**Σ**_{2}|

(ii) MVAR(1) with |**Σ**_{1}| > |**Σ**_{2}|

II. MVAR(2) with coefficient matrices **B**_{1} and **B**_{2}

(iii) MVAR(2) with |**Σ**_{1}| < |**Σ**_{2}|

(iv) MVAR(2) with |**Σ**_{1}| > |**Σ**_{2}|

To illustrate the potential to detect multiple change-points, we can apply the procedure with the subsegment after the change-point. Each time we reject H_{0}, we re-apply the method to the subsamples splitted by the proposed change-point.

We considered MVAR(1) and MVAR(2) cases with the different covariance structure after the change-point at *θ* = 0.5, 0.7. Tables 1, 2 provide MVAR(1) cases in which our estimators have smaller MSE's and higher matching proportions than Aue estimators. Tables 3, 4 give the change-point estimation results in MVAR(2) cases. Our method performs better than the Aue method in Tables 3, 4. Aue method depends on the covariance estimator ${\widehat{\mathbf{\Sigma}}}_{T}$ in the statistic in (19), and the change-point estimation results depend on the covariance structure accordingly. But our change-point estimator gives stable simulation results in the performance based upon mean and MSE. Overall in this simulation the proposed estimator performs very well.

**Table 1**. Change-point estimation results of the proposed $\widehat{\theta}$ and Aue estimator ${\widehat{\theta}}_{A}$ in case (i) MVAR(1) with *T* = 200 in 1, 000 repetitions.

**Table 2**. Change-point estimation results of the proposed $\widehat{\theta}$ and Aue estimator ${\widehat{\theta}}_{A}$ in case (ii) MVAR(1) with *T* = 200 in 1, 000 repetitions.

**Table 3**. Change-point estimation results of the proposed $\widehat{\theta}$ and Aue estimator ${\widehat{\theta}}_{A}$ in case (iii) MVAR(2) with *T* = 200 in 1, 000 repetitions.

**Table 4**. Change-point estimation results of the proposed $\widehat{\theta}$ and Aue estimator ${\widehat{\theta}}_{A}$ (iv) MVAR(2) with *T* = 200 in 1, 000 repetitions.

## 4. Application to the Epilepsy fMRI Data

### 4.1. Epilepsy fMRI Data

The same subjects from our previous study are used in this study except for one patient who showed a severe movement artifact (Jeong et al., 2019). Patients who underwent unilateral MTL resection (MTLR) to treat medically intractable temporal lobe epilepsy at Seoul National University Hospital at least 1 year before recruitment and were between 19 and 50 years of age are retrospectively recruited. We only include subjects who showed at least a low average or a higher level (scores> 80) of memory capacity and general intelligence evaluated by the standard neuropsychological test. A total of 34 patients (16 left and 18 right; median age = 32.5 years) and 24 age- and education-year-matched healthy controls (HC, median age = 32 years) are included in the present study. All subjects provided informed consent. This study was approved by the Institutional Review Board of Seoul National University Hospital.

The MR images were acquired on a research-dedicated 3T Magnetom Trio Tim Syngo (Siemens, Erlangen, Germany) using a 32-channel head coil. Five minutes of resting-state functional data (eyes open, fixation cross) were acquired both before (pre-encoding baseline) and after performing the in-scanner memory encoding paradigm of words and figures (post-encoding) using a T2^{*}-weighted gradient echo planar imaging sequence (36 axial slices, slice thickness = 3.4 mm, no gap, TR = 2,000 ms, TE = 30 ms, FOV the change-point 220 × 220 mm, flip angel = 80°, voxel size = 3.4 × 3.4 × 3.4*mm*^{3}, and interleaved). Whole-brain high-resolution anatomic T1-weighted images were obtained with the 3D TFL sequence (TR = 1,670 ms, TE = 1.89 ms, FOV = 250 × 250H mm, flip angle = 9°, voxel size = 1.0 × 1.0 × 1.0*mm*^{3}). The resting-state functional data underwent a number of preprocessing steps including motion correction, slice time correction, co-registration, spatial normalization, and spatial smoothing (AFNI, version: 16.0.00, https://afni.nimh.nih.gov/afni/). More details about the memory encoding paradigm and data preprocessing have been described elsewhere (Jeong et al., 2018, 2019).

Since previous studies have shown that FC during the post-encoding resting-state is related to the memory performance (Tambini et al., 2010; Staresina et al., 2013; Tambini and Davachi, 2013), we thought that data that showed evident FC changes of memory related areas rather than entire data of post-encoding resting-state may possibly provide the additional information of the memory consolidation network. In order to find the point that showed evident FC changes, we used change-point analysis.

For change-point estimation analysis, we select 20 ROIs within DMN, which is a well-known network that subserves memory function (Andrews-Hanna et al., 2010; Jeong et al., 2015). Table 5 shows the details of selected ROIs. In both pre- and post-encoding rest data, BOLD signals were extracted within each selected ROI (8 mm radius). For patients, one ROI, left hippocampus for left MTLR (LMTLR) and right hippocampus for right MTLR (RMTLR), was excluded from BOLD time series extraction. The extracted time series are then used for change-point estimation of post-encoding rest data to detect changes of FC after memory encoding. Next, pair-wise Pearson correlation coefficients between each pair of ROIs are calculated and averaged for each ROI, separately for time series after change-point and pre-encoding. FC is additionally calculated for the whole-segment of post-encoding time series in order to evaluate the usefulness of the change-point estimation method in the detection of memory encoding-related FC changes.

### 4.2. Change Analysis for Epilepsy Data

Change analysis starts with an individual test for existence of change. Under the significance of change, change-point estimation can be accomplished for its identification. With the epilepsy data, we performed the block permutation method to obtain the critical values. For dependent data like fMRI data, block permutation is preferable and useful (Kirch and Steinebach, 2006; Kirch, 2007; Adolf et al., 2011). This block permutation procedure is done individually because the change-point should be estimated subject-wise. The block permutation critical values are obtained with the block size *K* = 5 in 1,000 repetitions. At the significance level *α* = 0.05 we reject the null hypothesis for each individual subject except hc22 and lt16 (significant at *α* = 0.10). Since we focus on detectable change-point estimation, we get and use the change-points for all subjects. The change-point is estimable and informative for change, so we include these two subjects for the analysis. Therefore, we identify one change-point for each subject.

Table 6 shows a summary of change-point estimation and Figure 1 shows the change-points of individual subjects. The experiments are conducted during the two resting states. Change-point estimation is done during the second state because the first resting state is considered the baseline. Change-points are significantly different among groups [*F*_{(2,55)} = 4.29, *p* < 0.05, ANOVA]. Bonferroni *post-hoc* tests reveal that the LMTLR group shows significantly higher change-point values than the HC (*p* < 0.05) and RMTLR groups (*p* < 0.05). There is no significant differences between the HC and RMTLR groups (*p* > 0.10).

**Figure 1**. Change process ${G}_{t}^{2}$ according *t* and change-point estimation for each subject during post-encoding rest in HC **(Upper)**, LMTLR **(Middle)**, and RMTLR **(Lower)** groups.

Changes of FC after change-point and pre-encoding baseline are compared within each subject group by using the Wilcoxon signed rank test (Figure 2). FC without considering change-point estimation (whole-segment of post-encoding period) are also compared with baseline FC. FC values in areas that survived are presented in Table 7. In the HC group, the FC of the right parahippocampal cortex (PHC) shows a tendency to increase after the change-point (*p* = 0.067) and significantly increased in whole-segment of post-encoding rest (*p* < 0.05). In RMTLR group, none of ROIs shows significant FC differences. In the LMTLR group, although not statistically significant, the FC of the right temporal pole (TempP) shows a tendency to decrease after change-point (*p* = 0.088). In this ROI, whole-segment comparisons shows significant difference (*p* < 0.05). Interestingly, the FC of the right posterior inferior parietal lobule (pIPL) is significantly decreased after the change-point (*p* < 0.05) but shows no differences with the whole-segment of post-encoding rest (*p* = 0.642). Moreover, the FC difference (after change-point – baseline) in the right pIPL is negatively correlated with verbal immediate memory scores of the standard neuropsychological Rey-Kim memory test in LMTLR group (*r* = −0.548, *p* = 0.028; Spearman's correlation; Figure 2. There are no other significant clinical correlations of FC, to memory scores in either after-change point or whole-segment comparisons.

**Figure 2**. Changes of functional connectivity during pre- and post-encoding (after change-point) resting. Red color indicates increased FC and blue indicates decreased FC during post-encoding (Post-encoding minus Pre-encoding baseline). Size of circle indicates relative magnitude of FC changes in all ROIs in HC (healthy controls), LMTLR (left medial temporal lobe), and RMTLR (right MTLR). We displayed functional connections with gray line only in ROIs that showed significant (*p* < 0.05) or trend (*p* < 0.1) of changes. **(A,C,E)** Top view. **(B,D,F)** Front view. **(G)** Relationships of FC differences to memory capacity. ^{*}Significant at 0.05 level.

The relationships between the behavioral accuracy (d-prime) of the in-scanner memory task and the FC values of post-encoding rest after change-point are also analyzed. Our previous study describes behavioral accuracy in detail (Jeong et al., 2019). Figure 3 shows the significant correlations. The FC of the ventral medial prefrontal cortex (vMPFC) is positively correlated with word accuracy in LMTLR (*r* = 0.563, *p* < 0.05). In RMTLR, word accuracy is positively correlated with the FC of the dorsal medial prefrontal cortex (dMPFC; *r* = 0.556), left and right posterior cingulate cortex (PCC; *r* = 0.684, *r* = 0.552), and left and right temporal parietal junction (TPJ; *r* = 0.577, *r* = 0.608) (*p* < 0.05 for all). Figure accuracy is positively correlated with the FC of the right PCC in the RMTLR group (*r* = 0.497, *p* < 0.05).

**Figure 3**. Relationships between the behavioral accuracy of in-scanner memory task (words, figures) and the FC values of post-encoding resting after change-point. ^{*}Significant at 0.05 level.

## 5. Discussion

In this paper we develop a new procedure that can be used to determine a change-point in a individual FC system, using a function of maximum eigenvalues of suitably scaled covariance matrices. The proposed method is applicable to either resting-state or task fMRI studies providing individualized change-point estimation. The method yields individualized change-point estimation results, but it permits to group-level estimation. Our change-point estimation method generally works well in simulations and shows improved performance over the previous approach. Our work avoids the use of a sliding window usually applied to correlation measures, and reinforces the notion that time-varying techniques are required to study FC. Our results stress the importance of taking the change-point into account to characterize dynamic FC in order to describe how “information” is integrated and changed over time.

In this study, we applied change-point estimation analysis to post-encoding rest fMRI data of MTLR patients in order to evaluate the applicability and usefulness of the proposed method. We found that both segment after change-point and whole-segment of resting-state data after memory encoding detect similar FC changes. In HC, the MTL area (PHC) shows increased FC after memory encoding which is consistent with previous findings of memory consolidation studies (Tambini et al., 2010; Staresina et al., 2013). In LMTLR patients, the FC of the right pIPL and TempP decreased after memory encoding. Importantly, right pIPL is only detected when we use post-encoding data after the change-point, but not detected with the whole-segment of data. Moreover, the magnitude of FC changes in the right pIPL is associated with verbal memory capacity of LMTLR patients. Patients who displayed more reduction in FC values have better verbal memory. Recent studies suggest that the posterior parietal cortex contributes to memory consolidation (Brodt et al., 2016, 2018). Therefore, we speculate that DMN areas other than damaged MTL, especially the right pIPL, may play a compensatory role for the verbal memory function of resected left MTL. In summary, by applying change-point analysis, we can reveal additional information on the memory consolidation network in an MTL resected brain. Further studies with other data sets are warranted to ensure the general use of this change-point estimation analysis of fMRI data.

The other notable contribution of our work is that we develop a generalized multivariate extension to use an eigenvalue system based on RMT. While previous approaches have largely targeted univariate approaches examining a single time series or a single measure of association, our method use full information of connectivity via the eigenvalue system. Another interesting aspect of our change-point detection is the quantitative comparison before and after change-point. It reveals some properties of FC after change occurs. Different analytic approaches provides different perspective information on fMRI data, which necessitates integrating knowledge gained through the variety of models for the understanding of dynamic FC.

One limitation of the method is that the simulation results depend on the underlying covariance structure. When there are more change-points, simultaneous multiple change-points estimation procedure should be used. Further study with this extension is expected when there are multiple change-points in FC.

## Data Availability Statement

The datasets generated for this article are not readily available because: it can be available only if requested to the third author. Requests to access the datasets should be directed to Chun Kee Chung, chungc@snu.ac.kr.

## Ethics Statement

The studies involving human participants were reviewed and approved by Institutional Review Board of Seoul National University Hospital. The patients/participants provided their written informed consent to participate in this study.

## Author Contributions

Statistical method was developed and simulation work was done by JK. Epilepsy data were collected by CC. Epilepsy data analyses and interpretations were conducted by WJ. All authors contributed to the article and approved the submitted version.

## Funding

This research was supported by Mid-career Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education, Science and Technology (2018R1A2B6001664).

## Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

## Acknowledgments

We thank Hyeongrae Lee for helping with preprocessing of resting-state functional data.

## References

Adolf, D., Baecke, S., Kahle, W., Bernarding, J., and Kropf, S. (2011). Applying multivariate techniques to high-dimensional temporally correlated fMRI data. *J. Stat. Plan. Inform.* 141, 3760–3770. doi: 10.1016/j.jspi.2011.06.012

Anderson, T. W. (1963). Asymptotic theory for principal component analysis. *Ann. Math. Stat.* 34, 122–148. doi: 10.1214/aoms/1177704248

Andrews-Hanna, J. R., Reidler, J. S., Sepulcre, J., Poulin, R., and Buckner, R. L. (2010). Functional-anatomic fractionation of the brain's default network. *Neuron* 65, 550–562. doi: 10.1016/j.neuron.2010.02.005

Aue, A., Hörmann, S., Horváth, L., and Reimherr, M. (2009). Break detection in the covariance structure of multivariate time series models. *Ann. Stat.* 37, 4046–4087. doi: 10.1214/09-AOS707

Aue, A., and Horváth, L. (2013). Structural breaks in time series. *J. Time Ser. Anal.* 34, 1–16. doi: 10.1111/j.1467-9892.2012.00819.x

Bai, Z., and Silverstein, J. (2010). *Spectral Analysis of Large Dimensional Random Matrices*. New York, NY: Springer.

Biswall, B., Yetkin, F., Haughton, V. M., and Hyde, J. S. (1995). Functional connectivity in the motor cortex of resting human brain using echo-planar MRI. *Magn. Reson. Med*. 34, 537–541. doi: 10.1002/mrm.1910340409

Brodt, S., Gais, S., Beck, J., Erb, M., Scheffler, K., and Schönauer, M. (2018). Fast track to the neocortex: a memory engram in the posterior parietal cortex. *Science* 362:1045. doi: 10.1126/science.aau2528

Brodt, S., Pöhlchen, D., Flanagin, V. L., Glasauer, S., Gais, S., and Schönauer, M. (2016). Rapid and independent memory formation in the parietal cortex. *Proc. Natl. Acad. Sci. U.S.A.* 113, 13251–13256. doi: 10.1073/pnas.1605719113

Calhoun, V., Adali, T., Pearlson, G., and Pekar, J. (2001). A method for making group inferences from functional MRI data using independent component analysis. *Hum. Brain Mapp*. 14, 140–151. doi: 10.1002/hbm.1048

Chang, C., and Glover, G.H. (2010). Time-frequency dynamics of resting-state brain connectivity measured with fMRI. *Neuroimage* 50, 81–98. doi: 10.1016/j.neuroimage.2009.12.011

Chenouri, S., Mozaffari, A., and Rice, G. (2020). Robust multivariate change point analysis based on data depth. *Can. J. Stat*. 48, 417–446. doi: 10.1002/cjs.11541

Cheon, S., and Kim, J. (2010). Multiple change-point detection of multivariate mean vectors with the Bayesian approach. *Comput. Stat. Data Anal*. 54, 406–415. doi: 10.1016/j.csda.2009.09.003

Cribben, I., Wager, T. D., and Lindquist, M. A. (2013). Detecting functional connectivity changepoints for single-subject fMRI data. *Front. Comput. Neurosci*. 7:143. doi: 10.3389/fncom.2013.00143

El Karoui, N. (2007). Tracy-Widom limit for the largest eigenvalue of a large class of complex sample covariance matrices. *Ann. Prob.* 35, 663–714. doi: 10.1214/009117906000000917

El Karoui, N. (2008). Spectrum estimation for large dimensional covariance matrices using random matrix theory. *Ann. Stat.* 36, 2757–2790. doi: 10.1214/07-AOS581

Friston, K. J. (2011). Functional and effective connectivity: a review. *Brain Connect.* 1, 13–36. doi: 10.1089/brain.2011.0008

Goodman, N. R. (1963). The distribution of the determinant of a complex Wishart distributed matrix. *Ann. Math. Stat.* 34, 178–180. doi: 10.1214/aoms/1177704251

Hellyer, P. J., Shanahan, M., Scott, G., Wise, R. J., Sharp, D. J., and Leech, R. (2014). The control of global brain dynamics: opposing actions of frontoparietal control and default mode networks on attention. *J. Neurosci.* 34, 451–461. doi: 10.1523/JNEUROSCI.1853-13.2014

Husková, M. (2004). Permutation principle and bootstrap in change point analysis. *Inst. Commun.* 44, 273–291. doi: 10.1090/fic/044/15

Hutchison, R. M., Womelsdorf, T., Allen, E. A., Bandettini, P. A., Calhoun, V. D., Corbetta, M., et al. (2013b). Dynamic functional connectivity: promise, issues, and interpretations. *Neuroimage* 80, 360–378. doi: 10.1016/j.neuroimage.2013.05.079

Hutchison, R. M., Womelsdorf, T., Gati, J. S., Everling, S., and Menon, R. S. (2013a). Resting state networks show dynamic functional connectivity in awake humans and anesthetized macaques. *Hum. Brain Mapp*. 34, 2154–2177. doi: 10.1002/hbm.22058

Inclán, C., and Tiao, G.C. (1994). Use of cumulative sums of squares for retrospective detection of changes of variance. *J. Am. Stat. Assoc*. 89, 913–923. doi: 10.1080/01621459.1994.10476824

Jafria, M. J., Pearlsona, G. D., Stevensa, M., and Calhoun, V. (2008). A method for functional network connectivity among spatially independent resting-state components in schizophrenia. *Neuroimage* 39, 1666–1681. doi: 10.1016/j.neuroimage.2007.11.001

Jeong, W., Chung, C. K., and Kim, J. S. (2015). Episodic memory in aspects of large-scale brain networks. *Front. Hum. Neurosci.* 9:454. doi: 10.3389/fnhum.2015.00454

Jeong, W., Lee, H., Kim, J. S., and Chung, C. K. (2018). Neural basis of episodic memory in the intermediate term after medial temporal lobe resection. *J. Neurosurg.* 131, 790–798. doi: 10.3171/2018.5.JNS18199

Jeong, W., Lee, H., Kim, J. S., and Chung, C. K. (2019). Characterization of brain network supporting episodic memory in the absence of one medial temporal lobe. *Hum. Brain Mapp*. 40:2188. doi: 10.1002/hbm.24516

Johnstone, I.M. (2009). Approximate null distribution of the largest root in multivariate analysis. *Ann. Appl. Stat*. 4, 1616–1633. doi: 10.1214/08-AOAS220

Johnstone, I. M. (2001). On the distribution of the largest eigenvalue in principal components analysis. *Ann. stat*. 29, 295–327. doi: 10.1214/aos/1009210544

Kim, J. (2019). “Change-point estimation and testing for brain functional connectivity networks,” in *11th International Conference on Knowledge and Smart Technology (KST)* (Phuket), 226–231.

Kim, J., and Cheon, S. (2010). Bayesian multiple change-point estimation with annealing stochastic approximation Monte Carlo. *Comput. Stat.* 25, 215–239. doi: 10.1007/s00180-009-0172-x

Kim, J., and Hart, J. D. (2011). A change-point estimator using local fourier series. *J. Nonparametr. Stat*. 23, 83–98. doi: 10.1080/10485251003721232

Kirch, C. (2007). Block permutation principles for the change analysis of dependent data. *J. Stat. Plan. Inform.* 137, 2453–2474. doi: 10.1016/j.jspi.2006.09.026

Kirch, C., and Steinebach, J. (2006). Permutation principles for the change analysis of stochastic processes under strong invariance. *J. Comp. Appl. Math*. 186, 64–88. doi: 10.1016/j.cam.2005.03.065

Koike, T., Kan, S., Misaki, M., and Miyauchi, S. (2011). Connectivity pattern changes in default-mode network with deep non-REM and REM sleep. *Neurosci. Res.* 69, 322–330. doi: 10.1016/j.neures.2010.12.018

Lindquist, M. A., Waugh, C., and Wager, T. D. (2007). Modeling state-related fMRI activity using change-point theory. *Neuroimage* 35, 1125–1141. doi: 10.1016/j.neuroimage.2007.01.004

Lindquist, M. A., Xu, Y., Nebel, M. B., and Caffo, B. S. (2014). Evaluating dynamic bivariate correlations in resting-state fMRI: a comparison study and a new approach. *Neuroimage* 101, 531–546. doi: 10.1016/j.neuroimage.2014.06.052

Luger, R. (2006). Exact permutation tests for non-nested non-linear regression models. *J. Econometr.* 133, 513–529. doi: 10.1016/j.jeconom.2005.06.005

Marčenko, V., and Pastur, L. (1967). Distribution for some sets of random matrices. *Math USSR-Sb* 1, 457–483. doi: 10.1070/SM1967v001n04ABEH001994

Mayberg, H. S., Liotti, M., Brannan, S. K., McGinnis, S., Mahurin, R. K., Jerabek, P. A., et al. (1999). Reciprocal limbic-cortical function and negative mood: converging PET findings in depression and normal sadness. *Am. J. Psychiatry* 156, 675–682.

McKeown, M. J., Makeig, S., Brown, G. G., Jung, T.-P., Kindermann, S. S., Bell, A. J., et al. (1998). Analysis of fMRI data by blind separation into independent spatial components. *Hum. Brain Mapp.* 6, 160–188.

Okamoto, M. (1973). Distinctness of the eigenvalues of a quadratic form in a multivariate sample. *Ann. Stat.* 1, 763–765.

Petersen, A., and Muller, H.-G. (2016). Frechet integration and adaptive metric selection for interpretable covariances of multivariate functional data. *Biometrika* 103, 103–120. doi: 10.1093/biomet/asv054

Raj Rao, N., Mingo, J. A., Speicher, R., and Edelman, A. (2008). Statisticsl eigen-inference from large Wishart matrices. *Ann. Stat.* 36, 2850–2885. doi: 10.1214/07-AOS583

Robinsona, L. F., Wager, T. D., and Lindquist, M. A. (2010). Change point estimation in multi-subject fMRI studies. *Neuroimage* 49, 1581–1592. doi: 10.1016/j.neuroimage.2009.08.061

Salvador, R., Suckling, J., Coleman, M. R., Pickard, J. D., Menon, D., and Bullmore, E. (2005). Neurophysiological architecture of functional magnetic resonance images of human brain. *Cereb. Cortex* 15, 1332–1342. doi: 10.1093/cercor/bhi016

Sato, J. R., Junior, E. A., Takahashi, D. Y., Felix, M. D. M., Brammer, M. J., and Morettina, P. A. (2006). A method to produce evolving functional connectivity maps during the course of an fMRI experiment using wavelet-based time-varying Granger causality. *Neuroimage* 31, 187–196. doi: 10.1016/j.neuroimage.2005.11.039

Shao, X., and Zhang, X. (2010). Testing for change points in time series. *J. Am. Stat. Assoc.* 105, 1228–1240. doi: 10.1198/jasa.2010.tm10103

Staresina, B. P., Alink, A., Kriegeskorte, N., and Henson, R. N. (2013). Awake reactivation predicts memory in humans. *Proc. Natl. Acad. Sci. U.S.A.* 110, 21159–21164. doi: 10.1073/pnas.1311989110

Strasser, H., and Weber, C. (1999). On the asymptotic theory of permutation statistics. *Math. Methods Stat.* 8, 220–250.

Sun, F. T., Miller, L. M., and D'Esposito, M. (2004). Measuring interregional functional connectivity using coherence and partial coherence analyses of fMRI data. *Neuroimage* 21, 647–658. doi: 10.1016/j.neuroimage.2003.09.056

Supekar, K., Menon, V., Rubin, D., Musen, M., and Greicius, M. D. (2008). Network analysis of intrinsic functional brain connectivity in Alzheimer's disease. *PLoS Comput. Biol*. 4:e1000100. doi: 10.1371/journal.pcbi.1000100

Tambini, A., and Davachi, L. (2013). Persistence of hippocampal multivoxel patterns into postencoding rest is related to memory. *Proc. Natl. Acad. Sci. U.S.A.* 110, 19591–19596. doi: 10.1073/pnas.1308499110

Tambini, A., Ketz, N., and Davachi, L. (2010). Enhanced brain correlations during rest are related to memory for recent experiences. *Neuron* 65, 280–290. doi: 10.1016/j.neuron.2010.01.001

Tracy, C. A., and Widom, H. (1996). On orthogonal and symplectic matrix ensembles. *Commun. Math. Phys.* 177, 727–754. doi: 10.1007/BF02099545

Wied, D., and Krömer, W. (2012). Testing for a change in correlation at an unknown point in time using an extended functional delta method. *Econometr. Theor.* 28, 570–589. doi: 10.1017/S0266466611000661

Xu, Y., and Lindquist, M. A. (2015). Dynamic connectivity detection: an algorithm for determining functional connectivity change points in fMRI data. *Front. Neurosci*. 9:285. doi: 10.3389/fnins.2015.00285

Keywords: Tracy-Widom distribution, random matrix theory, fMRI, epilepsy, eigenvalue, dynamic functional connectivity, covariance, change-point

Citation: Kim J, Jeong W and Chung CK (2021) Dynamic Functional Connectivity Change-Point Detection With Random Matrix Theory Inference. *Front. Neurosci.* 15:565029. doi: 10.3389/fnins.2021.565029

Received: 23 May 2020; Accepted: 29 March 2021;

Published: 04 May 2021.

Edited by:

Plamen Ch. Ivanov, Boston University, United StatesCopyright © 2021 Kim, Jeong and Chung. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Jaehee Kim, jaehee@duksung.ac.kr