包星星，趙璨，包書圣，饒家聲，3，楊朝陽，李曉光，

1.北京航空航天大學(xué)生物與醫(yī)學(xué)工程學(xué)院/生物材料與神經(jīng)再生北京市重點實驗室，北京100083；2.中國康復(fù)科學(xué)所康復(fù)工程研究所，北京100068；3.北京航空航天大學(xué)生物醫(yī)學(xué)工程高精尖創(chuàng)新中心，北京100083；4.首都醫(yī)科大學(xué)基礎(chǔ)醫(yī)學(xué)院/神經(jīng)生物學(xué)系，北京100069

前言

脊髓損傷是一種高致殘率的中樞神經(jīng)性疾病，嚴(yán)重影響患者的正常生活，并給家庭帶來沉重的經(jīng)濟(jì)負(fù)擔(dān)［1］。在脊髓損傷急性期階段，通過創(chuàng)造適宜的微環(huán)境可以實現(xiàn)內(nèi)源性神經(jīng)干細(xì)胞在脊髓損傷局部向神經(jīng)元的分化，從而實現(xiàn)神經(jīng)組織的修復(fù)［2-5］。相比于急性期的研究而言，在脊髓損傷慢性期開展修復(fù)的研究相對較少，其中最主要的原因是膠質(zhì)瘢痕的存在。膠質(zhì)瘢痕的形成和穩(wěn)定是脊髓損傷慢性期的病理生理學(xué)特征之一，在損傷局部所形成的致密膠質(zhì)瘢痕會分泌硫酸軟骨素蛋白聚糖，連同其所包繞而成的囊腔共同構(gòu)成物理和化學(xué)雙重屏障，嚴(yán)重阻礙了軸突的再生與髓鞘化，從而影響神經(jīng)功能的恢復(fù)［6-12］。

對脊髓損傷慢性期的相關(guān)研究而言，膠質(zhì)瘢痕及其所包繞的囊腔構(gòu)成的脊髓壞死區(qū)域是脊髓組織修復(fù)過程中無法繞開的部位。針對特定時間點的橫斷面研究可以使用免疫組化的方法，通過病理切片來評估脊髓壞死區(qū)域的位置、尺寸、邊界范圍等信息。然而，針對長時間的縱向連續(xù)觀察研究而言，病理切片獲取壞死區(qū)域信息的方法存在著明顯的局限性，單個受試對象無法既提供不同時間節(jié)點的觀察結(jié)果同時又給出相應(yīng)的病理切片。為了解決這一矛盾，一些研究采用磁共振成像（Magnetic ResonanceImaging,MRI）的方法來進(jìn)行脊髓損傷后壞死區(qū)域的檢測研究，保障了縱向?qū)嶒灥捻樌M(jìn)行［13-16］。

在上述研究中無論采用何種成像模態(tài)，MRI的脊髓壞死區(qū)域大多依賴于資深放射科醫(yī)生的定性觀察或手動識別［17］?？紤]到MRI的部分容積效應(yīng)，壞死組織在MRI圖像中通常具有邊緣模糊的特點，這使得壞死區(qū)域的準(zhǔn)確識別存在一定難度。此外，不同損傷位置、程度、范圍、方式等導(dǎo)致?lián)p傷區(qū)域各異，由此形成的脊髓壞死區(qū)域也形態(tài)多樣且不規(guī)則，更增添了手動識別的困難度［17］。因此，基于操作者手動勾勒得到的脊髓壞死區(qū)域往往不具有可重復(fù)性和可驗證性，阻礙了對膠質(zhì)瘢痕及囊腔區(qū)域的準(zhǔn)確評估。

目前已知病理切片是識別脊髓損傷慢性期膠質(zhì)瘢痕及其包繞的囊腔區(qū)域的“金標(biāo)準(zhǔn)”，因此，應(yīng)用病理切片結(jié)果來指導(dǎo)MRI圖像壞死區(qū)域的識別無疑是一個較好的思路。要實現(xiàn)上述工作，首先應(yīng)當(dāng)實現(xiàn)病理切片與MRI圖像的對應(yīng)。本研究基于病理切片與MRI圖像的壞死區(qū)域空間相對位置具有一致性這一特點，通過一系列圖像處理建立了兩者之間的映射關(guān)系，并驗證了映射后獲得的MRI圖像實際壞死區(qū)域與“金標(biāo)準(zhǔn)”的膠質(zhì)瘢痕和囊腔范圍之間的吻合度，為將來計算機(jī)自動識別與分割MRI壞死區(qū)域奠定了基礎(chǔ)。

1 材料與方法

1.1 動物模型制備

本研究制備了9 只成年雌性Wistar 大鼠的脊髓砸傷模型，砸傷位置為脊髓的T8 節(jié)段，使用NYU 脊髓打擊器（NewYork University,New York,USA）砸傷脊髓，其中撞擊棒重10 g，并從50 mm 的高度自由下落完成施力。術(shù)后3個月內(nèi)，實驗大鼠飼養(yǎng)于恒溫恒濕的環(huán)境中，確保充足的食物和水。術(shù)后連續(xù)5 d 注射青霉素（10 萬單元/d/次），每天人工護(hù)理排尿直到實驗動物恢復(fù)自主排尿功能。實驗方案經(jīng)過北京航空航天大學(xué)生物與醫(yī)學(xué)倫理委員會批準(zhǔn)（批準(zhǔn)號：BM20180047）。

1.2 數(shù)據(jù)采集

本研究在大鼠脊髓損傷慢性期（術(shù)后3 個月）采集了其脊髓矢狀面的T1加權(quán)（T1-weighted,T1W）和T2加權(quán)（T2-weighted,T2W）兩種在體MRI結(jié)構(gòu)像數(shù)據(jù)，所有數(shù)據(jù)均采自首都醫(yī)科大學(xué)7.0T 小動物磁共振掃描儀（Bruker Biospec, Karlsruhe, Germany）。壞死區(qū)域在T1W 圖像中主要表現(xiàn)為低、中等信號強(qiáng)度，在T2W 圖像中主要表現(xiàn)為高信號強(qiáng)度。不同模態(tài)下各異的影像學(xué)表現(xiàn)為脊髓壞死區(qū)域的定位提供了豐富的信息。

1.3 病理切片制備

脊髓病理切片是確定脊髓壞死區(qū)域尺寸和范圍的“金標(biāo)準(zhǔn)”，是建立脊髓病理切片到MRI圖像映射關(guān)系的基礎(chǔ)。本研究對實驗動物脊髓矢狀位的病理組織切片（10 μm）進(jìn)行了膠質(zhì)纖維酸性蛋白（Glial Fibrillary Acidic Protein, GFAP）染色，并通過普通光學(xué)顯微鏡獲取了GFAP染色后的病理切片圖像。圖1的病理切片展示了脊髓損傷慢性期神經(jīng)元壞死液化形成的空洞（即囊腔），囊腔壁由致密的膠質(zhì)瘢痕構(gòu)成。

1.4 方法

圖1 為本研究所提出的方法示意圖。為建立脊髓病理切片到MRI圖像的映射關(guān)系，本研究所構(gòu)建的方法包括以下幾部分工作：（1）預(yù)處理：雙邊濾波去噪以及T1W 和T2W 圖像的剛體配準(zhǔn)；（2）脊髓分離及拉直：利用Snake模型實現(xiàn)MRI圖像中脊髓的粗分割，剔除冗余信息，并利用雙三次插值法對脊髓進(jìn)行拉直；（3）圖像融合：基于雙樹復(fù)小波變換（Dual Tree-Complex Wavelet Transform,DT-CWT）對T1W 和T2W圖像進(jìn)行融合；（4）脊髓MRI圖像壞死區(qū)域的確定：基于病理切片與MRI圖像的脊髓壞死區(qū)域空間相對位置具有一致性這一特點，建立脊髓病理切片到MRI圖像的映射關(guān)系，以實現(xiàn)MRI圖像中脊髓壞死區(qū)域的標(biāo)注。

圖1 脊髓壞死區(qū)域鑒別框架示意圖Fig.1 Schematic diagram of spinal cord necrosis area differentiation framework

1.4.1 預(yù)處理MRI圖像在獲取和生成的過程中會混雜一定的噪聲，這無疑會給圖像的后續(xù)處理增加干擾信息。雙邊濾波方法在進(jìn)行噪聲平滑的同時能較好地保留圖像的邊緣信息［18］，因此本研究將該方法應(yīng)用于脊髓MRI多模態(tài)圖像的去噪，取得了較好的去噪效果。

T1W 和T2W 兩種序列的采集存在時間差，生理因素及機(jī)器震動等均有可能導(dǎo)致組織發(fā)生空間位移，從而引起圖像的不匹配。本文以T1W 圖像為參考圖像，對T2W 圖像進(jìn)行剛性變換以改變T2W 的空間位置，實現(xiàn)T1W 和T2W 圖像的配準(zhǔn)，以消除兩種序列的空間位置差異。

1.4.2 脊髓分離及拉直原始脊髓MRI圖像不僅存在脊髓，還包括椎體、椎間盤、肌肉以及皮膚等其他組織，為去除冗余信息，本文采用了Snake 模型從MRI圖像中提取脊髓。Snake 模型在最小化能量泛函的動力下驅(qū)使曲線不斷靠近目標(biāo)物體邊緣以獲取目標(biāo)物體輪廓［19］。初始輪廓曲線是通過將T1W 圖像進(jìn)行對數(shù)變換后手動獲取的，根據(jù)曲線上點v在圖像中的空間位置可確定Snake模型的能量函數(shù)：

其中，Etotal為曲線演化過程中的總能量；第一項為內(nèi)部能量，保持著曲線的連續(xù)性和平滑性；第二項則為外部能量，表示灰度值或梯度等圖像特征信息對曲線的能量控制，本研究取的是v點處圖像的梯度信息：

由于在實際制備脊髓病理切片的過程中可能會產(chǎn)生任意角度的彎曲，為準(zhǔn)確獲取病理切片與MRI圖像的空間對應(yīng)關(guān)系，本研究采用插值效果最佳的雙三次插值法對病理切片和MRI圖像中的脊髓均進(jìn)行了拉直處理。脊髓拉直圖像中某位置的強(qiáng)度值，對應(yīng)于距離原圖像同位置像素點最近的16個像素強(qiáng)度值的加權(quán)平均值，各權(quán)重的計算如下［20］：

其中，x表示原圖像中距離給定數(shù)據(jù)點最近的16個像素點到目標(biāo)圖像對應(yīng)點的距離，S(x)表示16 個像素點對應(yīng)權(quán)重。

1.4.3 圖像融合為盡可能全面地獲取脊髓壞死區(qū)域的影像信息，本研究充分考慮了脊髓T1W 和T2W 圖像的信息，利用DT-CWT 方法將T1W 和T2W 圖像轉(zhuǎn)換到了頻率域，在頻域?qū)崿F(xiàn)兩種圖像的融合。DTCWT 不僅保留了小波變換的時頻分析特性，而且具備平移不變性、方向分析能力良好、數(shù)據(jù)冗余有限、計算效率較高以及重構(gòu)效率好等特點，成為近年來應(yīng)用于圖像融合的熱門方法之一［21］。圖2 為DT-CWT 的分解示意圖，通過樹A 和樹B 并行獲取實部和虛部的小波變換系數(shù)，從而形成雙樹結(jié)構(gòu)。每級分解通過低、高通濾波及間隔采樣，可獲取兩個低頻子圖像以及±15°、±45°和±75° 這6個方向的高頻子圖像1, 2,…, 6)。本文首先將脊髓拉直的T1W 和T2W圖像進(jìn)行DT-CWT 五級分解，然后對同方向高頻系數(shù)的融合取小波系數(shù)模的最大值，以使融合圖像突出細(xì)節(jié)信息。低頻系數(shù)的融合取其小波系數(shù)的均值，再經(jīng)過逆變換可得到T1W 和T2W 融合圖像（圖3）。

圖2 DT-CWT變換Fig.2 DT-CWT

圖3 融合流程圖Fig.3 Fusion flow chart

1.4.4 基于病理“金標(biāo)準(zhǔn)”標(biāo)定脊髓MRI圖像壞死區(qū)域的方法本研究的目的旨在構(gòu)建病理切片與MRI圖像的映射關(guān)系，從而依據(jù)病理切片“金標(biāo)準(zhǔn)”來客觀地標(biāo)記脊髓MRI圖像中的壞死區(qū)域。圖4 為映射關(guān)系構(gòu)建的示意圖，圖4a 為拉直的脊髓病理圖像，包含壞死區(qū)域的脊髓段已用黑色框標(biāo)記并對其進(jìn)行裁剪。病理圖像中壞死組織的邊界清晰，可手動獲取壞死區(qū)域并進(jìn)行二值化，將損傷區(qū)域賦值為1，其余賦值為0，得到圖4b 所示二值圖像。針對一張MRI圖像（1 mm 層厚）對應(yīng)多張病理切片（10 μm 層厚）的情況，本研究分別對多張病理切片進(jìn)行二值化處理，隨后分別與MRI圖像建立映射關(guān)系，最后再求取總的映射關(guān)系式，以盡量減小瘢痕在MRI圖像上被標(biāo)記的誤差。映射關(guān)系的建立是依據(jù)被拉直的病理切片和MRI圖像中壞死組織邊界到脊髓左側(cè)邊緣的水平距離占整個脊髓寬度的比例相等這一特點來實現(xiàn)的（研究中所用的距離、寬度和高度都是指相應(yīng)方向的像素個數(shù)），具體公式如下：

其中，s為對應(yīng)的病理切片數(shù)。對于MRI圖像中脊髓的裁剪，首先根據(jù)脊髓的影像學(xué)表現(xiàn)，可確定圖4c脊髓寬度n以及大致的損傷中心位點，再計算二值化圖像的高度、寬度比例的均值可獲取應(yīng)裁剪的位置和高度m：

位置的標(biāo)記分兩種情況：若在病理圖像中，瘢痕區(qū)域的數(shù)量等于1，則選取損傷區(qū)域y 軸方向上1/8、1/4、1/2、3/4 和7/8 這5 個位置，每個位置一對像素點，共10個像素點（分別對應(yīng)同一y軸水平損傷邊緣的兩個位置）；若瘢痕的區(qū)域數(shù)量大于1，則增加3/8、5/8、1/3 和2/3 這4 個位置，共18 個像素點，最后連接這些像素點以勾勒整個瘢痕邊界。

圖4 病理圖像與MRI圖像映射關(guān)系構(gòu)建示意圖Fig.4 Schematic diagram of mapping relationship between pathological image and MRI image

1.4.5 評估方法傳統(tǒng)評估方法在本研究中難以適用，這主要有兩個方面的原因：（1）傳統(tǒng)評估方法不適用于圖像尺寸不一的配準(zhǔn)情況；（2）本研究使用了二值化圖像和相對位置關(guān)系兩個要素，因此難以提供額外的圖像紋理、灰度以及空間位置等信息供傳統(tǒng)評估方法所使用。根據(jù)本研究的具體情況，筆者建立了一種基于形狀參數(shù)的評估方法，包括圓度和偏心率兩個特征，計算公式如下：

其中，PR 和PE 分別表示病理圖像中膠質(zhì)瘢痕及囊腔區(qū)域的圓度和偏心率，IR 和IE 則分別為MRI圖像中壞死區(qū)域的圓度和偏心率。式中S為瘢痕面積，C為瘢痕周長，c為瘢痕區(qū)域焦點間的距離，a為長軸長度。

本研究通過比較病理切片區(qū)域和MRI圖像中兩個區(qū)域的圓度和偏心率特征來反映MRI圖像壞死區(qū)域標(biāo)記的準(zhǔn)確性，計算公式如下：

其中，RR、ER 分別為病理圖像中膠質(zhì)瘢痕及囊腔區(qū)域與MRI圖像中壞死區(qū)域的圓度、偏心率的比值。在計算時，將較小的值設(shè)置為分子，較大的值設(shè)為分母，MRI圖像中所標(biāo)記的實際壞死區(qū)域與病理切片所給出的“金標(biāo)準(zhǔn)”越相似，則RR 和ER 的值越接近1。當(dāng)兩者完全相同時，計算結(jié)果即為1。

2 實驗結(jié)果與分析

2.1 脊髓分離結(jié)果

Snake 模型能量函數(shù)的外部能量大小與曲線上點的局部梯度大小呈負(fù)相關(guān)，圖5a 中，脊髓正常組織與瘢痕組織的信號突變使得曲線的演化方向偏離脊髓邊緣（圖5b），導(dǎo)致最后獲取的脊髓不完整，這將影響后續(xù)病理圖像與MRI圖像映射關(guān)系的準(zhǔn)確建立。而圖5c 中，T1W 圖像通過對數(shù)變換后，降低了脊髓正常組織與瘢痕強(qiáng)度的對比度，最終可獲取完整的脊髓（圖5d）。

2.2 圖像融合結(jié)果

表1 對T1W 和T2W 圖像在不同DT-CWT 分解級數(shù)下融合效果進(jìn)行了比較，評估的參數(shù)包括互信息（Mutual Information,MI）、熵（Entropy）以及峰值信噪比（Peak Signal to Noise Ratio,PSNR）。綜合來看，在五級分解的情況下可得到最優(yōu)的融合效果。由此，本研究采取DT-CWT 五級分解再融合進(jìn)行后續(xù)的實驗。

2.3 MRI圖像中實際壞死區(qū)域標(biāo)注的結(jié)果

通過建立病理切片與MRI圖像的映射，本研究實現(xiàn)了在MRI圖像上準(zhǔn)確標(biāo)記壞死區(qū)域（圖6）。

另外，本研究對病理圖像中的實際膠質(zhì)瘢痕及其包繞的囊腔區(qū)域與MRI圖像中標(biāo)記的壞死區(qū)域在形狀上的相似程度進(jìn)行了量化分析，結(jié)果發(fā)現(xiàn)，這兩種區(qū)域的相似度評估結(jié)果在圓度上為0.93±0.03，在偏心率上達(dá)到了0.97±0.02，表明了MRI圖像中所標(biāo)記的區(qū)域與病理“金標(biāo)準(zhǔn)”具有較高的吻合度。

圖5 脊髓粗分割Fig.5 Coarse segmentation of spinal cord

表1 基于不同DT-CWT分解級數(shù)下融合效果統(tǒng)計Tab.1 Statistics offusion effect based on different DT-CWT decomposition levels

圖6 MRI圖像中脊髓壞死區(qū)域的標(biāo)記結(jié)果Fig.6 The results of MRI in the necrotic area of spinal cord

3 總結(jié)與討論

本研究提出一種應(yīng)用病理切片來確定脊髓損傷慢性期MRI圖像中實際壞死區(qū)域的方法，基于病理切片膠質(zhì)瘢痕及其所包繞的囊腔區(qū)域與MRI圖像的脊髓壞死區(qū)域具有空間相對位置一致性的特點，建立脊髓病理切片到MRI圖像的映射關(guān)系，據(jù)此實現(xiàn)MRI圖像中脊髓壞死區(qū)域的客觀、準(zhǔn)確標(biāo)記。利用形狀信息的評估方法顯示MRI標(biāo)記的壞死區(qū)域與病理學(xué)“金標(biāo)準(zhǔn)”的膠質(zhì)瘢痕和囊腔區(qū)域之間高度相似，為今后應(yīng)用圖像識別技術(shù)構(gòu)建計算機(jī)自動輔助分割脊髓損傷慢性期壞死區(qū)域的方法奠定基礎(chǔ)。

后續(xù)的研究可以在兩個方面繼續(xù)深入開展：（1）進(jìn)一步改進(jìn)映射方法，采用加權(quán)求和的策略，結(jié)合主觀經(jīng)驗，將那些與MRI圖像中壞死區(qū)域的大小、形狀及位置等較為相似的病理切片樣本數(shù)據(jù)賦予更大的權(quán)重值，以獲取帶有加權(quán)系數(shù)的映射關(guān)系；（2）在已建立映射關(guān)系的基礎(chǔ)上，評價MRI中已標(biāo)記壞死區(qū)域的三維結(jié)構(gòu)形態(tài)與病理學(xué)“金標(biāo)準(zhǔn)”三維結(jié)構(gòu)形態(tài)之間的差異，從而更全面地評估映射關(guān)系的準(zhǔn)確性。