劉迎劉小君 齊貝貝田會娟

生物組織的δ-P1近似漫反射光學(xué)模型*

劉迎1)劉小君2)齊貝貝1)田會娟3)

1)(天津大學(xué)理學(xué)院，天津300072)
2)(精密儀器與光電子工程學(xué)院，光電信息技術(shù)科學(xué)教育部重點(diǎn)實(shí)驗(yàn)室，天津300072)
3)(天津工業(yè)大學(xué)，信息與通訊工程學(xué)院，天津300160)
(2010年8月2日收到;2010年9月2日收到修改稿)

本文研究了生物組織的改進(jìn)的δ-P1近似漫反射光學(xué)模型，推導(dǎo)了含有等效光源一階矩的雙點(diǎn)源近似空間分辨漫反射解Rδ－P1(ρ)．研究表明，考慮等效光源一階矩的光學(xué)模型，較好地描述了具有強(qiáng)的前向散射特性和較大吸收系數(shù)的生物組織散射特性;與漫射近似下的漫反射率RSDA(ρ)相比，新的光學(xué)模型能較好地描述光源附近的漫輻射強(qiáng)度分布，并且由于解析表達(dá)式中含有散射相函數(shù)的二階參量γ，這對于在近光源區(qū)域由空間分辨漫反射測量數(shù)據(jù)來反演生物組織的光學(xué)參量具有重要的意義．

組織光學(xué)，δ-P1近似，等效光源，微區(qū)漫反射

PACS:42．25．Dd，42．62．Be，87．64．Cc

1.引言

對混濁介質(zhì)中光子輸運(yùn)問題的研究為漫射光學(xué)光譜技術(shù)提供了理論基礎(chǔ)．Boltzmann輸運(yùn)方程已經(jīng)廣泛地用于研究光在生物組織中的傳播問題，它的解描述了輻射通量(fluence)、反射率(reflectance)或輻射強(qiáng)度(radiance)等可測量物理量與生物組織的光學(xué)參量之間的關(guān)系［1］．在輸運(yùn)理論中，生物組織的光學(xué)性質(zhì)用吸收系數(shù)μa、散射系數(shù)μs和散射相函數(shù)p(μ)描述．在實(shí)際問題中，散射相函數(shù)p(μ)的復(fù)雜性，限制了對Boltzmann輸運(yùn)方程的精確求解［2］．為了得到解析的解，通常將Boltzmann輸運(yùn)方程以球諧函數(shù)的形式展開成為一組封閉的方程組，P1近似只取零階和一階項(xiàng)［3］．在P1近似基礎(chǔ)上，如果在輸運(yùn)方程球諧函數(shù)的形式展開式中，光源項(xiàng)只取到零階，這就是漫射近似．對于漫反射問題，通常要求光源與探測器之間的距離遠(yuǎn)大于一個(gè)輸運(yùn)平均自由程，對于生物組織，一般需要大于幾個(gè)毫米．這個(gè)條件限制了漫射近似的應(yīng)用，例如，當(dāng)被測組織的體積比較小或組織的吸收比較大時(shí)，需要測量光源附近組織的輻射強(qiáng)度分布［4］;當(dāng)對表層組織測量時(shí)，需要研究散射次數(shù)較少的非漫射近似光子行為［5］．

近些年對解析理論模型的研究主要是圍繞克服漫射近似的這些缺點(diǎn)進(jìn)行的．一些研究者采用不同方法來擴(kuò)大漫射近似的適用范圍，例如對漫射系數(shù)D進(jìn)行改進(jìn)［6，7］和采用雙點(diǎn)源近似［8，9］．2005年，Ripoll等［6］將一個(gè)適合于吸收的漫射系數(shù)D用在漫射近似中，并以實(shí)驗(yàn)證實(shí)改進(jìn)的漫射近似能夠描述光在強(qiáng)吸收介質(zhì)中的傳輸．Hull和Foster［7］在研究P3近似漫反射時(shí)提出了一種P3近似與漫射近似的混合模型，該理論模型將漫射近似中的漫射系數(shù)D用描述P3近似的漫射系數(shù)Dasym替代．在研究漫反射問題時(shí)，雙點(diǎn)源近似的空間分辨漫反射解被證明能較好地描述光源附近的輻射強(qiáng)度分布［8，9］．

漫射近似能對接近各向同性的光輻射做出精確的描述，而光源附近的光輻射分布是各向異性的．Bevilacqua等［10］用Monte Carlo方法對光源附近約1個(gè)輸運(yùn)平均自由程以內(nèi)的漫反射進(jìn)行了研究，指出光源附近組織的光輻射與相函數(shù)的前兩階矩有很大關(guān)系，并引入了與相函數(shù)二階矩有關(guān)的光學(xué)參量γ．Kienle等人［11］的研究表明，用Monte Carlo (MC)模型反演組織的光學(xué)參量時(shí)，僅考慮g值是不夠的，還應(yīng)該考慮二階參量γ，否則產(chǎn)生誤差會高達(dá)100%．2008年Seo等人［12］研究了δ-P1近似下的空間分辨漫反射率問題時(shí)，通過引入含有高階參量的反照率，對光源附近的光輻射分布給出了更好的數(shù)學(xué)描述．我們近些年的研究表明，二階參量γ的影響主要是通過光源帶來的，因此考慮光源近似的高階項(xiàng)的必要的．

本文研究了生物組織的δ-P1近似漫反射光學(xué)模型，推導(dǎo)了含有等效光源一階矩的雙點(diǎn)源近似空間分辨漫反射解Rδ-P1(ρ)．與漫射近似RSDA(ρ)僅僅是μa和μ's的函數(shù)相比，Rδ-P1(ρ)還與參量γ有關(guān)，因此新的光學(xué)模型能較好地描述光源附近的漫輻射強(qiáng)度分布．

2.光輻射輸運(yùn)方程的δ-P1近似

穩(wěn)態(tài)情況下，介質(zhì)中r點(diǎn)^s方向的輻射強(qiáng)度L(r，^s)滿足Boltzmann方程［1］

其中，μt=μa+μs為介質(zhì)的衰減系數(shù)，μa和μs分別為介質(zhì)的吸收系數(shù)和散射系數(shù)．p(^s，^s')為散射相函數(shù)，表示來自^s'方向的光被散射到^s方向的概率，通常采用的實(shí)際散射相函數(shù)是Henyey-Greenstein相函數(shù)，

其中g(shù)的取值從－1到1，稱為各向異性因子．g=0表示各向同性散射，g=－1表示完全后向散射，g= 1表示完全前向散射．生物組織具有前向散射的特性，一些組織的g值達(dá)到0.99—0.999［13］．如果用delta函數(shù)δ(^s－^s')描述光子的前向散射，散射相函數(shù)可表示為

輻射強(qiáng)度通常包括尚未被散射的約化強(qiáng)度Lc(r，^s)(來自光源)和漫射強(qiáng)度Ld(r，^s)兩部分，即L(r，^s)能夠表示為

將(3)式和(4)式代入(1)式后得到

其中，μ*t=μa+μ*s，μ*s=μs(1－f)．式中Sc(r，^s)，是由約化入射強(qiáng)度Lc(r，^s)引起的等效的源函數(shù)，又稱為有效光源項(xiàng)．

在δ-P1近似下，(5)式中Ld(r，^s)，p*(^s，^s')和Sc(r，^s)在其球諧函數(shù)展開式中只取前兩項(xiàng)，

式中η=cosθ=^s·^r，η0=cosγ=^s·^s'．將(6)式代入(5)式得到φ0(r)和φ1(r)滿足的方程組

在δ-P1近似下，相函數(shù)pδ(^s，^s')與pHG(^s，^s')應(yīng)該滿足二階相似關(guān)系．pHG(^s，^s')的Legendre矩，因此g1=g， g2=g2．根據(jù)(2)式，pδ(^s，^s')的一階和二階Legendre矩可表示為f+(1－f)g1*和f+(1－f)g2*．根據(jù)相似關(guān)系［7］得到

根據(jù)(6b)式，有g(shù)2*=0，因此從(8b)式得到f=g2．利用f=g2，我們得到

根據(jù)(9)式，(7)式中μ*tr=μa+μs(1－g)=μtr，因此φ0(r)和φ1(r)滿足的方程組最終可表示為

漫射近似假設(shè)了q1(r)=0，因此有

對于均勻無限介質(zhì)，(11a)式格林函數(shù)解為

式中D=1/3μtr是漫射系數(shù)．輻射通量φ0(r)∝ΦG(r)?q0(r)，符號“?”表示卷積;輻射強(qiáng)度

3.δ-P1近似下的有效光源項(xiàng)及其空間分辨漫反射解

從上面分析我們看到，如果不考慮有效光源的一次項(xiàng)q1(r)，我們并不能在輸運(yùn)方程的δ-P1近似理論框架下解決任何新的問題．下面我們在(10)式基礎(chǔ)上，考慮q1(r)≠0情況下φ0(r)和φ1(r)的解．當(dāng)q1(r)≠0時(shí)，由(10 b)式得到

與漫射近似的(11b)式相比，式中多了q1(r)/μtr項(xiàng)．利用，得到

由約化入射強(qiáng)度Lc(r，^s)引起的等效的源函數(shù)Sc(r，^s)是復(fù)雜的，在實(shí)際應(yīng)用中往往將Lc(r，^s)作點(diǎn)源近似．對無限細(xì)的準(zhǔn)直光束沿z軸入射到一半無限均勻介質(zhì)的情況，考慮雙點(diǎn)源近似［9］

式中α*=μ*s/μ*t，表示反照率;z01=2/μ*t=2 z0， z02=0，z03=－2 zb，z04=－2zb－z01，如圖1所示;外推長度zb=2 AD，其中A=(1+Reff)/(1－Reff)，有效反射系數(shù)Reff與界面的相對折射率nrel有關(guān)，可近似表示為［14］

在本文中，均采用nrel=1.4．根據(jù)(15)和(16)式，(13)式中q1(r)可以表示為

q1(r)=0時(shí)的空間分辨漫反射率為［7］

其中

根據(jù)前面分析，考慮光源一次項(xiàng)q1(r)的δ-P1空間分辨漫反射率可以表示為

4.結(jié)果與討論

為評估所得到的雙點(diǎn)源δ-P1近似模型，我們采用由Wang等［15］編寫的MC模擬程序得到的反射率RMC(ρ)作為參考標(biāo)準(zhǔn)，將Rδ-P1(ρ)與漫射近似漫反射解RSDA(ρ)進(jìn)行了比較．對于RSDA(ρ)的計(jì)算，(19)式中φ0(ρ)和j(ρ)為［16］

在圖2和圖3中，圖示的橫坐標(biāo)均為光源與探測器之間的距離;圖2中的縱坐標(biāo)表示空間分辨漫反射率R(ρ)，圖3中的縱坐標(biāo)表示漫反射率R(ρ)與檢驗(yàn)標(biāo)準(zhǔn)RMC(ρ)之間的誤差ΔR/R=|R(ρ)－RMC(ρ)|/RMC(ρ)，R(ρ)分別表示Rδ-P1(ρ)和RSDA(ρ)．此外，進(jìn)行MC模擬計(jì)算時(shí)，考慮準(zhǔn)直光束入射到半無限均勻介質(zhì)的情況，模擬的光子個(gè)數(shù)為106．

在圖2中，我們比較了RMC(ρ)，RSDA(ρ)和Rδ-P1(ρ)三種反射率曲線，介質(zhì)的參數(shù)取值為μa= 0.02 mm－1，μs=20.0 mm－1，g=0.9．在此參數(shù)取值下，該介質(zhì)的輸運(yùn)反照率α'=μ's/μtr=0.99，輸運(yùn)平均自由程l't=1/μtr=0.495 mm．為了說明(20)式中光源項(xiàng)q1(ρ)/μtr對反射率的影響，曲線*Rδ-P1(ρ)是(20)式中不含q1(ρ)/μtr項(xiàng)的反射率曲線．從圖中可以看出，與漫射近似模型相比，我們所得到的雙點(diǎn)源δ-P1近似漫反射解在一個(gè)輸運(yùn)平均自由程之內(nèi)(ρ＜0.4 mm)與MC模擬結(jié)果仍然保持很好的一致性，并且由于q1(ρ)/μtr項(xiàng)的存在，使其與MC模擬結(jié)果之間的誤差減?。?/p>

保持μs=20.0 mm－1，g=0.9不變，圖3考慮了吸收系數(shù)μa值分別為0.01 mm－1，0.02 mm-1，0.05 mm－1，0.1 mm－1，0.2 mm－1，0.5 mm－1六種情況，相應(yīng)的a'分別為0.995，0.990，0.976，0.952，0.909，0.8，而l't數(shù)值上是a'一半，因此橫坐標(biāo)ρ＜4 mm的范圍，覆蓋了大約8—10 lt．由圖3所示的誤差曲線可以看出，在ρ≥1.5 lt范圍內(nèi)，RSDA(ρ)曲線的誤差均大于10%，并且隨著μa的增大，誤差高達(dá)20%以上．而對于Rδ-P1(ρ)曲線，當(dāng)a'＜0.952時(shí)，其誤差均小于10%．

圖4中顯示了各向異性常數(shù)g取不同值時(shí)Rδ-P1(ρ)與Monte Carlo模擬計(jì)算結(jié)果RMC(ρ)之間的誤差．圖中曲線具有相同的μs=20.0 mm－1和μa=0.02 mm－1，g在0.7—0.99范圍內(nèi)取值(對于大多數(shù)生物組織，其g值在這個(gè)范圍內(nèi))．結(jié)果顯示在ρ≥0.5 lt范圍內(nèi)，誤差均小于10%．

研究表明，雙點(diǎn)源近似下的δ-P1近似理論模型的適用范圍優(yōu)于漫射近似理論模型，即能夠應(yīng)用于吸收系數(shù)較大的散射介質(zhì)，并能較好地描述光源附近的輻射強(qiáng)度分布．此外，由于該模型考慮了有效光源項(xiàng)對近場輻射強(qiáng)度分布的影響，使得漫反射解的表達(dá)式中含有表征生物組織的光學(xué)特性的高階光學(xué)參量γ，而漫射近似僅含有表征生物組織的光學(xué)特性的一階光學(xué)參量μa和μ's．這項(xiàng)工作將為生物組織光學(xué)特性參量的在體無損檢測提供重要的理論依據(jù)．

［1］Kim A D 2004 J．Opt．Soc．Am．A 21 820

［2］Wang Q H，Zhang Y Y，Lai J C，Li Z H，He A Z 2007 Acta Phys．Sin．56 1203(in Chinese)［王清華、張穎穎、來建成、張小娟、高宗慧、田會娟2007物理學(xué)報(bào)56 1203］

［3］Liu Y，Wang L J，Guo Y F，Zhang X J，Gao Z H，Tian H J 2007 Acta Phys．Sin．56 2119(in Chinese)［劉迎、王利軍、郭云峰、張小娟、高宗慧、田會娟2007物理學(xué)報(bào)56 2119］

［4］Xu T，Zhang C P，Chen G Y，Zhang X J，gao Z H，Tian H J 2005 Chin．Phys．14 1813

［5］Xu L Q，Li H，Xie S S 2008 Acta Phys．Sin．57 6024(in Chinese)［徐蘭青、李暉、謝樹森2008物理學(xué)報(bào)57 6024］

［6］Ripoll J，Yessayan D，Zacharakis G，Ntziachristos V 2005 J．Opt．Soc．Am．A 22 546

［7］Hull E L，F(xiàn)oster T H 2001 J．Opt．Soc．Am．A 18 584

［8］Tian H J，Liu Y，Wang L J，Zhang Z B，Xiao L F 2009 Acta Phys．Sin．58 243(in Chinese)［田會娟、劉迎、王利軍、張智卜、肖立峰2009物理學(xué)報(bào)58 243］

［9］Wang R，Liu Y 2010 Chinese J．Lasers．37 1147(in Chinese)［王銳、劉迎2010中國激光37 1147］

［10］Bevilacqua F，Depeursinge C 1999 J．Opt．Soc．Am．A 16 2935

［11］Kienle A，F(xiàn)orster F K，Hibst R 2001 Opt．Lett．26 1571

［12］InSeok Seo，Carole K．Hayakawa，Vasan Venugopalan 2008 Med．Phys．35 681

［13］Zhang Z X 2008 New Technology and Appliccations of Biomedical Photonics(Beijing:Science Press)(in Chinese)［張鎮(zhèn)西2008生物醫(yī)學(xué)光子學(xué)新技術(shù)及應(yīng)用(北京:科學(xué)出版社)］

［14］Groenhuis R A J，F(xiàn)erwerda H A，Bosch J J T 1983 I:Theory Appl Opt．22 2456

［15］Wang L H，Jacques S L，Zheng L Q 1995 Computer Methods and Programs in Biomedicine 47 131

［16］Kienle A，Patterson M S 1997 J．Opt．Soc．Am．A 14 246

PACS:42．25．Dd，42．62．Be，87．64．Cc

*Project supported by the National Natural Science Foundation of China(Grant No．60278004)．

E-mail:liuying@tju．edu．cn

δ-P1approximation model of biologicaltissues*

Liu Ying1)Liu Xiao-Jun2)Qi Bei-Bei1)Tian Hui-Juan3)
1)(College of Science，Tianjin University，Tianjin 300072，China)
2)(Key Laboratory of Opt-electronics Information Technical Science，EMC Tianjin University，Tianjin 300072，China)
3)(School of Information and Communication Engineering，Tianjin Polytechnic University，Tianjin 300160，China)
(Received 2 August 2010;revised manuscript received 2 September 2010)

A modified delta-P1approximation is studied in this paper．An analytic solution Rδ-P1(ρ)of this model for reflectance in double source approximation is developed by taking the source series expansion terms to the first order．It is shown that when the source in extended to first order，a better description for high absorption and high scattering media will be developed．The new model can give satisfactory results in describing irradiance distribution close to source compared with the RSDA(ρ)of diffusion approximation．Moreover，Rδ-P1(ρ)provides the prediction of the second-order parameterγ，which is significant for developing the inversion algorithm to recover optical parameters from spatially-resolved diffuse reflectance by using smaller source-detector separations．

tissue optics，δ-P1approximation，effective source，diffuse reflectance close to source

*國家自然科學(xué)基金(批準(zhǔn)號:60278004)資助的課題．

E-mail:liuying@tju．edu．cn