小野下光子特征线混合笔束模型肺部剂量算法

崔凤洁，王宁宇，顾少娴，杨 鹏，尹楚欧，张盛元，胡金有,2，蔡芸竹，吴章文，汪 俊，勾成俊

(1.四川大学原子核科学技术研究所 辐射物理及技术教育部重点实验室，成都 610064;2.四川省医学科学院·四川省人民医院，成都 610072)

随着现代影像诊断技术和仪器设备的发展，适形调强、图像引导和立体定向等精准放疗技术已广泛应用于肿瘤临床治疗.这些新技术采用小射野光子束照射靶区并调整靶区的剂量分布，从而提高肿瘤控制率，降低正常组织并发症的概率[1,2].射野半径小于次级电子最大射程的照射野为小射野，通常认为小射野的尺寸小于3 cm×3 cm.小射野放疗本就存在横向电子不平衡，当射线入射到低密度肺组织中时，次级电子的输运范围增大，更易引起横向电子不平衡，使剂量计算更复杂[3-5].精准的剂量算法是精准放疗的关键，ICRU 24号报告[6]建议：剂量算法的误差应控制在3%以内.

现代放疗计划系统采用的剂量算法主要有解析法和蒙特卡罗(MC)模拟两类.解析法简化了粒子与人体组织的相互作用过程，计算速度快，但在小射野、非均匀组织等复杂条件下，其计算精度会降低[7,8].例如，卷积叠加法(CS)、笔形束卷积法(PBC)、筒串卷积叠加法(CCC)、快速傅里叶卷积法(FFT)和各向异性解析法(AAA)等非均匀修正算法得到的肺组织剂量明显高于MC模拟值[9,10].Hasani等[9]对CCC、FFT和CS算法的精度进行了评估.在10 cm×10 cm射野下，CCC、FFT、CS算法与MC模拟得到肺中的平均剂量误差分别为4.81%、-1.55%和4.96%，最大误差分别为7.3%、5.26%和7.4%.Cristiano等[10]将MC模拟的深度剂量数据与PBC、AAA的计算结果进行比较.当射野尺寸为1 cm×1 cm、肺密度为0.3 g/cm3时，PBC、AAA与MC计算的肺部平均剂量差异分别为21.9%和5.8%，最大差异分别为24.3%和11.5%.上述算法均不能精确模拟小射野下肺模体中散射光子和次级电子的输运行为，不适用于小射野肺部放疗.而MC方法将概率统计原理与计算机模拟实验相结合，可准确计算出小射野下肺模体中的剂量分布，但计算时间太长，不能完全用于临床剂量计算[11].因此，有必要开展小射野下非均匀修正算法的研究.

本文提出了光子特征线混合笔束模型(Photon characteristic line hybrid pencil beam model, PCL-HPBM).该算法首先利用特征线方法[12,13]计算无限宽光子束在肺模体中的深度剂量分布；
然后对水模体笔束核进行非均匀修正，计算笔束光子在肺模体中的剂量分布；
再利用有限笔束算法[14-16]计算任意照射野条件下肺模体中的三维剂量分布；
最后与MC模拟结果进行比较，验证该算法的精度.

2.1 蒙特卡罗模拟

本文使用程序为MC程序EGSnrc/Dosxyznrc.该程序具有完善的带电粒子输运机制和低能作用截面，能准确模拟在keV至GeV范围内的光子和电子在模体中的输运行为，从而得到不同射野条件下不同肺模体(肺组织深度、密度不同)中的剂量分布[17].模拟参数[18,19]如下：电子截止能量为0.521 MeV，光子截止能量为0.001 MeV，电子输运步长算法为PRESTA-Ⅱ，边界处理算法为EXACT，粒子源的类型选择为0号(以平行的矩形光束垂直入射)，光子能谱数据文件为mohan6.spectrum.在MC模拟中，随机抽样的光子数越多，精确度越高.为确保整个模拟过程的统计误差均小于1%，每次模拟的历史粒子数应足够大.

本文模拟了四个肺模体，大小为40 cm×40 cm×40 cm，体元大小为0.25 cm×0.25 cm×0.25 cm.水、肺组织元素的组成及比例来源于EGSnrc的材料库文件521cru.pegs4dat.其中，肺模体A、B和C均由水(5 cm)、肺(15 cm)、水(20 cm)三层介质组成，肺密度分别为0.1、0.26 和0.4 g/cm3；
肺模体D由水(15 cm)、肺(15 cm)和水(10 cm)三层介质组成，肺密度为0.26 g/cm3.笔束大小为0.25 cm×0.25 cm.射野尺寸为0.5 cm×0.5 cm、1 cm×1 cm、2 cm×2 cm和不规则射野(如图1所示).

图1 不规则射野的尺寸和形状

2.2 光子特征线混合笔束模型

放疗能区内的光子与人体组织的相互作用主要有光电效应、康普顿效应和电子对效应三种形式.光电效应和电子对效应引起光子的吸收，而康普顿效应使光子偏离初始的入射方向.其相互作用机理与光子能量有关，对于临床上广泛使用的6 MV光子束，康普顿效应占优势[20].这三种相互作用的总微观截面表示为

σtot=σpe+σcomp+σpair

(1)

其中σpe、σcomp和σpair分别为光电效应、康普顿效应和电子对效应的微观截面.当能量为E0的无限宽光子束垂直入射到介质表面时，光子在均匀介质中输运的Boltzmann方程表示为[21]

Ω′→Ω)+S(z,E,μ)，S(z,E,μ)=

δ(z)δ(E-E0)δ(1-μ)/2π

(2)

其中Ø(z,E,μ)表示光子的微分注量分布；
Ø(z,E,μ)dEdμ表示在深度z处，能量为[E,E+dE]、方向余弦为[μ,μ+dμ]的光子的注量；
K(z,E′→E,Ω′→Ω)表示深度z处，能量为E′、方向为Ω′的光子经过康普顿散射后，散射光子能量为[E,E+dE]、方向为[Ω,Ω+dΩ]的几率；
S(z,E,μ)表示源项；
∑(z,E)表示深度z处，能量为E的光子的线衰减系数，其表达式为

∑(z,E)=n(z)·σtot(z,E)

(3)

其中n(z)为介质的原子数密度，表达式为

(4)

其中NA为阿伏伽德罗常数；
A(z)为深度z处的物质的原子质量；
ρ(z)为深度z处的介质密度.

在光子的特征线算法中，光子注量Ø表示为

Ø(z,E,μ)=Ø0(z,E,μ)+Ø1(z,E,μ)+

Øm(z,E,μ)

(5)

其中Ø0(z,E,μ)、Ø1(z,E,μ)、Øm(z,E,μ)分别为初级光子注量、一次散射光子注量和多次散射光子注量.初级光子的输运方程表示为[12,13]

(6)

一次散射光子的输运方程表示为[12,13]

S1(z,E,μ)

(7)

多次散射光子的输运方程表示为[12，13]

(8)

(7)式、(8)式中λ表达式为

(9)

Ø0(z,E,μ)、Ø1(z,E,μ)有明确的解析解，Øm(z,E,μ)只能由数值方法得到数值解.在计算Øm(z,E,μ)的数值解时，可以将方向余弦μ按照特征线方向进行离散化.则(7)式解的表达式为[12,13]

(10)

(8)式解的表达式为[12,13]

Ø1(z,E,μ)=

(11)

(9)式解的表达式为[12,13]

Øm(z,E,μ)=

(12)

其中d为介质的厚度.在得到光子的注量分布后，可以得到次级电子源和次级正电子源的注量分布，表达式为

(13)

其中dz为计算单元在深度方向的长度.设电子的微分注量为N(z,E,μ)，则电子的输运方程可表示为

σe(z,E′,E′-E,Ω′·Ω)

(14)

其中∑e(z)为深度z处的电子在单位路径上的作用次数；
σe为电子总的作用截面，包括轫致辐射、弹性相互作用、非弹性相互作用.利用特征线方法和电子输运双亲模型，可以获得电子的微分注量分布.根据微分注量分布就可以获得能量为E0的无限宽光子束的深度剂量分布，表达式为

E,T)TdT

(15)

D(E0,z)=De(E0,z)+Dp(E0,z)

(16)

模体中的剂量分布可由入射光子束的能量注量分布和笔束核卷积得到.笔束核是指无限窄的射束垂直入射到半无限大平面介质中沉积能量分布.该算法的关键是精准表述笔束核.对于体积较大的肺组织，只要照射条件(如射线能量、源皮距等)和体元大小不变，就可以认为笔束核具有空间不变性，可重复使用[20].为了利用有限笔束模型计算任意照射野条件下肺部的剂量分布，需要得到笔束光子在肺部的三维剂量分布.小射野下的肺组织可近似为层状非均匀介质，入射能量为E的笔束光子在肺部的三维剂量分布为

Dp(E,x,y,z)=Dp(E,0,0,z)·f(E,x,y,z)

(17)

其中Dp(E,x,y,z)表示笔束核，是笔束光子沿z轴正方向入射到原点时(x,y,z)点的剂量，为方便计算，将其进行最大值归一化处理；
Dp(E,0,0,z)表示中心轴深度剂量分布；
f(E,x,y,z)表示散射因子.本文采用Dosxyznrc程序模拟笔束光子在均匀水模体中的三维剂量分布，则f(E,x,y,z)的表达式为

f(E,x,y,z)=Dp(E,x,y,z)/Dp(E,0,0,z)

(18)

本研究考虑了散射因子在深度方向和离轴方向上的非均匀修正.首先，对散射因子在深度方向上进行等效深度修正，即肺模体中的散射因子近似用水中等效深度处的散射因子表示，公式如下.

f′(E,x,y,z)=f(E,x,y,zeff)，

(19)

其中ρ(h)、ρw分别表示h深度处的密度和水的密度；
zeff为肺模体z深度处的等效深度.其次，在等效深度修正的基础上，对散射因子在离轴方向上进行加权密度修正，公式如下.

f″(E,x,y,z)=f′(E,xeff,y,z)，

xeff=β·ρ(z)·x，

β=α+(1-α)·ρ(z)

(20)

其中α为散射权重系数；
β为密度修正系数.当ρ(z)=1 g/cm3时，β=1；
当ρ(z)=0 g/cm3时，β=α.中间的密度修正系数可由线性插值得到.xeff为肺模体中离轴距离x处的等效距离.经上述修正后，根据(1)式，得到肺模体中的笔束核Dp′(E,x,y,z).其表达式为

Dp′(E,x,y,z)=Dp(E,0,0,z)·f″(E,x,y,z)

(21)

无限宽束的中心轴剂量分布可叠加每个笔束光子的剂量贡献得到.函数表达式为[14]

(22)

实际上是对x、y在有限的范围内积分，积分区间与光子和次级电子在非均匀介质中的散射程度有关.本文采用20 cm×20 cm的宽束光子作为参考射野，通过光子特征线方法得到参考射野下肺模体的中心轴深度剂量分布，表示为Dref(E,0,0,z).同时根据(22)式，得到参考射野下的中心轴深度剂量分布.其表达式为

(23)

其中a、b为10 cm.再利用Dref(E,0,0,z)、Dref′(E,0,0,z)重新定义笔束光子在肺模体中的中心轴深度剂量，表达式为

Dp″(E,0,0,z)=Dref(E,0,0,z)·

Dp(E,0,0,z)/Dref′(E,0,0,z)

(24)

根据(17)式，得到非均匀修正后肺组织中精确的笔束核.其表达式为

Dp″(E,x,y,z)=Dp″(E,0,0,z)·

f″(E,x,y,z)

(25)

因本研究采用平行光子束垂直照射肺模体，故不考虑锥形束校正.则根据有限笔束算法[14-16]可以计算出任意射野下人体组织中的剂量分布.其表达式为

Dp″(E,x-x′,y-y′,z)dx′dy′

(26)

其中I(x′,y′)是在入射点(x,y,0)处入射光子束的能量注量，表示射野形状;Dp″(E,x-x′,y-y′,z)是入射到点(x′,y′,0)的光子束对计算点(x,y,z)的剂量贡献.(26)式可离散化表示为

D(E,i,j,k)=∑i′∑j′I(i′,j′)·

Dp″(E,i-i′,j-j′,k)

(27)

其中I(i′,j′)是射野形状的离散化强度矩阵；
Dp″(E,i-i′,j-j′,k)是入射到点(i′,j′,0)的光子束对计算点(i,j,k)的剂量贡献.临床中使用的X射线通常具有连续的能谱分布，故本研究采用6 MV光子能谱，利用Dosxyznrc程序获得了具有能谱信息的笔束核.最后，我们将该算法得到的小射野下肺模体的中心轴深度剂量分布与MC模拟结果进行了比较，计算点剂量相对误差以验证算法的精度.剂量相对误差(Relative error)的表达式为

(28)

其中DMC表示MC模拟计算点的剂量；
DPCL-HPBM表示由PCL-HPBM算法得到同一计算点的剂量；
DMC,max表示肺模体中心轴上的最大剂量.

图2为在0.5 cm×0.5 cm、1 cm×1 cm、2 cm×2 cm和不规则射野下，PCL-HPBM算法和MC方法得到6 MV光子束穿过肺密度为0.1 g/cm3的肺模体A时的中心轴深度剂量的对比.横坐标为肺模体深度(Depth).纵坐标分别为以剂量最大值为单位进行归一化的百分深度剂量(Percentage Depth Dose，PDD)和两组剂量数据的相对误差.在不同射野下，两种方法得到的PDD基本符合，但在介质交界处存在明显的剂量差异.在重点关注的肺部，PCL-HPBM算法得到的深度剂量在水-肺交界处(图2中深度为5 cm处)以及肺组织前部低于MC模拟结果, 而在肺-水交界处(图2中深度为20 cm处)高于MC模拟结果.但在肺组织中部及后部，两种方法得到的剂量符合得很好.在0.5 cm×0.5 cm、1 cm×1 cm、2 cm×2 cm和不规则射野下，最大误差出现在水-肺交界处，最大相对误差分别为36.9%、35.5%、32.2%和33.2%；
在肺组织前部，相对误差大于3%剂量点的最大深度分别为6.0、6.75、7.5 和7.25 cm.

图2 不同射野下肺模体A中的PDD曲线和相对误差曲线.射野尺寸：(a)0.5 cm×0.5 cm;(b)1 cm×1 cm;(c)2 cm×2 cm;(d)不规则射野

图3为在0.5 cm×0.5 cm、1 cm×1 cm、2 cm×2 cm和不规则射野下，PCL-HPBM算法和MC方法得到肺密度为0.26 g/cm3的肺模体B中的PDD及相对误差.两种方法得到的PDD基本一致.在肺组织中部和后部，PCL-HPBM算法计算值和MC模拟值的平均相对误差小于1%，但在肺前部及介质交界处的相对误差大于3%.在0.5 cm×0.5 cm、1 cm×1 cm、2 cm×2 cm和不规则射野下，出现在水-肺交界处的最大相对误差分别为22.7%、21.1%、14.7%和16.7%；
在肺组织前部，相对误差大于3%的剂量点的最大深度分别为5.75、6.0、6.75和6.25 cm.

图3 不同射野下肺模体B中的PDD曲线和相对误差曲线.射野尺寸：(a)0.5 cm×0.5 cm;(b)1 cm×1 cm;(c)2 cm×2 cm;(d)不规则射野

图4为在0.5 cm×0.5 cm、1 cm×1 cm、2 cm×2 cm和不规则射野下，由PCL-HPBM算法和MC方法得到肺密度为0.4 g/cm3的肺模体C中的PDD及相对误差.两种方法得到的PDD基本一致，在0.5 cm×0.5 cm、1 cm×1 cm、2 cm×2 cm和不规则射野下，最大误差出现在第一水肺交界处，最大相对误差分别为15.4%、12.9%、7.8%和8.7%，肺前部相对误差大于3%剂量点的最大深度分别为5.5、5.75、6.0和5.75 cm.

图4 不同射野下肺模体C中的PDD曲线和相对误差曲线.射野尺寸：(a)0.5 cm×0.5 cm;(b)1 cm×1 cm;(c)2 cm×2 cm;(d)不规则射野

图5为在0.5 cm×0.5 cm、1 cm×1 cm、2 cm×2 cm和不规则射野下，PCL-HPBM算法和MC方法得到肺密度为0.26 g/cm3的肺模体D中的PDD及相对误差.PCL-HPBM算法得到的深度剂量在水-肺交界处(图5中深度为15 cm处)及肺组织前部小于MC模拟结果，而在肺-水交界处(图5中深度为30 cm处)大于MC模拟结果.在0.5 cm×0.5 cm、1 cm×1 cm、2 cm×2 cm和不规则射野下，出现在水-肺交界处的最大相对误差分别为14.1%、13.3%、8.6%和10.8%；
在肺组织前部，相对误差大于3%剂量点的最大深度分别为15.5、15.75、15.75和15.75 cm.而在其他深度处，剂量点的相对误差小于3%.

图5 不同射野下肺模体D中的PDD曲线和相对误差曲线.射野尺寸：(a)0.5 cm×0.5 cm;(b)1 cm×1 cm;(c)2 cm×2 cm;(d)不规则射野

从图2～图5可以看出，在不同条件下，PCL-HPBM算法得到的PDD与MC模拟结果在介质交界处和肺前部存在较大差异，最大相对误差在7.8%～36.9%之间.在水-肺交界处，最大相对误差随着射野尺寸和肺密度的增加而减小.同时，PCL-HPBM算法得到的深度剂量迅速跌落，而MC模拟结果变化得相对缓慢.这是因为介质密度的减小增大了次级电子的不平衡性，重建电子平衡需要一定的剂量沉积深度，使介质交界处的深度剂量变化缓慢，而采用PCL-HPBM算法进行非均匀修正时忽略了这一过程.并且随着射野尺寸和肺密度的增大，由横向电子不平衡导致的剂量跌落幅度减小.在肺组织前部，相对误差大于3%的剂量点的深度范围随着射野尺寸增大和肺密度的减小而增大，在距水-肺交界面2.5 cm处的剂量点的相对误差均能满足小于3%的要求.而在肺组织中部与后部，两种方法得到的深度剂量符合得较好，平均相对误差均小于1%，并且随着深度的增大，深度剂量逐渐减小.

表1列出了在不同射野下不同肺模体中，PCL-HPBM算法的散射权重系数和密度校正系数.从表1可以看出，这些系数随着射野尺寸及肺密度的增加而减小.

表1 在不同射野不同肺模体条件下，PCL-HPBM算法的散射权重系数和密度校正系数

本文通过构建光子特征线混合笔束模型来进行小射野光子束剂量计算.在对散射因子进行非均匀修正时，我们考虑了电子不平衡性对能量沉积的影响，并重新生成了笔束卷积核.在分析不同小射野照射不同肺模体(肺组织密度、深度不同)时，该算法得到的中心轴深度剂量分布与MC模拟结果基本相符，大部分剂量点的相对误差小于3%；
在重点关注的肺组织中部和后部，平均剂量相对误差小于1%，从而验证了该算法的准确性.该算法同样适用于不规则射野的剂量计算.当忽略在介质交界处和肺组织前部相对误差较大的几个深度时，射野尺寸越小，该算法的计算精度越高.

为了提高肿瘤与正常组织的剂量比值，现代放疗技术越来越多地使用具有高度可调性的小射野.这推动了小射野精准剂量算法的发展.然而，前面提到的PBC、AAA和CCC等现有的非均匀修正算法均高估肺组织中的剂量，平均误差范围约为6%～22%，算法精度低.通过上述分析，在计算小射野下非均匀模体的深度剂量分布时，PCL-HPBM算法的精度均高于这些修正算法.该算法在小射野肺部剂量计算中有潜在的研究价值，但在介质交界处和肺组织前部，误差较大.因此我们需要进一步提高该算法的性能，来获得更精确的剂量分布，以便应用于更复杂的条件.