CN114863165A - 一种基于影像组学与深度学习特征融合的椎体骨密度分类方法 - Google Patents

Abstract

本发明涉及医学图像分割与图像分类技术领域，具体涉及一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，将CT图像下椎体分为骨质疏松、低骨量与正常组，包括S1:建立基于CRF和注意力引导的椎体分割网络，获取L1、L2椎体松质骨掩膜；S2:通过GCAM‑Net对L1和L2融合后的特征图进行深度学习特征提取，并利用L1和L2的CT图像及掩膜进行影像组学特征提取；S3：运用差分进化算法在深度学习特征中提取最优特征集，并将提取后影像组学特征通过SVM‑RFE法进行特征筛选，最后将影像组学特征和深度学习特征通过最大相关性融合算法进行特征融合，并利用神经网络分类。本发明技术方案结合深度学习特征和影像组学特征，有效提高骨密度分类的准确性。

Description

一种基于影像组学与深度学习特征融合的椎体骨密度分类 方法

技术领域

本发明涉及医学图像分割与图像分类技术领域，具体涉及一种基于影像组学与深度学习特征融合的椎体骨密度分类方法。

背景技术

骨质疏松症是一种骨量减少、骨组织微结构破坏、骨骼脆性增加和易发生骨折的全身性疾病，其发病率随着年龄的增长而明显上升，骨量减少是导致骨质疏松的直接原因，对骨量的检测是预防和治疗骨质疏松的关键，腰椎因其骨质含量大且分布均匀多被选为测量部位，多个椎体的综合骨密度值被用做最终的评价标准，通常为L1、L2椎体。

目前临床上用于评估骨质疏松的常用方法包括，双能X线骨密度测量法(DXA)，定量CT测量法(QCT)以及定量核磁共振成像(QMRI)等方法。使用DXA对骨密度进行测量，速度快且X射线的辐射剂量低，但是该方法的扫描区域无法将皮质骨和松质骨区分开。QMRI的方法可以对椎体中的骨小梁结构进行研究和评估，这种方法受限于成像速度缓慢，成像价格昂贵，受众面不高。QCT应用已知密度的体模和对应的测量分析软件对骨密度进行测量，准确性高，但是其辐射剂量大，操作复杂，价格也相对昂贵。

人工智能技术与医学影像的结合是目前研究的一大热点，研究人员尝试挖掘影学图像的内在特征作为分类标准，如CT值，虽然有实验证明使用CT值作为判断的准确度可以达到82％，但是CT值受管电压影响极大使得模型稳定性差。为了在低辐射剂量下获取高精度与高时效的骨密度分类方法，在过去的十年中，研究人员寻求一种不受辐射剂量影响的骨密度分类方法，并取得一定结果。章轶立等人依据组学特征，如灰度共生矩阵、灰度大小区域矩阵、灰度运行长度矩阵、相邻灰度差矩阵、灰度依赖矩阵等，建立了基于Group Lasso的Logistic回归的骨质疏松症风险评估；SHioji等验证了多种CNN算法模型，通过评估矿物质密度，判断骨量减少，更准确的预测了骨质疏松，结果显示CNN模型优于Logistic回归分析，为骨质疏松症的早判断，早干预提供理论依据；Mehta等利用最优随机森林和支持向量机算法模型，训练模型对于骨质疏松判断的准确率、敏感度和特异度，分别达到了98.0％、77.8％和100.0％，结果显示CNN模型可用于对偶发性骨质疏松诊断的辅助工具。

目前对于骨质疏松预测模型的建立主要从两个方面出发，一是通过影像组学特征建立机器学习模型，二是通过CT图像建立深度学习模型。影像组学特征需要手动分割，耗费大量时间，深度学习模型构建缺乏标准的公共数据集，所用的数据集规模较小；放射组学特征可以很好的对骨小粱的内部微观结构进行解码，而深度学习特征可以更好的识别骨小梁的边缘特征，保证机器学习的质量和效率。此外，现有的骨密度分类方法回避了综合使用多个椎体进行评估的标准，将单张椎体的预测结果作为最终结果，缺乏合理性，所提取的全部特征用于分类导致模型难以优化。

发明内容

针对以上问题，本发明提供了一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，通过结合深度学习特征和影像组学特征，有效提高骨密度分类的准确性。

为了实现上述目的，本发明采用的技术方案如下：

一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，将CT图像下椎体分为骨质疏松、低骨量与正常组，包括以下步骤：

S1：建立基于CRF和注意力引导的椎体分割网络，获取L1椎体、L2椎体的松质骨分割掩膜；

S2：通过GCAM-Net对L1椎体和L2椎体融合后的特征图进行深度学习特征提取，并利用L1椎体和L2椎体的CT图像及掩膜进行影像组学特征提取；

S3：运用差分进化算法在深度学习特征中提取最优特征集，并将提取后影像组学特征通过SVM-RFE法进行特征筛选，最后将影像组学特征和深度学习特征进行最大相关性融合算法进行特征融合，并利用神经网络分类。

优选地，在S1中，所述CRF和注意力引导的椎体分割网络包括：

特征提取模块，将图形特征嵌入到卷积神经网络中，用于学习卷积特征和图形特征；

通道特征融合模块，用于抑制背景中噪音、伪影对分割的干扰；

特征推理模块，用于划分皮质骨与松质骨区域，细化松质骨分割边界，填充分割掩膜内部的孔洞。

优选地，所述特征提取模块利用2D残差块和图卷积进行特征提取，使用注意门控制来自不同特征提取器中的信息流。

优选地，所述特征推理模块利用条件随机场建立椎体像素邻域间的关系，将分割问题转化为能量成本最小化问题，用平均场更新来近似代替复杂的后验概率分布，利用KL散度表示使用理论分布Q拟合真实分布P时产生的信息损耗，引入最小化自由能函数来优化KL散度的求解，求解公式如下所示：

式中D_KL(Q||P)为分布P与分布Q之间的KL散度，Q(h)为使用平均场拟合分布P的后验概率分布，lnZ为常数项，E(h)为分割时的能量函数，其中

为自由能，自由能的第一项为标记每个体素的成本，第二项表示分布Q的熵，h为隐变量，利用卷积操作代替随机场中的平均场更新，定义使用图卷积网络提取的特征X^G，卷积提取特征X^R，注意门的权重矩阵A^GR，图中以及卷积中待推理的特征分别为H^G以及H^R具体的操作步骤为：

(1)计算注意门权重矩阵的更新矩阵：

为A^GR的更新矩阵，r^GR为一个3×3的卷积核，X^G为使用图卷积网络提取的特征，X^R为卷积提取特征；

(2)使用sigmoid函数σ对更新矩阵做归一化处理：

式中的

为归一化后注意门权重矩阵的更新，σ为sigmoid函数；

(3)H^G的平均场更新：H^G＝r^GR*H^G，*表示卷积操作，X^G为图卷积网络提取的特征，等式左边的H^G为更新后的H^G，r^GR为一个3×3的卷积核；

(4)

⊙表示元素相乘，

为卷积的隐藏特征H^G的更新，

为归一化后注意门权重矩阵的更新，X^G为使用图卷积网络提取的特征；

(5)通过使用剩余连接添加一元势来更新H^R:

表示元素相加，

为卷积的隐藏特征H^G的更新，为了便于迭代更新，H^G的贡献归算在H^R中。

优选地，所述通道特征融合模块，聚合来自编码层的低级信息以及解码层中的高级特征，通过挤压和激发操作校准通道特征响应，得到的权重向量与低层特征相乘，并将重新加权的低层特征与高层特征相加得到最终结果，抑制背景中噪音、以及相似组织对分割的干扰。

优选地，在S2中，所述GCAM用于椎体特征提取卷积神经网络包括：

特征融合模块，通过对L1椎体和L2椎体的图像在进行卷积神经网络特征提取前，先一步进行特征融合，特征融合最后采用特征的通道数相融合的方式将L1椎体图像和L2椎体图像进行特征融合；

门控通道注意力模块，用于建模通道之间的关系，在归一化之前嵌入全局上下文并控制每个通道的权重，以及一个门控适应算子，根据归一化的输出逐通道地调整输入特征。

其中，GCAM采用一种标准化方法来创建渠道之间的竞争或合作关系。值得注意的是，规范化操作是无参数的。为了使GCAM可学习，设计了一个全局上下文嵌入算子，它在归一化之前嵌入全局上下文并控制每个通道的权重，以及一个门控适应算子，它根据归一化的输出逐通道地调整输入特征。按信道可训练的参数是轻量且有效的，并且使得GCAM便于广泛部署，同时占用少量的参数。此外，门控适应算子的参数对于解释GCAM的行为是容易和直观的。总之，精心设计了高度轻量级的，可解释的，但有效的GCAM架构，该架构基于建模渠道关系的规范化操作。假设

作为卷积网络得激活特征，其中得H和W作为图像得高度和宽度，C为图像的通道数，一般来说，GCAM执行以下变换：

其中的α，β和γ都是可训练的参数，嵌入权重α主要负责调整嵌入输出、门控权重γ和偏置权重β负责调整门控的激活，

表示第c个通道的特征图，

表示变换后的特征，F表示函数变换，x表示原特征。它们决定了GCAM在每个通道中的行为。相对于通道注意力模块SE的时间复杂度O(C²),GCAM模块的时间复杂度只有O(C)。

其中，大感受野信息可以避免小感受野信息(例如,卷积层)引起的局部模糊，因此，首先设计了一个全局上下文嵌入模块来聚合每个通道中的全局上下文信息、该模块可以利用卷积层的小感受域之外的全局上下文信息。给定嵌入权重α＝[α₁,α₂,...α_c],模块S_c的定义为：

其中的∈是一个小常数，主要用来避免零点的求导问题，S_c为小感受域模块，α表示嵌入权重集，c表示通道数，|| ||₂表示2范数，H表示图像高度，W表示图像宽度，∑表示求和，x_c表示通道特征。相对于SE来说，GCAM不使用全局平均池化层来汇总通道的环境。此外，GCAM使用可训练的参数α_c来控制每个通道的权重，因为不同的通道应该具有不同的重要性，特别是α_c接近0时，通道C将不参与通道归一化，总之，门控权重α使得GCAM能够学习一个通道对于其他通道是独立的情况。

其中，归一化方法可以在神经元(或通道)之间建立竞争关系，具有轻量级计算资源和稳定的训练性能。类似于LRN，使用l2归一化来跨信道操作，即信道归一化。S＝[S₁,S₂,...S_c]通道归一化的公式为：

其中c是一个小常数，标量

用于归一化S_C的比例，为了避免c太大的时候S_C的规模太小，∑表示求和，|| ||₂表示2范数，S_c为小感受域模块，

为通道归一化。

其中，采用门控机制来适应原始特征，通过引入门控机制，GCAM可以在训练过程中促进竞争和合作。设选通权重γ＝[γ₁,γ₂,...γ_c]和门控偏置β，设计以下门控函数：

其中，

为门控函数，γ_c为通道偏差，每个原始通道x_c的标度将由其对应的门来适配，即

由于通道归一化是无参数的，设计了可训练的权重和偏差γ和β，用于学习以通道方式控制门的激活。LRN只受益于神经元之间的竞争。然而，通过结合标准化方法和门控机制，GCAM能够对不同渠道之间的更多类型的关系(即竞争和合作)建模。当一个通道的选通权重(γ_c)被积极激活时，GCAM促进该通道与其他通道竞争，如在LRN。当门控权重被负面激活时，GCAM鼓励通道与其他通道合作。

此外，当选通权重和偏差为零时，该选通函数允许原始特征传递到下一层，即

其中，

表示传递后的特征，F表示函数变换，x表示原始特征。对身份映射进行建模的能力可以有效地提高深度网络中退化问题的鲁棒性。ResNet也从这个想法中受益。因此，建议在初始化GCAM层时将γ和β初始化为0。这样做，训练过程的初始步骤会更稳定，GCAM的最终表现会更好。

优选地，在S3中，利用差分进化算法、SVM-RFE算法和最大相关性计算对特征进行优化；

差分进化算法：将差分进化算法模块嵌入到神经网络提取特征的池化层后，全连接层前，对特征进行优化；

SVM-RFE算法：将SVM-RFE算法嵌入影像组学特征提取后，对批提取特征进行递归特征筛选，去除冗余信息；

最大相关性计算：将差分进化算法后的特征和SVM-RFE后的特征进行最大相关度计算，进一步对特征进行优化。

其中，为了对卷积神经网络提取到的特征进行特征优化，并且消除后续模型训练中过拟合的问题，使用差分进化算法进行特征选择。

其中，所述差分进化算法从某一随机产生的初始群体开始，利用从种群中随机选取的两个个体的差向量作为第三个个体的随机变化源，将差向量加权后按照一定的规则与第三个个体求和而产生变异个体，该操作称为变异；然后，变异个体与某个预先决定的目标个体进行参数混合，生成试验个体，这一过程称之为交叉；如果试验个体的适应度值优于目标个体的适应度值，则在下一代中试验个体取代目标个体，否则目标个体仍保存下来，该操作称为选择；在每一代的进化过程中，每一个体矢量作为目标个体一次，算法通过不断地迭代计算，保留优良个体，淘汰劣质个体，引导搜索过程向全局最优解逼近。

具体的流程如下：输入：分割后的特征向量，输出：优化后的特征向量，对种群进行初始化

其中

表示第j维的下界，

表示第j维的上界，NP表示向量数，D表示向量维度。然后通过以下表达式进行计算边界：

其中rand(0,1)表示区间[0,1]上的随机数，

表示第j维的下界，

表示第j维的上界，x_i,j(0)表示边界计算；使用精细的KNN算法计算适应度函数；执行变异，变异的表达式如下：v_i(n+1)＝λ_r1(n)+F(λ_r2(n)-λ_r3(n))，其中F代表范围为(0，2)的缩放变异因子，n表示第n代，v_i(n+1)表示经过变异得种群，λ_r1(n)表示未经过变异得种群，源向量是从总体中随机选择的λ_r1、λ_r2、λ_r3，并且λ_r1、λ_r2、λ_r3要互不相同；执行交叉，设置交叉比例，如果交叉比例大于0到1之间的随机数，那么C_i,j(n+1)＝v_i,j(n+1)，其他情况下，C_i,j(n+1)＝λ_i,j(n+1)，其中C_i,j(n+1)表示交叉后得种群，v_i,j(n+1)表示经过变异得种群，λ_i,j(n+1)表示未经过变异得种群。最后进行选择，在差分进化算法中使用贪婪选择的策略，就是选择较优的个体作为新的个体。对于每一个个体来说，得到地解要好于或者持平于个体通过变异、交叉、选择达到全部最优。

选择的公式如下：

其中，X_i(n+1)为选择后的个体，f表示当前种群适应度平均值，X_i(n)表示原始种群，C_i(n+1)表示交叉变异后的种群，n表示第n代。

优选地，为了去除影像组学提取地特征中的冗余特征，提取到关键信息，使用SVM-RFE方法对影像组学进行特征选择，筛选步骤如下：

步骤1：输入需要的数据样本集

v_i∈{1,2,...,l}，其中，l为类别数，v_i为类别数集合，N为总样本数，x_i为样本。并对数据进行标准化处理

μ为均值，σ为标准差，x_i为样本；

步骤2：输入特征集合S＝{1,2,3...D}并初始化，设置特征排序集R＝[ ]，其中D为特征集合子集合；

步骤3：生成

个样本，并且在训练的样本

中找出不同类别的两两组合得到最后的训练样本；

v_i＝1时，y_i＝1，v_i＝j+1时，y_i＝-1；…………

v_i＝l-1时，y_i＝1，v_i＝l时，y_i＝-1；其中v_i为类别数集合，x_i为样本，N为样本数，x_j为组合后训练样本。

步骤4：循环以下过程直至S＝[]，其中S为特征集合：

获得l个用于训练的子样本x_j(j＝1,2,...l(l-1)/2)，其中l为个数，x_j为样本；

并用此样本x_j训练SVM，分别得到w_j(j＝1,2...l)，其中w_j为特征权重，其中l为个数；

计算排序准则总分c_k＝∑_jw_jk ²(k＝1,2...,|S|)，其中w_jk为第k个子集第j个特征的特征权重，∑表示求和，c_k表示第k个子集所有特征的总分，|S|表示特征集合中子集合的个数；

找出排序准则分数最小的特征

其中p表示排序准则分数最小特征，arg min为最小值时的变量的取值，c_k表示第k个子集所有特征的总分；

更新特征集R＝{p}∪R，其中R为特征排序集，p表示排序准则分数最小特征，∪表示两个集合的并集；

在S中去除此特征S＝S/p，其中S为特征集合，p表示排序准则分数最小特征。

优选地，为了获得最优的融合特征向量，对经过优化的深度学习特征和影像组学特征进行最大相关度计算，最大相关性计算为，通过最大相关系数对R₁和R₂两个特征向量进行计算，计算地表达式如下：

其中的ρ(R₁,R₂)的值始终介于(-1，1)之间，表示相关性，1表示两个特征之间的强相关性，-1表示两个特征间的弱相关性，COV表示协方差，Var表示变量，R₁和R₂是两个特征向量。由于感兴趣的是最大相关性，因此通过下列表达式实现两个向量的最大相关：

CC(R₁,R₂)＝T(ρ(g₁(R₁),g₂(R₂)))

其中R₁和R₂是两个特征向量，g₁和g₂是Borel函数，T表示取Borel函数g₁和g₂的最小上界，位于(0，1)，CC(R₁,R₂)表示两个向量的最大相关性，ρ(g₁(R₁),g₂(R₂))的值始终介于(-1，1)之间，表示相关性；融合的步骤为：输入两个最佳特征向量R₁和R₂，执行均值填充使得两个向量维度之间相同，计算两个向量的最大相关性，如果两个向量的最大相关性接近1，则在融合向量中加入这两个特征，如果相关性接近0，那么丢弃这两个特征，最后输出一个融合的特征向量。

本发明有益效果：

本发明以CT图像为基础，引入条件随机场特征推理模块，建立椎体与背景像素邻域关系模型，优化对皮质骨边界的区分能力。利用差分进化对卷积神经网络提取的深度学习特征进行特征优化，并利用支持向量机对影像组学特征进行递归筛选，对深度学习筛选特征和影像组学筛选特征进行最大相关度计算，建立骨质疏松分类模型，有效提高骨密度分类的准确性。

附图说明

图1为本发明的整体流程图；

图2为本发明中分割网络的示意图；

图3为本发明中椎体L1和椎体L2的特征融合的示意图；

图4为本发明中GCAM-Net分类网络的示意图；

图5为本发明中GCAM模块示意图；

图6为本发明中分割网络的可视化效果图；

图7为本发明中分类效果ROC曲线图。

具体实施方式

下面结合附图将对本发明实施例中的技术方案进行清楚、完整地描述，以使本领域的技术人员能够更好的理解本发明的优点和特征，从而对本发明的保护范围做出更为清楚的界定。本发明所描述的实施例仅是本发明一部分实施例，而不是全部的实施例，基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。

参照图1-7，一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，将CT图像下椎体分为骨质疏松、低骨量与正常组，包括以下步骤：

1、特征提取模块

在使用图卷积神经网络提取特征的过程中，通过G＝{V,G,A}构建特征交互空间，其中V是交互图中的节点集合，G是节点之间边的集合，A是邻接矩阵，将原始输入的特征图投影到特征交互空间得到新的投影特征

其中新的交互节点数为K，新的边集为E，C'为交互空间的维度，通常新特征中节点之间通过线性组合进行连接，不同切片之间的血管形态差异较大，缺乏自适应的采样能力，模型难以优化，采用自适应采样策略，运用三线性插值采样器，使用随机梯度下降的方法为每个原始特征学习移动距离，从而对给定的V个节点的图中，对特征节点X_m及其周围的位移特征节点进行采样。对投影后得到的新特征，遵循图卷积神经网络的卷积操作，定义新的邻接矩阵为A^G,权重矩阵为W^G,则图卷积网络的运算公式为X^G＝σ(A^GX^PROW^G)，在网络模型的实际操作中，首先应用拉普拉斯平滑和更新邻接矩阵

从而将节点特征传播到整个图中，I是输入的一组SMA序列。在实际操作中采用了梯度流最大的剩余连接代替,1×1卷积实现

和W^G，推理完成后再将特征从交互空间投影回坐标空间。

2、通道特征融合模块

为了抑制背景噪音干扰，没有直接将拼接后的特征直接传送到上采样阶段，而是将来自跳跃连接层中的低层特征图与解码层中上采样的高级特征在通道维度进行拼接，基于挤压和激发(SE)操作来指导特征层之间的融合。将低阶特征

与高阶特征

在通道维度上进行拼接得到

C,H,W分别代表通道式，图像的高度和宽度。接着通过挤压和激励(SE)操作建模通道之间的相关性，通过计算一个权重向量对低层特征进行重新加权并抑制无关背景噪声的干扰。拼接后的特征F_concat送入一个1×1的卷积层中实现降维，紧着利用全局平均池化(GAP)得到加权向量

将Sigmoid函数生成的权重向量与低层特征相乘，然后将重新加权的低层特征与高层特征相加得到最终结果：

式中F^(c)为当前层的融合特征，

为该层上层的融合特征，σ表示sigmoid函数，C_1×1代表1×1的卷积操作，c代表模型的层数，F_L为低阶特征，F_H为高阶特征，GAP为全局平均池化操作，

和

分别表示元素相加与元素相乘，采用无参数双线性上采样策略，在保持分割性能的同时减少了参数量。

3、特征推导模块

为了避免以最大化概率来分配分割标签时，由于相邻体素共享特征表达相似而导致不正确的边界分割，在预重建网络中引入条件随机场特征推理模块，使用不同特征提取器——图卷积神经网络(图卷积神经网络)与卷积神经网络(CNN)提取椎体特征，按照所提特征对最终结果的贡献程度，通过注意门调节来自双方的信息流并将该特征输入条件随机场进行推理，对每个像素点进行精准分类，具体位置放置在网络下采样过程中的第4层以及第5层中。

卷积神经网络下采样提取的卷积特征X^R安置在每层2D残差块后，注意门A^GR控制两种潜在特征的信息流，接着使用条件随机场完成推理过程，全连接条件随机场是一种概率图模型，将图像中每个像素点都看作一个节点，像素与像素之间的关系看作连接节点的边，由此构成条件随机场。使用条件随机场的最终目标是实现能量函数的最小化，在本发明中用图卷积提取的特征和使用条件随机场推理的图卷积中的潜在特征以及卷积提取特征和使用条件随机场推理的卷积中潜在特征之间的关系看作一元势(两个一元势)，将使用条件随机场融合的最终特征与注意门之间关系看作二元势，使用高斯函数将推测所得的潜在特征H表征到观测特征X上，本发明中的可观测特为最终图卷积和卷积的融合特征X^F，一元势，二元势以及和能量函数的公式如下所示：

式中φ(H,X)指代一元势，h_n为隐藏特征，x_n为可观测特征；

式中

指代二元势，h_m ^g为参与交互的隐藏图卷积特征，h_n ^r为隐藏卷积特征，

为注意门权重矩阵，

是与隐藏特征映射相关联的内核势；

式中E(H^R,H^G,X^F|I)为使用概率图模型分割时的能量函数，I是输入的一组SMA序列，φ^G(H^G,X^G)为隐藏图卷积H^G与可观测性图卷积特征X^G之间的一元势，φ^R(H^R,X^R)为为隐藏卷积H^R与可观测性卷积特征X^R之间的一元势,ψ^GR(H^G,H^R,A^GR)为隐藏卷积特征，隐藏图卷积特征与注意门权重矩阵之间的二元势，本实施例利用高斯函数将预测的隐藏特征H驱动到可观测特征X，右上角的角标表示特征来源，表示成对势。

条件随机场中，图像的分割问题被转化为能量函数最小化问题，采用平均场更新来近似代替复杂的后验概率分布，使用具有拟合参数的Q作为后验代理，形成隐藏变量的后验分布：

式中P(H^R,H^G,X^F,|I)用于描述可观测的图像分布，Q(H^R,H^G,X^F)为拟合P的后验概率分布，隐藏变量

以及注意门权重矩阵

间是相互独立的且可分解的，q_n为可因数分解表达式，H为潜在特征，右上角的角标表示特征来源，H^R,H^G分别表示隐藏卷积特征和图特征，X为可观测特征，X^F为融合隐藏卷积特征和图形特征的最终可观测特征；

引入KL散度表示使用理论分布Q拟合真实分布P时产生的信息损耗，

由于包含式P(h|x)，KL散度依旧难以求解，参考马尔科夫随机场，引入最小化自由能函数来优化KL散度的求解，求解公式如下所示：

式中D_KL(Q||P)为分布P与分布Q之间的KL散度，Q(h)为使用平均场拟合分布P的后验概率分布，lnZ为常数项，E(h)为分割时的能量函数，其中

为自由能，自由能的第一项为标记每个体素的成本，第二项表示分布Q的熵,h为隐变量。将能量函数表达式代入KL散度的求解方程中，对隐变量求导并令导数为零有：

与

同理，式中

表示注意门权重矩阵

的更新，

表示分布Q对各隐变量的期望，

为注意门权重矩阵，

为卷积的潜在特征表达，

为图卷积的隐藏特征表达，

为对

求导后的二元势表达式，

为对

求导后的可观测卷积特征与隐藏特征之间的一元势，

描述可观测卷积特征；

利用相邻体素之间潜在特征表达的重新加权学习体素之间的共生关系，体素n的潜在卷积和图形特征之间的注意力平衡允许重新加权来自体素n的领域间的成对潜在信息。定义

代入上述公式中有特征更新如下：

式中

为隐藏图卷积特征的更新，

描述可观测图卷积特征，

描述可观测卷积特征，

是与隐藏特征映射相关联的内核势，

为参与交互的隐藏卷积特征，

是与隐藏特征映射相关联的内核势，

表示按照不同信息流的权重矩阵，将卷积特征归算到图卷积特征中；

式中

为隐藏卷积特征的更新，

描述可观测卷积特征，

是与隐藏特征映射相关联的内核势，

为参与特征交互的图卷积特征，

将参与交互的图卷积特征按照注意门权重归算到卷积特中；

权重矩阵

可以使用sigmoid函数对其进行更新：

式中，

为更新后的

为参与交互的注意门权重矩阵，

是与隐藏特征映射相关联的内核势，

为参与交互的隐藏图卷积特征；

在实际的操作过程中，将条件随机场的平均场更新替换为卷积操作，以适应任何卷积网络的端对端的训练。定义使用图卷积网络提取的特征X^G，卷积提取特征X^R，注意门的权重矩阵A^GR，图中以及卷积中待推理的特征分别为H^G以及H^R具体的操作步骤为：

(1)计算注意门权重矩阵的更新矩阵：

为A^GR的更新矩阵，r^GR为一个3×3的卷积核，

为使用图卷积网络提取的特征，X^R为卷积提取特征；

(2)使用sigmoid函数σ对更新矩阵做归一化处理：

式中的

为归一化后注意门权重矩阵的更新；

(3)H^G的平均场更新：H^G＝r^GR*H^G，*表示卷积操作，X^G为图卷积网络提取的特征，等式左边的H^G为更新后的H^G，r^GR为一个3×3的卷积核；

(4)

⊙表示元素相乘，

为卷积的隐藏特征H^G的更新，

为归一化后注意门权重矩阵的更新，X^G为使用图卷积网络提取的特征；

(5)通过使用剩余连接添加一元势来更新H^R:

表示元素相加，

为卷积的隐藏特征H^G的更新，为了便于迭代更新，H^G的贡献归算在H^R中。

4、门控通道注意模块

门控通道注意力模块，用于建模通道之间的关系，在归一化之前嵌入全局上下文并控制每个通道的权重，以及一个门控适应算子，根据归一化的输出逐通道地调整输入特征。

其中，GCAM采用一种标准化方法来创建渠道之间的竞争或合作关系。值得注意的是，规范化操作是无参数的。为了使GCAM可学习，设计了一个全局上下文嵌入算子，它在归一化之前嵌入全局上下文并控制每个通道的权重，以及一个门控适应算子，它根据归一化的输出逐通道地调整输入特征。按信道可训练的参数是轻量且有效的，并且使得GCAM便于广泛部署，同时占用少量的参数。此外，门控适应算子的参数对于解释GCAM的行为是容易和直观的。总之，精心设计了高度轻量级的，可解释的，但有效的GCAM架构，该架构基于建模渠道关系的规范化操作。假设

作为卷积网络得激活特征，其中得H和W作为图像得高度和宽度，C为图像的通道数，一般来说，GCAM执行以下变换：

其中的α，β和γ都是可训练的参数，嵌入权重α主要负责调整嵌入输出、门控权重γ和偏置权重β负责调整门控的激活，

表示第c个通道的特征图，

表示变换后的特征，F表示函数变换，x表示原特征。它们决定了GCAM在每个通道中的行为。相对于通道注意力模块SE的时间复杂度O(C²),GCAM模块的时间复杂度只有O(C)。

其中，大感受野信息可以避免小感受野信息(例如,卷积层)引起的局部模糊，因此，首先设计了一个全局上下文嵌入模块来聚合每个通道中的全局上下文信息、该模块可以利用卷积层的小感受域之外的全局上下文信息。给定嵌入权重α＝[α₁,α₂,…α_c],模块S_c的定义为：

其中的∈是一个小常数，主要用来避免零点的求导问题，S_c为小感受域模块，α表示嵌入权重集，c表示通道数，|| ||₂表示2范数，H表示图像高度，W表示图像宽度，∑表示求和，x_c表示通道特征。相对于SE来说，GCAM不使用全局平均池化层来汇总通道的环境。此外，GCAM使用可训练的参数α_c来控制每个通道的权重，因为不同的通道应该具有不同的重要性，特别是α_c接近0时，通道C将不参与通道归一化，总之，门控权重α使得GCAM能够学习一个通道对于其他通道是独立的情况。

其中，归一化方法可以在神经元(或通道)之间建立竞争关系，具有轻量级计算资源和稳定的训练性能。类似于LRN，使用l2归一化来跨信道操作，即信道归一化。S＝[S₁,S₂,...S_c]通道归一化的公式为：

其中c是一个小常数，标量

用于归一化S_C的比例，为了避免c太大的时候S_C的规模太小，∑表示求和，|| ||₂表示2范数，S_c为小感受域模块，

为通道归一化。

其中，采用门控机制来适应原始特征，通过引入门控机制，GCAM可以在训练过程中促进竞争和合作。设选通权重γ＝[γ₁,γ₂,...γ_c]和门控偏置β，设计以下门控函数：

其中，

为门控函数，每个原始通道x_c的标度将由其对应的门来适配，即

由于通道归一化是无参数的，设计了可训练的权重和偏差γ和β，用于学习以通道方式控制门的激活。LRN只受益于神经元之间的竞争。然而，通过结合标准化方法和门控机制，GCAM能够对不同渠道之间的更多类型的关系(即竞争和合作)建模。当一个通道的选通权重(γ_c)被积极激活时，GCAM促进该通道与其他通道竞争，如在LRN，当门控权重被负面激活时，GCAM鼓励通道与其他通道合作。

此外，当选通权重和偏差为零时，该选通函数允许原始特征传递到下一层，即

其中，

表示传递后的特征，F表示函数变换，x表示原特征。对身份映射进行建模的能力可以有效地提高深度网络中退化问题的鲁棒性。ResNet也从这个想法中受益。因此，建议在初始化GCAM层时将γ和β初始化为0。这样做，训练过程的初始步骤会更稳定，GCAM的最终表现会更好。

5、差分进化特征优化模块

为了对卷积神经网络提取到的特征进行特征优化，并且消除后续模型训练中过拟合的问题。使用差分进化算法进行特征选择。差分进化算法作为一种全局搜索优化问题的进化算法，首先它在搜索空间中产生一些初始值，然后对输入数据进行变异和交叉，接着选择过程成中产生的新的种群。具体的流程如下：输入：分割后的特征向量，输出：优化后的特征向量，对种群进行初始化

其中

表示第j维的下界，

表示第j维的上界，NP表示向量数，D表示向量维度。然后通过以下表达式进行计算边界：

其中rand(0,1)表示区间[0,1]上的随机数，

表示第j维的下界，

表示第j维的上界，x_i,j(0)表示边界计算；使用精细的KNN算法计算适应度函数；执行变异，变异的表达式如下：v_i(n+1)＝λ_r1(n)+F(λ_r2(n)-λ_r3(n))，其中F代表范围为(0，2)的缩放变异因子，n表示第n代，v_i(n+1)表示经过变异得种群，λ_r1(n)表示未经过变异得种群，源向量是从总体中随机选择的λ_r1、λ_r2、λ_r3，并且λ_r1、λ_r2、λ_r3要互不相同；执行交叉，设置交叉比例，如果交叉比例大于0到1之间的随机数，那么C_i,j(n+1)＝v_i,j(n+1)，其他情况下，C_i,j(n+1)＝λ_i,j(n+1)，其中C_i,j(n+1)表示交叉后得种群，v_i,j(n+1)表示经过变异得种群，λ_i,j(n+1)表示未经过变异得种群。最后进行选择，在差分进化算法中使用贪婪选择的策略，就是选择较优的个体作为新的个体。对于每一个个体来说，得到地解要好于或者持平于个体通过变异、交叉、选择达到全部最优。

选择的公式如下：

其中，X_i(n+1)为选择后的个体，f表示当前种群适应度平均值，X_i(n)表示原始种群，C_i(n+1)表示交叉变异后的种群，n表示第n代。

6、SVM-RFE特征筛选模块

为了去除影像组学提取地特征中的冗余特征，提取到关键信息，使用SVM-RFE方法对影像组学进行特征选择，筛选步骤如下：

步骤1：输入需要的数据样本集

v_i∈{1,2,...,l}，其中，l为类别数，v_i为类别数集合，N为总样本数，x_i为样本。并对数据进行标准化处理

μ为均值，σ为标准差，x_i为样本；

步骤2：输入特征集合S＝{1,2,3...D}并初始化，设置特征排序集R＝[]，其中D为特征集合子集合；

步骤3：生成

个样本，并且在训练的样本

中找出不同类别的两两组合得到最后的训练样本；

v_i＝1时，y_i＝1，v_i＝j+1时，y_i＝-1；…………

v_i＝l-1时，y_i＝1，v_i＝l时，y_i＝-1；其中v_i为类别数集合，x_i为样本，N为样本数，x_j为组合后训练样本。

步骤4：循环以下过程直至S＝[ ]，其中S为特征集合：

获得l个用于训练的子样本x_j(j＝1,2,...l(l-1)/2)；

并用此样本x_j训练SVM，分别得到w_j(j＝1,2...l)，其中w_j为特征权重；

计算排序准则总分c_k＝∑_jw_jk ²(k＝1,2...,|S|)，其中w_jk为第k个子集第j个特征的特征权重，∑表示求和，c_k表示第k个子集所有特征的总分，|S|表示特征集合中子集合的个数；

找出排序准则分数最小的特征

其中p表示排序准则分数最小特征，arg min为最小值时的变量的取值，c_k表示第k个子集所有特征的总分；

更新特征集R＝{p}∪R，其中R为特征排序集，p表示排序准则分数最小特征，∪表示两个集合的并集；

在S中去除此特征S＝S/p，其中S为特征集合，p表示排序准则分数最小特征。

6、最大相关度计算模块

最大相关性计算为，通过最大相关系数对R₁和R₂两个特征向量进行计算，计算地表达式如下：

其中的ρ(R₁,R₂)的值始终介于(-1，1)之间，表示相关性，1表示两个特征之间的强相关性，-1表示两个特征间的弱相关性，COV表示协方差，Var表示变量，R₁和R₂是两个特征向量。由于感兴趣的是最大相关性，因此通过下列表达式实现两个向量的最大相关：

CC(R₁,R₂)＝T(ρ(g₁(R₁),g₂(R₂)))

其中R₁和R₂是两个特征向量，g₁和g₂是Borel函数，T表示取Borel函数g₁和g₂的最小上界，位于(0，1)，CC(R₁,R₂)表示两个向量的最大相关性，ρ(g₁(R₁),g₂(R₂))的值始终介于(-1，1)之间，表示相关性；融合的步骤为：输入两个最佳特征向量R₁和R₂，执行均值填充使得两个向量维度之间相同，计算两个向量的最大相关性，如果两个向量的最大相关性接近1，则在融合向量中加入这两个特征，如果相关性接近0，那么丢弃这两个特征，最后输出一个融合的特征向量。

以上对本发明的实施例进行了详细说明，但所述内容仅为本发明的较佳实施例，不能被认为用于限定本发明的实施范围。凡依本发明申请范围所作的均等变化与改进等，均应仍归属于本发明的专利涵盖范围之内。

Claims (10)

1.一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，将CT图像下椎体分为骨质疏松、低骨量与正常组，其特征在于，包括以下步骤：

S1：建立基于CRF和注意力引导的椎体分割网络，获取L1椎体、L2椎体的松质骨分割掩膜；

S2：通过GCAM-Net对L1椎体和L2椎体融合后的特征图进行深度学习特征提取，并利用L1椎体和L2椎体的CT图像及掩膜进行影像组学特征提取；

S3：运用差分进化算法在深度学习特征中提取最优特征集，并将提取后影像组学特征通过SVM-RFE法进行特征筛选，最后将影像组学特征和深度学习特征进行最大相关性融合算法进行特征融合，并利用神经网络分类。

2.根据权利要求1所述的一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，其特征在于，在S1中，所述CRF和注意力引导的椎体分割网络包括：

特征提取模块，将图形特征嵌入到卷积神经网络中，用于学习卷积特征和图形特征；

通道特征融合模块，用于抑制背景中噪音、伪影对分割的干扰；

特征推理模块，用于划分皮质骨与松质骨区域，细化松质骨分割边界，填充分割掩膜内部的孔洞。

3.根据权利要求2所述的一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，其特征在于，所述特征提取模块利用2D残差块和图卷积进行特征提取，使用注意门控制来自不同特征提取器中的信息流。

4.根据权利要求2所述的一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，其特征在于，所述特征推理模块利用条件随机场建立椎体像素邻域间的关系，将分割问题转化为能量成本最小化问题，用平均场更新来近似代替复杂的后验概率分布，利用KL散度表示使用理论分布Q拟合真实分布P时产生的信息损耗，引入最小化自由能函数来优化KL散度的求解，求解公式如下所示：

式中，D_KL(Q||P)为分布P与分布Q之间的KL散度，Q(h)为使用平均场拟合分布P的后验概率分布，lnZ为常数项，E(h)为分割时的能量函数，h为隐变量，其中

为自由能，自由能的第一项为标记每个体素的成本，第二项表示分布Q的熵，利用卷积操作代替随机场中的平均场更新，定义使用图卷积网络提取的特征X^G，卷积提取特征X^C，注意门的权重矩阵A^GR，卷积中待推理的特征分别为H^G以及H^R具体的操作步骤为：

(1)计算注意门权重矩阵的更新矩阵：

为A^GR的更新矩阵，r^GR为一个3×3的卷积核，X^G为使用图卷积网络提取的特征，X^R为卷积提取特征；

(2)使用sigmoid函数σ对更新矩阵做归一化处理：

式中的

为归一化后注意门权重矩阵的更新；

(3)H^G的平均场更新：H^G＝r^GR*H^G，*表示卷积操作，X^G为图卷积网络提取的特征，等式左边的H^G为更新后的H^G，r^GR为一个3×3的卷积核；

(4)

⊙表示元素相乘，

为卷积的隐藏特征H^G的更新，

为归一化后注意门权重矩阵的更新，X^G为使用图卷积网络提取的特征；

(5)通过使用剩余连接添加一元势来更新

表示元素相加，

为卷积的隐藏特征H^G的更新，为了便于迭代更新，H^G的贡献归算在H^R中。

5.根据权利要求2所述的一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，其特征在于，所述通道特征融合模块，聚合来自编码层的低级信息以及解码层中的高级特征，通过挤压和激发操作校准通道特征响应，得到的权重向量与低层特征相乘，并将重新加权的低层特征与高层特征相加得到最终结果，抑制背景中噪音、以及相似组织对分割的干扰。

6.根据权利要求1所述的一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，其特征在于，在S2中，所述GCAM用于椎体特征提取卷积神经网络包括：

特征融合模块，通过对L1椎体和L2椎体的图像在进行卷积神经网络特征提取前，先一步进行特征融合，特征融合最后采用特征的通道数相融合的方式将L1椎体图像和L2椎体图像进行特征融合；

门控通道注意力模块，用于建模通道之间的关系，在归一化之前嵌入全局上下文并控制每个通道的权重，以及一个门控适应算子，根据归一化的输出逐通道地调整输入特征。

7.根据权利要求1所述的一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，其特征在于，在S3中，利用差分进化算法、SVM-RFE算法和最大相关性计算对特征进行优化；

差分进化算法：将差分进化算法模块嵌入到神经网络提取特征的池化层后，全连接层前，对特征进行优化；

SVM-RFE算法：将SVM-RFE算法嵌入影像组学特征提取后，对批提取特征进行递归特征筛选，去除冗余信息；

最大相关性计算：将差分进化算法后的特征和SVM-RFE后的特征进行最大相关度计算，进一步对特征进行优化。

8.根据权利要求7所述的一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，其特征在于，所述差分进化算法从某一随机产生的初始群体开始，利用从种群中随机选取的两个个体的差向量作为第三个个体的随机变化源，将差向量加权后按照一定的规则与第三个个体求和而产生变异个体，该操作称为变异；然后，变异个体与某个预先决定的目标个体进行参数混合，生成试验个体，这一过程称之为交叉；如果试验个体的适应度值优于目标个体的适应度值，则在下一代中试验个体取代目标个体，否则目标个体仍保存下来，该操作称为选择；在每一代的进化过程中，每一个体矢量作为目标个体一次，算法通过不断地迭代计算，保留优良个体，淘汰劣质个体，引导搜索过程向全局最优解逼近。

9.根据权利要求7所述的一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，其特征在于，所述SVM-RFE算法对影像组学特征进行筛选的步骤为：

步骤1：输入需要的数据样本集

其中，l为类别数，v_i为类别数集合，N为总样本数，x_i为样本；并对数据进行标准化处理

μ为均值，σ为标准差,x_i为样本；

步骤2：输入特征集合S＝{1,2,3...D}并初始化，设置特征排序集R＝[]，其中D为特征集合子集合；

步骤3：生成

个样本，并且在训练的样本

中找出不同类别的两两组合得到最后的训练样本；

v_i＝1时，y_i＝1，v_i＝j+1时，y_i＝-1；…………

v_i＝l-1时，y_i＝1，v_i＝l时，y_i＝-1；其中v_i为类别数集合，x_i为样本，N为样本数，x_j为组合后训练样本；

步骤4：循环以下过程直至S＝[]，其中S为特征集合：

获得l个用于训练的子样本x_j(j＝1,2,...l(l-1)/2)，其中l为个数，x_j为样本；

并用此样本x_j训练SVM，分别得到w_j(j＝1,2...l)，其中w_j为特征权重，l为个数；

计算排序准则总分c_k＝∑_jw_jk ²(k＝1,2...,|S|)，其中w_jk为第k个子集第j个特征的特征权重，∑表示求和，c_k表示第k个子集所有特征的总分，|S|表示特征集合中子集合的个数；

找出排序准则分数最小的特征

其中p表示排序准则分数最小特征，argmin为最小值时的变量的取值，c_k表示第k个子集所有特征的总分；

更新特征集R＝{p}∪R，其中R为特征排序集，p表示排序准则分数最小特征，∪表示两个集合的并集；

在S中去除此特征S＝S/p，其中S为特征集合，p表示排序准则分数最小特征。

10.根据权利要求7所述的一种基于影像组学与深度学习特征融合的椎体骨密度分类方法，其特征在于，所述最大相关性计算为，通过最大相关系数对R₁和R₂两个特征向量进行计算，计算地表达式如下：

其中的ρ(R₁,R₂)的值始终介于(-1，1)之间，表示相关性，1表示两个特征之间的强相关性，-1表示两个特征间的弱相关性，COV表示协方差，Var表示变量，R₁和R₂是两个特征向量。由于感兴趣的是最大相关性，因此通过下列表达式实现两个向量的最大相关：

CC(R₁,R₂)＝T(ρ(g₁(R₁),g₂(R₂)))

其中R₁和R₂是两个特征向量，g₁和g₂是Borel函数，T表示取Borel函数g₁和g₂的最小上界，位于(0，1)，CC(R₁,R₂)表示两个向量的最大相关性，ρ(g₁(R₁),g₂(R₂))的值始终介于(-1，1)之间，表示相关性；融合的步骤为：输入两个最佳特征向量R₁和R₁，执行均值填充使得两个向量维度之间相同，计算两个向量的最大相关性，如果两个向量的最大相关性接近1，则在融合向量中加入这两个特征，如果相关性接近0，那么丢弃这两个特征，最后输出一个融合的特征向量。

Images