光在生物组织内传输的过程中会发生吸收、散射、反射和透射，但是不同的生物组织结构、形态、成分差异很大，会表现出截然不同的光学特性，这种差异也是临床上实现光学诊断病变组织的依据。^[1]光学成像技术是光学方法在生物医学领域应用的重要方向，发展出了光学共焦显微成像、荧光成像、光声成像、光学相干层析成像(Optical Coherence Tomography, OCT)技术等等，与传统的成像技术相比，具有非侵入式、无损、精度高的优点。

OCT技术通过探测样品的后向散射光实现光学成像，可以通过提取OCT信号中的衰减情况实现对生物组织状态的检测。^[2-3]因此，研究生物组织中各参量的衰减情况可以为光学诊断判断生物组织状态提供参考和依据。

蒙特卡罗方法通过光子随机行走的方式模拟光在组织中的传输，已经成为组织光学领域验证其他理论或模型的一种非实验标准。^[4]

本文通过蒙特卡罗方法模拟无限窄光束辐照吸收系数、散射系数和衰减系数不同的生物组织，分析生物组织中吸收能量密度、光能流率、漫反射率、透射率、光通量的差异，并分析其产生原因。

蒙特卡罗模拟在依据光学特性参数建立生物组织模型的基础上，通过跟踪统计大量光子在介质中的光程，估计出相关物理量。

由计算机生成伪随机变量ξ，通过χ=f(ξ)将ξ∈(0, 1) 变换为χ∈(a, b)，再由$\int_a^\chi {p\left( \chi \right)} d\chi = \xi $生成为非均匀分布的概率密度函数p(χ)，其中χ∈(a，b)。光子每次移动前，通过生成χ产生模拟过程中所需步长s、偏转角θ和方位角φ等随机变量。

光子的步长定义为光子与组织发生吸收或散射等相互作用之前，在组织中行走的距离。因此由吸收和散射系数的定义可得：${s_1} = \frac{{-{\rm{ln}}\left( \xi \right)}}{{{\mu _t}}}$，其中：衰减系数μ_t=μ_a+μ_s，μ_a为吸收系数，μ_s为散射系数。

光在组织内部散射的过程是轴对称的，采用极坐标中偏转角θ和方位角φ确定具体方向，如图 1所示。偏转角θ∈[0, π]由方向余弦的几率密度函数p(cosθ)=$\frac{{1-{g^2}}}{{2{{(1 + {g^2}-2g{\rm{cos}}\theta )}^{3/2}}}}$确定，其抽样值: $\theta = {\rm{arccos}}(\frac{1}{{2g}}\left[{1 + {g^2}-{{\left( {\frac{{1-{g^2}}}{{1-g + 2g\xi }}} \right)}^2}} \right]){\rm{ }}\left( {g \ne 0} \right)$，g为各项异性因子。方位角φ∈[0, 2π]的几率密度函数为p(ψ)=$\frac{1}{{2\pi }}$，其抽样值ψ=2πξ。

图 1

图 1 移动球坐标系示意图

光子垂直组织表面入射后，在组织内部每移动一个步长s，其权重、坐标和方向余弦的变化如表 1所示。^[5]光子初始权重ω=1，当传输至ω < ω_th时，由“轮盘赌”判断是否继续传输。

表 1

表 1 光子移动一个步长，权重、坐标、方向余弦的变化情况 光子位置 权重 坐标 方向余弦 当前 ω (x，y，z) (μx, μy, μz) 步长s ω′=ω-ωμa/(μa+μs) x′=x+μxs μ′x=sinθ(μxμzcosψ-μysinψ)/(1-μz2)$\frac{1}{2}$+μxcosθ y′=y+μys μ′y=sinθ(μyμzcosψ-μxsinψ)/(1-μz2)$\frac{1}{2}$+μycosθ z′=z+μzs μ′z=-sinθcosψ(1-μz2)$\frac{1}{2}$+μzcosθ

表 1 光子移动一个步长，权重、坐标、方向余弦的变化情况

蒙特卡罗方法模拟组织内光分布的流程图如图 2所示。模拟初始值采用633 nm下皮肤组织的光学特性参数^[6]：μ_a=2.7 cm^-1，μ_s=187 cm^-1，g=0.81。另外, 模拟光子数为10⁵个，d=0.1 cm，dz=0.001 cm，dr=0.005 cm, ω_th=0.000 01。如表 2所示，分别改变μ_a、μ_s、μ_t，获得不同参数下无限窄准直光束入射皮肤组织，吸收能量密度、光能流率、漫反射率、透射率、光通量的分布(图 3至图 7)，分析光学特性参数对组织中光子传输的影响。图 3至图 6中的曲线1至曲线5分别对应表 2中相应序号的光学特性参数。序号2是将吸收系数改变为初始值的2倍，散射系数不变；序号3是将散射系数改变为初始值的1/2，吸收系数不变；序号4是保持衰减系数不变，将吸收系数改变为初始值的3倍；序号5是保持衰减系数不变，将吸收系数改变为初始值的1/3。上述参数的选择是在兼顾皮肤组织正常的光学特性参数范围和模拟结果差异较明显的前提下选择的。

表 2

表 2 模拟所采用的光学特性参数 序号 吸收系数μa(cm-1) 散射系数μs(cm-1) 衰减系数μt(cm-1) 1 2.7 187 189.7 2 5.4 187 192.4 3 2.7 93.5 96.2 4 8.1 181.6 189.7 5 0.9 188.8 189.7

表 2 模拟所采用的光学特性参数

图 2

图 2 模拟流程图

图 3为不同光学特性参数生物组织中吸收能量密度的分布情况。散射系数相同，吸收系数越大(μ_a2>μ_a1)，吸收能量密度的衰减梯度越大，在浅表层(r < 0.07 cm)，吸收能量密度越大。两曲线在r=0.07 cm附近交叉，在深层组织(r>0.07 cm)，吸收能量密度的分布与浅表层正好相反，这是因为浅表层吸收能量密度较大，更多的光子被吸收，因此组织深层吸收能量密度较小。

图 3

图 3 光学特性参数对吸收能量密度A(z)的影响

吸收系数相同，散射系数越大(μ_s1>μ_s3), 吸收能量密度的衰减梯度越大，光子的散射能力越强，更多的光子被吸收，因此曲线1在浅表层吸收能量密度更大，二者在r=0.04 cm附近交叉。

衰减系数相同，吸收系数越小(μ_a5 < μ_a1 < μ_a4)，吸收能量密度衰减梯度越小，在浅表层吸收能量密度越小。

图 4为不同光学特性参数生物组织中光能流率的分布。衰减系数越小(μ_t3 < μ_t1、μ_t4、μ_t5 < μ_t2)，光能流率的衰减梯度越小。衰减系数相同(μ_t1=μ_t4=μ_t5)，光能流率的衰减梯度相同，吸收系数越大(μ_a4>μ_a1>μ_a5), 光能流率越小。

图 4

图 4 光学特性参数对光能流率的影响

散射系数相同，吸收系数越小(μ_a1 < μ_a2)，光能流率越大。吸收系数相同，散射系数不同(μ_s3 < μ_s1), 两曲线在r=0.04 cm附近交叉，这主要是由于光能流率在数值上等于A(z)与μ_a的比值^[4]，二者吸收系数相同，因此光能流率的分布情况与图 3中相应曲线类似。

图 5为不同光学特性参数的生物组织中漫反射率的分布，漫反射率在径向r迅速衰减。散射系数相同，吸收系数越大(μ_a2>μ_a1)，径向漫反射率越小；衰减系数相同，吸收系数越大(μ_a4>μ_a1>μ_a5)，径向漫反射率越小。这主要是由于更多的光子被吸收，导致由于漫反射从上表面离开组织的光子变少。吸收系数相同，散射系数越小(μ_s3 < μ_s1), 发生散射的光子数越少，径向漫反射率越小，其衰减梯度也明显较小。

图 5

图 5 光学特性参数对径向漫反射率的影响

图 6为不同光学特性参数的生物组织中透射率的分布。散射系数相同，吸收系数越小(μ_a1 < μ_a2)，透射率越大；衰减系数相同，吸收系数越小(μ_a5 < μ_a1 < μ_a4)，透射率越大，这主要是由于更多光子未被组织吸收，可以通过组织内部的传输，从下表面离开组织。吸收系数相同，散射系数不同(μ_s3 < μ_s1), 二者曲线差异较大。曲线3由于散射系数较小，更多的光子在中心轴附近通过组织内部，从下表面离开组织，因此透射率近轴明显高于曲线1，且从中心轴开始迅速衰减；而曲线1散射系数较大，因此透射率在近轴并未迅速衰减，甚至出现先增大后减小的现象。

图 6

图 6 光学特性参数对透射率的影响

图 7为不同光学特性参数的生物组织中光通量的二维分布图，图中各序号对应表 2中各光学特性参数。散射系数相同，吸收系数越小(μ_a1 < μ_a2)，同一位置光通量越大，但是光子在径向r上分布得更远，这主要是由于吸收系数小的组织中更多的光子被散射。吸收系数相同，散射系数越小(μ_s3 < μ_s1), 光通量的衰减梯度较小，光子在径向上分布得更远。衰减系数相同，吸收系数越小(μ_a5 < μ_a1 < μ_a4)，散射系数越大，同一位置光通量越大，光通量的衰减梯度越大，光子在径向上分布得更远。

图 7

图 7 光学特性参数对光通量的影响

同一光源辐照单层生物组织，不同的光学特性参数对组织中光分布的各参量影响各不相同。

(1) 增大吸收系数和散射系数均会增大吸收能量密度及其梯度。

(2) 光能流率和光通量分布类似，吸收系数影响光能流率和光通量的大小，散射系数影响二者的衰减梯度。

(3) 漫反射率与透射率分布类似，吸收系数和散射系数均对漫反射率和透射率大小有影响，散射系数影响二者衰减梯度。