一种用于小儿骨龄预测的方法

1.本发明涉及计算机图像视觉和大数据医疗图像应用领域，特别涉及一种用于小儿骨龄预测的方法。


背景技术：

2.传统的骨龄评估通常是对受试者的手和手腕进行拍片，然后由医生根据拍片的情况进行解释。解释的方法包括简单的计数法、图谱法、评分法和计算机骨龄评分系统等。最常用的方法是g-p图谱法和tw3评分法；受生长发育的种族差异与长期趋势的影响，目前最适合中国当代儿童的骨龄标准为chn-05，在2006年，中华-05骨龄标准被编入《中华人民共和国行业标准》，成为中国目前唯一的行业骨龄标准。
3.最近，深度神经网络在医学领域得到了广泛的应用，因为使用深度学习方法自动估计骨龄比人工判断骨龄快得多，同时准确率也远远超过传统方法。hyunkwanglee等人使用 googlenet作为主干，通过显示g&p地图集的三到五个参考图像来确定最终骨龄。 toanducbui等人分别使用faster-rcnn和inception-v4网络进行感兴趣区域(roi)的检测和分类。chuanbinliu等人引入了一个注意代理和一个识别代理分别用于提出有区别的骨骼部分和特征学习和年龄评估。
4.尽管这些方法在baa中可以获得良好的性能，而在实际检测过程中，因为加上年代的发展变化，致使这些方法都不能客观地、准确地评价中国当代儿童的骨龄。基于上述问题，所以急需寻找一种方法，在基于中华-05的标准化样本进行研究，中华-05的标准由2005年中国城市汉族正常儿童组成，能够真实反映当代中国儿童的发育趋势，以中国儿童掌指骨边缘以及掌指骨长度作为参考标准，从而得到正确的骨龄。


技术实现要素：

5.针对现有技术存在的上述问题，本发明要解决的技术问题是：如何对小儿骨龄进行准确预测。
6.为解决上述技术问题，本发明采用如下技术方案：一种用于小儿骨龄预测的方法，包括如下步骤：
7.s100：选用公开数据集，公开数据集包括带有骨龄标签的手部x光片；随机选取部分手部x光片组成数据集a，并将a中所有图片的大小调整为560
×
560；
8.s200：建立tenet模型，tenet模型包括手部拓扑模块、边缘特征增强模块和深度学习网络改进模型d；其中，边缘特征增强模块为改进后的canny边缘检测算法，所述改进后的 canny边缘检测算法使用的是双边滤波去噪算法和大津算法；
9.s300：构建tenet模型的训练集：
10.s310：使用自动色彩均衡算法对a中的图片进行图片预处理，得到a中每张图对应的预处理图；
11.s320：使用手部拓扑模块对a中的图片进行特征提取处理，得到a中每张图对应的
手部拓扑图；
12.s330：使用边缘特征增强模块对a中图片对应的预处理图进行边缘强化，得到a中图片对应的边缘特征增强图；
13.s340：重复s310-s330，遍历a中所有手部x光片，得到a中每张图片对应的预处理图、 a中每张图片对应的手部拓扑图和a中每张图片对应的边缘特征增强图；
14.s350：将s340得到的所有预处理图、手部拓扑图和边缘特征增强图共同组成tenet模型的训练集，其中，每张图对应的预处理图、手部拓扑图和边缘特征增强图组成一个训练样本； s400：预设最大迭代次数，预设d的初始学习率为μ，对d进行训练：
15.s410：将训练集中的所有训练样本作为d的输入，使用adam优化器对d进行训练，实现对d的参数更新，该adam优化器中使用的是mae目标损失函数；
16.当达到最大迭代次数时停止训练，得到训练好的深度学习网络改进模型d’；
17.s420：将使用d’的tenet模型作为训练好的tenet模型；
18.s500：将一个需要预测未知骨龄的小儿手部x光片输入训练好的tenet模型，输出即为对小儿骨龄的预测值。
19.作为优选，所述s200中深度学习网络改进模型d的结构是采用深度学习网络作为主干网络，顺序依上而下，去除掉最顶层的神经网络层，然后在原主干网络的最底层后，依次加上一个3x3的下采样层、一个3x3的最大池化层和一个带有32个神经元的全连接层，并且将a 中的手部x光片对应的性别信息通过全连接层融入图像特征之中。
20.作为优选，所述s310中得到a中图片对应的预处理图的具体步骤如下：
21.s311：从a中选取第c张图：
22.s312：计算第c张图中的任意一个像素点p的自适应滤波中间变量rc(p)；
[0023][0024]
其中，p和q表示第c张图上的两个像素点，ic(p)和ic(q)分别表示p和q的灰度差，d(p,q) 表示p和q之间的距离度量函数，s
α
(
·
)表示亮度表达式函数；
[0025]
s313：计算第c张图的l(p)值，计算表达式如下：
[0026][0027]
其中，l(p)表示像素点p的预处理结果，像素p的映射区间为[0,255]，[minrc(p)，maxrc(p)] 表示l(p)的全定义域；
[0028]
s314：重复s312-s313，遍历第c张图上所有的像素点，将第c张图中所有像素点进行预处理，得到第c张图对应的预处理图；
[0029]
对图片进行预处理可以解决手部x光片背景区域和roi对比度低的问题。
[0030]
作为优选，所述s320中得到a中图片对应的手部拓扑图的具体步骤如下：
[0031]
s321：从a中选取第c张图：
[0032]
s322：使用水平方向矩阵与c上的第t个像素点进行平面卷积，得到第t个像素点横向上的亮度差分近似值g
xt
，具体表达式如下：
[0033][0034]
其中，g
xt
表示第t个像素点横向上的亮度差分近似值，c
t
表示第t个像素点；
[0035]
使用垂直方向矩阵与c上的第t个像素点进行平面卷积，得到第t个像素点纵向上的亮度差分近似值g
yt
，具体表达式如下：
[0036][0037]
其中，g
yt
表示第t个像素点纵向上的亮度差分近似值；
[0038]
s323：计算c上第t个像素点的梯度近似值g
t
，表达式如下：
[0039][0040]
计算c上第t个像素点的方向，表达式如下：
[0041][0042]
其中，α(x
t
,y
t
)表示c上第t个像素点的方向；
[0043]
s324：重复s322-s323，遍历c上所有的像素点，得到每个像素点的卷积数值，所有像素点的卷积数值组成一个矩阵b，该矩阵b即为梯度图c’；
[0044]
s325：对c’进行阈值二值化处理，得到c的手部骨干粗略前景边缘图。
[0045]
s326：使用mediapipe应用框架对c对应的手部骨干粗略前景边缘图上的关节点进行位置预测，得到c的手部拓扑图。
[0046]
利用手部骨干粗略前景边缘图得到的手部拓扑结构图，可以获得结构化的roi信息，可以良好的表达chn-05中所对应的语义特征，即手部掌指骨长度；可以更有利于后期对骨龄的评估。
[0047]
作为优选，所述s330中得到a中图片对应的边缘特征增强图的具体步骤如下：
[0048]
s331：从a中选取第c张图：
[0049]
s332：计算第c张图上每个像素点的空间邻近度和像素值相似度对第c张图进行双边滤波去噪处理，空间邻近度的表达式如下：
[0050][0051]
其中，i,j,k,l表示像素点a和b的坐标值，即a(i,j)和b(k,l)，δd为高斯标准差；
[0052]
像素值相似度的表达式如下：
[0053][0054]
其中，f(i,j)表示图像在a(i,j)点的像素值，f(k,l)表示图像在b(k,l)点的像素值，δr为高斯标准差，wr为计算近点a与中心点b的差异程度；
[0055]
设第c张图经过双边滤波去噪处理后得到的图为c”；
[0056]
s333：设定阈值t，将c”分为背景和目标；
[0057]
采用大津算法计算得到c”的像素最佳阈值，计算公式如下：
[0058]
p0+p1＝1；(9)
[0059][0060][0061]
其中，σ2表示类间方差，表示c”的平均灰度值，p0表示背景在图片上所占比例，p1表示目标在图片上所占比例，且表达式如下：
[0062][0063][0064]
其中，k＝{1,
…
,t}，n表示灰度值的类别，pn表示每类灰度值出现的频率，p0表示0～k之间每类灰度值出现概率的总和，l表示图像的像素级，p1表示k+1～l之间每类灰度值出现概率的总和；
[0065]
m0表示背景的平均灰度值，m1表示目标的平均灰度值，且表达式如下：
[0066][0067][0068]
s334：将公式(10)代入公式(11)并进行简化得到公式(16)，表达式如下：
[0069]
σ2＝p0p1(m
0-m1)2；(16)
[0070]
s335：遍历c”上的所有像素点，将c”上的相同类别灰度值的像素点依次通过公式(9)
‑ꢀ
公式(16)进行计算得到n个σ2；
[0071]
s336：将得到的n个σ2进行降序排列，选择值最大的σ2作为c”的像素最佳阈值；
[0072]
s337：将经过336处理后的c”作为第c张图对应的边缘特征增强图。
[0073]
在传统的canny边缘检测算法中，采用的线性高斯滤波进行去噪，其中使用含有加权系数的卷积对图像进行卷积去噪，实质上就是一种均值模糊，只是高斯模糊是按照加权平均的，距离越近的点权重越大，距离越远的点权重越小，因此计算均值会将图像中的边缘信息以及特征信息“模糊”掉，会丢失很多特征。而双边滤波是一种非线性滤波方法，是结合图像的空间邻近度和像素值相似度的一种折衷处理，同时考虑空域信息和灰度相似性，即通过引入空间域核和值域核达到目的，此种方法不光考虑了像素在空间中位置远近程度的影响，还考虑了像素亮度相近程度的影响。因此双边滤波器可以良好地保留图像边缘特征且过滤掉低频分量的噪音。
[0074]
在传统canny边缘检测算法中，采用的是双阈值算法来进行像素点的分类，然而固定的阈值对于每一张不同的图像不能达到普适性，因此会出现tc、laplacian、prewitt和roberts 算子无法识别骨化中心和掌骨末端边界的固有问题。使用根据不同图片计算适合其最佳阈值的大津算法代替手动选择阈值，使算法性能的不确定性最小化，并且使得最终结果更加符合 chn-05所属语义，即骨化中心清晰可见，为圆盘形，有平滑连续的缘且能观察到骺和骨干的融合程度。
[0075]
相对于现有技术，本发明至少具有如下优点：
[0076]
1.本发明模型使用手部拓扑图、边缘特征增强模块和深度学习网络改进模型构建了骨龄预测评估模型——tenet模型，在tenet模型中使用手部拓扑图加强了对手部x光片的信息特征提取，使得可利用的信息增多；利用改进后的canny边缘检测算法，即边缘特征
增强模块，提高了对手部x光片边缘部分的信息特征提取量和准确度，同时使用大津算法代替手动选择阈值使算法性能的不确定性最小化，确保最终结果更加符合chn-05所属语义。
[0077]
2.在手部x光片中常存在脉冲噪声，传统采用的去噪方法如高斯滤波、中值滤波、方框滤波等，但这些方法不能够良好地保护高频信息，从而模糊图像的边缘细节；本发明使用的双边滤波去噪算法能够通过对图像的空间邻近度和像素相似度进行折衷处理，且考虑空域信息和灰度相似性，在去噪的同时能保留住图像的边缘特征，从而达到保边去噪的效果。
[0078]
3.使用根据不同图片计算适合其最佳阈值的大津算法代替手动选择阈值，使算法性能的不确定性最小化。
[0079]
4.在数据预处理中，使用自动色彩均衡算法解决手部x光片背景区域和roi对比度低的问题。
附图说明
[0080]
图1为tenet模型的整体框架图。
[0081]
图2为本发明所使用的网络结构。
[0082]
图3为改进的canny算法流程图。
[0083]
图4为不同算法处理手部x光片的效果图。
[0084]
图5为本发明方法所使用的数据集的处理结果图。(a)-(d)为使用私有数据库的样本图像， (e)-(h)为使用rsna数据库的样本图像。
具体实施方式
[0085]
下面对本发明作进一步详细说明。提出了一种新的自动骨龄评估模型tenet，通过对 chn-05中骨龄评价的语义描述进行特征抽象，并以算法实现多特征表达的水平融合目的。本发明首先设计了手部拓扑图模块，提供掌指骨的长度以及位置结构信息，以符合chn-05 四方位评估骨龄的语义；其次开发了强化边缘特征模块，识别手部结构，提取有效的边缘信息，满足chn-05对骨头边缘清晰的要求，也增强了骨龄评估的准确性。
[0086]
参见图1-5，一种用于小儿骨龄预测的方法，包括如下步骤：
[0087]
s100：选用公开数据集，公开数据集包括带有骨龄标签的手部x光片；随机选取部分手部x光片组成数据集a，并将a中所有图片的大小调整为560
×
560；
[0088]
s200：建立tenet模型，tenet模型包括手部拓扑模块、边缘特征增强模块和深度学习网络改进模型d，深度学习网络为现有技术，其中，边缘特征增强模块为改进后的canny边缘检测算法，canny边缘检测算法为现有技术，所述改进后的canny边缘检测算法使用的是双边滤波去噪算法和大津算法，双边滤波去噪算法和大津算法为现有技术；
[0089]
所述s200中深度学习网络改进模型d的结构是采用深度学习网络作为主干网络，顺序依上而下，去除掉最顶层的神经网络层，然后在原主干网络的最底层后，依次加上一个3x3的下采样层、一个3x3的最大池化层和一个带有32个神经元的全连接层，并且将a中的手部x 光片对应的性别信息通过全连接层融入图像特征之中。
[0090]
s300：构建tenet模型的训练集：
[0091]
s310：使用自动色彩均衡算法对a中的图片进行图片预处理，得到a中每张图对应
的预处理图；自动色彩均衡算法为现有技术，采用自动色彩均衡算法对原始图像进行区域自适应滤波，完成色差校正，得到空域重构图像，突出图像中的特征和有价值的信息，使得图像更加符合人眼的感知；
[0092]
所述s310中得到a中图片对应的预处理图的具体步骤如下：
[0093]
s311：从a中选取第c张图：
[0094]
s312：计算第c张图中的任意一个像素点p的自适应滤波中间变量rc(p)；
[0095][0096]
其中，p和q表示第c张图上的两个像素点，ic(p)和ic(q)分别表示p和q的灰度差，灰度差表示拟生物学上的侧抑制性，d(p,q)表示p和q之间的距离度量函数，距离度量函数映射滤波的区域适应性，s
α
(
·
)表示亮度表达式函数；
[0097]
s313：计算第c张图的l(p)值，计算表达式如下：
[0098][0099]
其中，l(p)表示像素点p的预处理结果，像素p的映射区间为[0,255]，[minrc(p)，maxrc(p)] 表示l(p)的全定义域，全定义域的目的是使图像实现全局白平衡；
[0100]
s314：重复s312-s313，遍历第c张图上所有的像素点，将第c张图中所有像素点进行预处理，得到第c张图对应的预处理图；
[0101]
s320：使用手部拓扑模块对a中的图片进行特征提取处理，得到a中每张图对应的手部拓扑图；
[0102]
所述s320中得到a中图片对应的手部拓扑图的具体步骤如下：
[0103]
s321：从a中选取第c张图：
[0104]
s322：使用水平方向矩阵与c上的第t个像素点进行平面卷积，得到第t个像素点横向上的亮度差分近似值g
xt
，具体表达式如下：
[0105][0106]
其中，g
xt
表示第t个像素点横向上的亮度差分近似值，c
t
表示第t个像素点；
[0107]
使用垂直方向矩阵与c上的第t个像素点进行平面卷积，得到第t个像素点纵向上的亮度差分近似值g
yt
，具体表达式如下：
[0108][0109]
其中，g
yt
表示第t个像素点纵向上的亮度差分近似值；
[0110]
s323：计算c上第t个像素点的梯度近似值g
t
，表达式如下：
[0111][0112]
计算c上第t个像素点的方向，表达式如下：
[0113][0114]
其中，α(x
t
,y
t
)表示c上第t个像素点的方向；
[0115]
s324：重复s322-s323，遍历c上所有的像素点，得到每个像素点的卷积数值，所有像素点的卷积数值组成一个矩阵b，该矩阵b即为梯度图c’；
[0116]
s325：对c’进行阈值二值化处理，阈值二值化处理为现有技术，得到c的手部骨干粗略前景边缘图。
[0117]
s326：使用mediapipe应用框架对c对应的手部骨干粗略前景边缘图上的关节点进行位置预测，mediapipe应用框架为现有技术，得到c的手部拓扑图。
[0118]
s330：使用边缘特征增强模块对a中图片对应的预处理图进行边缘强化，得到a中图片对应的边缘特征增强图；
[0119]
所述s330中得到a中图片对应的边缘特征增强图的具体步骤如下：
[0120]
s331：从a中选取第c张图：
[0121]
s332：计算第c张图上每个像素点的空间邻近度和像素值相似度对第c张图进行双边滤波去噪处理，空间邻近度的表达式如下：
[0122][0123]
其中，i,j,k,l表示像素点a和b的坐标值，即a(i,j)和b(k,l)，δd为高斯标准差；
[0124]
像素值相似度的表达式如下：
[0125][0126]
其中，f(i,j)表示图像在a(i,j)点的像素值，f(k,l)表示图像在b(k,l)点的像素值，δr为高斯标准差，wr为计算近点a与中心点b的差异程度；
[0127]
设第c张图经过双边滤波去噪处理后得到的图为c”；
[0128]
s333：设定阈值t，将c”分为背景和目标；
[0129]
采用大津算法计算得到c”的像素最佳阈值，计算公式如下：
[0130]
p0+p1＝1；(9)
[0131][0132][0133]
其中，σ2表示类间方差，表示c”的平均灰度值，p0表示背景在图片上所占比例，p1表示目标在图片上所占比例，且表达式如下：
[0134][0135][0136]
其中，k＝{1,
…
,t}，n表示灰度值的类别，pn表示每类灰度值出现的频率，p0表示0～k之间每类灰度值出现概率的总和，l表示图像的像素级，一般为255，p1表示k+1～l之间每类灰度值出现概率的总和；
[0137]
m0表示背景的平均灰度值，m1表示目标的平均灰度值，且表达式如下：
[0138]
[0139][0140]
s334：将公式(10)代入公式(11)并进行简化得到公式(16)，表达式如下：
[0141]
σ2＝p0p1(m
0-m1)2；(16)
[0142]
s335：遍历c”上的所有像素点，将c”上的相同类别灰度值的像素点依次通过公式(9)
‑ꢀ
公式(16)进行计算得到n个σ2，所述灰度级是指像素所属于的灰度级别；
[0143]
s336：将得到的n个σ2进行降序排列，选择值最大的σ2作为c”的像素最佳阈值，计算过程中是在灰度级范围为0-255内，循环遍历所有的像素点，从而得到使公式(16)最大化的灰度级作为该图片的像素最佳阈值；
[0144]
s337：将经过336处理后的c”作为第c张图对应的边缘特征增强图。
[0145]
s340：重复s310-s330，遍历a中所有手部x光片，得到a中每张图片对应的预处理图、 a中每张图片对应的手部拓扑图和a中每张图片对应的边缘特征增强图；
[0146]
s350：将s340得到的所有预处理图、手部拓扑图和边缘特征增强图共同组成tenet模型的训练集，其中，每张图对应的预处理图、手部拓扑图和边缘特征增强图组成一个训练样本；
[0147]
s400：预设最大迭代次数，预设d的初始学习率为μ，对d进行训练：
[0148]
s410：将训练集中的所有训练样本作为d的输入，使用adam优化器对d进行训练，实现对d的参数更新，adam优化器为现有技术，该adam优化器中使用的是mae目标损失函数；
[0149]
当达到最大迭代次数时停止训练，得到训练好的深度学习网络改进模型d’；
[0150]
s420：将使用d’的tenet模型作为训练好的tenet模型；
[0151]
s500：将一个需要预测未知骨龄的小儿手部x光片输入训练好的tenet模型，输出即为对小儿骨龄的预测值。
[0152]
实验验证
[0153]
小儿骨龄评估(baa)是一项重要的人体生理检查，它可以反映人类的生长潜力和性成熟趋势。在临床实践中，《中国青少年儿童手腕骨成熟度及评价方法》(chn-05)是中国放射科医生进行baa的一种广泛使用方法。chn-05以掌指骨长度(ml)以及掌指骨边缘(mm) 作为参考标准，从而估计骨龄。受到chn-05语义描述的启发，本发明提出了一种新的模型，称作tenet，用于自动骨龄评估。在tenet模型中设计了一个手部拓扑模块识别手部关键位置，用于提取结构化语义；设计了一个边缘特征增强模块用于在训练过程中提供精确的骨骼边缘信息。tenet模型可以检测边缘整体信息以及拓扑结构局部信息，从而达到多特征水平融合评估骨龄的目的。
[0154]
实验结果表明，tenet模型在公开数据集rsna上达到了最先进的模型性能5.35meanabsoluteerror(mae)，tenet模型的设计遵循chn-05语义逻辑，因此在临床上使用上来说，它具有非常好的可靠性和可解释性。
[0155]
本实验是在私有数据集和北美放射学会骨龄评估(rsna-baa)公共数据集上评估tenet 模型。私有数据集包含了4954张手部x光片在训练集中，各包含了275张手部x光片在验证集和测试集中；rsna-baa公开数据集包含了5611张手部x光片在训练集中，各包含了 275张手部x光片在验证集和测试集中。在通过tenet模型处理之前，所有手部x光片的大小都被调整为了560
×
560。其中也报告了基本真实骨龄与对应预测骨龄之间的平均绝对误差。
[0156]
表1本发明模型的消融研究
[0157]
exp.oiphiimcatl1mae1(私有数据集)√
ꢀꢀ
√7.152(私有数据集) √ √6.923(私有数据集) √√√5.764(公开数据集) √√√5.80
[0158]
(exp为experiment，表示实验；oi为originalimage，表示原图；phi为 pre-processedhandimage，表示手部预处理图；imcat为 improvedadaptivethresholdingofcannyedgedetection，表示手部边缘强化图；l1为loss1，表示 mae损失函数；mae为meanabsoluteerror,表示平均绝对损失。)
[0159]
实现细节
[0160]
通过tensorflow1.9实现tenet，并在nvidiatitanrtxgpu和32gram的系统上完成培训，耗时约8小时。我们用了100个epochs来训练tenet的全过程。批处理大小为32。初始学习率为3
×
10-3
,50个epoch后降为10-3
80个epoch后降为10倍。使用的优化器是adam。
[0161]
tenet模型消融实验
[0162]
首先研究了各模块在tenet模型中的作用，即原图、预处理手部图像(phi)和边缘特征增强(imcat)。实验结果如表1显示，没有这两个模块，网络退化到xception，结果mae得分为7.19个月；将预处理图代替原图应用于该模型时，所得mae分数为6.92个月；当边缘特征增强模块应用于该模型时，得到的mae分数为5.76个月，mae分数提高了1.16个月；通过实际的实验数据说明了获取边缘特征的必要性和imcat算法的有效性。
[0163]
结合这两个组成部分，tenet模型在公开数据集上实现了5.80个月的mae。综上所得，由私有数据集的结果优于公开数据集的结果可以佐证tenet模型设计由chn-05手腕骨标准启发而来的正确性。注：mae值越小说明效果越好。
[0164]
表2imcat检测的消融研究
[0165][0166][0167]
手部边缘特征消融实验
[0168]
tenet模型检测手部边缘特征，并将拓扑结构联合起来进行骨龄评估。通过观察，imcat 也得益于其模块的设计，利用峰值信噪比(psnr)、均方误差(mse)、结构相似度(ssim)、均方根误差(rmse)和平均绝对误差(mae)来衡量性能，结果如表2所示。将传统的canny (tc)算法作为基线，利用四个边缘运算器作为基线，本发明中使用的imcat算法明显优于其他算法，说明imcat有能力更好地捕捉手部的整体边缘细节。
[0169]
在rsna公共数据集和私人数据集的结果表明，tenet模型实现了不错的表现。此外， tenet模型由chn-05先验知识设计而来，因此对于临床实践它有良好的可解释性和可
靠性。此外，本发明还为在尝试把特征获取模块和评分回归模块进行融合，形成一个端到端的自动骨龄评估框架提供了可靠的实验基础。
[0170]
最后说明的是，以上实施例仅用以说明本发明的技术方案而非限制，尽管参照较佳实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，可以对本发明的技术方案进行修改或者等同替换，而不脱离本发明技术方案的宗旨和范围，其均应涵盖在本发明的权利要求范围当中。