專(zhuān)利名稱(chēng)：一種利用表面模型的大腦皮層溝回結(jié)構(gòu)生長(zhǎng)仿真方法
技術(shù)領(lǐng)域：
本發(fā)明涉及一種利用表面模型的大腦皮層生長(zhǎng)仿真方法，屬于計(jì)算神經(jīng)科學(xué)與仿真等領(lǐng)域。適用于對(duì)三維大腦皮層生長(zhǎng)形態(tài)的仿真。

背景技術(shù)：
人類(lèi)大腦皮層在大小、形狀以及結(jié)構(gòu)模式方面存在著巨大的個(gè)體差異。其中大腦皮層的溝回卷曲模式是對(duì)大腦功能的一個(gè)很好的預(yù)測(cè)，并引起了很多研究人員的注意。雖然人類(lèi)大腦都是從一個(gè)相似的簡(jiǎn)單結(jié)構(gòu)開(kāi)始發(fā)育，但是其主要的溝回結(jié)構(gòu)差異在懷孕八個(gè)月時(shí)就已經(jīng)產(chǎn)生了。由于大量的復(fù)雜的生物發(fā)育活動(dòng)存在這一過(guò)程中，目前對(duì)大腦皮層溝回結(jié)構(gòu)差異的產(chǎn)生原因還是不確定的。而通過(guò)實(shí)驗(yàn)或觀測(cè)的方法來(lái)對(duì)懷孕胎兒大腦的動(dòng)態(tài)研究方法還不成熟，用計(jì)算機(jī)仿真的方法來(lái)模擬這一過(guò)程將有利于生物和醫(yī)學(xué)界對(duì)這一問(wèn)題的理解。
近些年來(lái)大腦皮層溝回結(jié)構(gòu)發(fā)育的仿真方法已經(jīng)成為研究的熱點(diǎn)課題。各種各樣腦溝區(qū)域分割的提取方法已經(jīng)被提出來(lái)，但仍存在很多問(wèn)題。方法1基于二維連續(xù)流體模型的皮層仿真方法，該方法利用流體力學(xué)模型來(lái)對(duì)大腦和生長(zhǎng)過(guò)程建模，然后用準(zhǔn)靜態(tài)方法來(lái)最小化能量函數(shù)。其缺點(diǎn)在于，流體模型本身不能很好的表達(dá)腦組織的物理特性，并且該方法僅應(yīng)用于二維簡(jiǎn)單皮層模型而不是真實(shí)三維大腦模型；方法2基于計(jì)算形態(tài)學(xué)方法的皮層仿真模型，該方法利用二維有限元的方法來(lái)對(duì)皮層結(jié)構(gòu)以及膠質(zhì)纖維進(jìn)行建模，并通過(guò)擴(kuò)張有限元來(lái)模擬皮層生長(zhǎng)。該方法的缺點(diǎn)在于，其僅僅應(yīng)用于二位簡(jiǎn)單環(huán)狀模型，還無(wú)法適用于真實(shí)的三維大腦皮層模型；方法3基于三維有限元的皮層生長(zhǎng)仿真模型，該方法用羊胎兒的DTI數(shù)據(jù)來(lái)對(duì)大腦皮層建模，并通過(guò)對(duì)有限元進(jìn)行切向擴(kuò)張來(lái)模擬皮層生長(zhǎng)。缺點(diǎn)在于該方法僅能實(shí)現(xiàn)小范圍的皮層形變，而無(wú)法仿真大形變溝回結(jié)構(gòu)的產(chǎn)生。
目前已有的大腦皮層溝回結(jié)構(gòu)仿真方法具有以下兩個(gè)主要的缺陷其一、使用了二維的簡(jiǎn)單模型來(lái)代替大腦皮層結(jié)構(gòu)，而不是真實(shí)三維胎兒大腦皮層結(jié)構(gòu)。其二、目前的三維方法只是用于小形變的大腦皮層生長(zhǎng)仿真，而無(wú)法適應(yīng)存在大形變的大腦皮層溝回結(jié)構(gòu)。


發(fā)明內(nèi)容
要解決的技術(shù)問(wèn)題 為了避免現(xiàn)有方法的不足之處，本發(fā)明提出一種基于利用表面模型的大腦皮層生長(zhǎng)仿真方法。
技術(shù)方案 本發(fā)明的基本思想是我們用一個(gè)三角形化的表面來(lái)表達(dá)胎兒在發(fā)育中的大腦皮層，該皮層通常通過(guò)對(duì)胎兒大腦進(jìn)行MRI掃描得到。然后我們對(duì)每個(gè)三角形本身建立彈性塑性力學(xué)模型，以及三角形彎曲的彈性塑性力學(xué)模型。通過(guò)增加每個(gè)三角形的原始大小建立皮層的生長(zhǎng)模型。然后我們通過(guò)偏微分方程組來(lái)表示皮層褶皺的過(guò)程，并用迭代的方法求解這一過(guò)程。在這一求解的過(guò)程中，我們通過(guò)MRI數(shù)據(jù)對(duì)該模型設(shè)置邊界條件。
本發(fā)明的技術(shù)特征在于步驟如下 步驟1首先對(duì)胎兒三維大腦核磁共振圖像進(jìn)行組織分割去除非大腦組織，對(duì)大腦圖像進(jìn)行腦組織灰質(zhì)，丘腦，腦室和腦殼四個(gè)部分的標(biāo)識(shí)； 步驟2從上述組織分割后得到的大腦圖像中利用行進(jìn)立方體方法重構(gòu)一系列頂點(diǎn)和三角形表示的大腦皮層表面； 步驟3然后建立表面上每條邊的彈性塑性力學(xué)模型和其中

是邊ij的原始長(zhǎng)度，

是當(dāng)邊ij變形的長(zhǎng)度，τce∈[500，2000]是邊ij的范性形變系數(shù)，

是沿著i和j之間的力，Kc∈
是邊ij的彈性系數(shù)，

和

分別是兩個(gè)頂點(diǎn)的當(dāng)前坐標(biāo)； 步驟4建立表面上每條邊彎曲的彈性塑性力學(xué)模型和其中p和q兩個(gè)頂點(diǎn)是與i和j兩點(diǎn)分享同一個(gè)三角形的點(diǎn)，

是邊pq的原始長(zhǎng)度，

是頂點(diǎn)p和q之間的距離，τcb∈[500，2000]是邊pq的范性形變系數(shù)，

是沿著p和q之間的彎曲力，Kb∈
是邊pq的彎曲彈性系數(shù)，

和

分別是兩個(gè)頂點(diǎn)的當(dāng)前坐標(biāo)； 步驟5建立表面上每條邊的生長(zhǎng)模型lc0是邊的原始長(zhǎng)度，k是該條邊期望其生長(zhǎng)結(jié)束時(shí)的目標(biāo)長(zhǎng)度，m∈
為控制生長(zhǎng)速度的系數(shù)； 步驟6最后通過(guò)偏微分方程組來(lái)表示皮層褶皺的過(guò)程，并用迭代的方法求解這一過(guò)程，得到大腦生長(zhǎng)溝回形成的仿真結(jié)果；其中mi為頂點(diǎn)質(zhì)量，x是頂點(diǎn)坐標(biāo)，Δt∈
是迭代步長(zhǎng)；所述的頂點(diǎn)質(zhì)量mi為該頂點(diǎn)周?chē)切钨|(zhì)量的三分之一； 所述模型的邊界條件為提取的胎兒三維大腦核磁共振圖像數(shù)據(jù)中的腦殼、丘腦和腦室組織區(qū)域所在的空間區(qū)域通過(guò)體元圖記錄下來(lái)設(shè)置為模型的邊界條件，并限制頂點(diǎn)坐標(biāo)x不能變形到這些空間區(qū)域中。
所述步驟6中的迭代的方法通過(guò)隱式迭代方法 和求解偏微分方程組來(lái)表示皮層褶皺的過(guò)程，其中β∈
，γ∈
。
將彈性塑性力學(xué)模型建立在三角形的每一條邊上，包括該條邊的由于變形引起的彈性力以及由于變形引起的永久的范性形變，同時(shí)由該邊連接的兩個(gè)三角形之間角度(彎曲)的變化也用該相同的力學(xué)模型表達(dá)。而皮層的增長(zhǎng)則又一個(gè)隨時(shí)間變化的函數(shù)來(lái)表達(dá)，該函數(shù)決定了三角形每條邊在不受力情況下的長(zhǎng)度。通過(guò)對(duì)該長(zhǎng)度增長(zhǎng)的控制，我們可以控制皮層生長(zhǎng)的速度。
我們通過(guò)胎兒的MRI掃描數(shù)據(jù)同時(shí)還可以得到限制皮層生長(zhǎng)的生物物理?xiàng)l件，包括發(fā)育中的腦殼，腦室以及丘腦等。這些組織將通過(guò)手工標(biāo)注MRI圖像獲得，并保存在三維體元圖中，當(dāng)求解皮層溝回形成過(guò)程中，皮層上的點(diǎn)將會(huì)和這個(gè)體元圖作比較并保證腦皮層不會(huì)變形到這些組織的點(diǎn)中。
有益效果 本發(fā)明提出的利用表面模型的大腦皮層生長(zhǎng)仿真方法可行性體現(xiàn)在，首先，隨著磁共振成像設(shè)備的精度不斷提高和三維大腦磁共振圖像的預(yù)處理方法進(jìn)一步成熟，獲取胎兒準(zhǔn)確大腦皮層表面相對(duì)容易；同時(shí)，目前用表面模型仿真進(jìn)行仿真的方法已經(jīng)大量的用于衣服紙張等柔軟物質(zhì)的仿真中，所以用類(lèi)似方法對(duì)殼狀軟組織皮層的仿真是可行的。
本發(fā)明相對(duì)于其它方法具有以下優(yōu)點(diǎn)1、利用真實(shí)三維胎兒數(shù)據(jù)進(jìn)行仿真；2、將皮層生長(zhǎng)模型加入到傳統(tǒng)的表面模型中。



圖1本發(fā)明方法的基本流程圖 圖2胎兒大腦組織標(biāo)示與皮層重建示意圖 圖3三角形表面上的力學(xué)模型示意圖 圖4大腦皮層溝回結(jié)構(gòu)生長(zhǎng)仿真結(jié)果 (a)仿真發(fā)育過(guò)程在有腦殼限制的情況下的發(fā)育過(guò)程 (b)仿真發(fā)育過(guò)程在沒(méi)有腦殼限制的情況下的發(fā)育過(guò)程。

具體實(shí)施例方式 現(xiàn)結(jié)合附圖對(duì)本發(fā)明作進(jìn)一步描述 根據(jù)本發(fā)明提出的利用表面模型的大腦皮層生長(zhǎng)仿真方法，我們用C++語(yǔ)言實(shí)現(xiàn)了一個(gè)大腦皮層生長(zhǎng)仿真原型系統(tǒng)。圖像數(shù)據(jù)的來(lái)源是實(shí)際中正常胎兒的三維大腦核磁共振圖像。
本發(fā)明整個(gè)流程可以參考附圖1，具體的實(shí)施步驟如下 1.大腦皮層表面與體元圖重建 對(duì)三維大腦核磁共振圖像進(jìn)行去非大腦組織，大腦組織標(biāo)識(shí)和重建幾何結(jié)果精確，拓?fù)浣Y(jié)構(gòu)正確的大腦皮層表面，以及其他組織(丘腦，腦室和腦殼)的體元圖，如附圖2。
2.大腦皮層表面上的生長(zhǎng)模型計(jì)算 對(duì)于大腦皮層三角形上的每一條邊，我們記錄其原始長(zhǎng)度，并按照以下公式在每次迭代中增加其原始長(zhǎng)度 其中l(wèi)c0是邊的原始長(zhǎng)度，k是該條邊的目標(biāo)長(zhǎng)度，即期望其生長(zhǎng)結(jié)束時(shí)的長(zhǎng)度，m控制了生長(zhǎng)速度。
3.大腦皮層表面上的彈性范性力學(xué)模型計(jì)算 如圖3所示我們對(duì)三角型的每一條邊按以下公式計(jì)算其范性形變 其中

是邊ij的原始長(zhǎng)度，

是當(dāng)這條邊當(dāng)前的長(zhǎng)度，即頂點(diǎn)i和j之間的距離，τce是當(dāng)前邊的范性形變系數(shù)。
然后我們按照以下公式計(jì)算該條邊對(duì)兩個(gè)頂點(diǎn)的彈性力 其中

是沿著i和j之間的力，Kc是該邊的彈性系數(shù)，

是邊ij的原始長(zhǎng)度，

是當(dāng)這條邊當(dāng)前的長(zhǎng)度，

和

分別是兩個(gè)頂點(diǎn)的當(dāng)前坐標(biāo)。
同時(shí)考慮到腦殼有一定厚度，我們需要計(jì)算表面再?gòu)澢鷷r(shí)產(chǎn)生的范性形變 其中，p和q兩個(gè)頂點(diǎn)是與i和j兩點(diǎn)分享同一個(gè)三角形的點(diǎn)，如圖三所示。

是邊pq的原始長(zhǎng)度，

是即頂點(diǎn)p和q之間的距離，τce是當(dāng)前邊的范性形變系數(shù)。
然后我們按照以下公式計(jì)算其對(duì)頂點(diǎn)p和q的彈性力 其中

是沿著p和q之間的彎曲力，Kb是該邊的彎曲彈性系數(shù)，

是邊pq的原始長(zhǎng)度，

是當(dāng)這條邊當(dāng)前的長(zhǎng)度，

和

分別是兩個(gè)頂點(diǎn)的當(dāng)前坐標(biāo)。
然后對(duì)于每一個(gè)頂點(diǎn)，我們對(duì)其上所有彈性力求和得到該頂點(diǎn)的當(dāng)前受力fi。
4.大腦皮層的動(dòng)態(tài)演變模型計(jì)算 在得到受力后我們就可以計(jì)算每個(gè)頂點(diǎn)按時(shí)間變化的加速度，速度(初始為零)以及位移。對(duì)于任意頂點(diǎn)i，其加速度為 其中mi為頂點(diǎn)質(zhì)量，通常定義為該頂點(diǎn)周?chē)切钨|(zhì)量的三分之一。而速度則按以下公式計(jì)算 其中Δt是迭代時(shí)間間隔。而該頂點(diǎn)的位置則變?yōu)? 5.大腦皮層生長(zhǎng)的邊界條件計(jì)算 對(duì)于每個(gè)頂點(diǎn)的新的形變位置，我們?cè)谌S體元圖中找到對(duì)應(yīng)的體元，如果該體元是腦殼，丘腦或者腦室等組織，則該頂點(diǎn)不能變形到該位置。在這種情況下，我們將使該頂點(diǎn)退回到原來(lái)的位置。
6.迭代生長(zhǎng)過(guò)程計(jì)算 對(duì)于每一次使得皮層大小變化的生長(zhǎng)的過(guò)程(第二步)發(fā)生，我們都會(huì)迭代的進(jìn)行3～5步之至到達(dá)一定的迭代數(shù)量或者穩(wěn)定。然后我們返回第二步使得大腦皮層繼續(xù)生長(zhǎng)并重復(fù)這一過(guò)程直至大腦皮層形成穩(wěn)定的溝回結(jié)構(gòu)。
為了測(cè)試該大腦皮層溝回結(jié)構(gòu)生長(zhǎng)仿真方法的準(zhǔn)確度，我們將該方法用于一個(gè)真實(shí)胎兒的大腦核磁共振圖像所重建出的大腦皮層表面。圖4顯示了我們的方法成功的仿真了腦溝回結(jié)構(gòu)從無(wú)到有的動(dòng)態(tài)形成過(guò)程。其中圖4(a)和圖4(b)分別顯示了胎兒大腦在有和沒(méi)有腦殼限制的情況下仿真發(fā)育過(guò)程，其中參數(shù)設(shè)置為Kc＝1，τce＝1000，Kb＝5，τcb＝1000，m＝0.002，k＝3。從實(shí)驗(yàn)結(jié)果看出我們的大腦皮層溝回結(jié)構(gòu)生長(zhǎng)仿真方法具有很好的表現(xiàn)。
權(quán)利要求
1.一種利用表面模型的大腦皮層生長(zhǎng)仿真方法，其特征在于步驟如下
步驟1首先對(duì)胎兒三維大腦核磁共振圖像進(jìn)行組織分割去除非大腦組織，對(duì)大腦圖像進(jìn)行腦組織灰質(zhì)，丘腦，腦室和腦殼四個(gè)部分的標(biāo)識(shí)；
步驟2從上述組織分割后得到的大腦圖像中利用行進(jìn)立方體方法重構(gòu)一系列頂點(diǎn)和三角形表示的大腦皮層表面；
步驟3然后建立表面上每條邊的彈性塑性力學(xué)模型和其中
是邊ij的原始長(zhǎng)度，
是當(dāng)邊ij變形的長(zhǎng)度，τce∈[500，2000]是邊ij的范性形變系數(shù)，
是沿著i和j之間的力，Kc∈
是邊ij的彈性系數(shù)，
和
分別是兩個(gè)頂點(diǎn)的當(dāng)前坐標(biāo)；
步驟4建立表面上每條邊彎曲的彈性塑性力學(xué)模型和其中p和q兩個(gè)頂點(diǎn)是與i和j兩點(diǎn)分享同一個(gè)三角形的點(diǎn)，
是邊pq的原始長(zhǎng)度，
是頂點(diǎn)p和q之間的距離，τcb∈[500，2000]是邊pq的范性形變系數(shù)，
是沿著p和q之間的彎曲力，Kb∈
是邊pq的彎曲彈性系數(shù)，
和
分別是兩個(gè)頂點(diǎn)的當(dāng)前坐標(biāo)；
步驟5建立表面上每條邊的生長(zhǎng)模型lc0是邊的原始長(zhǎng)度，k是該條邊期望其生長(zhǎng)結(jié)束時(shí)的目標(biāo)長(zhǎng)度，m∈
為控制生長(zhǎng)速度的系數(shù)；
步驟6最后通過(guò)偏微分方程組來(lái)表示皮層褶皺的過(guò)程，并用迭代的方法求解這一過(guò)程，得到大腦生長(zhǎng)溝回形成的仿真結(jié)果；其中mi為頂點(diǎn)質(zhì)量，x是頂點(diǎn)坐標(biāo)，Δt∈
是迭代步長(zhǎng)；所述的頂點(diǎn)質(zhì)量mi為該頂點(diǎn)周?chē)切钨|(zhì)量的三分之一；
所述模型的邊界條件為提取的胎兒三維大腦核磁共振圖像數(shù)據(jù)中的腦殼、丘腦和腦室組織區(qū)域所在的空間區(qū)域通過(guò)體元圖記錄下來(lái)設(shè)置為模型的邊界條件，并限制頂點(diǎn)坐標(biāo)x不能變形到這些空間區(qū)域中。
2.根據(jù)權(quán)利要求1所述的利用表面模型的大腦皮層生長(zhǎng)仿真方法，其特征在于通過(guò)隱式迭代方法和求解偏微分方程組來(lái)表示皮層褶皺的過(guò)程，其中β∈
，γ∈
。
全文摘要
本發(fā)明涉及一種利用表面模型的大腦皮層生長(zhǎng)仿真方法，技術(shù)特征在于首先對(duì)胎兒三維大腦核磁共振圖像進(jìn)行預(yù)處理，從組織分割后的大腦圖像中重構(gòu)幾何結(jié)果精確，拓?fù)浣Y(jié)構(gòu)正確的大腦皮層表面，該皮層表面由一系列頂點(diǎn)和三角形表示。然后對(duì)每個(gè)三角形本身建立彈性塑性力學(xué)模型，以及三角形彎曲的彈性塑性力學(xué)模型。通過(guò)增加每個(gè)三角形的原始大小建立皮層的生長(zhǎng)模型。然后通過(guò)偏微分方程組來(lái)表示皮層褶皺的過(guò)程，并用迭代的方法求解這一過(guò)程。在這一求解的過(guò)程中，我們通過(guò)胎兒MRI數(shù)據(jù)對(duì)該模型設(shè)置邊界條件。本發(fā)明相對(duì)于其它方法具有以下優(yōu)點(diǎn)1.利用真實(shí)三維胎兒數(shù)據(jù)進(jìn)行仿真；2.將皮層生長(zhǎng)模型加入到傳統(tǒng)的表面模型中。
文檔編號(hào)G06F19/00GK101546358SQ20091002240
公開(kāi)日2009年9月30日 申請(qǐng)日期2009年5月7日 優(yōu)先權(quán)日2009年5月7日
發(fā)明者雷 郭, 聶晶鑫, 劉天明, 剛 李 申請(qǐng)人:西北工業(yè)大學(xué)