生物發(fā)光斷層成像的仿真研究

劉 勇，齊樹(shù)波，方勇軍，陳 鋒

0 引言

與傳統(tǒng)醫(yī)學(xué)成像模態(tài)，例如X線計(jì)算機(jī)斷層成像（X-ray computed tomography，XCT）不同，分子成像能在特異性分子探針的幫助下，在分子水平上反映生物體生理或病理變化，在分子生物學(xué)與臨床醫(yī)學(xué)之間架起了相互連接的橋梁[1]。根據(jù)成像原理不同，分子影像技術(shù)主要可分為單光子發(fā)射成像（single photon emission computed tomography，SPECT）、正電子發(fā)射成像（positron emission tomography，PET）、磁共振成像（magnetic resonance imaging，MRI）、超聲成像（ultrasound）以及光學(xué)成像（optical imaging，OI）等。相對(duì)于其他分子成像技術(shù)，光學(xué)分子成像技術(shù)具有無(wú)電離輻射、對(duì)生物體無(wú)害、靈敏度高且其成本相對(duì)較低等優(yōu)點(diǎn)[1]，目前已發(fā)展成為一種理想的成像方式。

作為光學(xué)分子成像的一個(gè)重要分支，近年來(lái)，生物發(fā)光斷層成像技術(shù)（bioluminescence tomography，BLT）得到了快速發(fā)展[2-5]。作為一種新的分子斷層成像模式，BLT能夠以非侵入的方式，在細(xì)胞和分子水平研究或監(jiān)測(cè)生理和病理過(guò)程，因而有望為疾病早期診斷、基因治療和藥物研發(fā)等研究提供有力工具。然而，BLT作為一種基于擴(kuò)散光的成像技術(shù)，對(duì)比XCT等傳統(tǒng)成像技術(shù)，具有較低的成像空間分辨率。主要原因?yàn)椋河捎诠庠谏锝M織中的強(qiáng)散射特性，一方面導(dǎo)致了光在生物組織傳播過(guò)程中的強(qiáng)烈衰減；另一方面，也導(dǎo)致了光在生物體內(nèi)的傳播不是沿直線方向進(jìn)行。上述2個(gè)原因致使BLT的求解極具病態(tài)性（ill-posed），進(jìn)一步影響了其生物醫(yī)學(xué)應(yīng)用。

在本文中，基于數(shù)值仿真實(shí)驗(yàn)，我們驗(yàn)證了BLT成像的可行性。同時(shí)，我們將稀疏理論應(yīng)用于重建過(guò)程，以克服BLT成像的病態(tài)性。實(shí)驗(yàn)結(jié)果表明，基于多角度、非接觸BLT成像系統(tǒng)（multi-view non-contactBLT imaging system），我們能夠解析自發(fā)光源在成像體內(nèi)的分布，定位誤差小于1mm，這將有望為疾病早期診斷及藥物研發(fā)等研究提供有力工具。

1 方法

1.1 BLT前向數(shù)學(xué)模型

對(duì)于BLT成像技術(shù)而言，準(zhǔn)確描述光在生物組織中的傳輸過(guò)程是成像的首要問(wèn)題。根據(jù)輸運(yùn)理論，光在生物組織中的傳播可以通過(guò)輻射傳輸方程（radiative transferequation，RTE）來(lái)描述。直接對(duì)輻射傳輸方程求解雖然可以較為準(zhǔn)確地描述光子在組織中的傳播過(guò)程，但計(jì)算量較大。為了在實(shí)際中更好地進(jìn)行模擬，需要對(duì)輻射傳輸方程進(jìn)行一定的簡(jiǎn)化。對(duì)于BLT成像而言，通常可采用連續(xù)波模式下的擴(kuò)散方程（diffusion equation，DE）[6]對(duì)其近似，如式（1）所示：

其中，Ω 描述了成像物體；S（r）為自發(fā)光光源；Φ（r）為光場(chǎng)分布；μa（r）與Φ（r）分別描述了組織的吸收和擴(kuò)散系數(shù)。

根據(jù)微分方程理論，擴(kuò)散方程（1）的求解需要有合適的邊界條件。對(duì)于光學(xué)成像而言，常用的邊界條件為外插邊界條件和Robin邊界條件。在本文中，我們使用Robin邊界條件[6]，如式（2）所示：

考慮到有限元方法（finiteelement method，F(xiàn)EM）[7]可以有效地處理復(fù)雜幾何形狀且各向異質(zhì)的成像問(wèn)題，因而，在本文中，我們使用有限元方法對(duì)式（1）進(jìn)行求解。當(dāng)重建區(qū)被離散成小的網(wǎng)格后，式（1）可被改寫為如下的線性方程：

式（3）中：

式中，φi（r）、φj（r）為相應(yīng)的試驗(yàn)函數(shù)，可從一個(gè)合適的試驗(yàn)空間中選取，且必須滿足平方可積條件。考慮到在式（3）中，矩陣A通常是正定的，于是，式（3）可進(jìn)一步改寫為：

式中，W為權(quán)重，用于映射成像體內(nèi)未知的生物發(fā)光源S到已知的成像表面測(cè)量值Φ。

1.2 BLT逆向問(wèn)題求解

BLT的逆向問(wèn)題是從式（3）中求解未知的生物發(fā)光源S。然而，如前所述，考慮到：（1）光子在生物組織中的傳播是高度散射的；（2）用于描述光子在成像對(duì)象體內(nèi)傳播的數(shù)學(xué)模型（前向模型）是近似的。因此，對(duì)比XCT或MRI成像，BLT的重建問(wèn)題有更大的病態(tài)性。此外，測(cè)量數(shù)據(jù)Φ通常包含噪聲，這進(jìn)一步復(fù)雜了BLT的成像問(wèn)題。于是，直接對(duì)式（3）進(jìn)行求解不太可行?？紤]到在大多數(shù)生物學(xué)應(yīng)用中，生物發(fā)光源在成像體內(nèi)的分布通常具有稀疏性或可稀疏性。為了提高成像的空間分辨率，在本文中我們擬基于稀疏理論進(jìn)行求解[8-10]：

式（7）中，ε用于描述測(cè)量數(shù)據(jù)包含的噪聲水平；S是待重建的納米熒光團(tuán)的密度。

2 仿真實(shí)驗(yàn)設(shè)置及實(shí)驗(yàn)結(jié)果

2.1 BLT成像系統(tǒng)

在本文中，我們所模擬的成像系統(tǒng)是多角度、非接觸的BLT成像系統(tǒng)（multi-view non-contact BLT imagingsystem）。該系統(tǒng)主要包括光子探測(cè)器（CCD）、旋轉(zhuǎn)臺(tái)及計(jì)算機(jī)等設(shè)備，如圖1所示。

圖1 多角度、非接觸BLT成像系統(tǒng)

2.2 仿真實(shí)驗(yàn)設(shè)置

在數(shù)字仿真實(shí)驗(yàn)中，一個(gè)直徑3.0 cm、高5.0 cm的圓柱體被使用作為仿體，如圖2所示。仿體內(nèi)的光學(xué)參數(shù)被配置為 μa=0.3 cm-1、μ′s=10.0 cm-1。這與小鼠組織的平均光學(xué)參數(shù)類似。此外，1個(gè)直徑0.4 cm、高0.6 cm的小圓柱體被浸入在仿體中，用以模擬生物發(fā)光源。

圖2 仿真實(shí)驗(yàn)參數(shù)設(shè)置示意圖

直徑3.0 cm、高5.0 cm的圓柱體被使用作為仿體；仿體內(nèi)的光學(xué)參數(shù)為 μa=0.3 cm-1、μs′=10.0 cm-1；直徑0.4 cm、高0.6 cm的圓柱體被使用作為生物發(fā)光源。

為了提高BLT成像的空間分辨率，在BLT成像過(guò)程中，我們總計(jì)獲取了6張自發(fā)光投影圖像，等間隔 60°。

2.3 實(shí)驗(yàn)結(jié)果

對(duì)于BLT重建，首先，圖2中的圓柱仿體被離散為5 675個(gè)節(jié)點(diǎn)及29 317個(gè)四面體單元，如圖3所示。其次，重建被執(zhí)行。在重建過(guò)程中，我們使用稀疏重建算法求解式（7）。對(duì)于每個(gè)投影，檢測(cè)點(diǎn)分布在圓柱仿體140°×4.6 cm的范圍內(nèi)。BLT的重建結(jié)果如圖4所示。

圖3 圓柱形仿體的有限元離散結(jié)果

圖4 BLT重建結(jié)果

為了定量評(píng)估BLT重建性能，我們計(jì)算了重建的定位誤差，見(jiàn)表1。這里，定位誤差被定義為：實(shí)際生物發(fā)光源的中心與BLT重建發(fā)光源的中心最大值之間的距離。

表1 BLT重建的定位誤差

實(shí)驗(yàn)結(jié)果顯示，基于多角度、非接觸的BLT成像系統(tǒng)，結(jié)合稀疏重建算法，我們能夠在較少（6個(gè)）的角度下，獲取較好的BLT重建結(jié)果，如圖4及表1所示。

3 結(jié)語(yǔ)

基于多角度、非接觸的BLT成像系統(tǒng)，結(jié)合稀疏重建算法，本文解析了生物發(fā)光源在成像體內(nèi)的三維分布信息。數(shù)值仿真實(shí)驗(yàn)結(jié)果顯示，定位誤差小于1mm。在今后的研究中，我們將重點(diǎn)搭建BLT成像系統(tǒng)，進(jìn)行物體仿體及在體小動(dòng)物實(shí)驗(yàn)。

[1] Ntziachristos V，Ripoll J，WANG Lyu，etal.Looking and listening to light：the evolution of whole-body photonic imaging[J].Nature Biotechnology，2005，23（3）：313-320.

[2] WANG Ge，CONG Wen-xiang，Durairaj K，et al.In vivo mouse studies with bioluminescence tomography [J].Optics Express，2006，14（17）：7 801-7 809.

[3] Lv Y，Tian J，Cong W，et al.A multilevel adaptive finite element algorithm for bioluminescence tomography[J].Optics Express，2006，14（18）： 8 211-8 223.

[4] Gu X，ZhangQ，Larcom L，etal.Three-dimensional bioluminescence tomography with model-based reconstruction[J].Optics Express，2004，12（17）： 3 996-4 000.

[5] Lv Y，Tian J，Cong W，et al.Spectrally resolved bioluminescence tomography with adaptive finite element analysis：methodology and simulation[J].Physics inmedicine and biology，2007，52（15）：4 497.

[6] Schweiger M，Arridge S R，Hiraoka M，et al.The finite element method for the propagation of light in scattering media：boundary and source conditions[J].Medical Physics，1995，22（11）：1 779-1 792.

[7] Joshi A，Bangerth W，Sevick-Muraca E.Adaptive finite element based tomography for fluorescence optical imaging in tissue[J].Optics Express，2004，12（22）：5 402-5 417.

[8] Candès E J，Wakin M B.An introduction to compressive sampling[J].IEEE Signal Processing Magazine，2008，25（2）：21-30.

[9] Elad M，Aharon M.Image denoising via sparse and redundant representationsover learned dictionaries[J].IEEE Transactionsan Image Processing，2006，15（12）：3 736-3 745.

[10] Tropp JA，Gilbert A C.Signal recovery from random measurements via orthogonal matching pursuit[J].IEEE Transactionsan Information Theory，2007，53（12）：4 655-4 666.