Sub-state space estimation method based on function connection and random process model

Abstract

A sub-state space estimation method based on function connection and a random process model comprises the following steps: S1, constructing a random process model which comprises an input module, a random process module and an output module; the input module receives function connection data, the random process module adopts an OU model, and the output module outputs a multi-dimensional OU sequence; s2, a function connection matrix is calculated through BOLD data to serve as input of an input module; s3, sending the function connection matrix into a random process module, and outputting a random process sequence through an output module; s4, calculating dynamic function connection through a random process sequence; and S5, based on the step S4, estimating a brain sub-state space to obtain a sub-state probability and a state transition matrix. Compared with an existing sub-state space calculation method completely based on data, the method has higher robustness and adaptability, deep connection between functional connection and the sub-state space can be revealed, and more efficient technical support is provided for research of the brain sub-state.

Description

Technical Field

[0001] The present invention belongs to the interdisciplinary technical field of biomedical engineering and brain science, and specifically relates to a sub-state space estimation method based on functional connectivity and random process models. Background Art

[0002] Blood-oxygen-level dependent fMRI (BOLD-fMRI) is an imaging technique that indirectly detects neural activity by measuring changes in blood oxygen levels in the brain. It mainly includes two paradigms: resting state and task state. Among them, resting state functional MRI (rs-fMRI) highlights the rich structure of brain activity in the absence of tasks or stimulation. In the past few decades, the spontaneous fluctuations of brain functional activity in the resting state have attracted the attention of many researchers. In recent years, it is generally believed that the spontaneous fluctuation signals detected by resting state functional MRI contain huge information about brain functional activities and have extremely important data information mining potential. Moreover, the analysis of spontaneous fluctuation information can get closer to the essence of brain functional network organization and operation. When the brain completes a specific cognitive task, it is not just a single brain area that is responsible, but multiple brain areas work together. In order to quantify the relationship between different brain areas, researchers have proposed the concept of functional connectivity (FC). It measures the synchronization of spontaneous activities of different brain areas in the resting state, usually obtained by calculating the Pearson correlation coefficient in a specific frequency band. Over the past two decades, much effort has been devoted to studying functional connectivity, i.e., the functional interactions between different brain regions, which for a long time were thought to be static.

[0003] Research based on static functional connectivity technology assumes that the brain functional network is static and stable over a period of time, but in fact the brain's activities and networks are dynamic and time-varying, and this dynamicity is reflected in the resting state and different cognitive tasks. To address this problem, researchers have proposed the dynamic functional connectivity method (Dynamic functional connectivity, dFC), which has been widely used to study the dynamic functional connectivity patterns between different brain regions. The dynamics can be characterized by the appearance and switching of probabilistic metastable substates (PMS). Past studies have found that these brain processes can be divided into several stable brain states by clustering dynamic functional connectivity or BOLD temporal features or constructing hidden Markov models. These states correspond to the working modes of the brain under different cognitive tasks, emotional states or physiological conditions. By calculating the probability of occurrence of substates and the state transition matrix, we can understand the differences in substates among different populations and understand the hierarchical structure of brain states at the spatial and temporal levels.

[0004] However, a series of methods for calculating sub-state space in the past all rely entirely on the original BOLD sequence. When the collected BOLD sequence is contaminated or of low quality, especially when its temporal resolution is not high or the number of time points is too short, it is often difficult to calculate the accurate sub-state space, which limits the application and comparison of sub-state space analysis across data sets, as well as its promotion in clinical practice. At the same time, the connection between resting-state functional connectivity and sub-state space is not clear, which affects researchers' joint analysis of data. Summary of the invention

[0005] In order to solve the above technical problems, the present invention is based on the basic mathematical characteristics of whole-brain functional connectivity and random processes, and at the same time uses numerical simulation methods to solve the shortcoming that the state of high-dimensional random processes is unanalyzable, and designs an estimation method for the sub-state space based on functional connectivity and random process models.

[0006] The object of the present invention is achieved through the following technical solution: A sub-state space estimation method based on functional connectivity and random process model comprises the following steps:

[0007] S1. Construct a random process model, including an input module, a random process module, and an output module; the input module receives functional connection data as the initial input of the model, representing the covariance matrix of the multidimensional OU process; the random process module adopts the OU model and performs numerical simulation through the Euler method to solve the numerical solution; the output module outputs the multidimensional OU sequence;

[0008] S2, calculate the functional connectivity matrix through BOLD data as the input of the input module; generate the functional connectivity matrix FC by calculating the Pearson correlation coefficient between the BOLD signals of different brain regions;

[0009] S3, sending the functional connection matrix obtained in step S2 to the random process module, and outputting the random process sequence through the output module;

[0010] The discretized form of the stochastic process model is expressed by the following difference equation:

[0011]

[0012] Where Z t is an N-dimensional standard normal random variable:

[0013]

[0014] X t represents the random process signal at time t, I N represents the identity matrix, θ is the temporal regression coefficient; ∑ is the covariance matrix, i.e., the functional connectivity matrix of the input module; Δt is the same as the sampling interval TR of the real BOLD data;

[0015] By numerically solving the stochastic process model using the Euler method, a multi-dimensional OU process is obtained as the output module output;

[0016] S4. Calculate dynamic functional connectivity through random process sequence: Based on the multidimensional OU process output by the output module, the dynamic functional connectivity is further calculated through phase difference method, sliding window method, and time derivative multiplication;

[0017] S5. Based on step S4, the brain sub-state space is estimated to obtain the sub-state probability and state transition matrix. For the calculation of brain sub-states, the dynamic functional connections at each moment are first reduced in dimension by the LEiDA method, and the three initial dynamic functional connections of each subject are reduced in dimension into two-dimensional data of N*T. This method is mainly achieved by extracting the maximum eigenvector:

[0018] dFC(t)=VΛV -1

[0019] Where dFC(t) is the initial matrix, i.e., the dynamic functional connectivity matrix at time t, i.e., dFC_pd(t,i,j), dFC_sw(t,i,j) or dFC_td(t,i,j); V is the eigenvector matrix, Λ is the eigenvalue vector, and LEiDA is the eigenvector with the largest eigenvalue;

[0020] For N*T two-dimensional data, the data at each moment is taken as a sample and clustered. Then, the Dunn coefficient is used as the clustering evaluation index, and the cluster number with the largest Dunn coefficient is selected. After the clustering operation is completed, the data with N features at each moment is summarized into one state according to the centroid of each cluster. Then, a state sequence with a length of T is obtained through the original N*T data.

[0021] Or a hidden Markov model is used to extract states from sequence features to obtain a state sequence of length T;

[0022] After the clustering algorithm or the hidden Markov algorithm estimates the state, a state sequence with a length of T is obtained, and the sub-state probability and the state transfer matrix are calculated through the state sequence.

[0023] Beneficial effects of the invention: The invention is based on the basic mathematical characteristics of whole-brain functional connectivity and random processes, and solves the shortcomings of the unanalyzable state of high-dimensional random processes through numerical simulation methods. It designs a sub-state space estimation method based on functional connectivity and random process models, which overcomes the shortcomings of the original completely sequence-based calculation method that is overly dependent on the time resolution and sampling length of the original data, and is helpful for the application and comparison of sub-state space analysis across data sets, as well as its promotion in clinical practice. Compared with the existing sub-state space calculation method that is completely based on data, it has stronger robustness and adaptability, helps to reveal the deep connection between functional connectivity and sub-state space, and provides more efficient technical support for the study of brain sub-states. BRIEF DESCRIPTION OF THE DRAWINGS

[0024] Figure 1 It is a flowchart of the estimation method of the sub-state space based on functional connectivity and random process model of the present invention.

[0025] Figure 2 Schematic diagram of the result of state occurrence probability estimation in an embodiment of the present invention.

[0026] Figure 3 Schematic diagram of the state transfer matrix in an embodiment of the present invention. DETAILED DESCRIPTION

[0027] The technical solution of the present invention is further described below in conjunction with the accompanying drawings.

[0028] like Figure 1 As shown, a sub-state space estimation method based on functional connectivity and random process model of the present invention comprises the following steps:

[0029] S1. Construct a random process model, including an input module, a random process module, and an output module; the random process model adopts a multidimensional OU process model, and the time series obtained by numerical simulation of the OU process is used as output;

[0030] The input module receives the functional connectivity data as the initial input of the model, which represents the covariance matrix of the multidimensional OU process; the random process module adopts the OU model and performs numerical simulation to obtain the numerical solution through the Euler method; the output module outputs the multidimensional OU sequence;

[0031] S2. Calculate the functional connectivity matrix through BOLD data as the input of the input module; generate the functional connectivity matrix FC by calculating the Pearson correlation coefficient between the BOLD signals of different brain regions. For the BOLD signal time series of brain regions i and j, the Pearson correlation coefficient r ij The calculation expression is as follows:

[0032]

[0033] Among them, T represents the overall time length, B i,t represents the BOLD signal value of brain area i at time t, represents the average value of BOLD signal in brain area i;

[0034] S3, sending the functional connection matrix obtained in step S2 to the random process module, and outputting the random process sequence through the output module;

[0035] The discretized form of the stochastic process model is expressed by the following difference equation:

[0036]

[0037] Where Z t is an N-dimensional standard normal random variable:

[0038]

[0039] X t represents the random process signal at time t, I N represents the identity matrix (indicating that each normal random variable is independent and identically distributed); θ is the time regression coefficient (how fast the reaction regresses to the mean); ∑ is the covariance matrix, i.e., the functional connectivity matrix of the input module; Δt is the same as the sampling interval TR of the real BOLD data;

[0040] By numerically solving the stochastic process model using the Euler method, a multi-dimensional OU process is obtained as the output module output;

[0041] S4. Calculate dynamic functional connectivity through random process sequence: Based on the multidimensional OU process output by the output module, further calculate dynamic functional connectivity through phase difference method, sliding window method, time derivative multiplication and other methods;

[0042] The dynamic functional connectivity under the phase difference definition is expressed as:

[0043] dFC_pd(t,i,j)=cos(angle o (t)-angle j (t))

[0044] Where dFC_pd(t,i,j) represents the dynamic functional connectivity strength between brain regions i and j at time t. The closer it is to 1, the more synchronized the two brain regions are at that moment. The closer it is to -1, the more reversely synchronized the two brain regions are at that moment. 0 represents asynchronous.

[0045] angle i (t) represents the phase of the brain region at that moment:

[0046]

[0047] Where x(t) represents the random process sequence output by S3; is the Hilbert transform of x(t), defined as:

[0048]

[0049] The sliding window method divides a time series into T moments by translating the observation window. The dynamic functional connection at a certain moment is expressed as the Pearson correlation coefficient between the sequences under a certain window:

[0050]

[0051] Where dFC_sw(t,i,j) represents the dynamic functional connectivity strength between brain regions i and j at time t, N represents the window size, and B i,t1 represents the BOLD signal value of brain area i at time t1 within the observation window, represents the average value of BOLD signal of brain area i within the observation window;

[0052] The time derivative multiplication method calculates the dynamic functional connection by calculating the product of the time derivatives at a certain moment and then performing a certain smoothing on the result:

[0053]

[0054] Where dFC_td(t,i,j) represents the dynamic functional connectivity strength between brain areas i and j at time t, w represents the smoothing coefficient, and MTD(t,i,j) represents the time derivative product:

[0055]

[0056] where dt it is the time series derivative, σ i is the standard deviation of dt in brain region i;

[0057] dt it =ts it -ts it-1

[0058] where ts it is the signal activity intensity of brain area i at time t;

[0059] S5. Based on step S4, the brain sub-state space is estimated to obtain the sub-state probability and state transition matrix. For the calculation of brain sub-states, the dynamic functional connectivity at each moment is first reduced in dimension by the LEiDA method, and the three initial dynamic functional connectivity strengths dFC of each subject (all three-dimensional data of N*N*T) are reduced in dimension to two-dimensional data of N*T. This method is mainly achieved by extracting the maximum eigenvector:

[0060] dFC(t)=VΛV -1

[0061] Where dFC(t) is the initial matrix, i.e., the dynamic functional connectivity matrix at time t, i.e., dFC_pd(t,i,j), dFC_sw(t,i,j) or dFC_td(t,i,j); V is the eigenvector matrix, Λ is the eigenvalue vector, and LEiDA is the eigenvector with the largest eigenvalue;

[0062] For N*T two-dimensional data, take the data at each moment as a sample, cluster it through kmean++ or deep learning methods, and then use the dunn coefficient as the clustering evaluation index, and select the cluster number with the largest dunn coefficient. The calculation method of dunn coefficient is:

[0063]

[0064] Where C is a set of clusters, including k clusters; dist(c i ,c j ) is the cluster c i and c j The minimum distance between i ) is cluster C i The maximum distance metric within the

[0065] After completing the clustering operation, the data with N features at each moment is summarized into one state according to the centroid of each cluster, so a state sequence with a length of T is obtained through the original N*T data.

[0066] Or a hidden Markov model is used to extract states from sequence features to obtain a state sequence of length T; the hidden Markov model is represented as follows:

[0067] λ=(A,B,II)

[0068] Where A represents the state transition probability matrix, B represents the observed state transition probability, and II represents the initial state probability;

[0069] After the clustering algorithm or the hidden Markov algorithm estimates the state, a state sequence with a length of T is obtained, and the sub-state probability and the state transfer matrix are calculated through the state sequence.

[0070] In this embodiment, the specific method for evaluating the estimation accuracy is as follows:

[0071] Estimation accuracy evaluation: Calculate the difference between the PMS space calculated based on this method and the PMS space based entirely on real data to evaluate the fitting accuracy, that is, use the KL distance as a quantitative indicator:

[0072]

[0073] Where X is the domain of the probability distribution (i.e. the space of possible events), and P(x) is the probability of event x under distribution P.

[0074] Figure 2 The results of the PMS space output based on this method are presented, as well as a comparison with the calculation results of real data.

[0075] The difference between the state transfer matrix calculated based on this method and the state transfer matrix based entirely on real data is calculated to evaluate the fitting accuracy, that is, the Markov entropy distance is used as a quantitative indicator:

[0076]

[0077] Where D MED (P1, P2) represents the Markov entropy distance between the two state transfer matrices, A represents the steady-state distribution probability, P ij represents the probability of going from state i to state j.

[0078] Figure 3 The results of the state transfer matrix output by this method are presented, and compared with the calculation results of real data.

[0079] In summary, the method of the present invention aims to solve the problem that the sub-state space is highly dependent on the quality of real data, and effectively estimates the sub-state space by only using functional connectivity and random process models. The present invention is based on the basic mathematical properties of whole-brain functional connectivity and random processes, and at the same time solves the shortcomings of the unanalyzable state of high-dimensional random processes through numerical simulation methods, overcoming the shortcomings of the original completely sequence-based calculation method that is overly dependent on the time resolution and sampling length of the original data, which is helpful for the application and comparison of sub-state space analysis across data sets, as well as its promotion in clinical practice. Compared with the existing sub-state space calculation method that is completely based on data, it has stronger robustness and adaptability, helps to reveal the deep connection between functional connectivity and sub-state space, and provides more efficient technical support for the study of brain sub-states.

[0080] Those skilled in the art will appreciate that the embodiments described herein are intended to help readers understand the principles of the present invention, and should be understood that the protection scope of the present invention is not limited to such specific statements and embodiments. Those skilled in the art can make various other specific variations and combinations that do not deviate from the essence of the present invention based on the technical revelations disclosed by the present invention, and these variations and combinations are still within the protection scope of the present invention.

Claims

1. A sub-state space estimation method based on functional connectivity and stochastic process models, characterized in that: The following steps are involved: S1. Construct a random process model, including an input module, a random process module, and an output module; The input module receives the functional connectivity data as the initial input of the model, which represents the covariance matrix of the multidimensional OU process; the random process module adopts the OU model and performs numerical simulation to obtain the numerical solution through the Euler method; The output module outputs a multi-dimensional OU sequence; S2. Calculate the functional connectivity matrix through BOLD data as the input of the input module; generate the functional connectivity matrix FC by calculating the Pearson correlation coefficient between the BOLD signals of different brain regions. For the BOLD signal time series of brain regions i and j, the Pearson correlation coefficient r ij The calculation expression is as follows: Among them, T represents the overall time length, B i,t represents the BOLD signal value of brain area i at time t, represents the average value of BOLD signal in brain area i; S3, sending the functional connection matrix obtained in step S2 to the random process module, and outputting the random process sequence through the output module; The discretized form of the stochastic process model is expressed by the following difference equation: Where Z t is an N-dimensional standard normal random variable: X t represents the random process signal at time t, I N represents the identity matrix, θ is the temporal regression coefficient; ∑ is the covariance matrix, i.e., the functional connectivity matrix of the input module; Δt is the same as the sampling interval TR of the real BOLD data; By numerically solving the stochastic process model using the Euler method, a multi-dimensional OU process is obtained as the output module output; S4. Calculate dynamic functional connectivity through random process sequence: Based on the multidimensional OU process output by the output module, the dynamic functional connectivity is further calculated through phase difference method, sliding window method, and time derivative multiplication; The dynamic functional connectivity under the phase difference definition is expressed as: dFC_pd(t,i,j)=cos(angle i (t)-angle j (t)) angle i (t) represents the phase of the brain region at that moment: Where x(t) represents the random process sequence output by S3, is the Hilbert transform of x(t); The sliding window method divides a time series into T moments by translating the observation window. The dynamic functional connection at a certain moment is expressed as the Pearson correlation coefficient between the sequences under a certain window: N represents the window size, B i,t1 represents the BOLD signal value of brain area i at time t1 within the observation window, represents the average value of BOLD signal of brain area i within the observation window; The time derivative multiplication method calculates the dynamic functional connection by calculating the product of the time derivatives at a certain moment and then performing a certain smoothing on the result: w represents the smoothing coefficient, MTD(t,i,j) represents the time derivative product; S5. Based on step S4, the brain sub-state space is estimated to obtain the sub-state probability and state transition matrix. For the calculation of brain sub-states, the dynamic functional connections at each moment are first reduced in dimension by the LEiDA method, and the three initial dynamic functional connections of each subject are reduced in dimension into two-dimensional data of N*T. This method is mainly achieved by extracting the maximum eigenvector: dFC(t)=VΛV -1 Where dFC(t) is the initial matrix, i.e., the dynamic functional connectivity matrix at time t, i.e., dFC_pd(t,i,j), dFC_sw(t,i,j) or dFC_td(t,i,j); V is the eigenvector matrix, A is the eigenvalue vector, and LEiDA is the eigenvector with the largest eigenvalue; For N*T two-dimensional data, the data at each moment is taken as a sample and clustered. Then, the Dunn coefficient is used as the clustering evaluation index, and the cluster number with the largest Dunn coefficient is selected. After the clustering operation is completed, the data with N features at each moment is summarized into one state according to the centroid of each cluster. Then, a state sequence with a length of T is obtained through the original N*T data. Or a hidden Markov model is used to extract states from sequence features to obtain a state sequence of length T; After the clustering algorithm or the hidden Markov algorithm estimates the state, a state sequence with a length of T is obtained, and the sub-state probability and the state transfer matrix are calculated through the state sequence.

Images