1 引言

由于生物组织的散射特性与组织的形态、亚显微结构有关,因此,近年来基于亚扩散散射的光谱测量和生物成像技术受到了很多关注^[1-8]。在单光纤反射探针技术中^[1],光源光信号的传输和散射光的收集由同一根光纤承担,因此能够用于临床中的微操作。但是,在这种情况下光纤只能收集到光源附近微区的散射光,这些散射光因在组织中经历了较少次数的散射而携带了更多的关于组织微观形态的信息,因此需要考虑与散射相函数相关的高阶光学参量。

在光传输理论中,用散射系数μ_s和散射相函数(SPF)P( s˙, s˙')(表示粒子将 s˙'方向的入射波散射到 s˙方向的概率)描述组织的散射特性。在扩散近似下,组织的散射特性用约化散射系数μ'_s表征,μ'_s=μ_s(1-g),g=< s˙, s˙'>是SPF的一阶Legendre矩。μ'_s是一个与组织的微观形态相关的光学参量,已被尝试用来标记恶性肿瘤^[6]。但是在亚扩散情况下,散射信号与SPF的更高阶Legendre矩有关。

Bevilacqu等^[9]利用Monte Carlo(MC)模拟研究了一个平均自由程内光源附近的漫反射,结果表明,反射率可以近似地由吸收系数μ_a、μ'_s和一个与散射相函数SPF相关的参量γ来表示。γ=(1-g₂)/(1-g₁),其中g₁=g,g₁和g₂分别是散射相函数SPF的前两阶矩。从玻尔兹曼方程得到亚扩散情况下的解析解是非常困难的,并且解的表达式是非常复杂的,例如P₃近似的表达式中不仅仅包含μ_a和μ'_s,还包含高阶参量γ和δ,其中δ=(1-g₃)/(1-g₂),其繁琐的表达形式限制了它的实际应用^[10]。

近些年来,研究人员提出了一些相对简单的模型,这些模型都是基于半经验公式和MC模拟结果得到的^[11-14]。根据朗伯比尔定律,Amelink等^[11]给出一个包含高阶参量γ的半经验公式,并利用这个公式计算出了一些生物组织的γ值。本课题组根据光子迁移理论,提出了一个包含高阶参量γ的半经验模型^[12]。尽管这些基于MC模拟结果得到的模型既具有解析形式,又具有MC模拟的优点,但是这些模型都是依据某些特定的SPF得到的。在实验研究中发现,不同的生物组织具有不同形式的SPF^[15-17]。本课题组之前用MC模拟研究了SPF对空间分辨漫反射的影响,并且证明了当将MC模拟结果用于反演算法时,需要考虑与SPF相关的高阶参量 ^[18]。

在之前的研究中,本课题组提出了一个适用于小孔径测量漫反射的半经验解析模型,确定了反射率与孔径、约化散射系数的函数关系,并对其适用范围进行了分析^[13]。而本文首先阐述了不同的SPF可以通过参量γ和δ的组合进行表征,之后研究了与相函数相关的两个高阶光学参量γ和δ对该模型的共同作用,并给出了简化的半经验模型及其适用条件。

2 理论和方法

图1是单光纤漫反射测量示意图,直径为d的光纤同时用作光源和探测器。假设无限细光束通过光纤垂直入射到界面,经历生物组织的随机吸收和散射后,有一部分背向散射光被光纤收集。生物组织的光学性质用吸收系数μ_a、散射系数μ_s和散射相函数P( s˙, s˙')表示。

图 1. 小孔径漫反射测量示意图

Fig. 1. Diagram of diffuse reflectance measurement by a small detector aperture

下载图片 查看所有图片

对于小的探测器孔径d,漫反射率R_d可以表示为^[13]

Rd=k1μ's1+k2μa,(1)

式中参数k₁和k₂均与生物组织的散射相函数P( s˙, s˙')有关,也与生物组织的折射率n有关。对于大部分的组织,折射率为1.335~1.62,如角质层的折射率为1.55,晶状体表面的折射率为1.386^[19],本文取折射率n=1.4,没有考虑折射率的影响。为了研究k₁、k₂与P( s˙, s˙')的关系,考虑组合相函数

ptissue(x,gHG,α)=αpHG(x,gHG)+38(1-α)(1+cos2θ),(2)

式中x=cosθ,θ为散射角,p_HG为Henyey-Greenstein(HG)相函数,g_HG为HG相函数的一阶矩,α为大粒子的相对密度。 等号右侧第一项中的p_HG(x,g_HG)描述了组织中大粒子的散射特性,第二项为描述小粒子散射的Rayleigh相函数。通过改变参量(α,g_HG)使该组合相函数发生变化。对(2)式进行Legendre展开,得到前三阶矩分别为

g1=αgHG,(3)g2=αgHG2+0.11-α,(4)g3=αgHG3。(5)

根据γ和δ的定义式可得到:

γ=1-g21-g1=1-g12/α+0.11-α1-g1,(6)δ=1-g31-g1=1-g13/α21-g1,(7)

式中g₁=g=αg_HG,g为各向异性因子。

由(3)~(6)式可以发现:1)一组参数(α,g_HG)唯一决定了g、γ和δ;2)由于各向异性因子g=αg_HG,因此只用g描述粒子的散射特性是不够的;3)当用一组参数(γ,g)表征单粒子的散射特性时,散射特性的差异反映在三阶参量δ上(图2)。图2显示了当g取不同值时,γ随δ的变化曲线。由图2可知:当γ很小,不同的SPF的δ值几乎相同,因此δ对漫反射率的影响可以忽略;不同的SPF的δ值的差异随着γ增大而增大,因此当γ值较大时需要考虑δ值的影响。

图 2. g取不同值时,γ随δ的变化

Fig. 2. Relation between γ and δ for different g

下载图片 查看所有图片

大部分组织的g值在0.6~0.99范围内^[15]。近年来的研究表明γ值在1.37~1.87范围内^[20](皮肤的γ值最低,肝脏的γ值最高),但也有一些较早的研究认为组织的γ值更低(γ=1.14)^[17,19,21]。在本研究中,在可以代表组织光学特性的范围内考虑了7个g值和8个γ值,g∈{0.65,0.70,0.75,0.80,0.85,0.90,0.95},γ∈{1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9},共得到47个不同的SPF。每对参数(γ,g)表示一个SPF,但根据(2)式可知其中有9组参数不能用来表示SPF,这9组参数为:(γ=1.7,g=0.65)、(γ=1.8,g=0.65)、(γ=1.9,g=0.65)、(γ=1.8,g=0.70)、(γ=1.9,g=0.70)、(γ=1.8,g=0.75)、(γ=1.9,g=0.75)、(γ=1.9,g=0.80)和(γ=1.9,g=0.85)。考虑了6个吸收系数和6个约化散射系数,参数μ_a的取值分别为0.001,0.00468,0.0214,0.1,0.468,2.14 mm^-1;参数μ'_s的取值分别为0.5,1.0,1.5,2.0,2.5,3.0 mm^-1;g的取值分别为0.65,0.70,0.75,0.80,0.85,0.90,0.95;γ的取值分别为1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9。光在介质中的传输过程用Wang等^[22]编写的Monte Carlo模拟(MCML)程序对μ_a、μ'_s、SPF取不同值时的漫反射过程进行模拟计算,发射并跟踪10⁶个光子,共得到47×6×6=1692组反射率数据,每一组反射率数据都记录了垂直入射到介质中的光子经随机散射后能从介质表面逸出的位置和光子数。光子逸出的位置用宽度为0.02 mm的环形网格标记,覆盖了直径为20 mm的圆形区域,不同探测孔径的反射率数据可以通过对覆盖环内的光子数求和得到。

3 结果与分析

图3(a)、(b)为SPF不同时介质的漫反射率随孔径的变化。在图2上共选取了7个SPF:A(γ=1.2,g=0.65),B(γ=1.4,g=0.65),C(γ=1.6,g=0.65),D(γ=1.5,g=0.65), E(γ=1.5,g=0.75),F(γ=1.5,g=0.85)和G(γ=1.5,g=0.95)。图3(a)中的SPF有相同的g值,即g=0.65;而图3(b)的SPF有相同的γ值,即γ=1.5。为了便于比较,其他光学参量取相同的数值,即μ_a=0.1 mm^-1,μ'_s=2 mm^-1。图3(c)、(d)分别为图3(a)、(b)中漫反射率的改变程度。图3(c)中漫反射率的改变以B(γ=1.4,g=0.65)为参考,图3(d)中漫反射率的改变则以E(γ=1.5,g=0.75)为参考。从图中不难发现,当d很小时,漫反射率对SPF很敏感,但敏感度随着d的增大而不断减小。当d很大时,如d>1.8 mm时, 漫反射率与SPF无关。当d较小时,对于不同的SPF,漫反射率的差别比较明显。漫反射率随g的增大而增大,这归因于SPF的后向散射分量的增加,因此g=0.60时的漫反射率大于g=0.80时的漫反射率。

图 3. (a) g=0.65时,漫反射率随孔径d的变化;(b) γ=1.5时,漫反射率随孔径d的变化;(c)以B点为参考,

Fig. 3. (a) Variation of diffuse reflectance with detector aperture when g=0.65; (b) variation of diffuse reflectance with detector aperture when γ=1.5; (c) relative change of diffuse reflectance in Fig. 3(a) with reference to point B; (d) relative change of diffuse reflectance in Fig. 3(c) with reference to point E

下载图片 查看所有图片

图3(a)中不同相函数的漫反射率随孔径d的变化;(d)以E点为参考,图3(b)中不同相函数的漫反射率随孔径d的变化

为了研究参数k₁(或k₂) 与相函数参量γ之间的关系,在图2中的曲线1上选取了5个SPF:g的取值相同而γ的取值不同,g=0.65,γ∈{1.2,1.3,1.4,1.5,1.6}。图4 (a)显示了当选取d=0.6 mm和μ_a=0.1 mm^-1时,不同SPF的漫反射率随μ'_s的变化曲线。对于吸收系数μ_a,也可以得到漫反射率随μ'_s的变化是线性的。总共获得5×6×6=180个反射率数据,代入(1)式进行拟合后可得到k₁、k₂随γ变化的曲线,如图4(b)所示。从图4(b)可以看出,k₁随γ的增大而增大,而k₂则随γ的增大而减小。在图2中的垂直虚线上选取了7个SPF。在这7个SPF中,γ取相同的值(γ=1.5),但g值不同,g∈{0.65,0.70,0.75,0.80,0.85,0.90,0.95},相应的δ∈{1.90,1.94,1.97,2.00,2.03,2.05,2.07}。按照与图4(a)、(b)一样的拟合过程可以得到图4(c)、(d)。由图4(c)可知,当d=0.6 mm,μ_a=0.1mm^-1时,对于7个不同的SPF,漫反射率均随μ'_s呈线性变化。图4(d)为k₁、k₂随δ变化的曲线,可以看出,k₁随δ的增大而减小,而k₂则随δ的增大而增大。

图 5. 当μa=0.1 mm-1,μ's=2 mm-1时,不同的g值下漫反射率随γ的变化

Fig. 5. Relation between diffuse reflectance and γ for various g when μa=0.1 mm-1 and μ's=2 mm-1

下载图片 查看所有图片

上述分析表明,与相函数相关的两个高阶光学参量γ和δ对漫反射的影响具有相反的效应。图5显示了γ和δ对反射率的共同作用,即当μ_a=0.1 mm^-1,μ'_s=2 mm^-1时,对于不同的g值,漫反射率随γ变化的曲线。对于其他的μ_a和μ'_s值,得到的R_d-γ曲线与图5是相似的。从图5中能够发现,漫反射率随γ的增大而增大,但随着g值的增加而增大的趋势较缓,这是由于γ增大的同时δ也增大了。

图5Fig. 5

当介质的g值较小时,随着γ增加,漫反射率的变化幅度较大,而δ对漫反射率的影响很弱, 因此可以只考虑漫反射率与γ之间的关系。对于图3(b)中g=0.65的情况,通过非线性拟合可以得到k₁=1.83γ²-4.33γ+5.97(相关度为0.998),k₂=-0.175γ²+0.328γ+0.102(相关度为 0.988),漫反射率的表达式具有如下形式:

Rd=1.83γ2-4.33γ+5.97)μ's1+(-0.175γ2+0.328γ+0.102)μa。(8)

值得注意的是k₂≈1[图4(b)],因此对于μ_a≪1的弱吸收介质,反射率R_d与γ呈二次函数关系。

当介质的g值在0.7~0.8之间时,R_d随γ的变化呈不同的变化形式,并且随机性增加。出现这种情况是由于γ随g变化的函数关系不是单调的,故而导致γ随g的变化有递增的,有递减的,还有先增后减的 ^[23]。

当介质的g值较大时,需要考虑γ和δ对R_d的共同效应。因此,对于给定的探测器孔径d,k₁和k₂是参量γ和δ的函数。考虑γ∈{1.4,1.5,1.6,1.7,1.8}和g∈{0.80,0.85,0.90,0.95}的情况,通过多变量的非线性拟合可以得到k₁=2.44+2.56γ-1.44δ(相关度为0.978),k₂=1.34-0.633γ+0.315δ(相关度为0.955),漫反射率可以表示为

Rd=(2.44+2.56γ-1.44δ)μ's1+(1.34-0.633γ+0.315δ)μa。(9)

根据(9)式,虽然R_d随g值增大而降低了,但散射相函数对R_d的影响增大了。另一方面,当组织SPF的形式确定以后,高阶参量γ和δ也就确定了,因此在特定孔径下可以根据反射率数据分析组织的光学特性,从而达到利用散射光监测生物组织形态的目的。

4 结论

在单光纤反射探针技术中,探测器只能收集到光源附近微区的散射光,这些散射光因为在组织中经历了较少的散射而携带了更多的关于组织微观形态的信息,因此需要考虑与散射相函数相关的高阶光学参量。在这种情况下,虽然表示可测量宏观物理量(例如反射率)与生物组织光学参量之间关系的数学模型具有非常复杂的形式,但这为获得组织微观形态信息提供了新的途径。给出了适用于各向异性因子较大和较小两种情况下的简化的半经验模型及其适用条件,这对促进漫反射光谱技术在生物医学领域中的应用具有非常重要的意义。