CN110363096B - 一种基于隐马尔可夫模型的大脑时间信号处理方法 - Google Patents

Abstract

本发明公开了一种基于隐马尔可夫模型的大脑时间信号处理方法，包括步骤：1)将采集所得功能磁共振图像进行预处理，并且要求所采集的功能磁共振图像有相同的回波时间；2)对预处理后每个被试脑区的时间序列进行中心化和标准化，使用健康被试分别训练每个脑区的隐马尔可夫模型；3)根据所求的隐马尔可夫模型参数求解每个被试每个脑区序列的似然值，并根据不同被试时间序列长度不同进行缩放，得到每个被试每个脑区的一个特征；4)使用步骤3)得到每个被试一个脑区的特征后，使用SVM‑RFE方法分类两组被试。基于本发明所提出的分析方法对孤独症谱系障碍的功能磁共振数据进行分类，准确率已经达到80.37％，相较于构建脑网络，深度学习等其他方法有很大的提升。

Description

一种基于隐马尔可夫模型的大脑时间信号处理方法

技术领域

本发明属于医学信号处理技术领域，特别涉及一种基于隐马尔可夫模型的大脑时间信号处理方法。

背景技术

功能磁共振技术是根据脱氧血红蛋白和氧合血红蛋白在高强磁场下所表现出的顺磁性和抗磁性的不同来测量神经元活动所引发的血液动力学变化。功能磁共振技术以其非倾入性，无损伤，较高的空间分辨率和较高的时间分辨率等优点，已经成为脑科学研究中的主流技术，其中静息态功能磁共振图像更是成为脑科学研究的基础信号。静息态功能磁共振信号(rest-fMRI)能反应在静息状态下大脑自发的神经活动产生的低频波动信号，对于解释大脑神经活动的机理有重要意义。对大脑静息态功能磁共振影像进行特征提取(大脑连接网络，深度学习等)已经在认知、记忆等大脑高级活动的探索和大脑神经类疾病的诊断方面表现出了巨大的潜力。同时科学界也在探索更新、更符合大脑神经活动的特征提取方法，以达到对大脑功能的探索。

目前最广泛使用的大脑分析方法是基于脑区的时间序列构建大脑全连接网络(使用的方法为：皮尔逊相关，偏相关)以及使用l₁范数正则化构建大脑连接的稀疏化网络，但目前的大脑分析方法主要存在以下不足：

1.基于皮尔逊相关或基于偏相关构建的全连接大脑网络会导致大脑出现较多的虚假连接，不能反应大脑神经活动中脑区稀疏连接的客观事实。

2.基于l₁范数正则化构建大脑连接的稀疏化网络从一定程度上改善了基于相关的大脑网络结构，但是由于l₁范数正则项的系数是根据数学准则(交叉验证，贝叶斯信息准则等)从每个个体的数据中计算得到，这导致个体脑网络的稀疏程度不一致，大脑网络结构不稳定。

3.现有的分析方法仅使用时间序列求解两脑区之间的关系，造成时间尺度上信息的浪费。

当使用现有分析方法所得到的特征对功能磁共振图像分类时，往往由于分类特征过多，导致分类器过拟合现象严重。以90个感兴趣区域(ROI)为例皮尔逊相关与偏相关构建的脑网络矩阵为实对称矩阵，因此会有

个特征，当使用l₁范数正则化构建的稀疏化脑网络，由于个体网络结构的不同，因此会有90*90＝8100个特征。

发明内容

为了解决在大脑时间信号处理过程中遇到的上述问题，本发明提供了一种基于隐马尔科夫模型(HMM)的大脑时间信号处理方法，该方法构建的隐马尔可夫模型参数使用所有健康被试同一脑区的时间序列训练模型参数，极大的利用了时间尺度的信息，并最终根据健康被试的训练模型，使用极大似然的方法求解出被试每个脑区的特征，从而使得所得到的特征能够明显的区分两类大脑时间信号，特征数量少而精，对提升分类的准确率，为大脑时间信号处理和临床神经类疾病的诊断提供新的思路和理论支持。

本发明采用如下技术方案来实现的：

一种基于隐马尔可夫模型的大脑时间信号处理方法，包括以下步骤：

1)将采集所得功能磁共振图像进行预处理，并且要求所采集的功能磁共振图像有相同的回波时间；

2)对预处理后每个被试脑区的时间序列进行中心化和标准化，使用健康被试分别训练每个脑区的隐马尔可夫模型；

3)根据所求的隐马尔可夫模型参数求解每个被试每个脑区序列的似然值，并根据不同被试时间序列长度不同进行缩放，得到每个被试每个脑区的一个特征；

4)使用步骤3)得到每个被试一个脑区的特征后，使用SVM-RFE方法分类两组被试。

本发明进一步的改进在于，步骤1)中，数据预处理方法步骤为：

步骤101对采集功能磁共振同时所采集的T1加权结构项数据进行去头皮，组织分割，将大脑组织分割为灰质、白质和脑脊液三部分；

步骤102：对功能磁共振图像进行大脑组织掩模处理，并将得到的掩模与每个时间点的功能磁共振图像矩阵对应位置相乘，去除脑组织之外的其他结构；

步骤103：对功能磁共振图像进行带通滤波，去除大脑低频活动范围之外的噪声信号；

步骤104：将功能磁共振图像分别与T1图像的灰质、白质和脑脊液部分配准，再将功能磁共振图像配准到MNI152标准模板上；

步骤105：使用步骤104配准后的功能磁共振图像提取大脑脑区的时间序列，提取时间序列的方法为在对应大脑区域计算功能磁共振图像灰度值的平均值。

本发明进一步的改进在于，步骤105中，采用的大脑脑区模板为DMN模板，AAL模板，或者CC200模板。

本发明进一步的改进在于，步骤2)中，模型训练方法步骤为：

步骤201：由于功能磁共振信号的涨落反应大脑神经活，信号实际值的大小并无绝对的物理含义，因此对每个被试每个脑区的时间序列进行中心化和标准化处理，X(t)＝{x₁(t),x₂(t),...,x_R(t)}表示一个被试提取所得的R个脑区的时间序列，对每个脑区的时间序列做中心化和标准化处理如下：

步骤202：

Y_P(t)表示所有N个健康被试第p个脑区的时间序列，对于隐马尔可夫模型来说则相当于是N个观测序列，学习目标是隐马尔可夫模型参数λ＝(A,B,π)，A为隐状态之间的概率转移矩阵，B为观测概率矩阵，π为初始状态概率向量，设定I为观测所对应的隐状态的序列，隐马尔可夫模型使用如下概率模型表示：

步骤203：使用EM算法求解上述概率模型。

本发明进一步的改进在于，步骤3)中，被试脑区的似然值计算步骤为：

步骤301：通过训练得到了R个脑区的隐马尔可夫模型λ¹，λ²，…λ^p，…λ^R，从而计算出第n个被试的第p个脑区的观测值在当前隐马尔可夫模型下出现的概率如下：

步骤302：由于不同被试之间功能磁共振扫描的时间长度不同，而可能造成似然值的大小不同，为了补偿由于扫描时间长度引起的差别，使用时间长度因子

对每个被试的似然值进行缩放，第n个被试的特征f_n如下：

本发明进一步的改进在于，步骤4)中，SVM-RFE分类方法步骤如下：

步骤401：输入M个被试的特征F＝[f₁，f₂，...f_i，...，，f_M]；

步骤402：对所有被试特征进行归一化处理；

步骤403：使用SVM分类器对两组特征进行分类，记录当前准确率；

步骤404：判断是否还有特征剩余，若还有特征剩余则去掉SVM分类器中权重最小的特征，重复步骤403，若没有特征剩余，则对记录的准确率进行排序；

步骤405：选取分类准确率最高的特征组合脑区作为分析当前大脑时间序列信号的重要脑区，以及用于区分不同大脑时间信号的类别。

本发明具有如下有益的技术效果：

在大脑时间序列信号(功能磁共振，脑电，PET等)分析中，使用最广泛的是基于皮尔逊相关构建大脑连接网络，但其忽略了大脑脑区连接稀疏性的客观事实，造成过多的虚假连接，导致信息冗余。并且极大地压缩了所采集数据的时间尺度信息，造成信息的浪费。

本发明提出一种基于隐马尔可夫模型的大脑时间信号的处理方法，采用大脑时间序列统一的特征提取方法，首先根据大脑相同位置脑区的时间序列建立该区域大脑活动的隐马尔可夫模型，不仅充分利用了时间尺度的信息，而且基于人群训练所得的模型更能从本质上反应该类人群的大脑活动特征，而不是仅从个体神经活动提取具有个体差异的大脑连接网络。其次隐马尔可夫模型是一种基于贝叶斯概率的统计模型，以大脑激活程度作为隐藏状态，以当前物理手段实际测得的数据作为观测值，两者结合能更全面的反映大脑活动的多种可能。最后本发明所提取的特征能够应用在更具有临床意义的大脑精神类疾病分类过程中，基于本发明所提出的分析方法对孤独症谱系障碍的功能磁共振数据进行分类，准确率已经达到80.37％，相较于构建脑网络，深度学习等其他方法有很大的提升。

附图说明

图1是本发明中功能磁共振数据的预处理流程。

图2是基于本发明的方法对73名孤独症谱系障碍患者和49名正常被试的数据进行分类的ROC曲线。

具体实施方式

以下结合附图和实施例对本发明做出进一步的说明。

本发明提供的一种基于隐马尔可夫模型的大脑时间信号处理方法(以功能磁共振图像为例)，包括以下步骤：

1)将采集所得功能磁共振图像进行预处理，并且要求所采集的功能磁共振图像有相同的回波时间(repetition time，TR)；数据预处理方法步骤为：

步骤101对采集功能磁共振同时所采集的T1加权结构项数据进行去头皮，组织分割，将大脑组织分割为灰质、白质和脑脊液三部分；

步骤102：对功能磁共振图像进行大脑组织掩模处理，并将得到的掩模与每个时间点的功能磁共振图像矩阵对应位置相乘，去除脑组织之外的其他结构；

步骤103：对功能磁共振图像进行带通滤波，去除大脑低频活动范围之外的噪声信号；

步骤104：将功能磁共振图像分别与T1图像的灰质、白质和脑脊液部分配准，再将功能磁共振图像配准到MNI152标准模板上；

步骤105：使用步骤104配准后的功能磁共振图像提取大脑脑区的时间序列，通常可使用的大脑脑区模板为默认模式网络(Default mode network,DMN)，AnatomicalAutomatic Labeling(AAL)分区模板，Craddock 200(CC200)模板，提取时间序列的方法为在对应大脑区域计算功能磁共振图像灰度值的平均值。

2)对预处理后每个被试脑区的时间序列进行中心化和标准化，使用健康被试分别训练每个脑区的隐马尔可夫模型。模型训练方法步骤为：

步骤201：由于功能磁共振信号的涨落反应大脑神经活，信号实际值的大小并无绝对的物理含义，因此对每个被试每个脑区的时间序列进行中心化和标准化处理，X(t)＝{x₁(t),x₂(t),...,x_R(t)}表示一个被试提取所得的R个脑区的时间序列，对每个脑区的时间序列做中心化和标准化处理如下：

步骤202：

Y_P(t)表示所有N个健康被试第p个脑区的时间序列，对于隐马尔可夫模型来说则相当于是N个观测序列，学习目标是隐马尔可夫模型参数λ＝(A,B,π)，A为隐状态之间的概率转移矩阵，B为观测概率矩阵，π为初始状态概率向量，设定I为观测所对应的隐状态的序列，隐马尔可夫模型使用如下概率模型表示：

步骤203：使用EM算法求解上述概率模型。

3)根据所求的隐马尔可夫模型参数求解每个被试每个脑区序列的似然值，并根据不同被试时间序列长度不同进行缩放，得到每个被试每个脑区的一个特征。被试脑区的似然值计算步骤为：

步骤301：通过训练得到了R个脑区的隐马尔可夫模型λ¹，λ²，…λ^p，…λ^R，从而计算出第n个被试的第p个脑区的观测值在当前隐马尔可夫模型下出现的概率如下：

步骤302：由于不同被试之间功能磁共振扫描的时间长度不同，而可能造成似然值的大小不同，为了补偿由于扫描时间长度引起的差别，使用时间长度因了

对每个被试的似然值进行缩放，第n个被试的特征f_n如下：

4)使用步骤3)得到每个被试一个脑区的特征后，使用SVM-RFE方法分类两组被试。SVM-RFE分类方法步骤如下：

步骤401：输入M个被试的特征F＝[f₁，f₂，...f_i，...，，f_M]；

步骤402：对所有被试特征进行归一化处理；

步骤403：使用SVM分类器对两组特征进行分类，记录当前准确率；

步骤404：判断是否还有特征剩余，若还有特征剩余则去掉SVM分类器中权重最小的特征，重复步骤403，若没有特征剩余，则对记录的准确率进行排序；

步骤405：选取分类准确率最高的特征组合脑区作为分析当前大脑时间序列信号的重要脑区，以及用于区分不同大脑时间信号的类别。

实施例

本实施例以功能磁共振数据为例，脑电等其他大脑时间信号处理方法相同。

步骤1：按照如图1所示的方式对功能磁共振数据进行预处理；

步骤2：使用AAL90模板提取大脑90个脑区的时间序列；

步骤3：每个被试的信号长度为T，对每个被试每个脑区的时间序列中心化和标准化处理，

x(t)表示任意被试任意脑区的时间序列；

步骤4：将所有健康被试统一脑区的时间序列作为一个集合，训练对应脑区的隐马尔科夫模型λ＝(A，B，π)，使用EM算法求解之后，模型参数具体计算方式如下：

A＝[a_ij]10x10表示大脑隐藏的激活状态转移矩阵，在此假设大脑的激活状态分为10个等级，那么A为10*10的方阵，B＝[b_j(k)]表示观测概率矩阵，π＝[π_i]表示初始时刻的状态概率向量。

在未开始训练之前先初始化λ⁰，接下来将所有被试的同一脑区的数据x(t)，带入下面公式，计算得到新的λ，再将新的λ应用于下一次迭代，直到λ最终收敛。

其中N＝49表示所有健康被试的数目，γ_t(i)表示模型在t时刻处于状态i的概率，ξ_t(i，j)表示模型在t时刻处于状态i，t+1时刻处于状态j的概率，以上参数可由向前向后算法计算得到，不再赘述。对每个脑区训练各自的模型，最终得到90个隐马尔可夫模型。

步骤4：根据每个脑区的模型，计算每个被试每个脑区的似然值，以一个被试一个脑区的时间序列x(t)为例，λ^r为当前脑区的隐马尔可夫模型。

表示在t＝1的时刻，隐状态为i时观测为x(1)的概率，按照这个方式递推

那么最终似然值由下面公式得到：

步骤5：利用每个被试计算所得的90个特征，利用线性SVM方法构建分类器，每次去掉一个权重最小特征，最后将分类准确率最高的特征排序，选择最优的分类组合。最终得到的分类准确率为80.37％，十折交叉验证ROC曲线下面平均为75.1％，ROC曲线如图2所示。

Claims (1)

1.一种基于隐马尔可夫模型的大脑时间信号处理方法，其特征在于，包括以下步骤：

1)将采集所得功能磁共振图像进行预处理，并且要求所采集的功能磁共振图像有相同的回波时间；数据预处理方法步骤为：

步骤101对采集功能磁共振同时所采集的T1加权结构项数据进行去头皮，组织分割，将大脑组织分割为灰质、白质和脑脊液三部分；

步骤102：对功能磁共振图像进行大脑组织掩模处理，并将得到的掩模与每个时间点的功能磁共振图像矩阵对应位置相乘，去除脑组织之外的其他结构；

步骤103：对功能磁共振图像进行带通滤波，去除大脑低频活动范围之外的噪声信号；

步骤104：将功能磁共振图像分别与T1图像的灰质、白质和脑脊液部分配准，再将功能磁共振图像配准到MNI152标准模板上；

步骤105：使用步骤104配准后的功能磁共振图像提取大脑脑区的时间序列，提取时间序列的方法为在对应大脑区域计算功能磁共振图像灰度值的平均值，采用的大脑脑区模板为DMN模板，AAL模板，或者CC200模板；

2)对预处理后每个被试脑区的时间序列进行中心化和标准化，将所有健康被试单个脑区的时间序列作为一个集合，使用健康被试分别训练每个脑区的隐马尔可夫模型；模型训练方法步骤为：

步骤201：由于功能磁共振信号的涨落反应大脑神经活，信号实际值的大小并无绝对的物理含义，因此对每个被试每个脑区的时间序列进行中心化和标准化处理，X(t)＝{x₁(t),x₂(t),...,x_R(t)}表示一个被试提取所得的R个脑区的时间序列，对每个脑区的时间序列做中心化和标准化处理如下：

步骤202：

Y_P(t)表示所有N个健康被试第p个脑区的时间序列，对于隐马尔可夫模型来说则相当于是N个观测序列，学习目标是隐马尔可夫模型参数λ＝(A,B,π)，A为隐状态之间的概率转移矩阵，B为观测概率矩阵，π为初始状态概率向量，设定I为观测所对应的隐状态的序列，隐马尔可夫模型使用如下概率模型表示：

步骤203：使用EM算法求解上述概率模型；

3)根据所求的隐马尔可夫模型参数求解每个被试每个脑区序列的似然值，并根据不同被试时间序列长度不同进行缩放，得到每个被试每个脑区的一个特征；被试脑区的似然值计算步骤为：

步骤301：通过训练得到了R个脑区的隐马尔可夫模型λ¹,λ²,···λ^p,···λ^R，从而计算出第n个被试的第p个脑区的观测值在当前隐马尔可夫模型下出现的概率如下：

步骤302：由于不同被试之间功能磁共振扫描的时间长度不同，而可能造成似然值的大小不同，为了补偿由于扫描时间长度引起的差别，使用时间长度因子

对每个被试的似然值进行缩放，第n个被试的特征f_n如下：

4)使用步骤3)得到每个被试一个脑区的特征后，使用SVM-RFE方法分类两组被试，SVM-RFE分类方法步骤如下：

步骤401：输入M个被试的特征F＝[f₁,f₂,…f_i,…,f_M]；

步骤402：对所有被试特征进行归一化处理；

步骤403：使用SVM分类器对两组特征进行分类，记录当前准确率；

步骤404：判断是否还有特征剩余，若还有特征剩余则去掉SVM分类器中权重最小的特征，重复步骤403，若没有特征剩余，则对记录的准确率进行排序；

步骤405：选取分类准确率最高的特征组合脑区作为分析当前大脑时间序列信号的重要脑区，以及用于区分不同大脑时间信号的类别。

Images