CN102509282B - Efficiency connection analysis method fused with structural connection for each brain area - Google Patents

Abstract

A method for analyzing the performance connectivity of each brain region that fuses structural connections. The steps are to extract and analyze the brain structure network in the region of interest, and convert the obtained structural information into the prior probability distribution space of the performance connection parameters; The performance connection model of the Bayesian framework; finally, the performance connection of each brain region is obtained through integrated learning and EM algorithm. Compared with other methods, the present invention has the following advantages: 1. Map the structural connection to the prior probability space of the performance connection parameters by converting the model, and optimize the model parameters in the subsequent integrated learning, so that the relationship between the structural connection and the performance connection is real 2. Combined with structural connection information, the results of brain activity analysis are more reliable, and it is convenient to discuss the actual situation of the individual.

Description

A Method for Efficient Connectivity Analysis of Each Brain Region Fused with Structural Connectivity

technical field

The invention relates to a calculation method of performance connection of fusion structure connection, which belongs to the fields of medical image methodology and medical image signal processing, and is applicable to the analysis and research of brain activity.

Background technique

Previous studies have found that mental illnesses including depression, schizophrenia, and epilepsy are all related to structural or functional impairments of the brain. The structural information (structural connectivity) and functional information (functional connectivity, functional connectivity) of the brain reflect the condition of the whole brain, and play an important role in the exploration and research of the pathological mechanisms of various mental diseases. Moreover, there is a close relationship between structural connectivity and functional (functional) functional connectivity. Functional (effective) functional connectivity reflects structural connectivity to a certain extent, but does not entirely depend on structural connectivity. Previous studies seldom combine them. If the structural connectivity and functional connectivity are integrated from a multimodal perspective, it will bring great help to the research of various mental diseases.

There have been some preliminary methods and explorations in the study of brain activity combining structural connectivity and functional connectivity, but there are some problems. Method 1: Using the respective advantages of DTI (diffusion tensor imaging) and MEG (magnetoencephalography), conduct research and analysis separately from different angles, and then conduct comprehensive analysis. However, the relationship between the two is not explored, let alone the fusion between them. Method 2: Use the analysis of structural connectivity to constrain and optimize the results of functional connectivity. Although this method uses and analyzes the relationship between the two, it still fails to really integrate the two. Method 3: Use the experimental data of structural connectivity to simulate the function of the brain on a large scale, or study what type of network structure will produce neural activity of specific functions. This method can integrate the two, but still only focuses on exploring the connection between the two.

Contents of the invention

In order to solve the deficiencies in the existing methods, the present invention proposes a new method for analyzing the performance connection of each brain region connected by the fusion structure, and the specific technical scheme is as follows:

The basic idea of the effective connection analysis method for each brain region that integrates structural connections is to extract and analyze the brain structural network in the region of interest, and convert the obtained structural information into the prior probability distribution space of the effective connection parameters; The performance connection model of the Bayesian framework; finally, the performance connection of each brain region is obtained through ensemble learning (Ensemble Learning) and EM algorithm (Expectation-Maximization Algorithm).

The specific steps of this method include:

1. A method for analyzing brain activity by combining structural connections and performance connections, characterized in that the steps include:

1) First use the DTI data of diffusion tensor magnetic resonance imaging to trace the nerve fibers of the whole brain, and establish the structural connection network of the whole brain; in addition, perform 3D source reconstruction on the acquired magnetoencephalogram MEG signals; the 3D source reconstruction of MEG The MSP (Multiple Sparse Priors) method is used, and this step is mainly processed by the SPM8 ( http://www.fil.ion.ucl.ac.uk/spm/ ) software. This method is a method of reconstructing distributed dipole sources using hierarchical Bayesian and empirical Bayesian. Its advantage is that it can automatically select multiple sparse cortical sources with simple prior knowledge.

2) Extract the backbone network in the connection network according to the region of interest that needs to be analyzed, and convert the backbone network into a graph;

Let any brain region of interest in the backbone network be i, and each region of interest ROI(i) in the backbone network is regarded as a node i, and the cortical area of the region of interest represented by node i is S(i);

和

其中，E_f为连接节点i和j的所有纤维，l_f为这些纤维的长度，N_f为神经纤维的数目；所以l为连接两区域的所有神经纤维的平均长度，w反映的是两区域的连接密度；The nerve fibers ROI(i) and ROI(j) connecting two brain regions i and j in the backbone network correspond to the edge E(i, j) connecting nodes i and j, and the length and weight of the edge are respectively

and

Among them, E _f is all fibers connecting nodes i and j, l _f is the length of these fibers, and N _f is the number of nerve fibers; so l is the average length of all nerve fibers connecting two regions, and w reflects the two regions the connection density;

3) Normalize the obtained structural information of each region of interest (the structural information here refers to the length and weight of the edge, corresponding to the length and density of the nerve fibers connecting the region of interest), and convert it to the performance connection The prior probability distribution space of , the conversion model is: ${\mathrm{\Σ}}_{\mathrm{ij}}=\frac{{\mathrm{\Σ}}_0}{1+{\mathrm{\Σ}}_0{e}^{a-{\mathrm{bs}}_{\mathrm{ij}}}}={\mathrm{\α}}_{\mathrm{ij}}^{-1}={\mathrm{\α}}_0+e^{a-{\mathrm{bs}}_{\mathrm{ij}}};$

The performance connection between any brain region i and j obeys the Gaussian distribution N(0, ∑ _ij ), s _ij is the normalized structural connection information, and ∑ ₀ , a, b are the adjustable parameters of the model;

4) The performance connection model based on the variational Bayesian framework is Y=XW+E. In the model Y, W corresponds to the performance connection parameter matrix representing the brain interval, and the matrix W is the autoregressive coefficient parameter matrix;

E is Gaussian noise with a mean value of zero and a precision matrix of Λ, and Λ ~ Γ(b, c); X, Y are the signal of the region of interest after 3D source reconstruction, for a given data set D={X, Y} has: $p\left(\mathrm{D.}\right|W,\mathrm{\Λ})={\left(2\mathrm{\π}\right)}^{-\mathrm{dN}/2}{\left|\mathrm{\Λ}\right|}^2{e}^{-\frac{1}{2}\mathrm{Tr}\left(\mathrm{\Λ}{\mathrm{E.}}_{\mathrm{D.}}\left(W\right)\right)};$ In order to facilitate the analysis of model Y, W is stretched into a vector w, and the distribution of w is as follows:

$p\left(w\right|\left\{{\mathrm{\α}}_k\right\})=\underset{k=1}{\overset{\mathrm{no}}{\mathrm{\Π}}}{\left(\frac{{\mathrm{\α}}_k}{2\mathrm{\π}}\right)}^{1/2}{e}^{-{\mathrm{\α}}_k{\mathrm{E.}}_k\left(w\right)};$ Among them, d is the number of brain area signals, N is the length of the brain area signal sequence, w is the vector stretched by matrix W, n is the number of performance connection parameters, ${\mathrm{E.}}_k\left(w\right)=\frac{1}{2}{w}^T{I}_{k}w;$

The parameters w, Λ and α obey the Gaussian distribution N(0, ∑ _ij );

5) Through the ensemble learning method and the maximum expectation EM algorithm, the efficiency connection Y of each brain region is obtained.

Beneficial effect

Compared with other methods, the present invention has the following advantages: 1. Map the structural connection to the prior probability space of the performance connection parameters by converting the model, and optimize the model parameters in the subsequent integrated learning, so that the relationship between the structural connection and the performance connection is real 2. Combined with structural connection information, the results of brain activity analysis are more reliable, and it is convenient to discuss the actual situation of the individual.

Description of drawings

Figure 1: Schematic diagram of the basic flow of the method;

Figure 2: Schematic diagram of the processing flow of structural connection information;

Figure 3: Schematic diagram of the structural information conversion model diagram;

Figure 4: Schematic diagram of performance connection results.

Detailed ways

Now in conjunction with accompanying drawing, the present invention will be further described:

Whole flow process of the present invention can refer to accompanying drawing 1, and concrete implementation steps are as follows:

1. Establishment of brain structural connection network and preprocessing and source reconstruction of magnetic brain signals

The following steps are required to generate a high-precision structural connection network from MRI diffusion image DTI to the whole brain: (1) Diffusion-weighted image processing, such as eddy current correction, head motion correction, diffusion tensor model fitting, diffusion tensor (2) Segmentation of gray and white matter; (3) Tracing imaging of white matter nerve fiber tracts; (4) Segmentation of cerebral cortex structure and selection of regions of interest; ( 5) Extract the backbone network according to the region of interest that needs to be analyzed, and convert it into a graph. Each ROI(i) in the network can be regarded as a node i, and the cortical area of the region it represents is S(i). The nerve fiber connecting ROI(i) and ROI(j) in the network corresponds to the edge E(i, j) connecting nodes i and j, and the length and weight of the edge are respectively $l=\frac{1}{N_f}{\mathrm{\Σ}}_{f\∈{\mathrm{E.}}_f}{l}_f,$ $w=\frac{2}{S\left(i\right)+S\left(j\right)}{\mathrm{\Σ}}_{f\∈{\mathrm{E.}}_f}\frac{1}{l_f}.$ Among them, E _f is all fibers connecting nodes i and j, l _f is the length of these fibers, N _f is the number of nerve fibers, so l refers to the average length of all nerve fibers connecting two regions, and w reflects the The connection density of the two regions. For the processing flow, refer to FIG. 2 .

After the preprocessing of the magnetic brain signal (including conversion, segmentation, filtering, removal of artifacts, and averaging), the 3D source reconstruction is performed. This processing flow uses the MSP (Multiple Sparse Priors) method. The steps of this method mainly require SPM8 ( http://www.fil.ion.ucl.ac.uk/spm/ ) software for processing. This method is a method of reconstructing distributed dipole sources using hierarchical Bayesian and empirical Bayesian. Its advantage is that it can automatically select multiple sparse cortical sources with simple prior knowledge.

2. Transformation of structural connection information

Normalize the structural information of each region of interest obtained above, and transform it into the prior probability distribution space of the performance connection. The transformation model is:

${\mathrm{<span\; class="notranslate">\; \&Sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; ij</span>}}<span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; \&Sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}}{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; \&Sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}{\mathrm{<span\; class="notranslate">\; e</span>}}^{\mathrm{<span\; class="notranslate">\; a</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; bs</span>}}_{\mathrm{<span\; class="notranslate">\; ij</span>}}}}<span\; class="notranslate">\; =</span>{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}_{\mathrm{<span\; class="notranslate">\; ij</span>}}^{<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}}<span\; class="notranslate">\; =</span>{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; e</span>}}^{\mathrm{<span\; class="notranslate">\; a</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; bs</span>}}_{\mathrm{<span\; class="notranslate">\; ij</span>}}}$

The performance connection between any brain region i and j obeys the Gaussian distribution N(0, ∑ _ij ), s _ij is the normalized structural connection information, ∑ ₀ , a, b are the adjustable parameters of the model, and in the subsequent They are gradually optimized in ensemble learning, so that the relationship between structural connectivity and performance connectivity can be truly reflected. The function diagram of the conversion model is shown in Figure 3.

3. Analysis of performance connection model based on variational Bayesian framework

The efficiency connection model based on the variational Bayesian framework is Y=XW+E. In the model, W is the matrix of autoregressive coefficients, which corresponds to the performance connection parameter matrix of the brain interval, and these parameters obey the Gaussian distribution N(0, ∑ _ij ); E is the Gaussian noise with zero mean and Λ precision matrix, and Λ～Γ (b, c). X, Y are the signal of the region of interest after source reconstruction. For a given data set D={X, Y}:

$\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; D.</span>}<span\; class="notranslate">\; |</span>\mathrm{<span\; class="notranslate">\; W</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; \&Lambda;</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; \&pi;</span>}<span\; class="notranslate">\; )</span>}^{<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; dN</span>}<span\; class="notranslate">\; /</span>\mathrm{<span\; class="notranslate">\; 2</span>}}{<span\; class="notranslate">\; |</span>\mathrm{<span\; class="notranslate">\; \&Lambda;</span>}<span\; class="notranslate">\; |</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}{\mathrm{<span\; class="notranslate">\; e</span>}}^{<span\; class="notranslate">\; -</span>\frac{\mathrm{<span\; class="notranslate">\; 1</span>}}{\mathrm{<span\; class="notranslate">\; 2</span>}}\mathrm{<span\; class="notranslate">\; Tr</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&Lambda;</span>}{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; D.</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; W</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; )</span>}$

In order to facilitate the analysis of the model, the coefficient parameter matrix W is elongated into a vector w, and its distribution is as follows:

$\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; |</span><span\; class="notranslate">\; \{</span>{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}_{\mathrm{<span\; class="notranslate">\; k</span>}}<span\; class="notranslate">\; \}</span><span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; k</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 1</span>}}{\overset{\mathrm{<span\; class="notranslate">\; no</span>}}{\mathrm{<span\; class="notranslate">\; \&Pi;</span>}}}{<span\; class="notranslate">\; (</span>\frac{{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}_{\mathrm{<span\; class="notranslate">\; k</span>}}}{\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; \&pi;</span>}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; /</span>\mathrm{<span\; class="notranslate">\; 2</span>}}{\mathrm{<span\; class="notranslate">\; e</span>}}^{<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}_{\mathrm{<span\; class="notranslate">\; k</span>}}{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; k</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; )</span>}$

Among them, the number of brain area signals in d, N is the length of the signal sequence, w is the vector stretched by the performance connection parameter matrix W, and n is the number of performance connection parameters,

4. To solve the model, an improved EM algorithm is given in this example:

For the above given data set D and parameters θ={w, α, Λ}, the logarithmic argument of model m is

log p(D|m)＝F(θ)+KL(q(θ|D)||p(θ|D, m))

F(θ) is the negative free energy of the model. When q(θ|D)=p(θ|D, m), that is, the approximate posterior probability distribution of the model parameters is equal to the real posterior probability distribution, the model is obtained from the data The lower limit F(θ), the parameters at this time are also the parameters required by the model. However, it is very difficult to solve the parameters directly through the equation. F(θ) is further decomposed to get F(θ)=∫q(θ|D)log p(D |θ, m)dθ-KL(q(θ|D )||p(θ|m)). Given the initial value of the parameter θ and fixing α and Λ, it can be obtained through simplification that when q(w|D)=e ^I(w) , F(θ) takes the maximum value, where I(w)=∫∫q (Λ|D)q(α|D)log(p(D|w,Λ)p(w|α))dαdΛ. Update parameter w and fix w and Λ, similarly, α can be solved, then update parameter α and fix α and w to solve Λ. Through this improved EM (Expectation-Maximization) algorithm, iteratively continues until convergence, and the parameters of the model and the performance connection between each brain area are obtained (results are shown in Figure 4).

Claims (3)

1. A method for analyzing the efficacy linkage between brain regions connected by a fusion structure is characterized by comprising the following steps:

1) firstly, carrying out nerve fiber tracking of the whole brain by using diffusion tensor magnetic resonance imaging DTI data, and establishing a structural connection network of the whole brain; in addition, 3D source reconstruction is carried out on the acquired magnetoencephalogram MEG signal;

2) extracting a backbone network in the connection network according to the region of interest to be analyzed, and converting the backbone network into a graph;

setting any interesting brain area in a backbone network as i, regarding each interesting region ROI (i) in the backbone network as a node i, and setting the cortex area of the interesting region represented by the node i as S (i);

the nerve fibers ROI (i) and ROI (j) connecting the two brain areas i and j in the trunk network correspond to edges E (i, j) connecting the nodes i and j, and the length and the weight of the edges are respectively

And

wherein E is_fTo connect all fibres of nodes i and j, l_fThe length of these fibers, N_fIs the number of nerve fibers; so l is the average length of all nerve fibers connecting the two regions, and u reflects the connection density of the two regions; f is a nerve fiber connecting nodes i and j of the brain;

3) normalizing the obtained structural information of each region of interest, and converting the structural information into a prior probability distribution space connected with efficiency, wherein the conversion model is as follows:

the structural information refers to the length and weight of the edge, corresponding to the length and density of the nerve fibers connecting the regions of interest;

the efficiency connection between any brain regions i and j obeys a Gaussian distribution N (0, Sigma)_ij)，s_ijFor normalized structural connection information, sigma₀A and b are adjustable parameters of the model; wherein,

α₀the method is characterized in that the method is a constant parameter, parameters a and b control the degree of constraint of structural connection on prior probability distribution variance, wherein a is an inflection point parameter, and b is a slope parameter;

4) the efficiency connection model based on the variational Bayesian framework is Y = XW + E, in the model Y, W corresponds to an efficiency connection parameter matrix representing a brain interval, and the matrix W is an autoregressive coefficient parameter matrix; x is a time sequence matrix;

e is Gaussian noise with mean value of zero and precision matrix of Lambda, and the precision matrix Lambda corresponds to element Lambda_ijObey the following gaussian distribution: lambda_ij～N(0,Σ_ij) (ii) a X, Y are the region of interest signals after 3D source reconstruction, for a given data set D = { X, Y } there is: <math> <mrow> <mi>p</mi> <mrow> <mo>(</mo> <mi>D</mi> <mo>|</mo> <mi>W</mi> <mo>,</mo> <mi>&Lambda;</mi> <mo>)</mo> </mrow> <mo>=</mo> <msup> <mrow> <mo>(</mo> <mn>2</mn> <mi>&pi;</mi> <mo>)</mo> </mrow> <mrow> <mo>-</mo> <mi>dN</mi> <mo>/</mo> <mn>2</mn> </mrow> </msup> <msup> <mrow> <mo>|</mo> <mi>&Lambda;</mi> <mo>|</mo> </mrow> <mn>2</mn> </msup> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <mn>1</mn> <mn>2</mn> </mfrac> <mi>Tr</mi> <mrow> <mo>(</mo> <msub> <mi>&Lambda;E</mi> <mi>D</mi> </msub> <mrow> <mo>(</mo> <mi>W</mi> <mo>)</mo> </mrow> <mo>)</mo> </mrow> </mrow> </msup> <mo>,</mo> </mrow> </math> E_D(Ｗ)＝(Y-XW)^T(Y-XW) is the unnormalized error covariance matrix, | Λ purple²Is the square of the matrix Λ determinant; to facilitate the analysis of model Y, W is elongated into a vector W whose distribution is as follows:

wherein d is the number of brain region signals, N is the length of brain region signal sequence, W is the vector of matrix W stretching, N is the number of efficacy connection parameters,

α_kthe parameter precision of the kth efficiency connection parameter; i is_kIs the k unit matrix;

the parameter w_ij,λ_ij,α_ijObeying Gaussian distribution N (0, Sigma)_ij)；

5) And (4) solving the efficiency connection Y between each brain region through an integrated learning method and a maximum expectation EM algorithm.

2. The method of claim 1, wherein in step 1), the step of establishing a network of structural connections throughout the brain comprises:

101) the DTI image is processed, and the DTI image processing method comprises the following steps:

eddy current correction and cephalotaxis correction, fitting of a diffusion tensor model, calculation of a diffusion tensor and an anisotropic value and conversion of a standard space;

102) dividing grey white matter;

103) white matter nerve fiber tract tracking imaging;

104) segmentation of cerebral cortical structures and selection of regions of interest.

3. The method as claimed in claim 1, wherein in the step 5), for the data set D and the parameters θ = { w, α, Λ }, the logarithm of the model Y is: logp (D | m) = F (θ) + KL (q (θ | D) | p (θ | D, m)), where F (θ) is the negative free energy of the model Y, and when q (θ | D) = p (θ | D, m), i.e., the approximate posterior probability distribution of the model parameters is equivalent to the true posterior probability distribution, the model Y takes a lower limit F (θ) on the theoretical data, the parameter θ is exactly the parameter required by the model Y; where m is the model that produces the data set D; KL (q (θ | D) | p (θ | D, m)) represents the KL Divergence of q (θ | D) and p (θ | D, m), which is an abbreviation of Kullback-Leibler difference (Kullback-Leibler Divergence), also called Relative Entropy (Relative Entropy), and has the physical meaning: in the same event space, the event space of the probability distribution q (theta | D) is encoded by the probability distribution p (theta | D, m), and the encoding length of each basic event symbol is increased by a certain number of bits on average; is a representation of distance;

further decomposing F (theta) to obtain:

F(θ)=∫q(θ|D)logp(D|θ,m)dθ-KL(q(θ|D)||p(θ|m))；

given the initial value of the parameter θ and fixing α and Λ therein, we obtain by simplification: when q (w | D) = e^I(w)When F (θ) takes a maximum value, where i (w) =; updating the parameter w, fixing w and Λ, solving alpha in the same way, then updating the parameter alpha, fixing alpha and w, and solving Λ;

through the algorithm, iteration is continuously carried out until convergence, and the efficiency connection between the parameters of the model and each brain region is solved.

Images