擴散磁共振圖像的神經(jīng)纖維追蹤算法研究綜述

葉偉紅，王遠軍

(上海理工大學(xué) 醫(yī)學(xué)影像工程研究所，上海 200093)

1 引 言

擴散磁共振成像(Diffusion Magnetic Resonance Imaging，dMRI)根據(jù)水分子的各向異性來進行神經(jīng)纖維微結(jié)構(gòu)成像，它是目前唯一不侵入活體研究纖維結(jié)構(gòu)特征的技術(shù).反映水分子擴散狀況可以由總體平均擴散密度函數(shù)(Ensemble Average Propagator，EAP)與方向分布函數(shù)(Orientation Distribution Function，ODF)等概率密度函數(shù)[1]來表示，擴散張量成像[2](Diffusion Tensor Imaging，DTI)是目前最廣泛使用的估計ODF的方法.但是DTI的張量模型不能全面檢測出纖維彎曲、分叉、交叉等復(fù)雜的纖維結(jié)構(gòu)，所以在DTI的基礎(chǔ)上又提出了高角分辨率成像[3](High Angular Resolution Diffusion Imaging，HARDI)以得到更精確的纖維方向分布函數(shù).神經(jīng)纖維追蹤技術(shù)分為兩部分：利用ODF描繪出局部纖維結(jié)構(gòu)方向信息是神經(jīng)纖維追蹤的第1步，使用不同的神經(jīng)纖維追蹤算法整合纖維局部信息勾勒出完整的全腦白質(zhì)纖維是第2步，也是關(guān)鍵的一步.

神經(jīng)纖維追蹤算法分為兩大類，即：局部型追蹤方法和全局型追蹤方法:1)局部型追蹤方法又可以分為確定型追蹤方法與概率型追蹤方法.確定型追蹤算法的主要思想是從種子點出發(fā)，根據(jù)各體素的纖維方向分布信息給出的纖維主方向，追蹤出一條纖維路徑.此算法易于追蹤操作，但由于噪聲、局部模型估計ODF誤差、追蹤方式簡明單一等因素會導(dǎo)致追蹤出的部分纖維軌跡有所偏差以及全腦纖維結(jié)構(gòu)不全面.概率型追蹤算法的主要思想是從每個種子點出發(fā)多次執(zhí)行確定型追蹤算法生成大量的可能軌跡的集合，再從這個集合中根據(jù)遍歷每個體素的軌跡數(shù)量的高低來選取高概率可能性的纖維路徑.此算法可以提高全腦纖維結(jié)構(gòu)的完整度，但也產(chǎn)生了大量的偽纖維；2)全局型追蹤方法通常是基于圖形模型或能量最小化的方法，它的主要思想是創(chuàng)建一組纖維配置，通過最優(yōu)化代價函數(shù)來找到最符合測量數(shù)據(jù)集的纖維結(jié)構(gòu)來重建所有纖維路徑.此算法能夠避免局部型追蹤算法所產(chǎn)生的誤差累積，但是全局型追蹤算法依賴于隨機優(yōu)化的過程，所以不能保證收斂到全局最優(yōu)解，并且因為追蹤時間長而難以真正地運用于實際.

纖維追蹤是一項復(fù)雜的工作，它涉及許多步驟：局部纖維方向建模、張量插值、種子點的選取、纖維傳播以及終止準則等.每一個步驟都會產(chǎn)生潛在的錯誤，研究者們在這些步驟上做出了許多優(yōu)化來減少誤差的產(chǎn)生.雖然纖維追蹤技術(shù)在近幾年來有了較大的發(fā)展，但是目前仍沒有一種追蹤算法能夠追蹤出完美的纖維體結(jié)構(gòu)，在使用纖維軌跡量化連接時不可避免的會產(chǎn)生假陽體.

總體而言，本文主要綜述了兩大類的神經(jīng)纖維追蹤算法的研究進展，闡明這幾者追蹤算法的優(yōu)點以及存在的局限性，

圖1 神經(jīng)纖維追蹤算法總體框架Fig.1 General framework of neural fiber tracking

然后在此基礎(chǔ)上介紹了在神經(jīng)纖維追蹤過程中作出優(yōu)化的具體方面，具體框架見圖1，最后對神經(jīng)纖維追蹤算法的未來發(fā)展趨勢提出建議.

2 神經(jīng)纖維追蹤算法

2.1 局部型追蹤算法分析

在局部型追蹤算法中，根據(jù)原始成像體素的信號進行纖維局部結(jié)構(gòu)方向估計，通過整合局部纖維結(jié)構(gòu)信息追蹤出全腦長纖維結(jié)構(gòu)，依據(jù)局部纖維結(jié)構(gòu)信息類型和纖維傳播的方向可分為確定型追蹤算法和概率型追蹤算法.

最早的確定型追蹤算法是由Basser提出的流線型追蹤[4](Streamlines Tracking，STT)，他認為在三維空間中，腦白質(zhì)可以看作是光滑的曲線集，那么單條纖維追蹤便是一個微分方程問題，由公式(1)所示：

(1)

其中s(t)表示纖維曲線，t表示體素內(nèi)纖維曲線上的任意點，v(t)表示局部纖維方向信息模型.次年，Mori等提出纖維聯(lián)絡(luò)連續(xù)跟蹤[5](Fiber Assignment by Continuous Tracking，F(xiàn)ACT)，他根據(jù)各個體素的方向信息利用數(shù)值積分把纖維段連接成一條相連折線.這兩者對于纖維追蹤的思想成為后續(xù)研究纖維追蹤確定型算法的奠基石，STT算法也是目前最流行、應(yīng)用最廣泛的FT技術(shù).具體使用確定型追蹤算法追蹤的流程框架由表1所示.而后，張量線算法[6](Tensorlines)提出依據(jù)上一步傳播方向與當前體素纖維結(jié)構(gòu)信息結(jié)合綜合給出這一步的纖維主方向、張量偏轉(zhuǎn)算法[7](Tensor Deflection，TEND)在張量線的基礎(chǔ)上提出對纖維傳播方向均給出偏轉(zhuǎn)、向量選擇跟蹤算法[8](Vector Criterion Tracking，VCT)提出在纖維跟蹤過程中排除角度突變的體元有效避免部分容積效應(yīng)、多纖維聯(lián)絡(luò)連續(xù)跟蹤[9](Multiple Fiber Assignment by Continuous Tracking，mFACT)提出FACT與HARDI結(jié)合考慮多方向的纖維取向、在FACT中引入彌散熵[10]使得算法有更好的魯棒性、結(jié)合拓撲學(xué)去除錯誤纖維連接[11]等.確定型追蹤算法簡單、易于實現(xiàn)，它所追蹤出來的單條纖維準確度也較高，但由于確定型追蹤算法給出的是單個主特征的纖維方向，無法表示出纖維的分叉、交叉等復(fù)雜的結(jié)構(gòu)狀態(tài)[12]，而人類大腦70%～90%的白質(zhì)是由兩種或兩種以上的纖維群在同一區(qū)域內(nèi)交叉構(gòu)成的[13].

表1 基于確定型纖維追蹤的流程框架Table 1 Process framework based on deterministic fiber tracking

概率型追蹤算法為了改進這一缺陷便不只采取纖維主特征方向，而是在每個種子點給出多個可能的纖維走向，用概率密度函數(shù)(Probability Density Functions，PDF)評估各體素內(nèi)最可能的纖維束走行方向，追蹤方向從連續(xù)的ODF里采樣獲得.多數(shù)纖維追蹤算法只計算先驗分布來獲取當前體素包含的纖維結(jié)構(gòu)信息而導(dǎo)致追蹤過早的停止在白質(zhì)邊界.于是Friman提出一種基于貝葉斯理論的概率型追蹤算法[14]，它根據(jù)纖維傳播的上一步方向信息引入貝葉斯概率并考慮噪聲的影響計算當前體素的后驗分布，從后驗分布中采樣獲取當前纖維的傳播方向，計算每一條纖維的連接概率且繪制成概率圖的方式產(chǎn)生大量的纖維.這樣便能較好的解決噪聲影響和追蹤過早停止等問題，同時也出現(xiàn)了在每一傳播點大量采樣耗時長的問題.為了改善它，許多學(xué)者在貝葉斯概率追蹤算法的理論框架的基礎(chǔ)上應(yīng)用了粒子濾波[15](Particle Filtering Tractography，PFT)、無跡卡曼濾波[16](Unscented Kalman Filter，UKF)、設(shè)定概率追蹤閾值等[17-19]方法進行纖維追蹤.概率型追蹤算法使得從一個給定的種子點發(fā)出的軌跡的分布完全取決于局部纖維取向估計的統(tǒng)計不確定性.有些概率型追蹤算法假定離散光纖方向受噪聲和模型誤差影響(Diffusion Orientation Distribution Function，dODF)，從而導(dǎo)致纖維方向的分布具有容積效應(yīng).有些概率型追蹤算法假設(shè)每個體素包含連續(xù)分布的纖維方向(Fiber Orientation Distribution Function，fODF)，通常忽略噪聲和模型不確定性.所以在軌跡傳播過程中，從fODF中抽取隨機樣本構(gòu)建的是一個統(tǒng)計不確定性的纖維分布圖，因此也產(chǎn)生大量的偽纖維.

2.2 全局型追蹤算法分析

全局型追蹤算法不同于利用局部微結(jié)構(gòu)信息一步一步連結(jié)出長纖維的思想，它在纖維群里設(shè)置一個測試集，企圖通過最優(yōu)化代價函數(shù)來找到同時重構(gòu)所有纖維的方法.在2008年，Kreher等[20]提出吉布斯跟蹤算法成為全局型追蹤的主流思想，其將每一段小纖維結(jié)構(gòu)看成是一個小直圓柱體，它的長度、位置以及方向可以任意改變，如圖2所示.它首先創(chuàng)建一組試驗纖維配置并給予相互作用能使這組柱面配置形成有低分支概率的長聚合鏈，而后使用模擬退火方法[21]并計算相應(yīng)的DW信號，在這個過程中調(diào)整試驗纖維配置到實驗測量的信號以使差異最小化，找到最優(yōu)的最終配置并應(yīng)用于所有纖維重建.具體使用全局型追蹤算法追蹤的流程框架由表2所示.計算性能低、纖維細節(jié)分布不夠理想是全局型追蹤算法需解決的問題，通過不同更簡化的模型的應(yīng)用使得全局型追蹤算法計算時間有了一些縮短[22-24].馮遠靜等[25]提出利用蟻群算法思想使用少量粒子采樣來獲取全局纖維；而后又有學(xué)者[26]提出將一個基于球面卷積的多殼多組織模型整合到全局型追蹤算法框架里，并有效處理了部分容積效應(yīng)；Gutierrez等[27]提出一個通用的數(shù)據(jù)驅(qū)動框架，整體優(yōu)化追蹤算法的參數(shù)，采用多準則決策分析選擇一組標準參數(shù)應(yīng)用于全局型優(yōu)化追蹤大腦纖維.近來，深度學(xué)習(xí)等[28]方法加入到全局型纖維追蹤里，研究者可以直接從擴散加權(quán)成像或者擴散張量成像數(shù)據(jù)中學(xué)習(xí)纖維流線的生成過程，這使得全局型追蹤算法更加的多元化.相對于局部型追蹤算法，全局型追蹤算法對噪聲和成像偽影有更好的魯棒性.但是全局方法的主要問題是使用代價函數(shù)時計算時間長卻不能保證收斂到全局最優(yōu)解.進一步的問題是定義先驗知識的任意性，對先驗知識的應(yīng)用過于強烈可能會導(dǎo)致與實際數(shù)據(jù)不一致而追蹤到不存在的纖維軌跡，另一方面，如果對先驗知識的應(yīng)用較少，雖然纖維軌跡與基礎(chǔ)數(shù)據(jù)完美匹配，但這在解剖學(xué)上來說沒有什么實用性意義.

圖2 吉布斯跟蹤算法纖維結(jié)構(gòu)鏈Fig.2 Gibbs tracking algorithm for fiber structure chain

表2 基于全局型纖維追蹤的流程框架Table 2 Process framework based on global fiber tracking

3 優(yōu)化纖維束追蹤算法的主要方面

神經(jīng)纖維追蹤算法面臨著許多的挑戰(zhàn)[29]，研究者們利用低分辨率的張量圖像對神經(jīng)纖維進行重建，做出了許多努力.本文概述了在局部型纖維追蹤作出的5個主要優(yōu)化方向.全局型纖維追蹤算法風(fēng)格多樣，優(yōu)化方向不能一概而論，但是在種子點的選取和終止準則這兩方面的優(yōu)化與局部型纖維追蹤有著共通性.

3.1 局部纖維方向建模

局部纖維方向建模是纖維束追蹤的第1步，它的目的是獲得體素內(nèi)纖維結(jié)構(gòu)的走向.DTI是最常用的dMRI建模方法，它假設(shè)水分子的擴散運動呈現(xiàn)高斯分布，擴散信號E(b)可以由公式(2)表示：

(2)

其中E(b)是高斯信號，S(0)是未施加敏感梯度的圖像強度，S(b)是施加敏感梯度方向上的圖像強度，b是擴散敏感系數(shù)，D是擴散張量.在此建模方法下采樣方向少能快速成像.但是在一般的擴散過程中，當纖維不在單一方向擴散時EAP便不是高斯分布，對于纖維交叉等復(fù)雜結(jié)構(gòu)無法很好的描述.所以Tuch等[30]又提出了高角分辨率成像，采用更高的角度分辨率進行采樣：混合張量模型、高階張量模型[31](High Order Tensor，HOT)以及擴散譜成像[32](Diffusion Spectrum Imaging，DSI)把DTI中的張量模型推廣到非高斯情況，它們可以檢測出纖維交叉.但是混合張量模型中的張量數(shù)需要人為經(jīng)驗選擇，高階張量模型中的SH基階數(shù)過大時不能準確表示纖維方向，DSI在q空間中進行笛卡爾采樣，雖然使用很多DWI圖像和大范圍b值變化使獲得的ODF峰銳利清晰，但是采樣時間過長不能應(yīng)用于臨床；混合擴散成像[33](Hybrid Diffusion Imaging，HYDI)采用q空間多球殼稀疏采樣減少了采樣時間，試圖使用不同模型處理不同b值的數(shù)據(jù)，但是受人為插值和重采樣誤差的影響很大；Q球成像[34](Q-ball Imaging，QBI)容易實現(xiàn)利用單球殼數(shù)據(jù)來估計ODF的極大值來推斷纖維方向，對多球殼采樣應(yīng)用困難.擴散密度成像[35](Diffusion Propagator Imaging，DPI)是QBI的推廣.為了處理多球殼數(shù)據(jù)，假設(shè)把擴散信號表示成拉普拉斯方程的解的形式，這便導(dǎo)致產(chǎn)生了許多建模誤差；球面反卷積[36](Spherical Deconvolution，SD)使用連續(xù)加權(quán)的方法得到更尖銳的ODF，選擇合適的核函數(shù)是模型依賴的關(guān)鍵.約束球面反卷積[37](Constrained Spherical Deconvolution，CSD)使用球諧基函數(shù)表示ODF，球級傅里葉成像[38](Spherical Polar Fourier Imaging，SPFI)將擴散信號表示成球級傅里葉基的線性組合，這兩者都易受噪聲和部分容積效應(yīng)的影響；模型依賴的建模方法都試圖使用現(xiàn)有特定的函數(shù)并加以假設(shè)限制性條件表示擴散信號變化，沒有模型依賴的建模方法雖然保留了原始數(shù)據(jù)的特征，但通過人為后處理：重采樣、插值、正則化項改進等方法也很難不引入誤差.總之，局部纖維方向建模是局部型追蹤算法的前提，特別是在灰質(zhì)與白質(zhì)交界處的ODF關(guān)乎后續(xù)纖維連接的準確性.

3.2 張量插值

擴散張量成像是臨床應(yīng)用最廣泛的磁共振成像技術(shù)，故對張量圖像進行插值來提高圖像分辨率是非常有必要的.擴散張量是一個二階對稱正定矩陣，可以由公式(3)表示：

(3)

基于歐幾里得空間的張量圖像插值方法是經(jīng)典的插值方法之一[39]，將歐幾里得空間中的擴散張量分解為6個擴散加權(quán)標量圖像的矩陣分量，分別對其線性插值.基于矩陣分量進行插值的方法其結(jié)果是連續(xù)的，但在插值過程中會產(chǎn)生負定張量、在某些點周圍產(chǎn)生部分擴散張量各向異性趨于零的擠壓效應(yīng).為了保證擴散張量的對稱正定性，Zhukov等[40]提出使用三分段線性插值重建圖像中的連續(xù)張量場，Wang等[41]提出使用Cholesky因子對張量進行正定約束.又有學(xué)者提出把張量空間變換為黎曼流形[42]，從線性空間變?yōu)樵谝?guī)則的流形上進行插值，由于其需要使用矩陣的對數(shù)等數(shù)據(jù)來處理張量導(dǎo)致計算量很大.為了提高DTI數(shù)據(jù)的處理效率，許多學(xué)者把矩陣張量用特征值和特征向量來表示，通過對特征值和特征向量進行系統(tǒng)地插值[43，44]來提高計算速度.Merino-Caviedes等[45]提出由插值點的特征值和主特征向量構(gòu)成相對于笛卡兒坐標系的旋轉(zhuǎn)角度，而后進行二維平面域上的二階擴散張量插值.基于特征的插值方法可以充分保留各向異性特征，但是在插值過程中會出現(xiàn)方向上的不連續(xù)導(dǎo)致追蹤方向異常改變.擠壓效應(yīng)與方向的不連續(xù)問題分別是基于矩陣的插值方法和基于特征的插值方法需要攻克的難題.

3.3 種子點選取

當劃定特定的感興趣區(qū)域以后，纖維追蹤從感興趣區(qū)域內(nèi)的種子點為起始點開始執(zhí)行.一般感興趣區(qū)的選擇是根據(jù)研究者的需求進行自主劃定，當然這一項任務(wù)需要解剖學(xué)知識介入.除了選取特定區(qū)域的纖維，也可以進行全腦的神經(jīng)纖維束追蹤，Ben等人[46]認為種子點可以從所有白質(zhì)體素中播撒，也可以僅從白質(zhì)與灰質(zhì)交界面附近的體素進行播撒.在全腦追蹤時，以全腦白質(zhì)體素播種為追蹤起點可以保證被更多的種子點探測到由此可以追蹤到所有白質(zhì).有研究者在應(yīng)用中[47]指出使用傳統(tǒng)全局型追蹤算法對纖維進行重建時希望在腦白質(zhì)區(qū)域種子點數(shù)量較大，如此便可創(chuàng)建更多的纖維段連接成最佳的纖維束.但是這些白質(zhì)束中可能會存在過度長或直徑過度大的錯誤纖維，同時帶來了追蹤耗時長的問題；在白質(zhì)與灰質(zhì)交界面進行播種減少了許多種子點的過度交互影響，可以防止超出解剖學(xué)范圍內(nèi)過度長的纖維束，也導(dǎo)致追蹤時可能提早終止使得纖維軌跡不能覆蓋全腦白質(zhì).

3.4 傳播方向

傳統(tǒng)的追蹤方法使用的傳播方向是基于纖維取向分布中提取的對稱峰值方向，這些峰值方向假設(shè)體素中的纖維軌跡直接近似于纖維取向分布的主方向.而在復(fù)雜的纖維結(jié)構(gòu)區(qū)域被證明是不正確的，彎曲、交叉等纖維其實不完全對稱.為了準確估計對稱殘差，研究者們嘗試使用局部纖維方向建模來估計纖維傳播方向，也提出了非對稱纖維取向分布函數(shù)來解決纖維形狀不對稱問題[48，49].當纖維追蹤到纖維交叉、分支等含有多個纖維峰值方向的體素時，設(shè)定余弦相似度閾值，找到與上一步纖維方向最接近的當前纖維方向作為前進方向.Lei Ye等[50]便提出通過分析體素之間的擴散方向角度和擴散張量相似性進行傳播方向判定.如何精準的判斷纖維的傳播方向，非對稱性纖維取向模型的應(yīng)用將是一個大趨勢.

3.5 終止準則

目前關(guān)于纖維束追蹤停止條件最廣泛使用的是Smith[51]提出的解剖學(xué)約束纖維追蹤成像終止標準(Anatomically-constrained tractography，ACT).ACT設(shè)定皮層灰質(zhì)、皮層下灰質(zhì)、白質(zhì)和腦脊液的體積分數(shù)閾值，允許纖維停止在達到相應(yīng)體積分數(shù)閾值的地方.擴散各向異性分數(shù)(fractional anisotropy，F(xiàn)A)閾值以及纖維彎曲的曲率閾值等也將限制纖維的解剖學(xué)結(jié)構(gòu)合理性.隨后，Gabriel等[52]提出將解剖學(xué)先驗信息融入整個纖維追蹤過程減少流線分布中的偏差，并使用粒子濾波追蹤方法為停止在皮層下灰質(zhì)的纖維尋找另一種可能替代的有效路徑.近來研究者們發(fā)現(xiàn)很多纖維追蹤方法仍然使用基于體素的指標，比如FA、廣義各向異性(Generalized fractional anisotropy，GFA)等[53，54].對于同一體素中纖維方向都具有相同的各向異性值，所以FA等會容易受部分容積效應(yīng)的影響.故Yeh[55]提出使用定量各向異性(Quantitative Anisotropy，QA)設(shè)定閾值判定纖維終止效果更好.

在全局型追蹤方法里，終止準則不能保證一些非常短的纖維段在合適的位置產(chǎn)生與終止，所以在最終的纖維束圖中會有短小錯誤的纖維存在.雖然擬合經(jīng)驗數(shù)據(jù)非常好，但實際上沒有提供可靠的宏觀連接信息.同樣的，局部型追蹤方法也面臨著這樣的問題.纖維終止準則目前使用比較統(tǒng)一，相關(guān)改進的研究較少，但是纖維停止的準確率還需要提高，尤其是遇到皮質(zhì)褶皺內(nèi)的纖維如何準確停止仍是一個難題.

4 總結(jié)與展望

神經(jīng)纖維追蹤算法對于臨床醫(yī)學(xué)分析來說有著很高的研究價值，如早期的阿爾茲海默癥可以依靠神經(jīng)纖維微結(jié)構(gòu)變化情況進行輔助診斷[56，57].本文分析了纖維追蹤算法里局部型追蹤算法與全局型追蹤算法的聯(lián)系與區(qū)別，又在此基礎(chǔ)上介紹了關(guān)于追蹤算法目前在做的具體優(yōu)化方向，具體文獻整理分類由表3所示.在現(xiàn)有算法的基礎(chǔ)上，仍有許多問題亟待解決，未來的工作可以參考如下幾個方面：

表3 文獻整理分類表Table 3 Literature classification table

1)在局部方向建模中，可以對張量圖像進行基于矩陣分量以及特征值與特征向量的聯(lián)合插值預(yù)處理.在沒有既能快速采樣成像又擁有足夠高分辨率的情況下，如何準確插值是需要攻克的難題.不管是在確定型追蹤和概率型追蹤時，約束型球面反卷積模型都表現(xiàn)出良好的穩(wěn)定性和較高的計算效益使得它有很好的臨床應(yīng)用前景，但是由于實際的體素纖維分布通常是不對稱的，所以非對稱的纖維取向分布與約束型球面反卷積模型相結(jié)合可能會使追蹤算法更加準確.

2)在種子區(qū)域的選擇更加青睞于播撒在白質(zhì)與灰質(zhì)的交界處，防止神經(jīng)纖維在白質(zhì)中產(chǎn)生和終止.終止準則是唯一判定纖維停止的系統(tǒng)性規(guī)則，但是其中許多限制描述不一定適用于所有纖維.終止準則中還需要加入更多的解剖學(xué)先驗知識，完善纖維連接體的結(jié)構(gòu)信息，特別是皮層下溝回之間白質(zhì)的連接和停止情況.

3)雖然目前的纖維追蹤效果與實際還有著差距，許多算法在實際應(yīng)用上又很有難度，但是近幾年隨著像深度學(xué)習(xí)等新穎的方法引入全局型纖維追蹤領(lǐng)域，預(yù)示著實現(xiàn)纖維準確重構(gòu)還有著很大的進步空間.全局型追蹤算法需要將個別纖維連接配置應(yīng)用到全腦纖維連接，這不僅需要找到纖維連接的共性還要克服不同區(qū)域纖維走向的差異，更要保證計算效率.