CN110289075A - 一种基于模糊熵的直接子野优化方法及系统 - Google Patents
一种基于模糊熵的直接子野优化方法及系统 Download PDFInfo
- Publication number
- CN110289075A CN110289075A CN201910352358.3A CN201910352358A CN110289075A CN 110289075 A CN110289075 A CN 110289075A CN 201910352358 A CN201910352358 A CN 201910352358A CN 110289075 A CN110289075 A CN 110289075A
- Authority
- CN
- China
- Prior art keywords
- ziye
- gradient
- fuzzy entropy
- fuzzy
- matrix
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Classifications
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61N—ELECTROTHERAPY; MAGNETOTHERAPY; RADIATION THERAPY; ULTRASOUND THERAPY
- A61N5/00—Radiation therapy
- A61N5/10—X-ray therapy; Gamma-ray therapy; Particle-irradiation therapy
- A61N5/103—Treatment planning systems
- A61N5/1031—Treatment planning systems using a specific method of dose optimization
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61N—ELECTROTHERAPY; MAGNETOTHERAPY; RADIATION THERAPY; ULTRASOUND THERAPY
- A61N5/00—Radiation therapy
- A61N5/10—X-ray therapy; Gamma-ray therapy; Particle-irradiation therapy
- A61N5/103—Treatment planning systems
- A61N5/1039—Treatment planning systems using functional images, e.g. PET or MRI
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
- G16H20/00—ICT specially adapted for therapies or health-improving plans, e.g. for handling prescriptions, for steering therapy or for monitoring patient compliance
- G16H20/40—ICT specially adapted for therapies or health-improving plans, e.g. for handling prescriptions, for steering therapy or for monitoring patient compliance relating to mechanical, radiation or invasive therapies, e.g. surgery, laser therapy, dialysis or acupuncture
Landscapes
- Health & Medical Sciences (AREA)
- Engineering & Computer Science (AREA)
- Biomedical Technology (AREA)
- General Health & Medical Sciences (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Public Health (AREA)
- Radiology & Medical Imaging (AREA)
- Pathology (AREA)
- Life Sciences & Earth Sciences (AREA)
- Animal Behavior & Ethology (AREA)
- Veterinary Medicine (AREA)
- Medical Informatics (AREA)
- Primary Health Care (AREA)
- Epidemiology (AREA)
- Urology & Nephrology (AREA)
- Surgery (AREA)
- Radiation-Therapy Devices (AREA)
Abstract
本发明针对利用局部梯度信息生成子野则易导致局部最优结果、利用全局梯度信息生成子野则计算量大等问题，提供了一种能够更加快速、准确地生成可交付使用子野的基于模糊熵的直接子野优化方法及其系统：首先采用基于模糊熵分割原理生成子野形状，然后进行子野权重优化，以得到满足临床要求的治疗方案,本发明可以快速、准确地利用全局梯度信息得到子野形状，计算量减少的同时提高了计划质量，优化所得子野数目较少，提高了优化效率；利用全局梯度信息生成子野形状，有利于提高靶区剂量适形度、降低靶区周围危及器官以及正常组织接受的放射剂量，从而达到提高治疗增益比、提高肿瘤放射治疗精确度以及降低正常组织由放射引起并发症的概率的临床治疗要求。
Description
技术领域
本发明一种基于模糊熵的直接子野优化方法及系统，属于智能计算方法技术领域。
背景技术
调强放射治疗广泛应用于临床癌症治疗。调强放射治疗方案优化通过调整照射野的参数使 高剂量区的形状与靶区形状达到基本一致，同时尽可能减少靶区周围危及器官以及正常组织接 受的放射剂量，从而达到提高肿瘤放射治疗精确度、降低正常组织由放射引起并发症的概率的 临床治疗要求。
当前应用最为广泛的是MLC调强技术。MLC调强技术主要有静态调强和动态调强两种模式。 实现静态调强主要有两步法和直接子野优化法，应用于系统对应称之为两步法系统和直接子野 优化系统。
传统调强技术采用两步法及其系统实现，首先生成符合临床要求的注量图，即求解注量图 优化问题；然后执行序列分割步骤以获得可交付使用的子野以及对应的机器跳数。直接子野优 化方法及其系统将MLC的物理约束加入到优化过程中，直接对子野的形状和强度进行优化，避 免了序列分割步骤。两步法及其系统的缺陷是在于生成计划过程中未考虑MLC的物理约束，同 时序列分割会导致计划质量退化、治疗时间延长。而相较于两步法及其系统，直接子野优化方 法及其系统能够直接生成满足MLC物理约束的子野，显著提高调强放射治疗计划质量，因此， 直接子野优化方法及其系统成为近几年众多学者研究的热点。
直接子野优化方法及其系统的每一次迭代分为两步：子野形状优化和子野权重优化。在子 野形状优化中，生成能够最大限度改善目标函数的可交付子野并加入到治疗计划中；在子野权 重优化中，所有已生成子野的权重被重新优化。子野形状优化在直接子野优化方法及其系统中 至关重要，并且影响调强放射治疗计划的质量。
应用于直接子野优化方法及其系统的算法，主要分为三类：随机搜索法、基于局部梯度法 和column-generation法。随机搜索法(具体参见Shepard D M,Earl M A,Li X A,et al. Direct aperture optimization:A turnkey solution for step-and-shootIMRT.Medical Physics,2002,29(6):1007-1018)应用于子野形状优化，需要靶区在对应射野方向上投影所 得到初始子野作为搜索起点，叶片随机左右移动对初始子野进行修正。随机搜索方法的优点是 能够以一定概率跳出局部最小值，缺点是由于搜索时所依据的是局部梯度信息，因此该方法不 能保证所得优化结果为全局最优。基于局部梯度的算法应用于子野形状优化，需要通过两步法 得到初始子野，由于只能根据局部梯度信息对子野进行生成和优化，而且在优化过程中极容易 陷入局部最优，因此同样不能保证所得优化结果为全局最优。而column-generation法(具体 参见Romeijn H，Ahuja R，Dempsey J，Arvind K.Acolumn generation approach to radiation therapy treatment planning usingaperture modulation.SIAM J Optim，2005， 15：838-862)具有严格的理论推导，是采用遍历的搜索策略，能够保证所得优化结果是全局 最优，是三种方法中精度最高的算法，但是这种方法的计算量是非常大的，增加了方案优化所 需的时间。此外，如果在MLC物理约束复杂时，column-generation法需要用到图论的方法， 通过求解最短路径问题来确定子野的形状，而图论的方法很难综合地考虑放射治疗设备物理约 束条件与目标函数值充分下降之间的关系，容易造成频繁增加\删除子野的问题。
发明内容
本发明针对利用局部梯度信息生成子野则易导致局部最优结果、利用全局梯度信息生成子 野则计算量大等问题，提供了一种能够更加快速、准确地生成可交付使用子野的基于模糊熵的 直接子野优化方法及其系统。
为了解决上述技术问题，本发明采用的技术方案为：一种基于模糊熵的直接子野优化方法， 按下述步骤实施：
第一步：输入直接子野优化所需基本信息，并根据输入基本信息中的目标函数等信息计算 子野的梯度矩阵；
将目标函数设置如下：
则子目标函数表达式为：
fgEUD(D)＝H(gEUD(D)-gEUD0)·(gEUD(D)-gEUD0)；
fTCP(D)＝H(TCP0-ln(TCP(D)))·(TCP0-ln(TCP(D)))；
fNTCP(D)＝H(ln(1-NTCP0)-ln(1-NTCP(D)))·(ln(1-NTCP0)-ln(1-NTCP(D)))；
其中：
D(x)为剂量分布，且D(x)＝Wx，即剂量分布是关于剂量沉积矩阵W和注量矩阵x的线性 函数；
fl(D(x))为第l个子目标函数；
ξl为表示对应子目标函数重要性的权重系数；
L为子目标函数个数；
Dmin为靶区处方剂量；
H(x)为阶跃函数；
N为器官内的所有体素点；
Di为器官第i个体素的剂量值；
D1，D2为两个剂量约束点处剂量值；
Dmax为耐受剂量；
Dmean为靶区处方剂量或正常组织的平均耐受剂量；
其中：
α为剂量体积效应因子；
其中：
表示在第i个体素单元内的克隆源性细胞在受到Di剂量的 照射后的存活率；
λ为内在放射敏感性参数，表示不可修复的放射损伤；
γ为两次分次照射之间可以修复的损伤；
其中：
为标准正态累积分布函数；
D50为全部体积和部分体积V受照射时，由放射损伤引起的NTCP值为50％时所需的剂量；
n为体积效应因子；
m为控制NTCP剂量效应曲线斜率的参数；
第二步：根据子野梯度矩阵，计算该子野梯度矩阵的模糊熵；
首先将子野梯度矩阵转化为具有有限个灰度阶的梯度图Ga(m，n)，再利用子野梯度图的直 方图以及隶属函数μ(g)，计算子野梯度图的模糊熵，模糊熵的计算如下：
其中：
梯度图大小为M×N；
有限个灰度级范围为[Ll，Lu]；
μ(g)为定义在Lu-Ll级灰度上的隶属函数；
Sn(x)＝-x ln x-(1-x)ln(1-x)为Shannon函数；
f(g)为梯度图中灰度级为g的元素的个数；
由模糊熵分割原理可得，梯度图模糊熵E(G)的计算取决于不同隶属函数μ(g)，若将整个梯 度图看作一个模糊集，则一般选取标准S函数，隶属函数形式如下：
或
若按照梯度图中梯度元素的正负分为两个模糊集，则两个模糊集的隶属函数μ1(g)和μ2(g)，可 定义如下：
式中：
交叉点
[a，c]为模糊区域；
Δb＝c-a为模糊区域带宽；
模糊熵E(G)由带宽Δb和交叉点b决定；
第三步：利用模糊熵原理对所得梯度矩阵进行分割并修正得到可交付使用的子野，加入到 子野集合中；
通过计算子野梯度矩阵的模糊熵，确定使模糊熵取得极值的分割阈值，即分割阈值为使得 模糊熵取得最大值的梯度值；梯度元素大于该阈值，对应子射束关闭，小于该阈值，对应子射 束打开，最终形成一个初始子野形状；然后，根据MLC的物理约束，对初始子野形状进行修正， 最终得到可交付使用的子野形状；
第四步：对已有的全部子野的权重进行优化；输出优化结果。
子野权重优化采用的梯度类算法，以L-BFGS-B算法为例，其算法原理具体如下：
设定初始点x0和决定有限内存修正存储次数的整数m，定义初始有限内存矩阵并令k：＝0；
(1)如果投影梯度满足收敛测试||P(xk-gk，l，u)-xk||∞＜10-5则停止；
(2)利用CP算法计算广义柯西点；
(3)通过直接法计算搜索方向dk；
(4)在可行域内，沿着dk方向执行限定最大步长的线搜索来计算步长λk，令xk+1＝xk+ λkdk。线搜索始于单位步长，满足并且尝试满足 其中，α＝10-4，β＝0.9；
(5)计算
(6)如果yk满足eps＝2.2×10-16，将sk和yk添加到Sk和Yk，如果超过m 次更新被存储，则从Sk和Yk中删除最早存储的列；
(7)更新Lk和Rk，并且设置
(8)令k：＝k+1，转(1)；
其中：
满足充分下降条件；
是子空间上带约束的二次型最小化的近似解；
fk是函数值；
gk是梯度值；
{sk，yk}是算法的修正对；
sk＝xk+1-xk，yk＝gk+1-gk，Yk＝[yk-m，…，yk-1]，Sk＝[sk-m，…，sk-1]，P(xk- gk，l，u)-xk是投影梯度；
Lk和Rk是m×m的矩阵；
优选地，若确定使模糊熵取得极值的分割阈值后，将子野梯度图中，大于分割阈值的梯度 置为0，小于分割阈值的梯度置为1，梯度图中所有为1的元素构成子野的初始形状，即梯度 图中为0的位置所对应的子射束关闭，梯度图中为1的位置所对应的子射束打开，最终形成一 个初始子野形状；对获得子野的初始形状结合MLC的物理约束进行修正，以得到可交付使用的 子野，加入到治疗计划中。
本发明涉及一种基于模糊熵的直接子野优化系统，其特征在于，包括：
信息输入模块，用于输入优化系统所需基本信息；
信息预处理模块，根据输入的基本信息，获得子野的梯度矩阵；
子野形状生成模块，利用模糊熵原理对子野梯度矩阵进行分割并修正，以得到可交付使用 子野加入到子野集合中；
子野权重优化模块，对已有的全部子野的权重进行优化；
优化结果输出模块，输出优化结果信息。
优选地，所述子野形状生成模块，先将子野的梯度矩阵转换为子野的梯度图，利用子野梯 度图的直方图以及隶属函数，计算子野梯度图的模糊熵，设定分割阈值为使得模糊熵取得最大 值的梯度值，将子野梯度图中，大于分割阈值的梯度置为0，小于分割阈值的梯度置为1，梯度 图中所有为1的元素构成子野的初始形状，结合MLC的物理约束，对获得子野的初始形状进行 修正，以得到可交付使用的子野，加入到治疗计划中。
本发明与现有技术相比具有的有益效果是：本发明首先采用基于模糊熵分割原理生成子野 形状，然后进行子野权重优化，以得到满足临床要求的治疗方案。相比一般方法中采用局部梯 度信息或者遍历方式生成子野形状，本发明的优化方法可以快速、准确地利用全局梯度信息得 到子野形状，计算量减少的同时提高了计划质量，优化所得子野数目较少，提高了优化效率； 此外，本发明优化方法利用全局梯度信息生成子野形状，有利于提高靶区剂量适形度、降低靶 区周围危及器官以及正常组织接受的放射剂量，从而达到提高治疗增益比、提高肿瘤放射治疗 精确度以及降低正常组织由放射引起并发症的概率的临床治疗要求。
附图说明
下面结合附图对本发明做进一步的说明。
图1为本发明基于模糊熵的直接子野优化方法具体流程图。
图2为本发明中子野形状的三种类别示意图。
图3为本发明基于模糊熵分割生成子野形状方法示意图
图4为本发明中隶属函数示意图。
图5为本发明基于模糊熵的直接子野优化系统的组成框图。
具体实施方式
本发明一种基于模糊熵的直接子野优化方法，按下述步骤实施：
第一步：输入直接子野优化所需基本信息，并根据输入基本信息中的目标函数等信息计算 子野的梯度矩阵；
将目标函数设置如下：
则子目标函数表达式为：
fgEUD(D)＝H(gEUD(D)-gEUD0)·(gEUD(D)-gEUD0)；
fTCP(D)＝H(TCP0-ln(TCP(D)))·(TCP0-ln(TCP(D)))；
fNTCP(D)＝H(ln(1-NTCP0)-ln(1-NTCP(D)))·(ln(1-NTCP0)-ln(1-NTCP(D)))； 其中：
D(x)为剂量分布，且D(x)＝Wx，即剂量分布是关于剂量沉积矩阵W和注量矩阵x的线性 函数；
fl(D(x))为第l个子目标函数；
ξl为表示对应子目标函数重要性的权重系数；
L为子目标函数个数；
Dmin为靶区处方剂量；
D1，D2为两个剂量约束点处剂量值；
H(x)为阶跃函数；
N为器官内的所有体素点；
Di为器官第i个体素的剂量值；
Dmax为耐受剂量；
Dmean为靶区处方剂量或正常组织的平均耐受剂量；
其中：
α为剂量体积效应因子；
其中：
表示在第i个体素单元内的克隆源性细胞在受到Di剂量的 照射后的存活率；
λ为内在放射敏感性参数，表示不可修复的放射损伤；
γ为两次分次照射之间可以修复的损伤；
其中：
为标准正态累积分布函数；
D50为全部体积和部分体积V受照射时，由放射损伤引起的NTCP值为50％时所需的剂量；
n为体积效应因子；
m为控制NTCP剂量效应曲线斜率的参数；
第二步：根据子野梯度矩阵，计算该子野梯度矩阵的模糊熵；
首先将子野梯度矩阵转化为具有有限个灰度阶的梯度图Ga(m，n)，再利用子野梯度图的直 方图以及隶属函数μ(g)，计算子野梯度图的模糊熵，模糊熵的计算如下：
其中：
梯度图大小为M×N；
有限个灰度级范围为[Ll，Lu]；
μ(g)为定义在Lu-Ll级灰度上的隶属函数；
Sn(x)＝--x ln x-(1-x)ln(1-x)为Shannon函数；
f(g)为梯度图中灰度级为g的元素的个数；
由模糊熵分割原理可得，梯度图模糊熵E(G)的计算取决于不同隶属函数μ(g)，若将整个梯 度图看作一个模糊集，则一般选取标准S函数，隶属函数形式如下：
或
若按照梯度图中梯度元素的正负分为两个模糊集，则两个模糊集的隶属函数μ1(g)和μ2(g)，可 定义如下：
式中：
交叉点
[a，c]为模糊区域；
Δb＝c-a为模糊区域带宽；
模糊熵E(G)由带宽Δb和交叉点b决定；
第三步：利用模糊熵原理对所得梯度矩阵进行分割并修正得到可交付使用的子野，加入到 子野集合中；
通过计算子野梯度矩阵的模糊熵，确定使模糊熵取得极值的分割阈值，即分割阈值为使得 模糊熵取得最大值的梯度值；梯度元素大于该阈值，对应子射束关闭，小于该阈值，对应子射 束打开，最终形成一个初始子野形状；然后，根据MLC的物理约束，对初始子野形状进行修正， 最终得到可交付使用的子野形状；
第四步：对已有的全部子野的权重进行优化；输出优化结果。
子野权重优化采用的梯度类算法，以L-BFGS-B算法为例，其算法原理具体如下：
设定初始点x0和决定有限内存修正存储次数的整数m，定义初始有限内存矩阵并令k：＝0；
(1)如果投影梯度满足收敛测试||P(xk-gk，l，u)-xk||∞＜10-5则停止；
(2)利用CP算法计算广义柯西点；
(3)通过直接法计算搜索方向dk；
(4)在可行域内，沿着dk方向执行限定最大步长的线搜索来计算步长λk，令xk+1＝xk+ λkdk。线搜索始于单位步长，满足并且尝试满足 其中，α＝10-4，β＝0.9；
(5)计算
(6)如果yk满足eps＝2.2×10-16，将sk和yk添加到Sk和Yk，如果超过m 次更新被存储，则从Sk和Yk中删除最早存储的列；
(7)更新Lk和Rk，并且设置
(8)令k：＝k+1，转(1)；
其中：
满足充分下降条件；
是子空间上带约束的二次型最小化的近似解；
fk是函数值；
gk是梯度值；
{sk，yk}是算法的修正对；
sk＝xk+1-xk，yk＝gk+1-gk，Yk＝[yk-m，…，yk-1]，Sk＝[sk-m，…，sk-1]，P(xk- gk，l，u)-xk是投影梯度；
Lk和Rk是m×m的矩阵；
确定使模糊熵取得极值的分割阈值后，将子野梯度图中，大于分割阈值的梯度置为0，小 于分割阈值的梯度置为1，梯度图中所有为1的元素构成子野的初始形状，即梯度图中为0的 位置所对应的子射束关闭，梯度图中为1的位置所对应的子射束打开，最终形成一个初始子野 形状；对获得子野的初始形状结合MLC的物理约束进行修正，以得到可交付使用的子野，加入 到治疗计划中。
本发明一种基于模糊熵的直接子野优化系统，包括：
信息输入模块，用于输入优化系统所需基本信息；
信息预处理模块，根据输入的基本信息，获得子野的梯度矩阵；
子野形状生成模块，利用模糊熵原理对子野梯度矩阵进行分割并修正，以得到可交付使用 子野加入到子野集合中；
子野权重优化模块，对已有的全部子野的权重进行优化；
优化结果输出模块，输出优化结果信息。
所述子野形状生成模块，先将子野的梯度矩阵转换为子野的梯度图，利用子野梯度图的直 方图以及隶属函数，计算子野梯度图的模糊熵，设定分割阈值为使得模糊熵取得最大值的梯度 值，将子野梯度图中，大于分割阈值的梯度置为0，小于分割阈值的梯度置为1，梯度图中所有 为1的元素构成子野的初始形状，结合MLC的物理约束，对获得子野的初始形状进行修正，以 得到可交付使用的子野，加入到治疗计划中。
为使本发明实施例的目的、特点和优点能够更为明显易懂，下面结合附图对本发明的具体 实施作详细的说明。
本发明提供了用于调强放射治疗方案优化的基于模糊熵的直接子野优化方法，图1为该方 法的具体实施流程图。
一种基于模糊熵的直接子野优化方法，输入方案优化所需基本信息，所述方案优化所需基 本信息包括病人CT(computed tomography)数据、组织勾画信息、放射源照射信息、MLC的物 理约束、目标函数信息以及在目标函数中所使用的剂量约束参数；根据计划所需的基本信息， 计算病例的剂量沉积矩阵，根据剂量沉积矩阵以及目标函数信息计算子野的梯度矩阵；生成子 野形状，利用模糊熵原理对所得梯度矩阵进行分割、并根据MLC的物理约束对子野形状进行修 正，得到可交付使用的子野形状；优化子野权重，对已有的所有子野的权重进行优化；输出优 化结果，所述优化结果包括剂量分布、DVH(dose-volumehistogram)曲线、子野数目、子野形 状、子野权重、优化所需时间。
进一步地，所述病人CT数据通过计算机成像获得；组织勾画信息是通过在病人CT数据上 对正常组织及靶区进行勾画而获得；所述放射源照射信息包括放射源个数、放射源能量、放射 源位置以及等中心点位置；所述MLC的物理约束根据所允许交付使用的子野类型分为三类：无 约束(允许交付图2(a)(b)(c)所示子野)、不允许交错(允许交付图2(a)(b)所示子野)、子野形 状连通且不允许交错(允许交付图2(a)所示子野)；所述目标函数采用子目标函数加权和形式构 成，其中子目标函数包括最大剂量子目标函数、最小剂量子目标函数、均匀剂量子目标函数、 DVH子目标函数、等效均一剂量(generalized equivalentuniform dose，gEUD)子目标函数、 肿瘤控制率(tumor control probability，TCP)子目标函数及正常组织并发症概率(normal tissue complication probability，NTCP)子目标函数；所述目标函数所使用的剂量约束参数 包括物理约束和生物约束，物理约束主要为DV(Dose-Volume)约束，生物约束主要为gEUD、 TCP和NTCP。其中，所述病人CT数据优选螺旋CT数据；所述勾画方式通过物理师手动勾画或 者勾画软件自动勾画。
更进一步地，所述目标函数设置如下：
子目标函数表达式为：
fgEUD(D)＝H(gEUD(D)-gEUD0)·(gEUD(D)-gEUDo)；
fTCP(D)＝H(TCP0-ln(TCP(D)))·(TCP0-ln(TCP(D)))；
fNTCP(D)＝H(ln(1-NTCP0)-ln(1-NTCP(D)))·(ln(1-NTCP0)-ln(1-NTCP(D)))；
其中，剂量分布D(x)＝Wx，即剂量分布是关于剂量沉积矩阵W和注量矩阵x的线性函数， fl(D(x))为第l个子目标函数，ξl为表示对应子目标函数重要性的权重系数，L为子目标函数个数， Dmin为靶区处方剂量，H(x)为阶跃函数，N为器官内的所有体素点，Di为器官第i个体素的剂量 值，D1，D2为两个剂量约束点处剂量值；Dmax为耐受剂量，Dmean为靶区处方剂量或正常组织 的平均耐受剂量，α为剂量体积效应因子， 表示在第i个体素单元内的克隆源性 细胞在受到Di剂量的照射后的存活率，λ为内在放射敏感性参数，表示不可修复的放射损伤，γ 为两次分次照射之间可以修复的损伤， 为标准正态累积分布函数，D50为全部体积和部分体积V受 照射时，由放射损伤引起的NTCP值为50％时所需的剂量，n为体积效应因子，m为控制NTCP剂 量效应曲线斜率的参数。
进一步地，所述剂量沉积矩阵的计算，首先通过放射源个数、放射源位置及等中心点位置 信息确定射野方向的数目以及每个射野方向上子射束的数目，然后根据病人的CT数据、组织 勾画信息在剂量计算平台上采用不同计算方法计算每个射野方向上的剂量沉积矩阵，最后根据 剂量沉积矩阵和目标函数信息计算子野的梯度矩阵。
更进一步地，所述的剂量计算平台所采用的算法可以是笔形束算法、点核算法、蒙特卡洛 算法以及确定性方法。
进一步地，所述子野形状生成，如图3所示，首先将目标函数的梯度矩阵Gb(m，n)转化为 具有有限个灰度阶的梯度图Ga(m，n)，表达式如下：
其中，Gb(m，n)表示梯度矩阵中第m行、n列处梯度元素；ceil(x)函数是向上取整，它返回的是 大于或等于函数参数，并且与之最接近的整数；floor(x)函数是向下取整，它返回的是小于或等 于函数参数，并且与之最接近的整数；ln(x)函数是对数函数。
进一步地，所述采用隶属函数计算子野梯度矩阵的模糊熵曲线，是转化后的具有有限个灰 度阶的梯度图Ga(m，n)利用隶属函数μ(g)来计算模糊熵曲线，模糊熵的计算如下：
其中，梯度图大小为M×N，有限个灰度级范围为[Ll，Lu]，μ(g)为定义在Lu-Ll级灰度上 的隶属函数，Sn(x)＝--x ln x-(1-x)ln(1-x)为Shannon函数，f(g)为梯度图中灰度级为g 的元素的个数。通过计算子野梯度矩阵的模糊熵，确定使模糊熵取得极值的分割阈值。在确定 使模糊熵取得极值的分割阈值后，将子野梯度图中，大于分割阈值的梯度置为0，小于分割阈 值的梯度置为1，梯度图中所有为1的元素构成子野的初始形状，即梯度图中为0的位置所对 应的子射束关闭，梯度图中为1的位置所对应的子射束打开，最终形成一个初始子野形状。然 后，根据MLC的物理约束，对初始子野形状进行修正，最终得到可交付使用的子野形状。
更进一步地，利用模糊熵分割原理，梯度图模糊熵E(G)的计算取决于不同隶属函数μ(g)， 若将整个梯度图看作一个模糊集，则一般选取标准S函数，如图4左侧部分所示，隶属函数形 式如下：
或
若按照梯度图中梯度元素的正负分为两个模糊集，则两个模糊集的隶属函数μ1(g)和μ2(g) 如图4右侧部分所示，可定义如下：
式中，交叉点[a，c]为模糊区域，Δb＝c-a为模糊区域带宽。模糊熵E(G)由带宽 Δb和交叉点b决定。
进一步地，所述子野权重优化，是当每个射束方向上的第一个子野生成后，只赋予该子野 初始子野权重，暂不进行权重优化；在所有射束方向的第一个子野生成完成后，每再生成一个 新的子野，对已生成的所有子野进行权重优化，直至所得优化结果满足优化停止条件则优化结 束。所述优化停止条件包括：迭代次数达到上限、优化结果达到预期目标等。满足任一停止条 件，方案优化结束。
更进一步地，所述子野权重优化采用的梯度类算法，以L-BFGS-B算法为例，其算法原理具 体如下：
设定初始点x0和决定有限内存修正存储次数的整数m，定义初始有限内存矩阵并令k：＝0。
(1)如果投影梯度满足收敛测试||P(xk-gk，l，u)-xk||∞＜10-5则停止。
(2)利用CP算法计算广义柯西点。
(3)通过直接法计算搜索方向dk。
(4)在可行域内，沿着dk方向执行限定最大步长的线搜索来计算步长λk，令xk+1＝xk+ λkdk。线搜索始于单位步长，满足并且尝试满足 其中，α＝10-4，β＝0.9。
(5)计算
(6)如果yk满足eps＝2.2×10-16，将sk和yk添加到Sk和Yk，如果超过m 次更新被存储，则从Sk和Yk中删除最早存储的列。
(7)更新Lk和Rk，并且设置
(8)令k：＝k+1，转(1)。
其中，满足充分下降条件，是子空间上带约束的二次型最小化的近似解，fk是函数值，gk是梯度值，{sk，yk}是算法的修正对，sk＝xk+1-xk，yk＝gk+1-gk，Yk＝[yk-m，…，yk-1]，Sk＝[sk-m，…，sk-1]，P(xk-gk，l，u)-xk是投影梯度，Lk和Rk是m×m的矩 阵。
如图5所示，本发明还提供一种包含上述方法的基于模糊熵的直接子野优化系统，包含信 息输入模块，用以输入病人的CT数据、组织勾画信息、放射源照射信息、MLC的物理约束、目 标函数信息以及在目标函数中所使用的剂量约束参数；信息预处理模块，在使用当前放射源照 射的情况下，通过剂量计算平台，计算不同射野方向上的剂量沉积矩阵；子野形状生成模块， 利用模糊熵原理生成新的可交付使用的子野形状；子野权重优化模块，新子野生成成功，采用 梯度类算法对已有的全部子野的权重进行优化；优化结果输出模块，优化结束后，输出剂量分 布、DVH曲线、子野数目、子野形状、子野权重、优化所需时间等信息。其中，子野形状和子 野权重分别用来控制MLC叶片位置以及对应MLC叶片位置下的射线源照射时长。这两个模块所 产生的子野形状、子野权重以及其他射线源参数信息构成放射治疗计划被输入放射治疗设备中， 按照放射治疗计划参数，控制放射治疗设备，对病人进行放射治疗。
进一步地，所述信息输入模块，病人CT数据通过计算机成像获得；组织勾画信息是通过 在病人CT数据上对正常组织及靶区进行勾画而获得；所述放射源照射信息包括放射源个数、 放射源能量、放射源位置以及等中心点位置；所述MLC的物理约束分为三类：无约束、不允许 交错、子野形状连通且不允许交错；所述目标函数采用子目标函数加权和形式构成，其中子目 标函数包括最大剂量子目标函数、最小剂量子目标函数、均匀剂量子目标函数、DVH子目标函 数、gEUD子目标函数、TCP子目标函数及NTCP子目标函数；所述目标函数所使用的剂量约束 参数包括物理约束和生物约束，物理约束主要为DV约束，生物约束主要为gEUD、TCP和NTCP。
进一步地，所述信息预处理模块，首先通过放射源个数、放射源位置及等中心点位置信息 确定射野方向的数目以及每个射野方向上子射束的数目，然后根据病人的CT数据、组织勾画 信息在剂量计算平台上采用不同计算方法计算每个射野方向上的剂量沉积矩阵，最后根据剂量 沉积矩阵和目标函数信息计算子野的梯度矩阵。
进一步地，所述子野形状生成模块，利用模糊熵分割原理，通过隶属函数计算子野梯度矩 阵的模糊熵曲线，以确定使模糊熵取得极值的分割阈值，梯度元素大于该阈值，对应子射束关 闭，小于该阈值，对应子射束打开，最终形成一个初始子野形状。然后，根据MLC的物理约束， 对初始子野形状进行修正，最终得到可交付使用的子野形状。
进一步地，所述子野权重优化模块，是当每个射束方向上的第一个子野生成后，只赋予该 子野初始子野权重，不进行权重优化；在所有射束方向的第一个子野生成完成后，每生成一个 新的子野，对已生成的所有子野进行权重优化，判断当前优化结果是否满足停止条件，满足则 停止优化、进入优化结果输出模块，否则转至信息预处理模块。所述优化停止条件包括：迭代 次数达到上限、优化结果达到预期目标等。满足任一停止条件，方案优化结束。
上面结合附图对本发明的实施例作了详细说明，但是本发明并不限于上述实施例，在本领 域普通技术人员所具备的知识范围内，还可以在不脱离本发明宗旨的前提下作出各种变化。
Claims (4)
1.一种基于模糊熵的直接子野优化方法，其特征在于，按下述步骤实施：
第一步：输入直接子野优化所需基本信息，并根据输入基本信息中的目标函数信息计算子野的梯度矩阵；
将目标函数设置如下：
则子目标函数表达式为：
fgEUD(D)＝H(gEUD(D)-gEUD0)·(gEUD(D)-gEUD0)；
fTCP(D)＝H(TCP0-ln(TCP(D)))·(TCP0-ln(TCP(D)))；
fNTCP(D)＝H(ln(1-NTCP0)-ln(1-NTCP(D)))·(ln(1-NTCP0)-ln(1-NTCP(D)))；
其中：
D(x)为剂量分布，且D(x)＝Wx，即剂量分布是关于剂量沉积矩阵W和注量矩阵x的线性函数；
fl(D(x))为第l个子目标函数；
ξl为表示对应子目标函数重要性的权重系数；
L为子目标函数个数；
Dmin为靶区处方剂量；
H(x)为阶跃函数；
N为器官内的所有体素点；
Di为器官第i个体素的剂量值；
D1，D2为两个剂量约束点处剂量值；
Dmax为耐受剂量；
Dmean为靶区处方剂量或正常组织的平均耐受剂量；
其中：
α为剂量体积效应因子；
其中：
表示在第i个体素单元内的克隆源性细胞在受到Di剂量的照射后的存活率；
λ为内在放射敏感性参数，表示不可修复的放射损伤；
γ为两次分次照射之间可以修复的损伤；
其中：
为标准正态累积分布函数；
D50为全部体积和部分体积V受照射时，由放射损伤引起的NTCP值为50％时所需的剂量；
n为体积效应因子；
m为控制NTCP剂量效应曲线斜率的参数；
第二步：根据子野梯度矩阵，计算该子野梯度矩阵的模糊熵；
首先将子野梯度矩阵转化为具有有限个灰度阶的梯度图Ga(m,n)，再利用子野梯度图的直方图以及隶属函数μ(g)，计算子野梯度图的模糊熵，模糊熵的计算如下：
其中：
梯度图大小为M×N；
有限个灰度级范围为[Ll,Lu]；
μ(g)为定义在Lu-Ll级灰度上的隶属函数；
Sn(x)＝-xlnx-(1-x)ln(1-x)为Shannon函数；
f(g)为梯度图中灰度级为g的元素的个数；
由模糊熵分割原理可得，梯度图模糊熵E(G)的计算取决于不同隶属函数μ(g)，若将整个梯度图看作一个模糊集，则一般选取标准S函数，隶属函数形式如下：
或
若按照梯度图中梯度元素的正负分为两个模糊集，则两个模糊集的隶属函数μ1(g)和μ2(g)，可定义如下：
式中：
交叉点
[a,c]为模糊区域；
Δb＝c-a为模糊区域带宽；
模糊熵E(G)由带宽Δb和交叉点b决定；
第三步：利用模糊熵原理对所得梯度矩阵进行分割并修正得到可交付使用的子野，加入到子野集合中；
通过计算子野梯度矩阵的模糊熵，确定使模糊熵取得极值的分割阈值，即分割阈值为使得模糊熵取得最大值的梯度值；梯度元素大于该阈值，对应子射束关闭，小于该阈值，对应子射束打开，最终形成一个初始子野形状；然后，根据多叶准直器(multileafcollimator,MLC)的物理约束，对初始子野形状进行修正，最终得到可交付使用的子野形状；
第四步：对已有的全部子野的权重进行优化；输出优化结果。
子野权重优化采用的梯度类算法，以L-BFGS-B算法为例，其算法原理具体如下：
设定初始点x0和决定有限内存修正存储次数的整数m，定义初始有限内存矩阵并令k:＝0；
(1)如果投影梯度满足收敛测试‖P(xk-gk,l,u)-xk‖∞<10-5则停止；
(2)利用CP算法计算广义柯西点；
(3)通过直接法计算搜索方向dk；
(4)在可行域内，沿着dk方向执行限定最大步长的线搜索来计算步长λk，令xk+1＝xk+λkdk。线搜索始于单位步长，满足并且尝试满足 其中，α＝10-4，β＝0.9；
(5)计算
(6)如果yk满足eps＝2.2×10-16，将sk和yk添加到Sk和Yk，如果超过m次更新被存储，则从Sk和Yk中删除最早存储的列；
(7)更新Lk和Rk，并且设置
(8)令k:＝k+1，转(1)；
其中：
满足充分下降条件；
是子空间上带约束的二次型最小化的近似解；
fk是函数值；
gk是梯度值；
{sk,yk}是算法的修正对；
sk＝xk+1-xk，yk＝gk+1-gk，Yk＝[yk-m,…,yk-1]，Sk＝[sk-m,…,sk-1]，P(xk-gk,l,u)-xk是投影梯度；
Lk和Rk是m×m的矩阵；
2.根据权利要求1所述的一种基于模糊熵的直接子野优化方法，其特征在于，确定使模糊熵取得极值的分割阈值后，将子野梯度图中，大于分割阈值的梯度置为0，小于分割阈值的梯度置为1，梯度图中所有为1的元素构成子野的初始形状，即梯度图中为0的位置所对应的子射束关闭，梯度图中为1的位置所对应的子射束打开，最终形成一个初始子野形状；对获得子野的初始形状结合MLC的物理约束进行修正，以得到可交付使用的子野，加入到治疗计划中。
3.一种根据权利要求1所述的基于模糊熵的直接子野优化系统，其特征在于，包括：
信息输入模块，用于输入优化系统所需基本信息；
信息预处理模块，根据输入的基本信息，获得子野的梯度矩阵；
子野形状生成模块，利用模糊熵原理对子野梯度矩阵进行分割并修正，以得到可交付使用子野加入到子野集合中；
子野权重优化模块，对已有的全部子野的权重进行优化；
优化结果输出模块，输出优化结果信息。
4.根据权利要求3所述的一种基于模糊熵的直接子野优化系统，其特征在于，所述子野形状生成模块，先将子野的梯度矩阵转换为子野的梯度图，利用子野梯度图的直方图以及隶属函数，计算子野梯度图的模糊熵，设定分割阈值为使得模糊熵取得最大值的梯度值，将子野梯度图中，大于分割阈值的梯度置为0，小于分割阈值的梯度置为1，梯度图中所有为1的元素构成子野的初始形状，结合MLC的物理约束，对获得子野的初始形状进行修正，以得到可交付使用的子野，加入到治疗计划中。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910352358.3A CN110289075B (zh) | 2019-04-29 | 2019-04-29 | 一种基于模糊熵的直接子野优化方法及系统 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910352358.3A CN110289075B (zh) | 2019-04-29 | 2019-04-29 | 一种基于模糊熵的直接子野优化方法及系统 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110289075A true CN110289075A (zh) | 2019-09-27 |
| CN110289075B CN110289075B (zh) | 2022-04-22 |
Family
ID=68001897
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910352358.3A Active CN110289075B (zh) | 2019-04-29 | 2019-04-29 | 一种基于模糊熵的直接子野优化方法及系统 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110289075B (zh) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111986778A (zh) * | 2020-07-31 | 2020-11-24 | 上海联影医疗科技股份有限公司 | 一种调强计划优化系统、装置及存储介质 |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20160089549A1 (en) * | 2013-05-06 | 2016-03-31 | Koninklijke Philips N.V. | An interactive dose gradient based optimization technique to control imrt delivery complexity |
| CN105709341A (zh) * | 2016-01-15 | 2016-06-29 | 中国科学院合肥物质科学研究院 | 一种基于梯度法和漫水填充法的调强子野优化方法 |
| CN107823806A (zh) * | 2017-09-15 | 2018-03-23 | 中北大学 | 一种用于调强放射治疗直接子野优化的方法及系统 |
-
2019
- 2019-04-29 CN CN201910352358.3A patent/CN110289075B/zh active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20160089549A1 (en) * | 2013-05-06 | 2016-03-31 | Koninklijke Philips N.V. | An interactive dose gradient based optimization technique to control imrt delivery complexity |
| CN105709341A (zh) * | 2016-01-15 | 2016-06-29 | 中国科学院合肥物质科学研究院 | 一种基于梯度法和漫水填充法的调强子野优化方法 |
| CN107823806A (zh) * | 2017-09-15 | 2018-03-23 | 中北大学 | 一种用于调强放射治疗直接子野优化的方法及系统 |
Non-Patent Citations (3)
| Title |
|---|
| PENGCHENG ZHANG等: "The Aperture Shape Optimization Based on Fuzzy Enhancement", 《IEEE ACCESS》 * |
| 杨婕等: "一种基于梯度信息的直接子野优化算法", 《生物医学工程学杂志》 * |
| 王捷等: "一种快速调强放射治疗直接子野优化方法", 《中国医学物理学杂志》 * |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111986778A (zh) * | 2020-07-31 | 2020-11-24 | 上海联影医疗科技股份有限公司 | 一种调强计划优化系统、装置及存储介质 |
| CN111986778B (zh) * | 2020-07-31 | 2024-02-20 | 上海联影医疗科技股份有限公司 | 一种调强计划优化系统、装置及存储介质 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN110289075B (zh) | 2022-04-22 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| US10850122B2 (en) | Optimization methods for radiation therapy planning | |
| JP3775993B2 (ja) | 放射線治療計画を作成するシステム | |
| CN110124214B (zh) | 预测剂量分布引导的调强放疗计划优化系统、介质及设备 | |
| Das et al. | Beam orientation selection for intensity-modulated radiation therapy based on target equivalent uniform dose maximization | |
| US7734010B2 (en) | Method and apparatus for planning and delivering radiation treatment | |
| CN107823806B (zh) | 一种用于调强放射治疗直接子野优化的方法及系统 | |
| JP2006142043A (ja) | 放射線治療線量を自動的に計画する方法 | |
| CN110327554B (zh) | 基于预测剂量分布引导的调强放疗计划优化方法及应用 | |
| Bedford et al. | Constrained segment shapes in direct‐aperture optimization for step‐and‐shoot IMRT | |
| CN108711447A (zh) | 基于体素权重因子的自动调强多目标优化方法及其应用 | |
| Sir et al. | Stochastic programming for off-line adaptive radiotherapy | |
| CN112566692A (zh) | 用于确定弧形疗法的弧形剂量的系统和方法 | |
| CN109499012A (zh) | 一种优化剂量引导的tps自动迭代优化算法 | |
| Zhang et al. | A two-stage sequential linear programming approach to IMRT dose optimization | |
| CN110289075A (zh) | 一种基于模糊熵的直接子野优化方法及系统 | |
| Wu et al. | IMRT optimization based on the generalized equivalent uniform dose (EUD) | |
| Xing et al. | Inverse planning in the age of digital LINACs: station parameter optimized radiation therapy (SPORT) | |
| CN113870976A (zh) | 一种基于剂量预评估的自适应放疗剂量调强优化计算方法 | |
| Haas et al. | On improving physical selectivity in the treatment of cancer: a systems modelling and optimisation approach | |
| Andersson | Mathematical optimization of radiation therapy goal fulfillment | |
| KR101996268B1 (ko) | 방사선의 전달 정도 관리 방법 및 이를 이용한 방사선의 전달 정도 관리 장치 | |
| Ghanbarzadeh et al. | The scatter search based algorithm for beam angle optimization in intensity-modulated radiation therapy | |
| Bos et al. | Comparison between manual and automatic segment generation in step‐and‐shoot IMRT of prostate cancer | |
| Bogner et al. | Fast direct Monte Carlo optimization using the inverse kernel approach | |
| EP4019085B1 (en) | Radiation therapy planning |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |