專(zhuān)利名稱(chēng)：基于拓?fù)浔３值目焖傩羞M(jìn)纖維跟蹤方法
技術(shù)領(lǐng)域：
本發(fā)明屬于磁共振擴(kuò)散張量成像領(lǐng)域，涉及一種拓?fù)浔３值目焖傩羞M(jìn)纖維跟蹤方 法。
背景技術(shù)：
擴(kuò)散張量成像(DTI)是在磁共振成像(MRI)基礎(chǔ)上發(fā)展起來(lái)的一種新的無(wú)創(chuàng)性成 像發(fā)法，它利用組織中水分子自由熱運(yùn)動(dòng)各向異性的原理探測(cè)組織的微觀結(jié)構(gòu)，達(dá)到研究 人體功能的目的。目前擴(kuò)散張量成像是唯一可在活體顯示腦白質(zhì)纖維束的無(wú)創(chuàng)性成像方 法，它使評(píng)估腦白質(zhì)纖維束的組織結(jié)構(gòu)及其連通性成為可能，具有一般MRI檢查無(wú)法比擬 的優(yōu)越性，并為神經(jīng)影像學(xué)開(kāi)辟了新的廣闊前景。近年來(lái)，DTI研究的重要進(jìn)展是將其應(yīng)用于腦白質(zhì)纖維的可視化，其中纖維跟蹤或 者稱(chēng)為腦白質(zhì)纖維成像是DTI可視化中研究的一個(gè)熱點(diǎn)。在腦白質(zhì)纖維中，沿著其長(zhǎng)軸方 向擴(kuò)散的水分子由于阻力較小而擴(kuò)散較快，而垂直于其長(zhǎng)軸方向擴(kuò)散的水分子由于阻力較 大而擴(kuò)散較慢。纖維跟蹤使用連續(xù)的曲線表示纖維的走向和分布，通過(guò)纖維跟蹤方法可以 得到三維連續(xù)的腦白質(zhì)組織結(jié)構(gòu)，可以顯示腦白質(zhì)纖維細(xì)節(jié)。目前對(duì)腦白質(zhì)神經(jīng)纖維的重 建是磁共振擴(kuò)散成像領(lǐng)域內(nèi)的研究重點(diǎn)，對(duì)腦白質(zhì)纖維走向和聯(lián)系進(jìn)行研究有助于獲得大 腦結(jié)構(gòu)和功能連接的信息，并可以為由于纖維缺失或結(jié)構(gòu)異常造成的疾病診斷提供有效的 信息。這不僅有助于深入的了解人腦纖維的結(jié)構(gòu)，而且在臨床上有很大的價(jià)值。在DTI中，擴(kuò)散張量是一個(gè)3*3的對(duì)稱(chēng)矩陣，其中有6個(gè)分量。將這個(gè)矩陣對(duì)角 化可以得到3個(gè)特征值及其所對(duì)應(yīng)的特征向量。其中最大特征值所對(duì)應(yīng)的特征向量(主 向量)的方向就是水分子擴(kuò)散的主要方向，通常也被認(rèn)為是該張量所在體素(將代表一定 厚度的三維空間的體積單元稱(chēng)為體素)內(nèi)纖維的方向。目前纖維跟蹤方法有很多種，它們 大概可以分成兩種基于張量域和基于全局能量最小化方法。基于張量域的纖維跟蹤算法 主要是利用局部張量信息進(jìn)行纖維跟蹤，DTI能產(chǎn)生每個(gè)體素的優(yōu)選擴(kuò)散方向，空間上每個(gè) 點(diǎn)張量的排列稱(chēng)為張量域。最初進(jìn)行纖維跟蹤是由一條纖維上的某點(diǎn)開(kāi)始，計(jì)算該點(diǎn)前進(jìn) 方向，沿該向量的方向跟蹤一段距離后，再以軌跡上新的一點(diǎn)作為開(kāi)始點(diǎn)，將這些點(diǎn)連接起 來(lái)，就可以顯示被跟蹤的纖維?；趶埩坑虻睦w維跟蹤算法關(guān)鍵在于當(dāng)前點(diǎn)擴(kuò)散方向的確 定，但其存在一個(gè)缺點(diǎn)，不能處理纖維分叉和交叉，算法往往在分叉區(qū)和交叉去停止。而基 于全局能量最小化的方法可以處理分叉和交叉纖維，該方法增加了纖維不確定性和隨機(jī)性 的考慮，通過(guò)最小化代價(jià)函數(shù)，尋求最優(yōu)路徑，具有最小代價(jià)的路徑與真實(shí)路徑相對(duì)應(yīng)?；?于水平集(一種用于界面追蹤和形狀建模的數(shù)值技術(shù))的快速行進(jìn)纖維跟蹤算法是一種基 于全局能量最小化的方法，它通過(guò)定義一個(gè)速度函數(shù)來(lái)控制曲面演化，每次演化迭代過(guò)程 中通過(guò)檢測(cè)與周?chē)従狱c(diǎn)的所有可能路徑來(lái)尋求最優(yōu)路徑。但是演化過(guò)程中引入了太多的 假性陽(yáng)支，并且對(duì)于纖維走向描述不明顯。纖維跟蹤提供了一種方法來(lái)研究腦白質(zhì)組織結(jié)構(gòu)和聯(lián)通性，但是纖維跟蹤還是有 一定的局限性。對(duì)活體纖維跟蹤結(jié)果的評(píng)價(jià)尚缺乏金標(biāo)準(zhǔn)。DTI是活體顯示神經(jīng)纖維束軌跡的唯一方法，因?yàn)榻M織標(biāo)本在進(jìn)行解剖、冷凍、脫水、固定、切片和溶解等處理過(guò)程中，其 微觀結(jié)構(gòu)必然發(fā)生變化，進(jìn)而產(chǎn)生幾何變形，應(yīng)用組織學(xué)方法在體外驗(yàn)證活體跟蹤結(jié)果有 很大難度。目前還沒(méi)有一種算法能得到大家的普遍認(rèn)可，因此提出一種可靠的、有效的、快 速的纖維跟蹤算法是當(dāng)前研究領(lǐng)域的一個(gè)研究重點(diǎn)
發(fā)明內(nèi)容
針對(duì)基于全局能量最小化的快速行進(jìn)(Fast Marching, FM)的算法存在較多假性 陽(yáng)支，纖維走向不明顯等問(wèn)題，本發(fā)明提出了一種拓?fù)浔３值目焖傩羞M(jìn)纖維跟蹤方法。為 此，本發(fā)明采用如下的技術(shù)方案—種拓?fù)浔３值目焖傩羞M(jìn)纖維跟蹤方法，其特征在于，包括下列步驟第一步DTI數(shù)據(jù)處理；第二步種子點(diǎn)選取及初始化(4)指定初始的活動(dòng)點(diǎn)，即種子點(diǎn)r'，作為曲面演化的起始位置； (5)確定速度函數(shù)“~ _ JTmin((| ；(r) · φ) I), (j 約(r'). n(r) |), (\ J1 (r) ·；(r') |))‘其(6)確定窄帶點(diǎn)和遠(yuǎn)離點(diǎn)，并根據(jù)下列公式確定窄帶點(diǎn)的到達(dá)時(shí)間T，遠(yuǎn)離點(diǎn)的到 達(dá)時(shí)間初始為T(mén) = TIME_MAXE，并建立窄帶點(diǎn)的最小排序堆T{r) = T(r')+廠(廠)其中T(r')為？的到達(dá)時(shí)間； 第三步采用曲率加權(quán)的速度函數(shù)進(jìn)行拓?fù)浔３值那嫜莼?4)從最小排序堆里輸出到達(dá)時(shí)間T最小的點(diǎn)r'，將其標(biāo)記為活動(dòng)點(diǎn)，并從最小 排序堆中刪除；(5)考察活動(dòng)點(diǎn)r'的鄰接點(diǎn)r，對(duì)于每個(gè)鄰接點(diǎn)r，如果是活動(dòng)點(diǎn)，則不變；如果是 遠(yuǎn)離點(diǎn)則修改為窄帶點(diǎn)，計(jì)算其到達(dá)時(shí)間T (r)，記錄其父節(jié)點(diǎn)r ‘，并將其放入最小排序堆 中；如果是窄帶點(diǎn)，則重新計(jì)算其到達(dá)時(shí)間T(r)，如果到達(dá)時(shí)間T(r)小于上次迭代的到達(dá) 時(shí)間，則更新該點(diǎn)的到達(dá)時(shí)間及父節(jié)點(diǎn)r'，調(diào)整最小排序堆；(6)如果曲面演化超過(guò)體數(shù)據(jù)容積范圍，則停止迭代；設(shè)定一個(gè)時(shí)間閾值，當(dāng)點(diǎn)的 到達(dá)時(shí)間超過(guò)該時(shí)間閾值時(shí)，停止迭代；或者預(yù)先定義迭代次數(shù)，達(dá)到迭代次數(shù)停止，則循 環(huán)結(jié)束，否則轉(zhuǎn)到(1)，繼續(xù)執(zhí)行；第三步采用時(shí)間梯度下降法確定所有纖維路徑；第四步利用聯(lián)通矩陣挑選真實(shí)纖維路徑。本發(fā)明的快速行進(jìn)纖維跟蹤方法，在演化過(guò)程中加入拓?fù)浔３帜Ｐ?，并將曲率?入到速度函數(shù)，將彎曲能量考慮在全局能量范圍內(nèi)，從而更好的控制演化過(guò)程，避免了原有 的快速行進(jìn)算法使用局部相似性作為唯一測(cè)度的問(wèn)題。本發(fā)明提供的方法能夠很好的反應(yīng)纖維的分叉信息，符合真實(shí)纖維的模型，保持了良好的拓?fù)浣Y(jié)構(gòu)，相對(duì)于現(xiàn)有技術(shù)而言，具 有纖維跟蹤結(jié)果精確，假性陽(yáng)支減少，對(duì)噪聲魯棒性良好等優(yōu)點(diǎn)。
圖1是本發(fā)明流程總圖；圖2本發(fā)明采用的快速行進(jìn)算法的流程圖；圖3是本發(fā)明采用的快速行進(jìn)算法示意圖；圖4本發(fā)明采用的拓?fù)浔３值难莼Ｐ褪疽鈭D；圖5采用本發(fā)明的方法處理投射纖維放射冠得到的計(jì)算機(jī)屏幕輸出結(jié)果；圖6計(jì)算機(jī)屏幕輸出的圖5的局部放大圖。
具體實(shí)施例方式參見(jiàn)圖1本發(fā)明鑒于快速行進(jìn)纖維跟蹤方法存在的問(wèn)題，提供了一種拓?fù)浔３值?快速行進(jìn)纖維跟蹤方法，包括以下幾個(gè)方面A. DTI數(shù)據(jù)處理B.種子點(diǎn)選取及初始化C.采用曲率加權(quán)的速度函數(shù)進(jìn)行拓?fù)浔３值那嫜莼疍.采用時(shí)間梯度下降法確定所有路徑E.利用聯(lián)通矩陣挑選真實(shí)纖維路徑參見(jiàn)圖2，下面詳細(xì)介紹本發(fā)明拓?fù)浔３值目焖傩羞M(jìn)纖維跟蹤方法和如何利用它 來(lái)實(shí)現(xiàn)纖維跟蹤?！?DTI數(shù)據(jù)處理讀取圖片并組織成三維體數(shù)據(jù)，在成像的過(guò)程中由于電子干擾和外界環(huán)境的影 響，醫(yī)學(xué)圖像經(jīng)常含有噪聲，因此使用濾波器進(jìn)行去噪。·種子點(diǎn)選取及初始化用戶(hù)指定種子點(diǎn)，作為曲面演化的起始位置(即纖維跟蹤的起始位置)，然后進(jìn)行 初始化1)活動(dòng)點(diǎn)活動(dòng)點(diǎn)就是網(wǎng)格中到達(dá)時(shí)間已知的點(diǎn)，初始的時(shí)候也就是用戶(hù)指定的 種子點(diǎn)，演化曲面經(jīng)過(guò)該點(diǎn)的時(shí)間T = 0的點(diǎn)。2)窄帶點(diǎn)所有活動(dòng)點(diǎn)的鄰域，其到達(dá)時(shí)間時(shí)間T = 1/F，并將窄帶點(diǎn)放入一個(gè)最 小排序堆里面，其中F是曲率加權(quán)的速度函數(shù)，是控制曲面演化方向和速度的重要函數(shù)，在 本發(fā)明中按公式(2)計(jì)算。3)遠(yuǎn)離點(diǎn)到達(dá)時(shí)間為無(wú)窮遠(yuǎn)的點(diǎn)，T = MAX_TIME,除了活動(dòng)點(diǎn)和窄帶點(diǎn)之外的都 是遠(yuǎn)離點(diǎn)?！げ捎们始訖?quán)的速度函數(shù)進(jìn)行拓?fù)浔３值那嫜莼?)從最小排序堆里輸出到達(dá)時(shí)間T最小的點(diǎn)r'，將其標(biāo)記為活動(dòng)點(diǎn)，并從最小排 序堆中刪除。2)考察r'的鄰接點(diǎn)即該點(diǎn)在空間上的鄰居點(diǎn)，如圖3所示，對(duì)于每個(gè)鄰接點(diǎn)r，如 果是活動(dòng)點(diǎn)，則不變；如果是遠(yuǎn)離點(diǎn)則修改為窄帶點(diǎn)，并根據(jù)公式(1)計(jì)算其到達(dá)時(shí)間，記錄其父節(jié)點(diǎn)r'，并將其放入最小排序堆中。如果是窄帶點(diǎn)，則根據(jù)公式(1)重新計(jì)算其到 達(dá)時(shí)間，如果到達(dá)時(shí)間T小于上次迭代的到達(dá)時(shí)間，則更新該點(diǎn)的到達(dá)時(shí)間及父節(jié)點(diǎn)r'， 調(diào)整最小排序堆。在每一次迭代過(guò)程中，曲面就會(huì)從一點(diǎn)r'演化到r，i 是r'的鄰居節(jié)點(diǎn)，并 且在窄帶中具有最小的到達(dá)時(shí)間，同理，r'是由r"演化而來(lái)，那么在演化模型中認(rèn)為 r" -r' — r為演化過(guò)程中的拓?fù)湫畔?，并且r〃 一 r' — r在纖維路徑上，如圖4所示。 其中T(r')為r'的到達(dá)時(shí)間，F(xiàn)(r)速度函數(shù)，保證了時(shí)間T(r)從種子點(diǎn)，沿著 n(r)方向，以速度F(r)前進(jìn)。
F(r) =_----⑵ 其中r'為活動(dòng)點(diǎn)，n(r')表示曲面到r'的演化方向，n(r)表示曲面到r的演 化方向，ei(r)表示r的主向量方向，ei(r')表示r'的主向量方向。演化方向n(r)使用 n(r) = |r-r' |近似，保留的拓?fù)鋘(r' ) = |r' ~r" |是纖維曲率。當(dāng)C(r) < 0時(shí)， F(r) =0，以保證曲面演化總是擴(kuò)張，不會(huì)收縮。曲率加權(quán)速度函數(shù)基于以下四個(gè)方面的考慮a.從同一活動(dòng)點(diǎn)出發(fā)到達(dá)不同鄰域(不同窄帶點(diǎn))速度相等時(shí)，選擇曲率最小的 方向前進(jìn)。b.從不同活動(dòng)點(diǎn)出發(fā)到達(dá)相同鄰域(相同窄帶點(diǎn))到達(dá)時(shí)間相等時(shí)，選擇曲率較 小的活動(dòng)點(diǎn)為父節(jié)點(diǎn)。c.當(dāng)纖維走向變得不明顯時(shí)，向量場(chǎng)方向不一致時(shí)，曲率偏轉(zhuǎn)較小的方向優(yōu)先，保 持纖維的連續(xù)和光滑性。d.保證曲面演化沿著纖維延伸的方向，避免出現(xiàn)X和T字形狀纖維路徑。曲率加權(quán)的速度函數(shù)把彎曲能量考慮在內(nèi)，當(dāng)演化曲面達(dá)到纖維邊界，因?yàn)榍?的限制，不規(guī)則的向量場(chǎng)將被排除在外，因?yàn)樗斐蓮澢芰孔兇蟆Ｍ瑫r(shí)在纖維束內(nèi)部，曲 率加權(quán)的速度函數(shù)通過(guò)最小化彎曲能量，是的纖維跟蹤過(guò)程能夠沿著張量場(chǎng)的方向前進(jìn)， 減少了假性陽(yáng)支。3)如果曲面演化超過(guò)體數(shù)據(jù)容積范圍，則停止迭代；定義一個(gè)時(shí)間閾值，當(dāng)點(diǎn)的 到達(dá)時(shí)間超過(guò)該閾值時(shí)，認(rèn)為該點(diǎn)與種子點(diǎn)關(guān)聯(lián)程度較低，與種子點(diǎn)之間不存在腦白質(zhì)纖 維，停止迭代；或者預(yù)先定義迭代次數(shù)，達(dá)到迭代次數(shù)停止，則循環(huán)結(jié)束，否則轉(zhuǎn)到1)，繼續(xù) 執(zhí)行。 采用時(shí)間梯度下降法確定所有路徑假設(shè)存在一條纖維、，并假設(shè)長(zhǎng)度是L，T是Y上的一點(diǎn)，時(shí)間T是在速度F的 作用下的|v:r|的累計(jì)結(jié)果?？焖傩羞M(jìn)算法保證從種子點(diǎn)a到r的最小路徑發(fā)生的演化過(guò)程 中，并且T(r)就是最小代價(jià)。如果A到r存在一條路徑則滿(mǎn)足
那么從迭代結(jié)束之后，從波前位置通過(guò)時(shí)間梯度下降法，計(jì)算返回到種子點(diǎn)的最 短路徑就是一條纖維路徑。 利用聯(lián)通矩陣挑選真實(shí)纖維路徑 快速行進(jìn)算法迭代結(jié)束之后從任何到達(dá)時(shí)間為T(mén)的點(diǎn)，通過(guò)最小路徑算法都可以 與種子點(diǎn)相連，真實(shí)路徑的判斷就要通過(guò)聯(lián)通矩陣 ，設(shè)定一個(gè) 閾值，通過(guò) 函數(shù)映射， 與種子點(diǎn)具有較強(qiáng)連通性的最可能路徑被選擇出來(lái)，而那些不滿(mǎn)足條件路徑被認(rèn)為是不可 能存在的纖維，則被刪除。采用本發(fā)明提出的拓?fù)浔３值目焖傩羞M(jìn)法對(duì)于腦白質(zhì)纖維進(jìn)行重建，重建算法從 種子點(diǎn)開(kāi)始。重建數(shù)據(jù)來(lái)自GE MEDICAL SYSTEMS，采用EP/SE序列。種子點(diǎn)選擇胼胝體和 矢狀圖的交界處，圖5為選擇投射纖維放射冠，使用本發(fā)明的拓?fù)浔３值目焖傩羞M(jìn)法得到 的結(jié)果，圖中很好的區(qū)分不同纖維走向。從圖6可以看出纖維路徑呈現(xiàn)樹(shù)形結(jié)構(gòu)，很好的反 映了纖維的分叉信息，符合真實(shí)纖維的模型，結(jié)果保持了很好的拓?fù)浣Y(jié)構(gòu)。
權(quán)利要求
一種拓?fù)浔３值目焖傩羞M(jìn)纖維跟蹤方法，其特征在于，包括下列步驟第一步DTI數(shù)據(jù)處理；第二步種子點(diǎn)選取及初始化(1)指定初始的活動(dòng)點(diǎn)，即種子點(diǎn)r′，作為曲面演化的起始位置；(2)確定速度函數(shù)其中，r′為活動(dòng)點(diǎn)，n(r′)表示曲面到r′的演化方向，n(r)表示曲面到r的演化方向，e1(r)表示r的主向量方向，e1(r′)表示r′的主向量方向，演化方向n(r)使用n(r)＝|r-r′|近似，纖維曲率n(r′)＝|r′-r″|；(3)確定窄帶點(diǎn)和遠(yuǎn)離點(diǎn)，并根據(jù)下列公式確定窄帶點(diǎn)的到達(dá)時(shí)間T，遠(yuǎn)離點(diǎn)的到達(dá)時(shí)間初始為T(mén)＝TIME_MAXE，并建立窄帶點(diǎn)的最小排序堆 <mrow><mi>T</mi><mrow> <mo>(</mo> <mi>r</mi> <mo>)</mo></mrow><mo>=</mo><mi>T</mi><mrow> <mo>(</mo> <msup><mi>r</mi><mo>′</mo> </msup> <mo>)</mo></mrow><mo>+</mo><mfrac> <mrow><mo>|</mo><mi>r</mi><mo>-</mo><msup> <mi>r</mi> <mo>′</mo></msup><mo>|</mo> </mrow> <mrow><mi>F</mi><mrow> <mo>(</mo> <mi>r</mi> <mo>)</mo></mrow> </mrow></mfrac><mo>,</mo> </mrow>其中T(r′)為r′的到達(dá)時(shí)間；第三步采用曲率加權(quán)的速度函數(shù)進(jìn)行拓?fù)浔３值那嫜莼?1)從最小排序堆里輸出到達(dá)時(shí)間T最小的點(diǎn)r′，將其標(biāo)記為活動(dòng)點(diǎn)，并從最小排序堆中刪除；(2)考察活動(dòng)點(diǎn)r′的鄰接點(diǎn)r，對(duì)于每個(gè)鄰接點(diǎn)r，如果是活動(dòng)點(diǎn)，則不變；如果是遠(yuǎn)離點(diǎn)則修改為窄帶點(diǎn)，計(jì)算其到達(dá)時(shí)間T(r)，記錄其父節(jié)點(diǎn)r′，并將其放入最小排序堆中；如果是窄帶點(diǎn)，則重新計(jì)算其到達(dá)時(shí)間T(r)，如果到達(dá)時(shí)間T(r)小于上次迭代的到達(dá)時(shí)間，則更新該點(diǎn)的到達(dá)時(shí)間及父節(jié)點(diǎn)r′，調(diào)整最小排序堆；(3)如果曲面演化超過(guò)體數(shù)據(jù)容積范圍，則停止迭代；設(shè)定一個(gè)時(shí)間閾值，當(dāng)點(diǎn)的到達(dá)時(shí)間超過(guò)該時(shí)間閾值時(shí)，停止迭代；或者預(yù)先定義迭代次數(shù)，達(dá)到迭代次數(shù)停止，則循環(huán)結(jié)束，否則轉(zhuǎn)到(1)，繼續(xù)執(zhí)行；第三步采用時(shí)間梯度下降法確定所有纖維路徑；第四步利用聯(lián)通矩陣挑選真實(shí)纖維路徑。FSA00000103553200011.tif,FSA00000103553200012.tif
全文摘要
本發(fā)明屬于共振擴(kuò)散成像領(lǐng)域，涉及一種基于拓?fù)浔３值目焖傩羞M(jìn)纖維跟蹤方法，包括讀取DTI數(shù)據(jù)；人工選取種子點(diǎn)并初始化；利用快速行進(jìn)法采用曲率加權(quán)的速度函數(shù)從種子點(diǎn)向周?chē)従狱c(diǎn)演化，在演化的過(guò)程中記錄拓?fù)湫畔?，即每個(gè)演化點(diǎn)的源節(jié)點(diǎn)；采用時(shí)間梯度下降法計(jì)算出所有路徑；通過(guò)聯(lián)通矩陣挑選真實(shí)的纖維路徑。本發(fā)明具有纖維跟蹤更符合真實(shí)的纖維路徑，減小假性陽(yáng)支的優(yōu)點(diǎn)，得到的纖維光滑，很好的反應(yīng)了纖維走向。
文檔編號(hào)G06F19/00GK101872385SQ20101016049
公開(kāi)日2010年10月27日 申請(qǐng)日期2010年4月30日 優(yōu)先權(quán)日2010年4月30日
發(fā)明者張加萬(wàn), 張怡, 張勝平, 米博會(huì) 申請(qǐng)人:天津大學(xué)