基于透射解析模型的時(shí)域擴(kuò)散熒光層析原理與實(shí)驗(yàn)

高 峰，李 嬌，張麗敏，趙會(huì)娟

(天津大學(xué)精密儀器與光電子工程學(xué)院，天津 300072)

熒光分子層析(fluorescence molecular tomography，F(xiàn)MT)[1-3]是結(jié)合近紅外擴(kuò)散熒光層析(fluorescence diffuse optical tomography，F(xiàn)DOT)理論和近紅外熒光探針技術(shù)以實(shí)現(xiàn)生物活體內(nèi)特異大分子生化過(guò)程的無(wú)損三維定量觀測(cè)方法，它可提供其他分子成像和現(xiàn)有平面熒光強(qiáng)度成像技術(shù)所不具備的超靈敏度和三維層析能力．FDOT作為 FMT的成像基礎(chǔ)技術(shù)已在連續(xù)光和頻域調(diào)制測(cè)量模式下獲得了原理性成功和初步應(yīng)用，而時(shí)域FDOT模式研究則處于基本測(cè)量技術(shù)和理論體系建立階段．由于時(shí)域 FDOT擁有多組分熒光產(chǎn)率和壽命圖像的同時(shí)重建能力和高信噪比完整測(cè)量信息的獲取能力，從而在小動(dòng)物模型應(yīng)用中獲得了廣泛的重視．相對(duì)于熒光產(chǎn)率，熒光壽命客觀地提供了特異分子的生化微環(huán)境信息，對(duì)生物醫(yī)學(xué)基礎(chǔ)研究和臨床診斷均有重要的意義．

本文提出了基于透射解析模型的時(shí)域擴(kuò)散熒光層析圖像重建方法．算法應(yīng)用廣義脈沖譜技術(shù)和歸一化玻恩比逆模型，其中，廣義脈沖譜技術(shù)通過(guò)拉普拉斯變換將耦合擴(kuò)散方程從時(shí)域轉(zhuǎn)換到復(fù)頻域，可實(shí)現(xiàn)多參數(shù)、多組份熒光圖像的快速重建，而歸一化玻恩比逆模型表示則可免除系統(tǒng)耦合系數(shù)的標(biāo)定及消除系統(tǒng)響應(yīng)的影響．本文采用外推條件下無(wú)限平板透射光學(xué)擴(kuò)散方程解析解作為逆模型權(quán)重矩陣的計(jì)算基礎(chǔ)，由于求解熒光參數(shù)線(xiàn)性方程的過(guò)程的欠定性，將基于行操作的代數(shù)重建技術(shù)用于線(xiàn)性求逆過(guò)程，其應(yīng)用優(yōu)勢(shì)在于每次迭代只需單行運(yùn)算與存儲(chǔ)，這有助于采用高密度空間剖分實(shí)現(xiàn)高分辨率的圖像重建．在算法驗(yàn)證研究中，首先進(jìn)行基于有限元正向模型數(shù)據(jù)的數(shù)值模擬研究，對(duì)算法的可行性及空間分辨率和魯棒性等關(guān)鍵性能方面做初步的驗(yàn)證與評(píng)估；然后采用所提算法和多通道時(shí)間相關(guān)單光子計(jì)數(shù)(time-correlated single photon counting，TCSPC)測(cè)量系統(tǒng)對(duì)自制仿體進(jìn)行圖像重建實(shí)驗(yàn)，進(jìn)一步驗(yàn)證了重建算法的有效性和準(zhǔn)確性．

1 算法原理

采用擴(kuò)散方程作為光在組織體中的傳輸模型，熒光檢測(cè)量的數(shù)學(xué)表達(dá)式可由激發(fā)光和出射熒光之間耦合擴(kuò)散方程得到[4]．另外采用了在解決以波動(dòng)和擴(kuò)散為主導(dǎo)現(xiàn)象的逆向問(wèn)題中，已被證明的一種有效方法——廣義脈沖譜技術(shù)(generalized pulse spectrum technique，GPST)，通過(guò)拉氏變換將時(shí)域擴(kuò)散方程轉(zhuǎn)化到復(fù)頻域，從而使原四維全時(shí)空問(wèn)題降為與時(shí)間相關(guān)的三維問(wèn)題[5-6]．

時(shí)域擴(kuò)散熒光的逆向問(wèn)題即時(shí)域 FDOT圖像重建，其目的是在背景光學(xué)參數(shù)分布已知的情況下，實(shí)現(xiàn)對(duì)熒光產(chǎn)率和熒光壽命空間分布的同時(shí)重建．即基于已得光子傳輸模型的正向算子，來(lái)求解所需的熒光參數(shù)分布，從而實(shí)現(xiàn)擴(kuò)散熒光層析成像的圖像重建．具體過(guò)程如下：源在sr處接收在dr處的熒光光子密度可看作是所有熒光體元 dV在整個(gè)體積上的積分，從而熒光密度可表示為[7]

式中：Φm( rd, rs, q)為源位于 rs處在 rd處探測(cè)到的熒光光子密度的拉普拉斯變換值；Φx( r, rs, q)為源位于 rs處接收在r處的激發(fā)光光子密度的拉普拉斯變換值；Gm(rd, r, q)為源位于r處接收在 rd處的擴(kuò)散方程的格林函數(shù)的拉普拉斯變換值[8]．

在檢測(cè)量處理方面，為了重建熒光參數(shù)，采用了歸一化玻恩比[7]，即用在檢測(cè)器處探測(cè)到的熒光光流量除以相應(yīng)的激發(fā)光的光流量，這樣可以有效地消除不同的源和探測(cè)器對(duì)計(jì)算引起的不便．由式(1)推導(dǎo)及檢測(cè)量的處理后得

在無(wú)限平板透射模型中， Φx( ri)為外推邊界條件下的拉普拉斯變換后的擴(kuò)散方程的解析解[8]，即

式中：r±,m為根據(jù)外推邊界近似中對(duì)稱(chēng)的正負(fù)源位置；zb為外推邊界的坐標(biāo)；z0= 1 /μs′為單一的近似各向同性源入射光線(xiàn)位置；為有效衰減系數(shù)．

由式(3)中的離散過(guò)程導(dǎo)致了所需求解的體元上的熒光參數(shù)個(gè)數(shù)遠(yuǎn)遠(yuǎn)多于測(cè)量數(shù)據(jù)個(gè)數(shù)，使得對(duì)體元熒光參數(shù)的求解過(guò)程變?yōu)榍范▎?wèn)題，從而意味著測(cè)量數(shù)據(jù)上的微小變化可能引起重建圖像的完全變異．且線(xiàn)性方程組的解易受到噪聲干擾，因此很難直接用矩陣求逆的方法得到．在此情形下，只能按照某種近似原則將原問(wèn)題的直接求解轉(zhuǎn)化為求原問(wèn)題的穩(wěn)定的、合理的近似解，該過(guò)程稱(chēng)為正規(guī)化過(guò)程．本文采用了代數(shù)重建技術(shù)(algebraic reconstruction technique，ART)[9]，它是根據(jù)求解線(xiàn)性方程組思想，以解決離散化數(shù)字圖像重建為目的，通過(guò)對(duì)矩陣進(jìn)行逐行計(jì)算，從而有效快速地對(duì)方程進(jìn)行求解[10]．為了同時(shí)求出熒光產(chǎn)率af()ημr及熒光壽命()τr，本文還引入了一對(duì)實(shí)變換因子對(duì)：．其中分別為背景在激發(fā)光和熒光波長(zhǎng)下的吸收系數(shù)為背景的熒光壽命．由式(1)，可以得出其熒光參數(shù)為

2 模擬研究的重建結(jié)果討論及分析

本文在模擬研究中應(yīng)用了同課題組的有限元方法對(duì)非均勻目標(biāo)體產(chǎn)生正向問(wèn)題的數(shù)據(jù)，再將所得數(shù)據(jù)代入設(shè)計(jì)的算法中進(jìn)行圖像重建．圖 1(a)所示為設(shè)計(jì)的非均勻目標(biāo)體，通過(guò)在正向數(shù)據(jù)中混入白噪聲及改變目標(biāo)體的中心間距(center-to-center separation，CCS)來(lái)對(duì)算法進(jìn)行噪聲魯棒性和空間分辨率的研究．首先討論噪聲對(duì)圖像重建的影響，結(jié)果表明圖像的抗噪性受到q取值的影響：當(dāng)時(shí)，如圖 1(b)所示，重建圖像的熒光壽命可以承受 50,dB以上的噪聲，抗噪聲的能力較差；而產(chǎn)率幾乎不受噪聲影響．當(dāng)時(shí)，如圖1(c)所示，重建的熒光壽命可以承受 35,dB以上的噪聲，有較好的抗噪能力；而產(chǎn)率相對(duì)的抗噪聲性能變差．且熒光壽命和產(chǎn)率的量化值都變大，從而影響了重建算法的分辨率．由以上的重建結(jié)果可以看出，q值增大可以提高結(jié)果中熒光壽命的抗噪性，但也影響了成像的質(zhì)量，犧牲了圖像的分辨率．

圖1 模型的模擬結(jié)果Fig.1 Simulated results of model

對(duì)于算法空間分辨能力的討論結(jié)果，如圖 1(d)和(e)所示，隨著 CCS的減小，在熒光產(chǎn)率和熒光壽命的重建中，兩個(gè)目標(biāo)體之間可分度均逐漸下降．當(dāng)CCS大于 10,mm 時(shí)，熒光產(chǎn)率可以分辨；而當(dāng) CCS大于14,mm時(shí)，熒光壽命可以分辨．

由上述結(jié)果可以看出不同熒光壽命和熒光產(chǎn)率區(qū)域可以明顯區(qū)分，重建圖像質(zhì)量較好．同時(shí)值得注意的是，熒光產(chǎn)率的重建深度及量化精度要比熒光壽命的好些．對(duì)重建算法的空間分辨率進(jìn)行測(cè)試，結(jié)果表明熒光產(chǎn)率比壽命好，壽命在深度信息上有偏差．

3 實(shí)驗(yàn)及分析

3.1 實(shí)驗(yàn)材料及設(shè)備

3.1.1 測(cè)量系統(tǒng)

實(shí)驗(yàn)研究基于本實(shí)驗(yàn)室所發(fā)展的多通道TCSPC皮秒時(shí)間分辨測(cè)量系統(tǒng)，采用透射檢測(cè)模式．測(cè)量過(guò)程為：皮秒脈沖激光器發(fā)出波長(zhǎng)為660,nm的激發(fā)光，經(jīng)直徑為 62.5,μm、數(shù)值孔徑為 0.22的光纖導(dǎo)出，之后1∶16光開(kāi)關(guān)將激發(fā)光源依次導(dǎo)入16個(gè)源位置．對(duì)于每次源入射，4個(gè) 4∶1光開(kāi)關(guān)分4次、每次選擇 4個(gè)檢測(cè)點(diǎn)做并行TCSPC 檢測(cè)，其中4個(gè)并行檢測(cè)點(diǎn)被置于相鄰區(qū)域，探測(cè)光纖的直徑為 500,μm、數(shù)值孔徑為 0.37．隨后切換源位置重復(fù)上述測(cè)量，直至 16個(gè)源位置“掃描”完畢．其中，在測(cè)量熒光時(shí)，接收光纖需經(jīng)過(guò)帶通濾波器將激發(fā)光濾除．實(shí)驗(yàn)測(cè)量系統(tǒng)示意如圖2所示．

圖2 TCSPC測(cè)量系統(tǒng)示意Fig.2 Measurement system of TCSPC

3.1.2 實(shí)驗(yàn)仿體選用固體模型作為背景，其主要組成成分是聚甲醛(Polyformaldehyde)．如圖3所示，仿體為100,mm×25,mm×70,mm的長(zhǎng)方體，吸收系數(shù)μa(r)=0.003,8 mm-1，約化散射系數(shù)(實(shí)驗(yàn)中所有光學(xué)參數(shù)均基于反射模式的光學(xué)參數(shù)時(shí)域測(cè)量方法測(cè)得).根據(jù)經(jīng)驗(yàn)選擇其熒光差率 ημaf(r)＝0.000,01 mm-1，熒光壽命 τ(r)＝10,000,ps．在本實(shí)驗(yàn)中，設(shè)計(jì)了2個(gè)固態(tài)仿體：仿體 1為背景中含有一個(gè)目標(biāo)體，仿體 2背景中則包含了兩個(gè)相同幾何尺寸的目標(biāo)體且其中心距為 20,mm．目標(biāo)體是直徑為 5,mm、高為60,mm的圓柱形孔，由1% Intralipid 溶液和Cy5.5熒光染料(激發(fā)和發(fā)射峰值分別在 670,nm 和 710,nm)的混合溶液填充．其吸收系數(shù)μa(r)=0.001,7,mm-1，約化散射系數(shù)目標(biāo)體的熒光產(chǎn)率和壽命為需重建的參數(shù)．在重建過(guò)程中只需在視場(chǎng)范圍(field of view，F(xiàn)OV)進(jìn)行重建．由于正向模型計(jì)算采用了外推邊界條件，因此重建區(qū)域?yàn)?64,mm×29,mm×60,mm的長(zhǎng)方體．光源的位置在外推邊界條件下，分布于 Z＝2,mm 的 XY平面，視場(chǎng)范圍為21,mm×21,mm，光源間隔為 7,mm×7,mm，共有 16個(gè)光源．16個(gè)探測(cè)器的位置分布于Z＝27,mm的XY平面，與光源位置相對(duì)．目標(biāo)體中心坐標(biāo)(X＝32,mm，Z＝15,mm)．離散化后體積元為V＝L×W×H＝1,mm×1,mm×1,mm的立方體．整個(gè)計(jì)算模型可離散成111,360個(gè)立方體單元．圖 3(a)為仿體、計(jì)算模型及坐標(biāo)系示意．在 ART計(jì)算過(guò)程中迭代次數(shù)為15，松弛因子選為0.5．

圖3 計(jì)算和實(shí)驗(yàn)?zāi)Ｐ褪疽釬ig.3 Calculating and experimental model

3.2 實(shí)驗(yàn)測(cè)量及結(jié)果分析

在實(shí)驗(yàn)測(cè)量時(shí)，TCSPC各部分的設(shè)定如下：TAC范圍為70.01×10-9；CFD最低限定為-15.69；PMT增益為 80%，共有 4,096個(gè)采樣點(diǎn)，采樣時(shí)間間隔為17.1,ps．在整個(gè)測(cè)量過(guò)程中，整個(gè)系統(tǒng)在黑暗狀態(tài)下，以防止外部光源影響．激發(fā)光和熒光測(cè)量過(guò)程中積分時(shí)間均固定在 10,s，改變激發(fā)光的光強(qiáng)．圖 4(a)顯示了第 1個(gè)光源激發(fā)時(shí)，16個(gè)探測(cè)通道測(cè)得的激發(fā)光(EXTdata)和熒光(EMSdata)的歸一化時(shí)間擴(kuò)展曲線(xiàn)．可以看出，熒光信號(hào)較激發(fā)光有一定的時(shí)間延遲，且有一定程度的展寬．

圖4(b)和(c)分別為仿體 1和仿體 2的熒光產(chǎn)率和壽命空間分布的重建結(jié)果．由重建圖像的XY平面的剖面圖看出算法對(duì)熒光產(chǎn)率圖像進(jìn)行了很好的重建，重建目標(biāo)的形狀、位置、大小均與實(shí)際目標(biāo)相符．重建目標(biāo)的熒光壽命形狀和位置與實(shí)際情況較符合，但體積要比實(shí)際目標(biāo)大，與模擬實(shí)驗(yàn)結(jié)果也是相符的．特別在仿體 2的重建圖像中可以很好分辨兩個(gè)目標(biāo)體，從而體現(xiàn)了算法良好的分辨力．兩個(gè)仿體的成功重建效果說(shuō)明了算法具有可靠性且適用于不同目標(biāo)體的成像．但在重建圖像中目標(biāo)邊界處出現(xiàn)了偽逆，尤其熒光壽命圖像較明顯，我們認(rèn)為此現(xiàn)象是由于目標(biāo)與背景交接處的光學(xué)特性(主要是熒光特性)變換較大，使得熒光擴(kuò)散光方程產(chǎn)生了誤差．因此，由重建圖像及實(shí)驗(yàn)討論可知，實(shí)驗(yàn)結(jié)果很好地驗(yàn)證了本文提出的透射解析模型時(shí)域擴(kuò)散熒光層析重建算法可行性，且與模擬結(jié)果相似．

圖4 測(cè)量和重建結(jié)果Fig.4 Measurement and reconstruction results

4 結(jié) 語(yǔ)

本文提出了基于透射解析模型時(shí)域擴(kuò)散熒光層析重建方法．應(yīng)用有限元方法提供正向問(wèn)題的數(shù)據(jù)，在數(shù)值模擬的基礎(chǔ)上評(píng)價(jià)了算法的空間分辨率及魯棒性，驗(yàn)證了算法具有強(qiáng)度測(cè)量成像模式無(wú)可比擬的空間分辨率和多參數(shù)成像功能．另外采用此重建方法對(duì)自行制作非均勻仿體進(jìn)行圖像重建，結(jié)果中對(duì)于仿體的熒光產(chǎn)率和熒光壽命的良好重建，驗(yàn)證了本重建算法的可行性和實(shí)用前景．

［1］ Ntziachristos V，Tung C H，Bremer C，et al. Fluorescence molecular tomography resolves protease activity in vivo [J]. Nature Medicine，2002，8(7)：757-760.

［2］ 徐可欣，高 峰，趙會(huì)娟．生物醫(yī)學(xué)光子學(xué)[M].北京：科學(xué)出版社，2007.Xu Kexin，Gao Feng，Zhao Huijuan. Biomedical Photonics [M]. Beijing：Science Press，2007(in Chinese).

［3］ Cherry S R. In vivo molecular and genomic imaging：New challenges for imaging physics[J]. Physics in Medicine and Biology，2004，49：13-48.

［4］ O’Leary M A，Boas D A，Li X D. Fluorescence lifetime imaging in turbid media[J]. Optics Letters，1996，21(2)：158-160.

［5］ Gao Feng，Zhao Huijuan，Yamada Y，et al. Timeresolved diffuse optical tomography using a modified generalized pulse spectrum technique [J]. IEICE Transactions on Information and Systems，2002，E85-D (1)：133-142.

［6］ Gao Feng，Zhao Huijuan，Yamada Y．A linear，featured-data scheme for image reconstruction in timedomain fluorescence molecular tomography[J]. Optics Express，2006，14(16)：7109-7124.

［7］ Ntziachristos V，Weissleder R. Experimental threedimensional fluorescence reconstruction of diffuse media by use of a normalized Born approximation[J]. Optics Letters，2001，26(12)：893-895.

［8］ Andrea C D，Spinelli L，Comelli D，et al. Localization and quantification of fluorescent inclusions embedded in a turbid medium[J]. Physics in Medicine and Biology，2005，50：2313-2327.

［9］ 赫爾曼 G T. 由投影重建圖像：CT的理論基礎(chǔ) [M].北京：科學(xué)出版社，1985：25-38.Herman G T. Image Reconstruction from Projections：The Fundamentals of Computerized Tomography [M].Beijing：Science Press，1985：25-38(in Chinese).

［10］ Kienle A，Patterson M S. Improved solutions of the steady-state and the time-resolved diffusion equations for reflectance from a semi-infinite turbid medium [J].Journal of the Optical Society of America，1997，A14：246-254.

［11］ Kumar A T N，Skoch J，Bacskai J．Fluorescent-lifetimebased tomography for turbid media[J]. Optics Letters，2005，30(24)：3347-3349.