基于脑功能网络的帕金森患者DBS术后效果预测系统及方法与流程

基于脑功能网络的帕金森患者dbs术后效果预测系统及方法
技术领域
1.本发明涉及一种基于脑功能网络的帕金森患者dbs(脑深部电刺激术，deep brain stimulation，简称dbs)术后效果预测系统及预测方法，使用计算机通过帕金森患者fmri影像信息融合基于图论的指标构建脑功能网络，实现对帕金森患者dbs术后效果的预测，属于计算机应用技术领域。


背景技术：

2.目前，随着人口老龄化趋势发展，中老年人罹患帕金森病的数量逐年增加。作为神经系统退行性疾病之一的帕金森病，严重威胁着老年人生活质量。
3.在日常生活中，帕金森病常见的运动症状包含运动迟缓、冻结步态和四肢静止性震颤等，以及伴有认知功能障碍、睡眠障碍、抑郁或嗅觉变差等常见症状，给患者及其护理者带来了严重负担。手术治疗作为缓解帕金森病的有效手段之一，在提升患者生活质量中发挥了重要作用。
4.帕金森病作为典型的神经退行性疾病，其发展过程伴随着患者神经元与髓鞘的丧失。脑深部电刺激术(deep brain stimulation，dbs)通过在大脑特定区域插入电极，并使电极发射电信号对周围组织进行刺激，以达到改善人体神经信号传输模式的目的。因其具有可逆性及可调节性的特点，被认为是治疗药物控制不佳帕金森病的理想疗法。但是，该类手术需在人体大脑内放置电极及体外调控设备，具有较高的经济成本，且属于有创手术，会面临出血、感染等相关风险，给患者带来不可避免的经济压力身体压力。
5.因此，开发一种帕金森dbs术后效果辅助评估系统，有效预测帕金森患者dbs术后的恢复效果，作为患者是否进行该手术的参考之一，十分必要。
6.目前，现有的帕金森dbs术后恢复效果预测系统所使用的技术手段，大多是基于对临床指标的统计分析，未结合帕金森病作为神经系统退行性疾病的特点，缺乏有效的生物标志物，预测效果不佳。尽管有一些基于fmri的研究表明刺激区域到其他脑区的连接对患者dbs术后恢复效果有预测作用，但其未考虑全脑功能网络作为一个整体的协调作用，导致系统的预测效果较差。


技术实现要素：

7.本发明的目的是为了克服现有的帕金森患者dbs术后恢复效果预测系统存在的不足，提出一种基于脑功能网络的帕金森患者dbs术后效果预测系统及方法，通过建立两种新型动态脑功能网络，融合基于图论的指标，能够更好捕捉帕金森患者的脑内功能连接信息，提高dbs术后效果预测的准确率。
8.本发明的创新点在于：提出一种新的预测系统，将基于图论的全脑功能网络与机器学习算法进行融合，实现在术前基于脑功能网络的帕金森患者dbs术后恢复效果预测。通过分析系统预测恢复效果与真实恢复效果之间的误差评估指标，大幅提升预测效果。
9.本发明设计的两种新型动态脑功能网络，具有以下优点：
10.采用滑动窗口策略构建动态脑功能网络，具有低时间粒度，考虑动态性因素，通过不同窗口构建的功能矩阵捕捉大脑在一段时间内动态的变化。
11.随机动态脑功能网络使用大小不同的窗口，具有多尺度的优势，能更为全面地提取不同小时间段内的动态特征。
12.基于图论的分析方法，依据脑图谱对大脑进行roi划分，计算任意两个roi间的最短路径长度时序序列或全局效率时序序列的皮尔逊相关系数作为最终脑网络中这两个roi的脑连接强度，引入了各个roi节点信息传输能力的相关关系。
13.本发明采用的技术方案如下：
14.一种基于脑功能网络的帕金森患者dbs术后效果预测系统，包括数据预处理模块、动态脑功能网络构建模块、特征筛选模块、预测模型训练模块和预测模块。
15.其中，数据预处理模块，负责对预处理后的患者mri图像进行脑区roi划分，以及生成roi节点的bold时序信号(txt格式)。该模块的输入为患者的fmri图像。
16.动态脑功能网络构建模块，负责根据bold时序信号构建固定动态脑功能网络和随机动态脑功能网络。该模块的输入为数据预处理模块生成的bold时序信号。
17.特征筛选模块，负责使用随机森林算法对输入的脑功能网络进行特征筛选。该模块的输入为动态脑功能网络构建模块生成的脑功能网络。
18.预测模型训练模块，负责进行回归模型训练。该模块的输入为特征筛选模块得到的脑功能网络的特征，以及相对应患者的术前udprs-iii与术后udprs-iii评分。
19.预测模块，负责利用回归模型进行患者dbs术后效果预测。该模块的输入为待预测患者的fmri图像。
20.上述模块之间的连接关系为：
21.数据预处理模块的输出端与动态脑功能网络构建模块的输入端相连；
22.动态脑功能网络构建模块的输出端与特征筛选模块的输入端相连；
23.特征筛选模块的输出端与预测模型训练模块、预测模块的输入端分别相连；
24.预测模型训练模块的输出端与预测模块的输入端相连。
25.本系统的工作过程如下：
26.首先，向数据预处理模块输入多个病人的fmri图像数据。数据预处理模块在进行数据预处理后，得到病人bold时间序列，发送给脑功能网络构建模块。
27.然后，脑功能网络构建模块根据病人bold时间序列，构建固定动态脑功能网络和随机动态脑功能网络，并输入至特征筛选模块。
28.之后，特征筛选模块利用随机森林算法，对病人的固定动态脑功能网络和随机动态脑功能网络进行特征筛选，得出病人重要的脑网络特征。
29.然后，预测模型训练模块根据病人的脑网络特征与该病人对应的术前udprs-iii和术后udprs-iii评分，进行回归模型训练，并将训练好的回归模型输入预测模块，用于后续预测使用。
30.当进行新患者术后效果预测时，向数据预处理模块输入患者的fmri图像数据，然后，依次进行数据预处理、构建脑功能网络、特征筛选，在预测模块中利用训练好的回归模型，对患者的术后恢复情况进行预测。
31.进一步地，基于上述系统，本发明提出了一种基于脑功能网络的帕金森患者dbs术
后效果预测方法，包括以下步骤：
32.步骤1：对原始mri图像数据进行图像预处理(例如，进行图像时间校正、图像头动校正、图像空间标准化、图像空间平滑、图像线性趋势消除、图像降噪和脑图谱选择，等等)。对预处理完后的mri图像进行脑区roi的划分，并生成roi节点的bold信号，可以使用gretna软件等实现。
33.步骤2：构建动态脑功能网络，包括构建固定动态脑功能网络和随机动态脑功能网络。
34.具体地，包括以下步骤：
35.步骤2.1：构建固定动态脑功能网络。包括以下步骤：
36.步骤2.1.1：选取长度为l的滑动窗口，按照固定长度为s的步长，将每个roi的bold时序信号划分成n个有重叠的子序列，每个子序列即为一个窗口。
37.步骤2.1.2：分别计算每个窗口中每个roi的最短路径长度。将每个窗口作为有k个节点的加权网络g，对于两个节点i和j，在它们的所有通路中，权重累加最小的路径被称为两节点间的最短路径，记为对于共有k节点的加权网络g中的某节点i，其最短路径长度定义如下：
[0038][0039]
步骤2.1.3：计算所有roi中任意两个roi x和roi y之间的最短路径长度序列的皮尔逊相关系数r，得出最终的固定动态脑功能网络：
[0040][0041]
其中，xa表示向量x的均值，ya表示向量y的均值。
[0042]
步骤2.2：构建随机动态脑功能网络。包括以下步骤：
[0043]
步骤2.2.1：选取长度范围在l_min～l_max的滑动窗口，按照长度范围在s_min～s_max的步长，将bold时序信号划分成n个有重叠的、不同的窗口间无时序，窗口内有时序的子序列，每个子序列即为一个窗口。
[0044]
步骤2.2.2：分别计算每个窗口中每个roi的节点全局效率：
[0045][0046]
其中，k表示总节点数，g表示该网络，表示节点i的全局效率，表示节点i和j之间的加权最短距离。
[0047]
步骤2.2.3：计算所有roi中任意两个roi全局效率序列的皮尔逊相关系数，得出最终的随机动态脑网络。
[0048]
步骤3：采用随机森林算法，对所构建的脑网络进行特征筛选。
[0049]
根据随机森林算法计算出每个脑连接特征的重要性，将重要性映射到该连接的左右脑区后进行归一。对每个脑连接特征重复上述操作，得到该网络下每个脑区的重要性权
重。
[0050]
步骤4：训练dbs术后预测模型。
[0051]
采用常见的线性回归模型与非线性回归模型进行dbs术后效果预测，并对预测效果进行比较。具体地，包括以下步骤：
[0052]
步骤4.1：进行内层交叉验证。
[0053]
在每个机器学习模型的训练过程中，采用网格搜索的方式对模型进行优化，即，不断循环遍历设置的所有候选参数，并比对获得的每一种参数组合及其对应的预测值，当预测值与真实值最接近时，则此时的模型即为最优模型。
[0054]
步骤4.2：进行外层交叉验证。将源数据进行分割，通过k折交叉验证，保证每次使用k-1份数据进行训练、1份数据进行验证。最后，将k次的验证结果取平均值，得到最终的预测效果。
[0055]
步骤4.3：分别训练最小二乘法线性回归、岭回归、套索回归、支持向量机回归和极端梯度回归模型。
[0056]
在嵌套交叉验证内层循环基于网格搜索算法，确定每个模型的最优参数，在外层循环进行交叉验证，获得每个受试者的术后udprs-iii提升率的预测值，并计算其与真实值的皮尔逊相关性系数、平均绝对误差(mean absolute error,mae)、均方根误差(root mean square error,rmse)，作为预测效果评估指标。
[0057][0058]
其中，updrs-iii
ratio
表示udprs-iii提升率，updrs-iii
preoperation
表示患者术前的udprs-iii值，updrs-iii
postoperation
表示患者术后的udprs-iii值。
[0059]
从中获得最优的机器学习预测模型，同时对该最优参数机器学习回归模型进行存储，以备预测模块使用。
[0060]
步骤5：利用预测模型，进行患者术后效果预测。
[0061]
输入患者的fmri图像数据，然后，按照步骤1至步骤3所述方法，依次进行数据预处理、构建脑功能网络和特征筛选，最后利用步骤4训练好的预测模型，对患者的术后恢复情况进行预测。
[0062]
有益效果
[0063]
本发明方法，对比现有技术，具有以下优点：
[0064]
(1)本发明系统，采用两种脑功能网络，均考虑动态性因素且具有低时间粒度，采用滑动窗口策略构建，通过不同窗口构建的功能矩阵捕捉大脑在一段时间内动态的变化。
[0065]
(2)本发明系统采用的随机动态脑功能网络具有多尺度的优势，能够更加全面地提取不同小时间段内的动态特征。
[0066]
(3)本发明系统采用了固定动态脑功能网络和随机动态脑功能网络均融合图论的分析方法，结合动态脑活动、roi节点的最短路径长度或全局效率、roi之间的相关关系三种信息，形成的脑网络具有稳定且有效的脑连接特征。
附图说明
[0067]
图1为本发明提出的基于脑功能网络的帕金森患者dbs术后效果预测系统组成示意图。
[0068]
图2为本发明提出的用于帕金森dbs术后恢复效果预测的固定动态脑网络构建方法的流程图。
[0069]
图3为本发明提出的用于帕金森dbs术后恢复效果预测的随机动态脑网络构建方法的流程图。
[0070]
图4为固定动态脑功能网络构建流程示意图。
[0071]
图5为随机动态脑功能网络构建流程示意图。
[0072]
图6为嵌套交叉验证示意图。
[0073]
图7为基于静态脑功能网络的术后udprs-iii提升率预测值与真实值拟合图。
[0074]
图8为基于均方根的动态脑功能网络的术后udprs-iii提升率预测值与真实值拟合图。
[0075]
图9为基于锥形滑动窗口的动态脑功能网络的术后udprs-iii提升率预测值与真实值拟合图。
[0076]
图10为基于固定动态脑功能网络的术后udprs-iii提升率预测值与真实值拟合图。
[0077]
图11为基于随机动态脑功能网络的术后udprs-iii提升率预测值与真实值拟合图。
具体实施方式
[0078]
下面结合附图和实施例对本发明做进一步详细说明。
[0079]
一种基于脑功能网络的帕金森患者dbs术后效果预测系统，如图1所示，包括数据预处理模块、动态脑功能网络构建模块、特征筛选模块、预测模型训练模块和预测模块。
[0080]
首先，向数据预处理模块输入多个病人的fmri图像数据。数据预处理模块在进行数据预处理后，得到病人bold时间序列，发送给脑功能网络构建模块。
[0081]
然后，脑功能网络构建模块根据病人bold时间序列，构建固定动态脑功能网络和随机动态脑功能网络，并输入至特征筛选模块。
[0082]
之后，特征筛选模块利用随机森林算法，对病人的固定动态脑功能网络和随机动态脑功能网络进行特征筛选，得出病人重要的脑网络特征。
[0083]
然后，预测模型训练模块根据病人的脑网络特征与该病人对应的术前udprs-iii和术后udprs-iii评分，进行回归模型训练，并将训练好的回归模型输入预测模块，用于后续预测使用。
[0084]
当进行新患者术后效果预测时，向数据预处理模块输入患者的fmri图像数据，然后，依次进行数据预处理、构建脑功能网络、特征筛选，在预测模块中利用训练好的回归模型，对患者的术后恢复情况进行预测。
[0085]
一种基于脑功能网络的帕金森患者dbs术后效果预测方法，如图2所示，包括以下步骤：
[0086]
步骤1：数据预处理。
[0087]
步骤1.1：因mri图像由按层扫描的方式获得，本方法选择大脑中间层作为标准，通过gretna中的时间校正步骤保证一个扫描周期内各层扫描时间的统一。
[0088]
步骤1.2：本方法评估了受试者的头动状况，将每个受试者图像序列的第一个图像作为基准，并把剩余图像向基准配准，调整因受试者头动造成的不同时刻的图像错位。
[0089]
步骤1.3：本方法考虑了标准化精度与时间成本，使用统一分割算法，借助受试者的t1结构像，通过将受试者功能图像配准到标准脑模板的方式来对齐所有受试者功能图像。
[0090]
步骤1.4：本方法中，受试者的fmri图像经过一个4mm
×
4mm
×
4mm的半高宽(即纵坐标峰值一半所对应的横坐标波长宽度)高斯滤波器，使得图像的信噪比有了显著的提升。
[0091]
步骤1.5：由于被扫描者的生理变化或仪器不稳定，rsfmri数据可能会出现系统性的信号增加或减少。线性趋势消除可以消除仪器在进行扫描时所产生的偏移对后续数据处理的影响。本方法在线性趋势消除步骤使用线性拟合策略，减小了在获取数据时因仪器本身原因产生的误差。
[0092]
步骤1.6：因自发性脑功能活动在低频段最为活跃，本实验使0.01hz～0.1hz的盒状滤波器进行滤波，以消除因机器或受试者心跳等产生的其他频率噪声。
[0093]
步骤1.7：bna由中国科学院脑网络组研究中心提出，该图谱将人脑分为246个区域(36个皮层下核团亚区，210个脑皮层区)，同时该246个脑区可被分为24个脑回或7个脑叶。该图谱的数据来源于40名健康成人的脑影像，且其绘制融合了多模态mri数据，使得分区更为精细。因为实验数据集为54人的静息态脑功能影像，且希望基于精细的脑区划分构建脑功能网络，所以本方法选择bna模板将246个脑区作为roi构建脑功能网络。
[0094]
步骤2：构建脑功能网络。
[0095]
本实施例中，一共构建5种脑功能网络，其中，后两种(固定动态与随机动态脑网络)为本方法所提出用于帕金森dbs术后恢复效果预测的新型脑功能网络。前3种使用目前常见的网络构建方式，用于与本方法进行预测效果对比。
[0096]
步骤2.1：构建静态功能网络
[0097]
步骤2.2.1：计算246个roi中任意两个时序bold信号间的皮尔逊相关系数，可得出每个受试者的静态脑功能网络。
[0098]
步骤2.2：构建基于均方根的动态脑功能网络
[0099]
步骤2.2.1：选取长度为50的滑动窗口，按照固定长度为1的步长，将bold时序信号划分成191个有重叠的时序子序列，每个子序列即为一个窗口。
[0100]
步骤2.2.2：通过公式(5)循环计算第i个roi与第j个roi之间的均方根值rms，n为窗口个数，表示第n个滑动窗口中第i个roi与第j个roi间的连接强度，形成基于均方根的动态脑功能网络。
[0101][0102]
其中，dfc
ij
表示基于均方根的动态脑功能网络中第i行、第j列位置。
[0103]
步骤2.3：构建基于锥形滑动窗口的动态脑功能网络。
[0104]
步骤2.3.1：选取长度为50的锥形滑动窗口，其中锥形滑动窗口为锥形余弦窗口，窗口的参数相关系数r值为0.5。按照固定长度为1的步长，将时间序列bold信号划分成191个有重叠的时序子序列，每个子序列即为一个窗口。
[0105]
步骤2.3.2：计算每个锥形窗口内246个roi中任意两个roi间的皮尔逊相关系数，形成一个246
×
246的相关性矩阵。根据公式(6)对该矩阵执行费舍尔变换，得到每个窗口变换后的脑连接矩阵。
[0106][0107]
其中，z表示费舍尔变换得到的结果。
[0108]
步骤2.3.3：分别计算第i个roi与第j个roi之间的连接在191个窗口中的方差，得到最终基于锥形滑动窗口的动态脑功能网络。
[0109]
步骤2.4：构建固定动态脑功能网络。
[0110]
步骤2.4.1：选取长度为50的滑动窗口，按照固定长度为1的步长，将每个roi的bold时序信号划分成191个有重叠的子序列，每个子序列即为一个窗口。
[0111]
步骤2.4.2：通过公式(1)分别计算246个roi在每个小窗口中的最短路径长度。
[0112]
步骤2.4.3：通过公式(2)计算246个roi中任意两个roi的最短路径长度序列的皮尔逊相关系数r，得出最终的固定动态脑功能网络。
[0113]
构建流程如图3所示。
[0114]
步骤2.5：构建随机动态脑功能网络
[0115]
步骤2.5.1：选取长度范围在25～75的滑动窗口，按照长度范围在1～90的步长，将bold时序信号划分成191个有重叠的、不同的窗口间无时序，窗口内有时序的子序列，每个子序列即为一个窗口。
[0116]
步骤2.5.2：通过公式分别计算246个roi在每个窗口中的节点全局效率
[0117]
步骤2.5.3：计算246个roi中任意两个roi全局效率序列的皮尔逊相关系数，得出最终的随机动态脑网络。
[0118]
随机动态脑功能网络构建流程如图4所示。
[0119]
步骤3：采用随机森林算法，对所构建的脑网络进行特征筛选。
[0120]
步骤3.1：根据随机森林算法计算出每个脑连接特征的重要性。
[0121]
步骤3.2：将重要性映射到该连接的左右脑区后进行归一。
[0122]
步骤3.3：对每个脑连接特征重复上述步骤，得到该网络下每个脑区的重要性权重。
[0123]
步骤4：dbs术后预测模型训练。
[0124]
应用常见线性回归模型与非线性回归模型进行dbs术后效果预测并对预测效果进行比较。
[0125]
步骤4.1：内层交叉验证。在每个机器学习模型的训练过程中，采用网格搜索的方式对模型进行优化，即不断循环遍历设置的所有候选参数，并通过比对获得的每一种参数组合及其对应的预测值，当预测值与真实值最为接近时，此时的模型即为最优模型。以ridge模型为例，其参数包括：alpha、solver等。假设alpha有5种可能，solve有3种可能，如果只对两参数进行优化，对其情况进行穷举。
[0126]
步骤4.2：外层交叉验证。本实验将54人数据进行分割，且通过十折交叉验证保证每次使用九份数据进行训练，一份数据进行验证。最后将十次的验证结果进行平均，得到最终的预测效果。嵌套交叉验证示意图如图5所示。
[0127]
步骤4.3：分别使用最小二乘法线性回归、岭回归、套索回归、支持向量机回归、极端梯度回归进行训练。在嵌套交叉验证的内层循环中使用网格搜索算法得到每个模型的最优参数，然后在外层循环中实行十折交叉验证，根据公式(4)获得每个受试者的术后udprs-iii提升率的预测值，并计算其与真实值的皮尔逊相关性系数、平均绝对误差(mae)、均方根误差(rmse)作为预测效果评估指标。从中获得最优的机器学习预测模型，存储最优参数机器学习回归模型。
[0128]
步骤5：利用预测模型，进行患者术后效果预测。
[0129]
输入患者的fmri图像数据，然后，按照步骤1至步骤3所述方法，依次进行数据预处理、构建脑功能网络和特征筛选，最后利用训练好的预测模型，对患者的术后恢复情况进行预测。
[0130]
当脑连接特征来源于静态脑功能网络时，各个回归模型的皮尔逊相关系数、mae和rmse如表1所示。各个回归模型术后udprs-iii提升率的预测值与实际值的拟合如图6所示。
[0131]
表1基于静态脑功能网络的机器学习模型预测结果
[0132][0133]
当脑连接特征来源于基于均方根的动态脑功能网络时，各个回归模型的皮尔逊相关系数、平均绝对误差和均方根误差如表2所示，每个模型术udprs-iii提升率的预测值与实际值的拟合如图7所示。
[0134]
表2基于均方根的动态脑功能网络的机器学习模型预测结果
[0135][0136]
当脑连接特征来源于基于锥形滑动窗口的动态脑功能网络时，各个回归模型的皮尔逊相关系数、mae和rmse如表3所示，每个模型术后udprs-iii提升率的预测值与实际值的拟合如图8所示。
[0137]
表3基于锥形滑动窗口的动态脑功能网络的机器学习模型预测结果
[0138][0139]
当脑连接特征来源于固定动态脑功能网络时，各个回归模型的皮尔逊相关系数、mae和rmse如表4所示，每个模型术后udprs-iii提升率的预测值与实际值的拟合如图9所示。
[0140]
表4基于固定动态脑功能网络的机器学习模型预测结果
[0141][0142]
当脑连接特征来源于随机动态脑功能网络时，各个回归模型的皮尔逊相关系数、mae和rmse如表5所示，每个模型术后udprs-iii提升率的预测值与实际值的拟合如图10所示。
[0143]
表5基于随机动态脑功能网络的机器学习模型预测结果
[0144][0145]
根据本实验所用的数据库，相比于本方法提出的两种动态网络，静态功能网络、基于均方根的动态脑功能网络与基于锥形滑动窗口的动态脑功能网络预测效果均比较差。
[0146]
基于均方根的动态脑功能网络与基于锥形滑动窗口脑动态网络引入了脑连接的动态特征，但由于其没有进行深层次的信息挖掘，创建的动态脑网络不适合本文预测帕金森患者dbs术后恢复效果的目标，使得最终预测效果不好。
[0147]
在固定动态脑功能网络与随机动态脑功能网络的构建过程中，借鉴前两种对比动态脑网络的优势，同样采取了滑动窗口策略来提取脑连接的动态变化。其引入了图论中的节点最短路径长度指标与全局效率指标，从节点入手，计算246个roi的最短路径长度或全局效率，从而获得了在某一个窗口下的每个roi的信息传输能力。然后，计算任意两个roi间的最短路径长度时序序列或全局效率时序序列的皮尔逊相关系数作为最终脑网络的脑连接强度，引入了各个roi节点信息传输能力的相关关系。所以该两种脑功能网络预测准确度与其他相比方法更高。