一种基于移动最小二乘算法的双能X射线骨密度检测方法与流程

一种基于移动最小二乘算法的双能x射线骨密度检测方法
技术领域
1.本发明属于及医疗器械技术领域，具体涉及应用双能x-射线骨密度测量领域。


背景技术：

2.骨质疏松疾病是世界六种多发疾病之一。但该疾病在早期无明显症状，并不会引起人们的重视，所以如何预防和诊断骨质疏松疾病十分重要。骨密度测量是目前用来评价骨质流失，诊断骨质疏松疾病的重要方法。人体的骨骼是由致密的皮质骨和骨髓腔内许多类似海绵状的的松质骨构成，其主要成分是羟基磷灰石，主要是以钙和磷元素存在。骨质疏松疾病就是钙、磷元素的流失，导致骨矿物质含量下降，医学上测量骨密度的主要就是羟基磷灰石的密度。
3.双能x射线吸收法是目前世界卫生组织认证最有效测量骨密度的方法，主要是由x射线球管通过开关脉冲法或k边缘过滤法产生两种不同能量的x射线，穿过等效人体骨骼和软组织的模体后发生衰减，再由探测器接收衰减后的光子能量，通过计算机及相关asic集成电路和算法得到骨密度与模体厚度之间对应关系再由屏幕显示roi区域的骨密度。双能x射线吸收法具有扫描速度快，精度高，辐射低，可以测量人体脊椎，股骨，腰椎等全身骨骼范围，具有较高的准确度和精确度。
4.双能x射线骨密度诊断系统通过探测器采集衰减的光子能量后，需要进一步用算法处理x射线衰减能量与骨骼、软组织的厚度数据，以进行后续的骨密度信息分析与诊断。传统的骨密度算法主要有查找表法、多项式拟合法。查找表法是将高低能衰减值与不同骨骼厚度一一对应，通过扫描骨骼与软组织模体得到的x射线衰减能量，根据对应表格找到相关厚度值。多项式拟合法是结合物质分解算法和x射线的衰减公式来进一步推导，运用泰勒公式展开，设置变化步长，对其进行一阶求导，找到目标函数的最小值。根据步长变化，得到新的目标函数，继续对其一阶求导，重复上述步骤，直至所有的离散数据点都遍历结束。即可得到高低能衰减值与物体厚度的关系。
5.双能x射线的骨密度诊断系统的常用传统算法，如查找表法是通过建立一个数据库，将采集到的x射线高低能衰减值与不同厚度骨骼，软组织相互对应，在测量骨密度时，根据x射线衰减值找到对应厚度，再进行厚度-骨密度的线性转换，即可获得对应部位的骨密度值。该方法的运行时间较长，需要较大数据内存。多项式拟合法通过二次或三次多项式进行近似表达，计算过程使用梯度法或牛顿迭代法求值，该方法应对噪声，探测器像素不均匀性具有不稳定性。多项式拟合计算量大，且计算过程繁琐，对拟合系数要求较高，否则拟合曲面难以达到最佳而容易增加误差。
6.基于光子计数探测器的双能x射线骨密度仪是采用碲锌镉晶体的光子计数探测器，光子计数探测器在理想情况下，能够充分利用x射线能谱信息，有效提高图像质量，获取物质成分信息，增强软组织对比度，但是由于制作工艺不成熟，光子计数探测器在模块边缘容易出现坏点，以及像素单元之间存在不均匀性，容易相互影响，探测器像素单元之间不均匀性对传统算法误差影响大。
7.移动最小二乘算法最常应用于力学、遥感地理中。基于传统的骨密度算法，为提高骨密度拟合精度，研究人员加强对系统校正和校正模体材料的研究。光子计数探测器的像素不均匀性对传统骨密度算法误差影响较大，移动最小二乘算法旨在解决光子计数探测器的边缘不均匀性并以其特有的紧支性实现更好的拟合以及校正。移动最小二乘法在二维散点拟合具有出色的表现，骨密度测量涉及三维数据，需在公式中设置两个自变量并确定权函数后使得数值计算过程的矩阵为非奇异矩阵，即矩阵可逆才能在三维空间成功拟合。


技术实现要素：

8.本发明旨在解决以上现有技术的问题。提出了一种基于移动最小二乘算法的双能x射线骨密度检测方法。本发明的技术方案如下：
9.一种基于移动最小二乘算法的双能x射线骨密度检测方法，其包括以下步骤：
10.获取模体在光子计数探测器的高能、低能x射线下的探测数据，并进行去噪在内的预处理；以及获取高能、低能x射线的入射能量；
11.根据经过预处理的探测数据和入射能量计算得到高低能衰减值；
12.在离散的高低能衰减值数据点上，按照模体厚度进行网格划分；
13.确定移动最小二乘算法的权函数后，在每个网格内，以目标节点为中心划分影响域半径，在空间域内，对每个高低能衰减值的数据点进行遍历，在影响域内设置变化步长；
14.采用移动最小二乘算法得到最终的拟合关系，当所有的节点遍历结束后，得到模体厚度-高低能衰减值拟合曲面；
15.确定需要测量的骨密度roi区域，光子计数探测器扫描待测部位得到x射线衰减值，计算得到对应的模体厚度，再进行一次线性变换，模体厚度与模体密度的乘积就是对应部位的骨密度值。
16.进一步的，所述根据经过预处理的探测数据和入射能量计算得到高低能衰减值，具体包括：
17.通过公式(1)，(2)分别得到高低能下的x射线衰减值，ih、i
l
分别代表高能、低能x射线在放置模体下的探测器采集数据，i
oh
、i
ol
分别代表高能、低能x射线的入射能量，yh、y
l
分别指代高能、低能条件下x射线的光子能量的衰减程度；
[0018][0019][0020]
进一步的，所述移动最小二乘法的权函数的类型有高斯权函数、三次样条权函数及双插值权函数，经实验证明，高斯权函数拟合效果最佳。
[0021]
进一步的，所述移动最小二乘算法的表达式如(3)所示，其中α(x)为系数向量矩阵，p(x)为基函数向量，i代表实际数据点，f(x)为拟合函数，对(3)式求解范数得到公式(4)，w(x)为权函数，x是目标节点，xi为其他离散数据点，j表示对(3)式取范数，yi为实际数据点，为求解得到其最小值，公式(4)对α(x)求偏导得到(5)；
[0022]
[0023][0024][0025]
将公式(5)中的矩阵a、b展开，(6)、(7)可以看出，两个矩阵是基函数与权函数构成的系数矩阵；
[0026]
a(x)＝p
t
w(x)p
ꢀꢀ
(6)
[0027]
b(x)＝p
t
w(x)
ꢀꢀ
(7)
[0028]yt
＝[y1,y2,...,yn]
ꢀꢀ
(8)
[0029]
w(x)为权函数、p代表基函数、yn代表实际x射线衰减值。
[0030]
基函数为k阶多项式，将其展开如(9)所示，权函数为对角矩阵，其表达式由(10)所示，将公式(5)进行变形，带入公式(3)最终得到移动最小二乘的标准式为(11)；
[0031][0032][0033]
f(x)＝p
t
a-1
(x)b(x)y
ꢀꢀ
(11)
[0034]
p为基函数、pm(xn)为m阶完备多项式下第n个x的值、m代表基函数p的阶数、n表示x的个数。
[0035]
由于数据为三维空间拟合，所以采用如公式(12)、(13)所示移动最小二乘方法，其中x，z分别代表模块al和聚碳酸酯高分子化合物的厚度，f
l
代表低能下的拟合值，fh为高能下的拟合值，带入y
l
、yh，可以分别得到模体厚度在高低能下的拟合关系
[0036]fl
(x,z)＝p
t
(x,z)a-1
(x,z)b(x,z)y
l
ꢀꢀ
(12)
[0037]fh
(x,z)＝p
t
(x,z)a-1
(x,z)b(x,z)yhꢀꢀ
(13)
[0038]
进一步的，所述获取模体在光子计数探测器的高能、低能x射线下的探测数据，并进行去噪在内的预处理；以及获取高能、低能x射线的入射能量，具体包括：
[0039]
搭建光学平台，用来采集模体的x射线高低能衰减值，主要包括以下部分：x射线光机、光子计数探测器、光学平台、校正模体、asic集成电路和中央计算机；
[0040]
x射线光机用于产生并发射x射线，光子计数探测器在骨密度测量系统中用于接收衰减后的x射线光子能量，光学平台用于放置x射线光机、校正模体和光子计数探测器，使其位于同一水平线上，具有准直功能，校正模体在系统中获得校正参数，建立骨密度与高低能衰减值的对应关系。asic集成电路将探测器采集到的光电信息转化成数字信息经串口电路传输至中央计算机由其实现roi图像对应骨密度值的屏幕显示；
[0041]
在室温22℃，湿度50％的环境下打开光子计数探测器，对光子计数探测器进行设
备校正，调试工作参数，光子计数探测器在测试前必须达到热稳定性，光子计数探测器在室温22℃下应该放在+5vdc电源下持续至少60分钟，并通过无源扫描、本底扫描在内的校正方法对光子计数探测器每一个像素功能性能进行检验，记录光子计数探测器的本底噪声；
[0042]
打开x射线光机，将x射线球管预热两分钟，保证光子能量达到稳定的输出，然后分别在70kev和40kev的x射线能量下，测量高低能的入射能量。
[0043]
进一步的，所述光子计数探测器是采用碲锌镉晶体的光子计数探测器。
[0044]
进一步的，所述在每个网格内，以目标节点为中心划分影响域，在影响域内，目标节点仅受相邻节点的不同权重影响，影响域外其他节点对目标节点没有影响作用。在空间域内，对每个高低能衰减值的数据点进行遍历，在影响域内设置变化步长；变化步长与权函数类型相关，选取高斯权函数，变化步长为影响域半径的一半。以目标节点为中心，确定一个圆形区域，影响域半径为目标节点到圆上的距离，需要保证矩阵可逆情况下设置影响域半径。
[0045]
8、进一步的，确定需要测量的骨密度roi区域，光子计数探测器扫描待测部位得到x射线衰减值，由公式(12)、(13)计算得到对应的模体厚度，再进行一次线性变换得到式(14)，模体厚度x与模体密度ρ的乘积就是对应roi部位的骨密度(bmd)值。
[0046]
roi
bmd
＝xρ
ꢀꢀ
(14)
[0047]
本发明的优点及有益效果如下：
[0048]
传统的骨密度检测系统中，常用光子计数探测器接收x射线衰减后的光子能量，光子计数探测器相较于闪烁体探测器具有以下优点：对环境温度不敏感，无潮解，化学物质稳定，无需光电倍增管，具更高的电荷收集特性和能量分辨率，但由于制备工艺不成熟和晶体生长缺陷，会导致像素之间的不均匀性和探测器边缘容易出现坏点，常用的传统算法受光子计数探测器边缘像素影响较大，会使得传统骨密度算法出现较大误差，影响骨密度值的诊断结果。本发明采用了一种基于移动最小二乘算法的双能x射线骨密度测量系统，移动最小二乘算法具有独特的紧支性，通过加强相邻节点的影响，使得拟合曲面与离散点的距离更小，有效提高拟合精度。本发明在基于光子计数探测器的双能x射线骨密度测量系统上使用移动最小二乘算法计算骨密度值，需在该算法公式中设置两个自变量并确定权函数后使得数值计算过程的矩阵为非奇异矩阵，即矩阵可逆才能在三维空间成功拟合，改进得到公式(12)、(13)，使其具有自适应性，局部拟合特点。在运用代码的过程中，调用不同权函数验证数据会出现运行错误，则需要进行代码调试。应注意在采集高低能衰减值时，每组厚度数据应当采集10份，为了得到探测器稳定状态的数据，最终应挑选其中数据值最为接近的6份数据平均值用来进行最终拟合。该算法选择高斯权函数，通过灵活设置影响域半径大小，使得该点的函数值只会受到附近影响域内节点的影响，不会被区域外的节点影响，改善探测器像素之间的不均匀性和边缘坏点的影响。移动最小二乘算法拟合x射线高低能投影值与厚度之间的数据关系，通过权函数控制样点影响权重以及样点影响半径，移动最小二乘法较好解决了拟合光滑性和局部化问题，在数据拟合方面具有很大优势，其具有对散点适应性强、具备局部拟合，优化探测器性能优点。
附图说明
[0049]
图1是本发明提供优选实施例双能x射线骨密度诊断系统示意图；
[0050]
图2是骨密度系统流程图；
[0051]
图3是移动最小二乘算法流程图；
[0052]
图4是移动最小二乘法空间曲面拟合示意图；
[0053]
图5是传统多项式与移动最小二乘拟合结果对比。
具体实施方式
[0054]
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、详细地描述。所描述的实施例仅仅是本发明的一部分实施例。
[0055]
本发明解决上述技术问题的技术方案是：
[0056]
本发明提出一种基于移动最小二乘算法的双能x射线骨密度检测系统，以解决光子计数探测器像素不均匀性对传统骨密度拟合算法的误差影响，本系统提高骨密度测量精度，有效校正探测器像素不均匀性，优化探测器性能。
[0057]
为实现上述目的，本发明实施提供如下技术方案。
[0058]
1、搭建如图(1)所示的光学平台，用来采集模体的x射线高低能衰减值。主要包括以下部分：x射线光机，光子计数探测器，光学平台，校正模体，asic集成电路和中央计算机等。
[0059]
2、在室温(22℃)，湿度50％的环境下打开探测器，对探测器进行设备校正，调试工作参数，探测器在测试前必须达到热稳定性，所以探测器在室温(22℃)下应该放在+5vdc电源下持续至少60分钟，并通过无源扫描、本底扫描等校正方法对探测器每一个像素功能性能进行检验，记录探测器的本底噪声。
[0060]
3、打开x射线光机，将x射线球管预热两分钟，保证光子能量达到稳定的输出，然后分别在70kev和40kev的x射线能量下，测量高低能的入射能量。
[0061]
4、在自然界中，金属al的相对原子质量与人体骨骼的相对原子质量相近，聚碳酸酯高分子化合物与人体软组织的相对原子质量相近，由于相对原子序数相近的物体对x射线的吸收率可视为等效，所以实验室常用al块和聚碳酸酯高分子化合物代替人体中的骨骼，软组织进行模体校正实验。将不同厚度的al块和聚碳酸酯高分子化合物由图(1)所示成梯形叠放，扫描一次，即可获得不同厚度组合下的x射线高低能衰减值。
[0062]
5、由于探测器直接采集到的光子能量数据较大，需要进行数据预处理，通过公式(1)，(2)分别得到高低能下的x射线衰减值，ih、i
l
分别代表高能、低能x射线在放置模体下的探测器采集数据。i
oh
、i
ol
分别代表高能、低能x射线的入射能量，yh、y
l
分别指代高能、低能条件下x射线的光子能量的衰减程度。
[0063][0064][0065]
6、移动最小二乘算法的表达式如(3)所示，其中α(x)为系数向量矩阵，p(x)为基函数向量，i代表实际数据点，f(x)为拟合函数，对(3)式求解范数得到公式(4)，w(x)为权函数，x是目标节点，xi为其他离散数据点，yi为实际数据点，为求解得到其最小值，公式(4)对α(x)求偏导得到(5)。
[0066][0067][0068][0069]
7、将公式(5)中的矩阵a、b展开，(6)、(7)可以看出，两个矩阵是基函数与权函数构成的系数矩阵。
[0070]
a(x)＝p
t
w(x)p
ꢀꢀ
(6)
[0071]
b(x)＝p
t
w(x)
ꢀꢀ
(7)
[0072]yt
＝[y1,y2,...,yn]
ꢀꢀ
(8)
[0073]
8、基函数为k阶多项式，将其展开如(9)所示，权函数为对角矩阵，其表达式由(10)所示,将公式(5)进行变形，带入公式(3)最终得到移动最小二乘的标准式为(11)。
[0074][0075][0076]
f(x)＝p
t
a-1
(x)b(x)y
ꢀꢀ
(11)
[0077]
9、由于本系统设计到的数据为三维空间拟合，所以在采用如公式(12)、(13)所示移动最小二乘方法，其中x，z分别代表模块al和聚碳酸酯高分子化合物的厚度，f
l
代表低能下的拟合值，fh为高能下的拟合值。带入y
l
、yh，可以分别得到模体厚度在高低能下的拟合关系
[0078]fl
(x,z)＝p
t
(x,z)a-1
(x,z)b(x,z)y
l
ꢀꢀ
(12)
[0079]fh
(x,z)＝p
t
(x,z)a-1
(x,z)b(x,z)yhꢀꢀ
(13)
[0080]
7、最后，扫描人体，选定roi区域，得到对应的x射线高低能衰减数据，计算出对应的骨骼与软组织的厚度，厚度与模块密度的乘积就是人体对应部位的骨密度。
[0081]
使用移动最小二乘算法拟合模体厚度-高低能衰减值的步骤如图(3)所示，首先，对探测器采集的数据进行相应的预处理，去除本底噪声和电子噪声，根据公式(1)、(2)计算得到高低能衰减值。在离散的数据点上，按照模体厚度的个数进行网格划分，移动最小二乘法常用的权函数的类型有高斯权函数，三次样条权函数，双插值权函数。确定权函数后，划分影响域半径。在空间域内，对每个数据点进行遍历。如图(4)所示，在每个网格内，以目标节点为中心划分影响域半径，在影响域内设置变化步长，越靠近目标节点的数据点对其影响权重越大，反之，随着距离越远则对目标节点影响越小。处在影响域外的数据节点则对目
标节点的作用效果为零，根据公式(12)、(13)得到最终的拟合关系。移动最小二乘算法可以灵活设置影响域半径，较好地解决了拟合光滑性和局部化问题，有效校正光子计数探测器像素不均匀特性。当所有的节点遍历结束后，可以得到模体厚度-高低能衰减值拟合曲面。最后，根据用户个人需求确定需要测量的骨密度roi区域，探测器扫描人体得到的x射线衰减值，计算得到对应的模体厚度，再进行一次线性变换，模体厚度与模体密度的乘积就是对应部位的骨密度值。
[0082]
本发明阐述了一种基于移动最小二乘算法的双能x射线骨密度诊断系统的结构以及移动最小二乘算法在拟合模体厚度-高低能衰减值的具体流程，其有益效果主要包括：
[0083]
从光子计数探测器的双能x射线骨密度系统上看，移动最小二乘算法拟合骨密度更加精确。在临床实践中，更方便医生诊断骨质疏松疾病，帮助病人提前预防和治疗骨骼疾病。移动最小二乘算法相较于传统的查找表法对内存的需求更小，运行时间更短，有效提高了骨密度系统的性能。对比传统的多项式拟合算法，移动最小二乘算法在一般的拟合算法基础上添加了权函数和影响域，使得目标节点在影响域内仅受相邻节点的影响，随距离越远影响权重越小，并影响域外数据点对目标节点没有权重影响，在拟合过程中更具有灵活性和解决局部化问题。
[0084]
在4
×
16通道的光子计数探测器接收x射线能量时，由图(5)可见，探测器边缘像素的不均匀性对传统算法的误差影响较大，由此可见，传统算法在拟合厚度-高低能衰减值时，无法校正探测器不均匀性对系统带来的误差影响。传统的多项式拟合方法无法兼顾局部化处理需求且对于大量数据会引起模型设定困难和计算不稳定问题。采取移动最小二乘算法的双能x射线骨密度系统在测量骨密度时具有明显的优越性。移动最小二乘算法在数学原理上比传统算法增加权函数和影响域，有效处理多维数据，而移动最小二乘法较好解决了拟合光滑性和局部化问题，在数据拟合方面具有很大优势，其具有对散点适应性强、具备局部拟合或插值特点以及精度高等优点，适用于像素单元信号均匀性较差、边缘坏点较明显的光子计数探测器的探测数据拟合。
[0085]
本发明将双能x射线骨密度系统与移动最小二乘算法相结合，对拟合骨密度精度有明显的提升，针对光子计数探测器像素不均匀性有明显提升。具有现实可行性，和广泛的应用价值。
[0086]
还需要说明的是，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、商品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、商品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个
……”
限定的要素，并不排除在包括所述要素的过程、方法、商品或者设备中还存在另外的相同要素。
[0087]
以上这些实施例应理解为仅用于说明本发明而不用于限制本发明的保护范围。在阅读了本发明的记载的内容之后，技术人员可以对本发明作各种改动或修改，这些等效变化和修饰同样落入本发明权利要求所限定的范围。