一种基于功能子网络特征的认知功能状态预测方法

本发明涉及生物医学信息领域，尤其涉及一种基于功能子网络特征的认知功能状态预测方法。

背景技术：

1、作为思维和意识器官，大脑的结构和功能都十分复杂；数以千亿计的神经元通过突触形成了错综复杂的连接性结构网络，而大脑执行功能无论简单与否总是依赖于多个脑区之间的广泛交互，形成了复杂的功能网络，因此从网络的角度来研究大脑给我们提供了一个新的视角。

2、影像学技术的不断成熟，图论以及复杂网络研究的快速发展，为脑科学研究提供了强有力的工具。人脑可以被划分为多个脑区，脑区功能的异常改变会导致认知功能状态的改变，同时脑区之间的连接关系也发生相应的改变。传统的预测认知功能状态的方法是从全局的大脑功能网络出发，忽略了大脑功能网络的模块化特征；或者是集中在研究单一层次上的模块性，没有提取多个层次的子网络对认知功能状态进行分析。

3、一些研究显示，在认知功能状态不同的脑网络中，并不是所有的脑区都具有明显的差异，认知功能状态的变化可能仅仅存在于某些局部模块中，而这些局部模块可能在全局网络中并不显著，忽略了脑功能网络的局部模块化结构会导致预测结果的不确定性。

技术实现思路

1、本发明的目的是为了提供一种基于功能子网络特征的认知功能状态预测方法，提高对认知功能状态预测的准确性。

2、为解决以上技术问题，本发明的技术方案为：

3、一种基于功能子网络特征的认知功能状态预测方法，包括：

4、步骤1：获取每位受试者的认知功能状态量表评分；采集每位受试者的功能磁共振图像并对图像进行预处理，基于预处理后的功能磁共振图像构建大脑功能网络；

5、步骤2：将大脑功能网络划分为多个子网络，通过递归聚类细分当前级别的每个子网络获得下一级的子网络，最终得到分层子网络；

6、步骤3：利用图论算法计算每个功能子网络拓扑属性参数的取值，并计算拓扑属性参数对应的受试者工作特征曲线下面积与认知功能状态量表之间的相关系数；

7、步骤4：使用统计学方法筛选出相关性较高的拓扑属性参数，其对应的受试者工作特征曲线下面积作为功能子网络的特征合集；

8、步骤5：寻找支持向量回归机的高斯径向基核函数中参数的最优值，利用特征合集对支持向量回归机进行优化；

9、步骤6：利用支持向量回归拟合功能子网络特征合集与认知功能状态量表评分，根据两者之间的映射关系预测认知功能状态。

10、进一步地，步骤1具体包括：

11、步骤1.1：对所有受试者进行蒙特利尔认知功能评估量表测试，评估他们的认知功能状态；

12、步骤1.2：采用梯度回波的单次激发平面回波成像序列扫描受试者大脑横断面成像，得到静息态功能磁共振图像；

13、步骤1.3：对功能磁共振图像进行预处理，预处理步骤包括：时间校正、头动校正、空间配准、空间标准化和平滑滤波操作；

14、步骤1.4：采用标准大脑分区模板划分脑区，提取每个脑区内部各体素在不同时间点上的血氧水平依赖强度，再将各脑区内所有体素的时间序列进行算术平均，得到每个脑区内所有体素的平均时间序列；

15、步骤1.5：采用皮尔逊相关分析方法，计算脑区之间平均时间序列的皮尔逊相关系数ρx,y：

16、

17、其中，n表示全脑分割的区域个数；xi表示脑区x在第i个时间点上的数值，表示脑区x在整个扫描期间内所有时间点的平均值；yi表示脑区y在第i个时间点上的数值，表示脑区y在整个扫描期间内所有时间点的平均值；

18、步骤1.6：根据相关系数得到平均时间序列的相关关联矩阵，并将相关关联矩阵取绝对值后得到大脑功能连接矩阵r：

19、

20、其中，ρx,y表示脑区x和脑区y之间平均时间序列的皮尔逊相关系数。

21、进一步地，步骤2具体包括：

22、步骤2.1：对于每位受试者的大脑功能网络，用相似性图g＝(v,e)来表示，其中v＝{v1,v2,…,vn}表示顶点，e表示边，图的邻接矩阵a＝(aij),i＝1,…n,j＝1,…n表示顶点间的相似性；设n维向量u＝{u1,u2,…,un}t为图g中n个节点的标签，令a·u＝λ·u，得到矩阵a的特征值λ和相应的特征向量u(v)(v＝1,2,…,h)，h为受试者的数量；

23、步骤2.2：对邻接矩阵进行归一化：首先，矩阵a乘以度矩阵d-1，进一步将其拆分为两个d-1/2，得到对称且归一化的矩阵a1(v)＝d-1/2ad-1/2，其中d是对角线上节点度数的对角矩阵(d1,d2,…,dn)；在传统的谱聚类算法中每次优化一个特征向量u(v)，k为聚类簇数，将特征向量u(v)正则化为公共质心u(*)(v＝1,2,…,h)，每个正则化项由反映顶点v的重要性的加权参数λv加权；优化后的目标函数为：

24、

25、为u(v)提供一个标准的光谱聚类，求解u(*)特定的特征向量：

26、

27、运用矩阵迹的循环性质，得到：

28、

29、步骤2.3：再进行对称优化，改进后的归一化邻接矩阵为m表示相似性图的总数；通过特征值分解计算特征向量u(v)，代入式(6)得到u(*)，不断更新迭代直到收敛，产生最终估计特征向量和特征值eu；

30、

31、步骤2.4：在进行分层聚类时，根据fiedler向量的正负将网络节点分层划分为当前级别下游的两个集群；若每次聚类中，组中所有受试者的fiedler特征值的中值为正，则重复上述步骤，在当前级别获得的两个集群被递归地送到步骤2.1-步骤2.3中，否则停止执行下一步；通过计算两类数据点之间的相似性，对所有数据点中最为相似的两个数据点进行组合，并反复迭代这一过程，形成自上而下的聚类树；

32、步骤2.5：连续迭代时，通过连续迭代之间的目标差异来监控收敛性，当目标函数中的数值差低于最小阈值时聚类停止，否则聚类继续；达到聚类停止条件时得到迭代的最优解，构成图v的分割，最终得到功能子网络的聚类结果；

33、在第一次划分之后，得到第一层划分结果，分为低阶感知网络和高阶认知网络；通过第二次划分后的得到第二层划分结果，其中低阶感知网络分为视觉网络和感觉运动网络，高阶认知网络分为默认模式网络和注意网络；这些功能专门化的子网络可以进一步划分为更细的子网络，最终得到分层子网络的树结构。

34、进一步地，步骤3具体包括：

35、步骤3.1：利用图论算法计算每个功能子网络拓扑属性参数的取值，并计算各个拓扑属性参数对应的受试者工作特征曲线下面积，其中，拓扑属性参数包括度、聚类系数、特征路径长度、小世界属性以及节点中心度；

36、步骤3.2：使用统计学方法，将每个功能子网络的拓扑属性参数对应的受试者工作特征曲线下面积与认知功能状态量表评分进行相关性分析，用皮尔逊方法计算相关系数r，公式如下：

37、

38、其中，g表示拓扑属性参数对应的受试者工作特征曲线下面积，w表示认知功能状态量表评分，σg、σw分别表示g、w的标准差。

39、进一步地，步骤4具体包括：

40、使用t检验方法在分布表中查找对应的临界值，公式如下：

41、

42、其中，t值是统计量，n表示特征合集中的元素个数；

43、根据t的数值查阅分布表，得到对应的p值，根据显著性水平，判断检验结果是否显著，如果p值小于显著性水平，则拒绝原假设，认为选定的拓扑属性参数对认知功能状态量表得分具有显著影响，最后将整体相关性较高的两个拓扑属性参数筛选出来，将这两个拓扑属性参数对应的受试者工作特性曲线下面积作为功能子网络的特征合集。

44、进一步地，步骤5具体包括：

45、步骤5.1：采用高斯径向基核函数的支持向量机模型，决策函数如下：

46、

47、其中，||xi-x||2表示模型输入与核函数中心欧式距离的平方；h为受试者的数量，b是决策函数中的偏执参数；g为径向基函数中的特性参数；其中，模型输入为步骤4中得到的功能子网络的特征合集；

48、步骤5.2：寻找最优的惩罚因子c以及特性参数g：首先初始化粒子群，计算每个粒子的适应度，其次在每一次的迭代中，粒子通过跟踪两个极值pbest和gbest来更新数据，当达到全局最优位置时得到全局最优解；

49、步骤5.3：基于特征合集对支持向量回归机进行训练：在算法合理范围内设置种群的规模以及迭代次数，对种群粒子的位置和速度进行初始化操作，每个个体就是一组(c,g)向量组；以支持向量机的网络均方误差作为适应度函数，计算粒子群算法中粒子群的适应度；对每次迭代过程中的每个粒子表示的参数值，带入到模型中训练得到的支持向量机模型，然后利用测试数据来评估模型的网络均方误差，即参数的好坏。

50、进一步地，步骤6具体包括：

51、步骤6.1：设功能子网络特征合集中包含n个元素，每个特征对应的特征向量表示为xi＝{xi1,xi2,…,xij}，其中j为特征维度，xij表示第i个受试者的第j个特征的取值；对应认知功能状态量表的目标值为yi，则支持向量机模型中的训练数据集表示为：

52、t＝{(x1,y1),(x2,y2),…,(xi,yi)},xi∈rh,yi∈{1,-1}(i＝1,2,…,h) (10)

53、得到支持向量到决策面的距离最小值公式：

54、

55、其中，ω是超平面的法向量，||ω||表示ω的模长，b是超平面的截距；

56、引入松弛变量ξi，ξi≥0，则支持向量机表示为yi(ωtxi+b)≥1-ξi，它需要最小化惩罚因子c>0；c越大，样本惩罚越大；c越小，样本惩罚越小；

57、步骤6.2：引入拉格朗日乘子法，得到：

58、

59、其中，ξ是指松弛变量ξi所有取值组成的向量，ai,μi为拉格朗日乘子，并且α＝(α1,α2,…,αh)t,μ＝(μ1,μ2,…,μh)t，α和μ是指拉格朗日乘子ai,μi所有取值组成的向量；令l(ω,b,ξ,α,μ)对ω,b,ξi的偏导为零，可以得出将结果代入上式可得：

60、

61、

62、上式是用于求解拉格朗日函数，目标是要使得划分超平面具有最大间隔，将最大间隔的凸二次规划问题转变为拉格朗日函数求解问题；其中，当时，0≤αi≤c；αj、yj、xj是映射到高维空间上对偶问题中的另一个变量；

63、求解得到拉格朗日乘子的极大值然后计算得到：

64、

65、引入非线性映射，并将其替换为核函数，得到非线性支持向量机模型决策函数，其预测值可表示为：

66、

67、其中，k(xi,x)是高斯径向基函数，表示为：

68、

69、步骤6.3：利用筛选出的功能子网络的特征，导入到参数优化后的支持向量回归机中，形成映射关系来预测认知功能状态。

70、本发明具有如下有益效果：

71、本发明从大脑功能网络中分层地将全局网络聚类为子网络，提供跨多个级别的分层网络拓扑结构；在子网络中对认知功能状态进行多角度分析，更具体地描述不同认知功能状态下大脑功能网络的差异性和相似性；根据大脑功能子网络变化判断不同的认知功能状态，利用优化后的支持向量回归机进行预测，提高对认知功能状态预测的准确性。