基于多信息融合和棧式自編碼網(wǎng)絡的擴散光學層析方法

孫志龍，張杰，劉宗洋，高峰，3，張麗敏，3

（1 天津大學 國際工程師學院，天津 300072）

（2 天津大學 精密儀器與光電子工程學院，天津 300072）

（3 天津市生物醫(yī)學檢測技術與儀器重點實驗室，天津 300072）

0 引言

擴散光學層析（Diffuse Optical Tomography，DOT）作為一種新型的功能光學成像技術，其利用波長范圍在600～900 nm 的近紅外光入射生物組織體，根據(jù)在生物組織邊界得到的光流量，從而重構組織內(nèi)部光學特性的三維分布［1］。該技術具有無創(chuàng)、無電離輻射以及可提供生物組織生理學信息等優(yōu)勢，在乳腺腫瘤成像和腦功能成像等領域極具應用潛力［2-4］。在DOT 成像過程中，近紅外光在組織中傳播具有強散射特性，且邊界測量數(shù)據(jù)遠小于需要重構的組織體光學參數(shù)個數(shù)［5］。因此，為提高DOT 圖像重建質量，國內(nèi)外研究人員發(fā)展了多種圖像重建優(yōu)化算法。其中，常用的策略為非線性迭代算法中添加不同類型的正則化項，基本思想是在每次迭代過程中添加額外的懲罰項進行正則化約束，以抑制測量噪聲和減小模型誤差所帶來的影響，從而獲得不適定問題的穩(wěn)定近似解。典型的正則化方法包括L1 范數(shù)正則化和L2 范數(shù)正則化［6］?；贚2 范數(shù)的正則化方法可有效的重建DOT 圖像，但L2 范數(shù)具有平滑效應，對圖像邊緣恢復較差。L1 范數(shù)則可以保持邊緣，但易出現(xiàn)階梯效應，導致圖像產(chǎn)生虛假邊緣。此外，國內(nèi)外學者提出將擴散光學層析成像技術與現(xiàn)有的醫(yī)學成像模態(tài)結合發(fā)展成為多模態(tài)成像技術。如BROOKSBY B A 和鄧勇等分別把磁共振成像（Magnetic Resonance Imaging，MRI）和微型計算機斷層掃描系統(tǒng)提供的解剖結構圖像作為結構先驗信息，以Laplace 算子的形式引入到正則化方法中以提高DOT 圖像重建質量［7-8］。實驗結果均表明，多模態(tài)方法可以提高重構光學參數(shù)圖像的空間分辨率和量化度。但是，多模態(tài)成像方法首先需要對結構圖像進行分割，且圖像分割的準確性取決于操作者的經(jīng)驗水平。此外，這些方法無法從根本上避免DOT 逆向問題的不適定性，導致圖像空間分辨率低和量化度低，以及由于多次迭代導致重建耗時長。

近年來，人工神經(jīng)網(wǎng)絡以其抽象特征提取和識別、較強的模型構建和泛化能力被廣泛應用于圖像重建領域。FENG Jinchao 等提出了基于反向傳播神經(jīng)網(wǎng)絡（Back Propagation Neural Network，BPNN）的DOT圖像獲取方法，通過直接建立邊界光信息與吸收系數(shù)的對應關系，克服逆問題求解的不適定性，但BPNN 訓練速度較慢并且模型泛化能力不足［9］。YEDDER H B 等和YOO J 等提出基于卷積神經(jīng)網(wǎng)絡的DOT 圖像重建方法。神經(jīng)網(wǎng)絡泛化能力進一步提高，但是網(wǎng)絡參數(shù)較多，計算量更大，對于計算機性能要求高［10-11］。在各種神經(jīng)網(wǎng)絡中，棧式自編碼（Stacked Auto-Encoder，SAE）具有網(wǎng)絡模型簡單，需要調整的網(wǎng)絡參數(shù)少，網(wǎng)絡訓練速度快的優(yōu)點。此外，SAE 的無監(jiān)督預訓練過程，可為多層神經(jīng)網(wǎng)絡訓練找到一個良好的初始權重和偏差，有利于抑制輸入數(shù)據(jù)中的噪聲干擾，并在基于誤差反向傳播算法的微調過程中減少訓練時間［12］。目前，基于人工神經(jīng)網(wǎng)絡的DOT 圖像重建研究多對規(guī)則的圓形異質體進行圖像重建，但針對形狀欠規(guī)則異質體的圖像重建還沒有相關研究。

為此，本文面向DOT 技術在乳腺腫瘤成像中的實際應用，發(fā)展了基于多信息融合和SAE 網(wǎng)絡的DOT圖像重建方法。該方法將MRI 獲取的解剖圖像與DOT 中測量得到的邊界光流量經(jīng)過預處理后作為輸入值，經(jīng)過SAE 網(wǎng)絡模型，最終輸出DOT 光學參數(shù)圖像。為驗證該方法的有效性，本文進行了基于圓形和橢圓形的單異質體以及雙異質體模擬。并與無結構圖像先驗信息時的基于SAE 的DOT 成像結果進行了對比。此外，設置了平均絕對值誤差、均方誤差和定量重建率三個指標對圖像重建質量進行評估。

1 研究方法

1.1 光在生物組織中的傳播

對乳腺組織等較厚組織體，近紅外光在其內(nèi)部的傳輸規(guī)律可由擴散方程，即輻射傳輸方程的一階球諧展開近似描述。其中，穩(wěn)態(tài)擴散方程描述為

式中，D(r)=c/(3μs'(r))為擴散系數(shù)，c為組織內(nèi)光速，μ's(r)為約化散射系數(shù)，μa(r)為吸收系數(shù)，Φ(r，rs)為光子密度，δ(r，rs)表示在r處的迷向光源。為求解上述方程常結合Robin 邊界條件

1.2 基于SAE 神經(jīng)網(wǎng)絡的圖像重建

基于神經(jīng)網(wǎng)絡的DOT 圖像重建策略摒棄了求解反向問題的復雜和耗時過程。本文只有吸收系數(shù)μa是需要重建，約化散射系數(shù)為已知常數(shù)。以未融合先驗信息時基于SAE 網(wǎng)絡的DOT 圖像重建為例，網(wǎng)絡的輸入數(shù)據(jù)為邊界光流量Φi(i=1，2，3，...，m)，輸出數(shù)據(jù)集為相應的吸收系數(shù)μa，j(j=1，...，N)。其中，m為測量數(shù)據(jù)個數(shù)，N為有限元節(jié)點數(shù)。

1.2.1 棧式自編碼器網(wǎng)絡模型

本文使用的棧式自編碼人工神經(jīng)網(wǎng)絡包含兩個自編碼器（Autoencoder，AE），SAE 的訓練過程包括無監(jiān)督的預訓練及有監(jiān)督的微調兩個過程［14-15］，如圖1 所示。

圖1 SAE 訓練過程Fig.1 SAE training process

在預訓練階段，第一個AE1 器將原始輸入數(shù)據(jù)Φ映射到SAE 網(wǎng)絡的第一個隱藏層，通過最小化重構誤差獲得網(wǎng)絡的權重W1和偏置b1。隨后第二個AE2 器將第一個隱藏層的輸出結果作為其輸入，繼續(xù)最小化重構誤差獲得網(wǎng)絡權重W2和偏置b2。在無監(jiān)督預訓練結束后，添加全連接層作為SAE 的輸出層并對網(wǎng)絡整體的權重及偏置進行微調。預訓練得到的權重和偏置Wk，bk（k=1，2），被用作各神經(jīng)層間權重和偏置的初始值。SAE 網(wǎng)絡第二個隱藏層和輸出層之間的權重W3和偏置b3被隨機初始化。最后網(wǎng)絡基于誤差反向傳播原理

通過計算預測輸出值和期望輸出值間的均方誤差來最小化其損失函數(shù)，以求對網(wǎng)絡模型中的權重和偏置進行微調得到最佳的參數(shù)。其損失函數(shù)定義為

式中，μa，j和分別為第j個樣本的實際值和神經(jīng)網(wǎng)絡預測值，N為總的數(shù)據(jù)集樣本數(shù)。

朝堂之上，七八位大臣對皇上說，你再不把秀容月明調回來，秀容兵就擁他為王了?；噬险f，我要聽你們的，就是中了胡人的離間之計。

1.2.2 自編碼AE 網(wǎng)絡

AE 是一種由編碼器和解碼器組成的無監(jiān)督學習神經(jīng)網(wǎng)絡模型，如圖2［16］所示。編碼器作用是將輸入數(shù)據(jù)的特征提取至隱藏層，而解碼器將輸入數(shù)據(jù)從隱藏層中重構出來。AE 網(wǎng)絡的訓練過程使得輸入和輸出大小相同，主要原理描述為

圖2 自編碼器網(wǎng)絡結構圖Fig.2 Autoencoder network diagram

對于編碼過程，對AE 網(wǎng)絡的輸入：m個輸入樣本Φi（i=1，2，3，...，m）進行編碼處理，通過式（4）的編碼函數(shù)f對原始輸入Φ編碼，輸出到隱藏層向量hk(k=1，...，n)，n為隱藏層神經(jīng)元個數(shù)

式中，W是一個n×m的權重矩陣，b1∈Rn是輸入層與隱藏層間的偏置。

對于解碼過程，通過式（5）的解碼函數(shù)將隱藏層輸出h映射到輸出層，完成對輸入樣本Φ的重構

式中，WT是一個m×n的權重矩陣，b2∈Rm是隱藏層與輸出層間的偏置，Φ'表示對原始輸入數(shù)據(jù)Φ進行重構后的信息。AE 網(wǎng)絡最終訓練的結果為Φ'=(f(Φ))≈Φ，保留輸入數(shù)據(jù)的主要特征，降低數(shù)據(jù)中的噪聲影響。

為了得到模型的最優(yōu)參數(shù)W、WT、b1、b2，通過計算原始樣本輸入Φ和重構后的Φ'間均方誤差來最小化其損失函數(shù)。

1.2.3 網(wǎng)絡有監(jiān)督微調階段

SAE 網(wǎng)絡有監(jiān)督微調階段主要包括：輸入層（L0）、全連接隱藏層（L1）和預測輸出層（L2）。人工神經(jīng)網(wǎng)絡進行微調原理圖如圖3 所示。

圖3 網(wǎng)絡微調原理圖Fig.3 Schematic of network fine-tuning

使用已經(jīng)初始化權重W[l]和偏差變量b[l]，根據(jù)迭代結果，計算輸出層預測值和真實值μa的均方誤差。將誤差分攤給各層所有單元，從而獲得各層單元的誤差，進而調整各層之間的權重和偏差。本研究使用Adam 優(yōu)化器基于學習率來控制神經(jīng)網(wǎng)絡的調整速度，激活函數(shù)使用的是Sigmoid 函數(shù)。

1.3 評價指標

為了定量評估兩種方法在DOT 重構中的性能，計算了平均絕對誤差（Mean Absolute Error，MAE）、均方誤差（Mean Square Error，MSE）和定量重建率（Quantitativeness Ratio，QR）三個度量指標［17-19］

式中，μa，true(i)、μa，pre(i)分別是有限元節(jié)點第i處真實值以及神經(jīng)網(wǎng)絡預測的吸收系數(shù)。μˉa，pre為預測吸收系數(shù)的平均值，N為有限元節(jié)點總數(shù)。此處，MAE 和MSE 均用于評估重建圖像與真實圖像之間的誤差。前者，MAE 計算了平均誤差，MSE 同時計算了平均誤差和標準偏差。QR 計算重建目標體吸收系數(shù)與真實吸收系數(shù)的比值，用于重建方法的定量性評估。理論上MAE 和MSE 越低，說明重建結果越好，QR 越接近100%，表示重建圖像的吸收系數(shù)值越接近真實值。

2 結果分析

2.1 數(shù)據(jù)集的生成

在半徑40 mm 的二維圓形仿體中嵌入異質體來模擬乳腺腫瘤組織，將仿體離散為8 762 個三角形和4 526 個有限元節(jié)點。為了生成用于神經(jīng)網(wǎng)絡訓練的數(shù)據(jù)集，嵌入大小、位置和吸收系數(shù)不同異質體，其中，圓形異質體半徑分別為6 mm、8 mm 和10 mm，橢圓異質體的半長軸變化范圍為6～10 mm，半短軸變化范圍為4～8 mm，步長均為2 mm。利用隨機函數(shù)生成120 個異質體中心坐標；背景吸收系數(shù)設為0.003 4 mm-1，約化散射系數(shù)設為1 mm-1。異質體-背景的吸收系數(shù)的對比度（異質體吸收系數(shù)與背景吸收系數(shù)比值）為1.5～3.0，步長為0.5，約化散射系數(shù)與背景一致。雙異質體與單異質體仿體生成基本一致，區(qū)別在于雙異質體仿體內(nèi)部有兩個位置隨機分布的異質體，且兩異質體不重合。

設置共16 個同軸源探測對均勻布配在仿體的邊界，如圖4 所示。當某一光源激發(fā)時，除自身以及兩側鄰近的位置外，共11 個位置探測，一組得到16×11（176）個測量數(shù)據(jù)用于圖像重建。為了模擬真實測量值，在計算測量值中加入40 dB 的隨機高斯噪聲。

圖4 二維圓形仿真模型Fig.4 2D circular simulation model

為了使得數(shù)據(jù)特征更加明顯，需要對測量數(shù)據(jù)進行預處理

式中，Φbac是在仿體的光學參數(shù)均為背景值時測得一組邊界光流量，Φ是加入異質體后測量得到的邊界光流量，Φave為單個樣本的平均值，σ為單個樣本的標準差，經(jīng)過數(shù)據(jù)預處理之后數(shù)據(jù)Φsca的均值為0，方差為1。

與光學仿體相對應模擬生成MRI 圖像，如圖5 所示。生成灰度圖像的大小為51×51，總共2 601 個像素點，在輸入網(wǎng)絡之前做最大值最小值歸一化預處理。將經(jīng)過預處理的測量值（176）和灰度圖像信息（2 601）進行融合，作為一組數(shù)據(jù)的輸入（2 777），網(wǎng)絡輸出為吸收系數(shù)（4 526），共獲得5 670 組樣本數(shù)據(jù)，劃分為訓練集3 780 組，驗證集1 620 組，測試集270 組。

圖5 模擬異質體分布的MRI 圖像和對應的光學圖像Fig.5 MRI image and corresponding optical image that simulate target distributions

對于SAE，根據(jù)經(jīng)驗將隱藏層1 和隱藏層2 的神經(jīng)元數(shù)量分別設置為100 和300。此外，神經(jīng)網(wǎng)絡的學習率均設置為0.001，最小批次分別設置為20 組，預訓練次數(shù)為2 000，微調次數(shù)為300。

2.2 數(shù)值模擬實驗

為了評估和比較融合先驗信息前后SAE 網(wǎng)絡重建圖像的性能，進行了一系列的數(shù)值模擬，并對實驗結果進行了定性和定量的比較。主要對四種情況進行了對比。

1）隨機分布的圓形異質體，目標體半徑固定10 mm，吸收對比度分布為1.5、2.0、2.5 和3.0。圖6 的第一行是無先驗信息方法的重建圖像（No prior），第二行是融合先驗信息方法的重建圖像（Prior），第三行是沿著穿過仿體和異質體中心線的吸收系數(shù)剖線，黑色圓圈代表真實異質體。

從圖6 可以看出當異質體半徑均為10 mm，在不同吸收對比度下，添加先驗信息與否均可重建出異質體的位置和尺寸，且隨著吸收對比度的提高，重建圖像的質量均有所提高。當吸收對比度為1.5 時，無先驗信息方法的重建圖像中異質體邊緣輪廓模糊，偽影較多，且吸收系數(shù)過高估明顯，而融合先驗信息方法對重建結果有明顯的改善。

圖6 不同吸收對比度下，基于無先驗信息和融合先驗信息的圓形單異質體重建及剖線Fig.6 Reconstructed images and profiles of a single circular target with different absorption contrasts based on no prior and prior methods

為進一步定量比較兩種方法，表1 列出了重建圖像的MAE、MSE 和QR 值，且為了容易對比，給出了表1對應的條形圖（圖7）。從表1 和圖7 可以發(fā)現(xiàn)，融合先驗信息方法的重建圖像的MAE 均明顯小于無先驗信息方法所得的值，而MSE 僅在吸收比度為2 時略微大于無先驗信息方法對應的值。其原因是吸收對比度為2 時，無先驗信息方法獲取的圖像較真實值偏小，導致其標準偏差值偏小。此外，融合先驗信息方法重建圖像的QR 值也更加接近100%。定量分析表明融合先驗信息方法圖像重建質量更高。

圖7 表1 對應的條形圖Fig.7 Bar charts corresponding to Table 1

表1 圓形單異質體不同目標-背景吸收系數(shù)對比度重建圖像評估指標計算結果Table 1 Calculation results of evaluation indexes for a single circular target reconstruction with different absorption contrasts

2）隨機分布的圓形異質體，吸收對比度固定為3.0，半徑分別為6 mm、8 mm 和10 mm。圖8 是采用融合先驗信息和無先驗信息兩種方法獲得的重建圖像以及對應的異質體吸收系數(shù)剖線?？梢杂^察到兩種方法均可較準確地重建圖像，但是，小半徑時，融合先驗信息方法重建圖像的輪廓更加清晰，且量化度更優(yōu)。

圖8 不同半徑下基于無先驗信息和融合先驗信息的圓形單異質體重建及剖線Fig.8 Reconstructed images and profiles of a single circular target with different radii based on no prior and prior methods

表2 和圖9 分別是重建圖像的MAE、MSE 和QR 值以及對應的條形圖。進一步定量比較可以發(fā)現(xiàn)融合先驗信息方法重建圖像的MAE 在不同半徑下均優(yōu)于無先驗信息方法的重建值，表明重建圖像與真實圖像差距較小。與第一種情況相同，融合先驗信息方法重建圖像的MSE 值沒有顯示出絕對優(yōu)勢。從QR 值可以看出，無先驗信息方法的QR 值在小半徑和大半徑時分布顯示了低估和過高估計，而融合先驗信息方法的QR 值具有一定程度的過高估。綜合對比發(fā)現(xiàn)，吸收對比度與異質體尺寸相比對重建結果影響較大，因此下文對橢圓形異質體重建圖像結果評估時，僅研究了不同吸收對比度的情況。

表2 圓形單異質體相同吸收對比度不同半徑重建圖像評估指標計算結果Table 2 Calculation results of evaluation indexes for a single circular target reconstruction with different radii and the same absorption contrasts

圖9 表2 對應的條形圖Fig.9 Bar charts corresponding to Table 2

3）隨機分布的單個橢圓形異質體，吸收對比度分別為1.5、2.0、2.5 和3.0。圖10 是兩種方法獲得的重建圖像以及對應的吸收系數(shù)剖線。當異質體為橢圓形時，融合先驗信息方法的圖像重建結果在異質體位置、尺寸、形狀以及量化度方面均優(yōu)于無先驗信息方法的結果。尤其當吸收系數(shù)對比度較小時，無先驗信息方法重建的異質體中心位置有明顯的偏差且形狀有明顯的膨脹。在其它吸收系數(shù)對比度下，無先驗信息方法重建圖像存在膨脹且邊緣輪廓模糊的現(xiàn)象。此外，從剖線圖可以看出融合先驗信息方法重建圖像的量化度更接近真實值，而無先驗信息方法重建結果在低吸收對比度和高吸收對比度時分別表現(xiàn)出明顯的過高估和低估。

圖10 不同吸收對比度下，基于無先驗信息和融合先驗信息的橢圓形單異質體重建及剖線Fig.10 Reconstructed images and profiles of a single elliptic target with different absorption contrasts based on no prior and prior methods

為進一步定量比較兩種方法，單個橢圓異質體重建圖像的MAE、MSE 和QR 值比較見表3 以及相應的條形圖見圖11。可以看到融合先驗信息方法的三個評估指標均表現(xiàn)出明顯的優(yōu)勢。以吸收系數(shù)對比度1.5為例，融合先驗信息方法比無先驗信息方法重建圖像的MAE 和MSE 分別減少了40.1%和46.8%，且其QR值（101%）比無先驗信息方法的QR 值（130%）更接近100%。其它組數(shù)據(jù)也可以觀察到類似結果，即融合先驗信息方法重建圖像的MAE 和MSE 顯著降低，并且QR 值也更接近100%。這些結果表明融合先驗信息方法重建圖像有更高的精度。

圖11 表3 對應的條形圖Fig.11 Bar charts corresponding to Table 3

表3 橢圓形單異質體不同目標-背景吸收系數(shù)對比度重建圖像評估指標計算結果Table 3 Calculation results of evaluation indexes for a single elliptic target reconstruction with different absorption contrasts

4）隨機分布的雙橢圓異質體，吸收系數(shù)對比度分別為1.5、2.0、2.5 和3.0。圖12 中，兩種方法均可有效預測異質體的位置。比較兩種方法獲取的重建圖像可以看到，無先驗信息方法重建圖像的背景偽影較多，異質體的邊緣輪廓也不夠清晰，且除吸收對比度為2.5 外，同一仿體中兩個相同的異質體的重建結果在形狀、尺寸和量化度上有明顯的差異。進一步觀察吸收系數(shù)剖線可以看到融合先驗信息方法的剖線更接近真實值，而無先驗信息方法獲取的吸收系數(shù)在低吸收對比度時量化度較差，隨著吸收對比度增大吸收系數(shù)有明顯的過高估現(xiàn)象，而融合先驗信息方法在不同對吸收對比度情況下均與真實值接近。

表4 為圖12 各組圖像的MAE、MSE 和QR 值，并繪制對應條形圖13。同樣以吸收系數(shù)對比度1.5 為例，融合先驗信息方法重建的圖像比無先驗信息方法重建圖像的MAE 和MSE 分別減少了34.9%和53.7%，并且QR 值分別是98%和56%，融合先驗信息方法的QR 值更接近100%，表4 中的其它數(shù)據(jù)也可以觀察到類似結果?？傮w來說，可以看出融合了先驗信息使得重建圖像MAE 和MSE 都顯著減少，并且QR 也更加接近100%，即融合先驗信息可以得到更好的重建結果。

表4 雙橢圓異質體不同目標-背景吸收系數(shù)對比度重建圖像評估指標計算結果Table 4 Calculation results of evaluation indexes for double elliptic targets reconstruction with different absorption contrasts

圖12 不同吸收對比度下，基于無先驗信息和融合先驗信息的雙橢圓形異質體重建及剖線Fig.12 Reconstructed images and profiles of double elliptic targets with different absorption contrasts based on no prior and prior methods

圖13 表4 對應的條形圖Fig.13 Bar charts corresponding to Table 4

3 結論

為了提高DOT 的重建圖像質量，本文發(fā)展了一種融合先驗信息SAE 神經(jīng)網(wǎng)絡方法的DOT 圖像重建方法。為驗證該方法的有效性，對融合先驗信息與無先驗信息兩種方法的實驗結果進行了比較和定量評估。其中，具有不同吸收對比度和不同尺寸的圓形單異質體模擬結果表明，融合先驗信息方法重建結果與真實圖像更接近，且當吸收對比度比較小時，融合先驗信息方法對于重建圖像質量改善較為明顯。在異質體吸收對比度為1.5 時，與無先驗信息方法相比，MAE 降低62%，MSE 降低11%，QR 由139%降低為107%，更接近100%。需要說明的是，當圓形異質體吸收對比度大于1.5 且半徑較大時，理論上更容易獲得異質體圖像，因此兩種方法均取得了較好的圖像重建質量。對于包含形狀欠規(guī)則的單個和兩個的異質體進行圖像重建，融合先驗信息方法能夠讓背景清晰、異質體尺寸和位置準確、異質體吸收系數(shù)更接近真實值。尤其當異質體吸收系數(shù)比較小，融合先驗信息可以減小預測誤差。定量分析表明，融合先驗信息的方法使得重建圖像的MAE、MSE 都顯著減小，與無先驗信息方法相比MAE 至少降低8%，MSE 至少降低5%，QR 更接近100%。綜合評估表明，融合先驗信息能有效提高成像精度和改善重建圖像的質量。