CN102592278A - 一种基于脑功能成像的脑内多区域协同竞争分析方法 - Google Patents

Abstract

本发明公开了一种基于脑功能成像的脑内多区域协同竞争分析方法，包括以下步骤：输入脑核磁共振图像实验数据序列；判断脑激活区；确定感兴趣区；感兴趣区协同竞争分析。本发明把脑图像中每个体素点分别在任务和控制时的信号作为两个样本，利用t检验，判断两个样本是否具有显著性差异，从而判断该体素点在实验任务中是否被激活。本发明划定了感兴趣区中t值最强点为中心点的簇，并把簇内的体素点信号的平均值，作为该感兴趣区的值，从而将感兴趣区量化。本发明根据感兴趣区的量化特征，两两计算不同感兴趣区的相关系数，利用相关系数的值确定感兴趣区间协同或竞争的关系。

Description

一种基于脑功能成像的脑内多区域协同竞争分析方法

技术领域

本发明涉及一种利用三维大脑磁共振图像发现各脑区间相互作用和协同竞争关系的分析方法，属于医学图像处理领域。

背景技术

随着功能性磁共振成像(functional Magnetic Resonance Imaging，fMRI)的不断成熟与发展，越来越多的学者运用fMRI研究大脑的高级神经活动。标准的fMRI实验设计之一称为组块设计(Blocked Designs)。实验有两种不同的条件，分别形成任务组块(Task blocks)和控制组块(Control block)，在任务组块中需要被试者完成某项设定的认知任务，在控制组块中被试者则处于休息状态，什么也不需要做，每个组块持续一定时间(一般为14至20秒)，并在连续的扫描中交替出现。在整个的实验过程中，每隔一段固定时间(其时间分辨率依机器的扫描速度而定，一般为2至3秒)，会采到一幅全脑加权像，它们构成一个离散的时间序列，可记为{f(n)}(n为正整数)。因此，组块持续时间可以使用扫描图像的数量代替，{f(n)}的长度l即总图像数由实验的周期数L与周期长度(每周期取得的图像数)所决定，即l＝L×(c+t)，其中，c和t分别为每个实验周期内所获取的完成控制组块的图像数和完成任务组块的的图像数。对于自变量n的每一个取值，{f(n)}都唯一地确定一幅全脑加权像，记为BOLD(n)。它既可能是属于控制组块的BOLDc(n)，又可能为属于任务组块的BOLDt(n)。在这类的实验中，通过对比脑扫描数据在任务和控制之间的差异而确定结果，得出许多有意义的结果，如语言区的精确定位、情绪区的新发现、记忆的脑神经网络等等。

但是当前对于从整体和动态角度考察认知过程中参与的脑区以及脑区间的反应模式和时空关系，并建立脑内信息加工的相关网络与模型，相对研究较少。由于认知在信息加工过程中，往往会有多个区域参与了活动，而且彼此之间存在着协同竞争关系。如果只通过统计分析得到不同脑区的激活图谱，很难解决在认知研究中激活脑区间的相互作用及活动规律(与任务相关的特定网络与模式)等问题。

发明内容

为解决现有技术存在的上述问题，本发明要设计一种基于脑功能成像的脑内多区域协同竞争分析方法，在获得感兴趣区和脑激活图的基础上，进行功能连接分析，发现脑区间的相互作用和协同竞争的关系。

为了实现上述目的，本发明的技术方案如下：一种基于脑功能成像的脑内多区域协同竞争分析方法，包括以下步骤：

A、输入脑核磁共振图像实验数据序列

脑fMRI图像实验数据序列，在整个的实验过程中有两种不同的条件，分别形成任务组块和控制组块，每个组块持续一定时间，并在连续的扫描中交替出现；每隔一段固定时间，机器会采到一幅全脑加权像，它们构成一个离散的时间序列，记为{f(n)}，n为正整数；因此，组块持续时间可以使用扫描图像的数量代替，{f(n)}的长度l即总图像数由实验的周期数L与周期长度所决定，即l＝L×(c+t)其中，c和t分别为每个实验周期内所获控制和任务的图像数；对于自变量n的每一个取值，{f(n)}都唯一地确定一幅全脑加权像，记为BOLD(n)；它既可能是属于控制组块的BOLDc(n)，又可能为属于任务组块的BOLDt(n)；所述的持续一定时间为14至20秒，所述的一段固定时间为2至3秒；

B、判断脑激活区

实验过程中的脑图像BOLD(n)可以看作一个四维数组，对于一幅脑图像的任何一个体素点，都可记为bold(x，y，z，n)，其中x、y和z是体素点的空间位置，n是这幅脑图像的序号；根据脑图像属于任务或控制组块的不同，还要细分为boldt(x，y，z，n_t)和boldc(x，y，z，n_c)；对于一个确定位置的单一体素点，x、y和z为一组常数，则boldt(x，y，z，n_t)和boldc(x，y，z，n_c)两个四维数组均成为一维数组，如表1所示；它们分别代表同一体素点在不同任务或控制下的信号；因此，可将其看作两个来自同一总体、接受不同处理的样本；

表1单一体素点t检验的比较对象

对两个样本进行t检验计算t值，其计算公式为：

$\mathrm{tValue}(x,y,z)=\frac{\stackrel{\‾}{\mathrm{boldt}(x,y,z,n_t)}-\stackrel{\‾}{\mathrm{boldc}(x,y,z,n_c)}}{{\mathrm{SD}}_c}---\left(1\right)$

式中，

和

分别表示两个样本的均值，SD_c为其合并标准差，它可按下式计算：

${\mathrm{SD}}_c=\sqrt{\frac{\underset{t=1}{\overset{T}{\mathrm{\Σ}}}{[\mathrm{boldt}(x,y,z,n_t)-\stackrel{\‾}{\mathrm{boldt}(x,y,z,n_t)}]}^2+\underset{c=1}{\overset{C}{\mathrm{\Σ}}}{[\mathrm{boldc}(x,y,z,n_c)-\stackrel{\‾}{\mathrm{boldc}(x,y,z,n_c)}]}^2}{T+C-2}(\frac{1}{T}+\frac{1}{C})}---\left(2\right)$

两组值显著性差异，不仅取决于它们的t值，而且与主观确定的显著性水平及样本自由度F的大小有关；所述的自由度F为T+C-2；

$\mathrm{tValue}(x,y,z)=\left\{\begin{array}{c}\mathrm{tValue}(x,y,z)=\mathrm{tValue}(x,y,z),\mathrm{iftValue}(x,y,z)\&GreaterEqual;t_{P/2}(T+C-2)\\ \mathrm{tValue}(x,y,z)=0,\mathrm{iftValue}(x,y,z)<t_{P/2}(T+C-2)\end{array}\right.---\left(3\right)$

如公式(3)所示，对于任意一个体素点，当它的t值落在拒绝域内，则认为其被激活，否则将不予考虑；

C、确定感兴趣区

在判断为激活的区域中，根据脑解剖结构预先圈定一些需要重点观测的区域，称为感兴趣区，在每个感兴趣区中，仍然有许多体素点存在，要找出在该感兴趣区中t值最强点，以及它所在的簇；对簇的定义为：设体素点(x_n，y_n，z_n)为t值最强点，则由27个体素点(x_n-1，y_n-1，z_n-1)(x_n-1，y_n-1，z_n)(x_n-1，y_n-1，z_n+1)(x_n-1，y_n，z_n-1)(x_n-1，y_n，z_n)(x_n-1，y_n，z_n+1)(x_n-1，y_n+1，z_n-1)(x_n-1，y_n+1，z_n)(x_n-1，y_n+1，z_n+1)(x_n，y_n-1，z_n-1)(x_n，y_n-1，z_n)(x_n，y_n-1，z_n+1)(x_n，y_n，z_n-1)(x_n，y_n，z_n)(x_n，y_n，z_n+1)(x_n，y_n+1，z_n-1)(x_n，y_n+1，z_n)(x_n，y_n+1，z_n+1)(x_n+1，y_n-1，z_n-1)(x_n+1，y_n-1，z_n)(x_n+1，y_n-1，z_n+1)(x_n+1，y_n，z_n-1)(x_n+1，y_n，z_n)(x_n+1，y_n，z_n+1)(x_n+1，y_n+1，z_n-1)(x_n+1，y_n+1，z_n)(x_n+1，y_n+1，z_n+1)组成的立方体为以体素点(x_n，y_n，z_n)为中心的簇；在以后根据体素点(x_n，y_n，z_n)坐标位置取信号值bold时，如果没有特殊说明，都取的是27个体素点信号的平均值；

D、感兴趣区协同竞争分析

对每一个感兴趣区内t值最强体素点的时间序列S＝S(x₀，y₀，z₀，n)进行正规化；对于控制组块的正规化如下式：

$\mathrm{Sc}(x_0,y_0,z_0,n)=\frac{[\mathrm{boldc}(x_0,y_0,z_0,n)-\stackrel{\&OverBar;}{\mathrm{boldc}(x_0,y_0,z_0,n)}]}{\stackrel{\&OverBar;}{\mathrm{boldc}(x_0,y_0,z_0,n)}]}---\left(4\right)$

而对于任务组块的正规化如下式：

$\mathrm{St}(x_0,y_0,z_0,n)=\frac{[\mathrm{boldt}(x_0,y_0,z_0,n)-\stackrel{\&OverBar;}{\mathrm{boldc}(x_0,y_0,z_0,n)}]}{\stackrel{\&OverBar;}{\mathrm{boldc}(x_0,y_0,z_0,n)}]}---\left(5\right)$

这里

是体素点(x₀，y₀，z₀)在控制组块中fMRI信号强度的均值，把该均值作为正规化的基线；

最后计算体素点任务与控制之间的相对变化差异，即任务时第一幅脑图像体素点(x₀，y₀，z₀)的St减去控制时第一幅脑图的Sc，并以此类推计算ΔS，如公式(6)：

ΔS(x₀，y₀，z₀，n)＝St(x₀，y₀，z₀，n)-Sc(x₀，y₀，z₀，n)        (6)

对于一位被试者，假设研究人员共为他划定了J个要观测的感兴趣区、在某一个感兴趣区j中，计算其t值最强点的ΔS_j，则可组成如t_j＝[ΔS_j(x₀，y₀，z₀，1)，ΔS_j(x₀，y₀，z₀，2)，......ΔS_j(x₀，y₀，z₀，T)]序列，在J个序列中，两两求相关系数，比如第i个感兴趣区和第j个感兴趣区的相关系数为：

$c_{\mathrm{ij}}=\mathrm{corr}(t_i,t_j)=\frac{\mathrm{Cov}(t_i,t_j)}{\sqrt{D\left(t_i\right)D\left(t_j\right)}}---\left(7\right)$

其中Cov(t_i，t_j)是两个不同周期数据的协方差，D(t_i)和D(t_j)分别为t_i和t_j的方差；

如果两个脑感兴趣区之间相关系数趋近于1或-1，这意味着它们之间存在协同或竞争的关系；如果两个脑区之间没有表现出与任务相关的激活变化，那么这两个脑感兴趣区的相关系数趋近于0，这意味着即使它们之间有解剖连接，也没有功能上的连接。

与现有技术相比，本发明具有以下有益效果：

1、本发明把脑图像中每个体素点分别在任务和控制时的信号作为两个样本，利用t检验，判断两个样本是否具有显著性差异，从而判断该体素点在实验任务中是否被激活。

2、本发明划定了感兴趣区中t值最强点为中心点的簇，并把簇内的体素点信号的平均值，作为该感兴趣区的值，从而将感兴趣区量化。

3、本发明根据感兴趣区的量化特征，两两计算不同感兴趣区的相关系数，利用相关系数的值确定感兴趣区间协同或竞争的关系。

附图说明

本发明共有附图5张，其中：

图1是本发明的基本流程图；

图2是脑fMRI成像数据序列示意图；

图3是组块设计实验示意图；

图4是确定感兴趣区内t值最强点簇及均值示意图；

图5是感兴趣区协同竞争关系举例示意图。

具体实施方式

下面结合附图对本发明进行进一步地描述，步骤方法如图1所示：

步骤一：脑fMRI图像实验数据序列，在整个的实验过程中有两种不同的条件，分别形成任务组块(Task blocks)和控制组块(Control block)，每个组块持续一定时间(一般为14至20秒)，并在连续的扫描中交替出现，如图2所示。每隔一段固定时间(其时间分辨率依机器的扫描速度而定，一般为2至3秒)，机器会采到一幅全脑加权像，如图3所示，它们构成一个离散的时间序列，可记为{f(n)}(n为正整数)。因此，组块持续时间可以使用扫描图像的数量代替，f(n)}的长度l即总图像数由实验的周期数L与周期长度(每周期取得的图像数)所决定，即l＝L×(c+t)其中，c和t分别为每个实验周期内所获控制和任务的图像数。

步骤二：脑激活区判断，通过比对脑区在任务和控制时信号的相对变化程度，可以判断该脑区是否激活。由于t检验是判定两组样本间差别是否显著常用的数理统计方法，因此使用t检验处理fMRI数据，判断脑区是否激活。

步骤三：确定感兴趣区，根据脑解剖结构，确定J个感兴趣区，在每个感兴趣区中，仍然有许多体素点存在，要找出在该感兴趣区中t值最强点，以及它所在的簇，取簇内体素点的平均值，作为感兴趣区中t值最强点的值。举例如图4，在某项实验中圈定了5个感兴趣区，在每个感兴趣区提取t值最强点，并取以该点为中心的簇，计算簇内体素点的平均值。

步骤四：感兴趣区协同竞争分析，将划分的J个感兴趣区中的t值最强点的序列，两两求相关系数，如果两个脑感兴趣区之间表现出协同或竞争，则相关系数趋近于1或-1，如果没有表现出与任务相关的激活变化，那么这两个脑区的相关系数趋近于0，举例如图5，左边圆圈1-5分别表示五个感兴趣区的左脑部分，右边圆圈1-5分别表示五个感兴趣区的右脑部分，圆圈之间连线的数字表示脑区间的相关系数。

从图5可以直观的看出感兴趣区之间的相关性。脑区构成的回路相关性强弱以及正负相关性，从整体上呈现了脑区之间的有效性连接。可以看到感兴趣区3和感兴趣区5的变化趋势非常一致，相关系数达到0.9以上，很好的说明了在处理任务时感兴趣区3和感兴趣区5的相关性高，即协同性好。由此可见。通过本发明可以直观的看出感兴趣区之间的协同或竞争的关系。

Claims (1)

1.一种基于脑功能成像的脑内多区域协同竞争分析方法，其特征在于：包括以下步骤：

A、输入脑核磁共振图像实验数据序列

脑fMRI图像实验数据序列，在整个的实验过程中有两种不同的条件，分别形成任务组块和控制组块，每个组块持续一定时间，并在连续的扫描中交替出现；每隔一段固定时间，机器会采到一幅全脑加权像，它们构成一个离散的时间序列，记为{f(n)}，n为正整数；因此，组块持续时间可以使用扫描图像的数量代替，{f(n)}的长度l即总图像数由实验的周期数L与周期长度所决定，即l＝L×(c+t)其中，c和t分别为每个实验周期内所获控制和任务的图像数；对于自变量n的每一个取值，{f(n)}都唯一地确定一幅全脑加权像，记为BOLD(n)；它既可能是属于控制组块的BOLDc(n)，又可能为属于任务组块的BOLDt(n)；所述的持续一定时间为14至20秒，所述的一段固定时间为2至3秒；

B、判断脑激活区

实验过程中的脑图像BOLD(n)可以看作一个四维数组，对于一幅脑图像的任何一个体素点，都可记为bold(x，y，z，n)，其中x、y和z是体素点的空间位置，n是这幅脑图像的序号；根据脑图像属于任务或控制组块的不同，还要细分为boldt(x，y，z，n_t)和boldc(x，y，z，n_c)；对于一个确定位置的单一体素点，x、y和z为一组常数，则boldt(x，y，z，n_t)和boldc(x，y，z，n_c)两个四维数组均成为一维数组，如表1所示；它们分别代表同一体素点在不同任务或控制下的信号；因此，可将其看作两个来自同一总体、接受不同处理的样本；

表1单一体素点t检验的比较对象

对两个样本进行t检验计算t值，其计算公式为：

$\mathrm{tValue}(x,y,z)=\frac{\stackrel{\&OverBar;}{\mathrm{boldt}(x,y,z,n_t)}-\stackrel{\&OverBar;}{\mathrm{boldc}(x,y,z,n_c)}}{{\mathrm{SD}}_c}---\left(1\right)$

式中，

和

分别表示两个样本的均值，SD_c为其合并标准差，它可按下式计算：

${\mathrm{SD}}_c=\sqrt{\frac{\underset{t=1}{\overset{T}{\mathrm{\&Sigma;}}}{[\mathrm{boldt}(x,y,z,n_t)-\stackrel{\&OverBar;}{\mathrm{boldt}(x,y,z,n_t)}]}^2+\underset{c=1}{\overset{C}{\mathrm{\&Sigma;}}}{[\mathrm{boldc}(x,y,z,n_c)-\stackrel{\&OverBar;}{\mathrm{boldc}(x,y,z,n_c)}]}^2}{T+C-2}(\frac{1}{T}+\frac{1}{C})}---\left(2\right)$

两组值显著性差异，不仅取决于它们的t值，而且与主观确定的显著性水平及样本自由度F的大小有关；所述的自由度F为T+C-2；

$\mathrm{tValue}(x,y,z)=\left\{\begin{array}{c}\mathrm{tValue}(x,y,z)=\mathrm{tValue}(x,y,z),\mathrm{iftValue}(x,y,z)\&GreaterEqual;t_{P/2}(T+C-2)\\ \mathrm{tValue}(x,y,z)=0,\mathrm{iftValue}(x,y,z)<t_{P/2}(T+C-2)\end{array}\right.---\left(3\right)$

如公式(3)所示，对于任意一个体素点，当它的t值落在拒绝域内，则认为其被激活，否则将不予考虑；

C、确定感兴趣区

在判断为激活的区域中，根据脑解剖结构预先圈定一些需要重点观测的区域，称为感兴趣区，在每个感兴趣区中，仍然有许多体素点存在，要找出在该感兴趣区中t值最强点，以及它所在的簇；对簇的定义为：设体素点(x_n，y_n，z_n)为t值最强点，则由27个体素点(x_n-1，y_n-1，z_n-1)(x_n-1，y_n-1，z_n)(x_n-1，y_n-1，z_n+1)(x_n-1，y_n，z_n-1)(x_n-1，y_n，z_n)(x_n-1，y_n，z_n+1)(x_n-1，y_n+1，z_n-1)(x_n-1，y_n+1，z_n)(x_n-1，y_n+1，z_n+1)(x_n，y_n-1，z_n-1)(x_n，y_n-1，z_n)(x_n，y_n-1，z_n+1)(x_n，y_n，z_n-1)(x_n，y_n，z_n)(x_n，y_n，z_n+1)(x_n，y_n+1，z_n-1)(x_n，y_n+1，z_n)(x_n，y_n+1，z_n+1)(x_n+1，y_n-1，z_n-1)(x_n+1，y_n-1，z_n)(x_n+1，y_n-1，z_n+1)(x_n+1，y_n，z_n-1)(x_n+1，y_n，z_n)(x_n+1，y_n，z_n+1)(x_n+1，y_n+1，z_n-1)(x_n+1，y_n+1，z_n)(x_n+1，y_n+1，z_n+1)组成的立方体为以体素点(x_n，y_n，z_n)为中心的簇；在以后根据体素点(x_n，y_n，z_n)坐标位置取信号值bold时，如果没有特殊说明，都取的是27个体素点信号的平均值；

D、感兴趣区协同竞争分析

对每一个感兴趣区内t值最强体素点的时间序列S＝S(x₀，y₀，z₀，n)进行正规化；对于控制组块的正规化如下式：

$\mathrm{Sc}(x_0,y_0,z_0,n)=\frac{[\mathrm{boldc}(x_0,y_0,z_0,n)-\stackrel{\&OverBar;}{\mathrm{boldc}(x_0,y_0,z_0,n)}]}{\stackrel{\&OverBar;}{\mathrm{boldc}(x_0,y_0,z_0,n)}]}---\left(4\right)$

而对于任务组块的正规化如下式：

$\mathrm{St}(x_0,y_0,z_0,n)=\frac{[\mathrm{boldt}(x_0,y_0,z_0,n)-\stackrel{\&OverBar;}{\mathrm{boldc}(x_0,y_0,z_0,n)}]}{\stackrel{\&OverBar;}{\mathrm{boldc}(x_0,y_0,z_0,n)}]}---\left(5\right)$

这里是体素点(x₀，y₀，z₀)在控制组块中fMRI信号强度的均值，把该均值作为正规化的基线；

最后计算体素点任务与控制之间的相对变化差异，即任务时第一幅脑图像体素点(x₀，y₀，z₀)的St减去控制时第一幅脑图的Sc，并以此类推计算ΔS，如公式(6)：

ΔS(x₀，y₀，z₀，n)＝St(x₀，y₀，z₀，n)-Sc(x₀，y₀，z₀，n)        (6)

对于一位被试者，假设研究人员共为他划定了J个要观测的感兴趣区、在某一个感兴趣区j中，计算其t值最强点的ΔS_j，则可组成如t_j＝[ΔS_j(x₀，y₀，z₀，1)，ΔS_j(x₀，y₀，z₀，2)，......ΔS_j(x₀，y₀，z₀，T)]序列，在J个序列中，两两求相关系数，比如第i个感兴趣区和第j个感兴趣区的相关系数为：

$c_{\mathrm{ij}}=\mathrm{corr}(t_i,t_j)=\frac{\mathrm{Cov}(t_i,t_j)}{\sqrt{D\left(t_i\right)D\left(t_j\right)}}---\left(7\right)$

其中Cov(t_i，t_j)是两个不同周期数据的协方差，D(t_i)和D(t_j)分别为t_i和t_j的方差；

如果两个脑感兴趣区之间相关系数趋近于1或-1，这意味着它们之间存在协同或竞争的关系；如果两个脑区之间没有表现出与任务相关的激活变化，那么这两个脑感兴趣区的相关系数趋近于0，这意味着即使它们之间有解剖连接，也没有功能上的连接。

Images