專利名稱：基于插值算法的柔性體點(diǎn)分布模型的形狀、形變及協(xié)方差矩陣計(jì)算方法
技術(shù)領(lǐng)域：
本發(fā)明涉及柔性體點(diǎn)分布模型分析技術(shù)，結(jié)合圖形學(xué)、矩陣?yán)碚摗⒕€性代數(shù)、 醫(yī)學(xué)圖像分析、計(jì)算機(jī)圖像處理技術(shù)和柔性體模型的參數(shù)計(jì)算方法。
背景技術(shù)：
柔性體相對(duì)于剛性體，在計(jì)算機(jī)中進(jìn)行建模和分析比較復(fù)雜，而在現(xiàn)實(shí)生活中 物體隨時(shí)間柔性變化的情況是非常普遍的。對(duì)于許多生物醫(yī)學(xué)組織，如人的心臟， 雙手，骨骼，腎臟、大腦、細(xì)胞等等，這些人們可以很好地描述其外部形狀特征 但卻很難給出準(zhǔn)確的剛性模型。點(diǎn)分布模型(PDM)是一種新近崛起的技術(shù)，它 能很好的對(duì)研究對(duì)象進(jìn)行描述，可以成功地應(yīng)用于非剛性模型的分解與合成處理。 近年來(lái)，柔性體的統(tǒng)計(jì)學(xué)建模及其與圖像解釋技術(shù)的結(jié)合催生了一些高級(jí)應(yīng)用軟 件，這些軟件不僅針對(duì)生物醫(yī)學(xué)，還涉及目標(biāo)跟蹤、識(shí)別、以至影視等領(lǐng)域。
PDM是一種基于對(duì)象樣本分析的統(tǒng)計(jì)模型。從Cootes等人開始，這種基于 訓(xùn)練樣本上標(biāo)記點(diǎn)的形體建模思想一直被貫徹下來(lái)。最典型的想法是，使標(biāo)記點(diǎn) 自動(dòng)對(duì)齊以確定圖形的位置和形變類型。對(duì)齊圖像后可以很容易的比較同組標(biāo)記 點(diǎn)在不同例圖中的位置變化，只需要計(jì)算它們的坐標(biāo)就可以了。因而，PDM模型 對(duì)于形狀的描述是通過(guò)輪廓點(diǎn)向量來(lái)完成的，它由中間形狀和形狀變化樣式組成。 通過(guò)對(duì)這種分布情況建模，我們可以得到與原始訓(xùn)練集相似的新形狀。
在三維空間，如果訓(xùn)練集圖像中的每個(gè)圖形都由一組數(shù)量為n的標(biāo)記點(diǎn)描 述，那么它可以用一個(gè)3n維的矢量表示(如果是二維對(duì)象即為2n維矢量)。訓(xùn) 練集中的標(biāo)記工作非常重要，通常是手動(dòng)完成的。另外，在研究中我們要收集數(shù) 百名例柔性體的個(gè)例以得到比較好的統(tǒng)計(jì)結(jié)果。訓(xùn)練集形成的概率模型維度可能 非常的大(在一個(gè)心臟模型中通?？蛇_(dá)到上萬(wàn)維)，所以建立PDM的過(guò)程是非?？菰镉仲M(fèi)時(shí)費(fèi)力的。然而，在一些應(yīng)用中，單獨(dú)的一個(gè)PDM模型或是有限個(gè)PDM 模型都不足以很好的描述對(duì)象變化。事實(shí)上，不論是連續(xù)PDM模型還是模型插 入方法都是為了滿足對(duì)象的性質(zhì)，即對(duì)象時(shí)間和空間上的變化。例如在心臟建模 時(shí)，需要不同時(shí)間段的模型才能很好地描述心室形狀在心臟跳動(dòng)不同階段的變 化。但是，在構(gòu)建整個(gè)心臟運(yùn)動(dòng)周期的連續(xù)模型時(shí)，我們只能由已知的模型序列 得到特定幾個(gè)時(shí)刻的模型，而其他任意時(shí)刻的模型就不得而知了，盡管在基于模 型的圖像分割領(lǐng)域PDM建模出現(xiàn)已有十來(lái)年，但對(duì)于模型插值及相關(guān)內(nèi)容的深 入研究還沒(méi)有開展過(guò)。

發(fā)明內(nèi)容
為了克服已有的柔性體點(diǎn)分布模型的不能形成實(shí)時(shí)分布的模型、模型精度差、 實(shí)用性差的不足，本發(fā)明提供一種能夠形成實(shí)時(shí)分布的中間過(guò)渡模型、模型精度
高、實(shí)用性良好的基于插值算法的柔性體點(diǎn)分布模型的形狀、形變及協(xié)方差 矩陣計(jì)算方法。
本發(fā)明解決其技術(shù)問(wèn)題所采用的技術(shù)方案是
一種基于插值算法的柔性體點(diǎn)分布模型的形狀、形變及協(xié)方差矩陣計(jì)算 方法，包括以下步驟
1) 、載入兩個(gè)相鄰時(shí)刻i;和Tb的柔性體分布模型Ma和Mb，其計(jì)算式為-
Ma = <", A,， Mb = <",①b, Xb>;
2) 、通過(guò)所述柔性體分布模型Ma和Mb，分別計(jì)算平均形狀"和"、特征矩陣Oa
和①b、以及特征值U和"，設(shè)定柔性體從la到Tb呈連續(xù)變化，
2.1)中間模型的平均形狀
用線性插值得到y(tǒng)&和ybi之間的中間模型形狀y"為
yti = ayai + bybi (12)
a和b由式(10)和(11)定義_ "。 Z) = 1 - a
(10) (11)
a和b反映了中間模型相對(duì)于Ma和Mb的距離；
插入的平均形狀為
■S" f=i S /=l S S ,=i
2.2)、中間模型的形變
中間模型的形變由下式算出
=少,,-y, =d +^)-(ay。 +6yfc)
(13);
(14)
y6f
2.3)、中間模型的協(xié)方差矩陣
中間模型的協(xié)方差矩陣由下式推導(dǎo)得出
— 1 ,=1 *i 一 1 ,=1
(15)
■y—If
令
因?yàn)?br> —丄,=i
盧
(16)
s。 =一1^。,' - y。)(y。, -y。)T =一1>戸《。
S 一 1 ，=i S — 1 ,=i
盧
(17)
(18)
顯然有
S,S:。 (19);
3)、根據(jù)得到中間模型的平均形狀、形變和協(xié)方差矩陣，得到t時(shí)刻的中間模型 Mt。
本發(fā)明的技術(shù)構(gòu)思為PDM方法是一種先驗(yàn)建模方法，它建立在目標(biāo)對(duì)象的
一組樣本(即訓(xùn)練集)已知的基礎(chǔ)上。這組樣本集可 提供對(duì)象的形狀及形變的統(tǒng)計(jì)學(xué)描述。在實(shí)際應(yīng)用中，這組對(duì)象樣本的形狀描述通常是使用輪廓來(lái)進(jìn)行的 (如圖片上的一組像素的坐標(biāo)值)。進(jìn)一步的，可以在輪廓上選取一些標(biāo)記點(diǎn)來(lái)描 述目標(biāo)物形狀，而這些標(biāo)記點(diǎn)往往與目標(biāo)物的某種特征相關(guān)聯(lián)。
不同樣本上的對(duì)應(yīng)標(biāo)記點(diǎn)在選取的時(shí)候遵循著"特征位置大致相同"的原則，但 往往這種定位是不夠精確的，存在一定偏差。這種偏差要?dú)w因于不同個(gè)體間的自 然形變。我們期望相對(duì)于整個(gè)形體的尺度而言，這種偏差是十分微小的。PDM方 法允許我們模仿這種微小的差異，并且指出哪些偏差的確是微小的，哪些相對(duì)而 言更大一些。
由于訓(xùn)練集中的形體是從不同的圖片集中獲得的，而不同圖片集的坐標(biāo)系并不 統(tǒng)一，因此有必要先將所有訓(xùn)練集形體進(jìn)行坐標(biāo)歸一化，即將它們對(duì)齊。對(duì)齊的 過(guò)程實(shí)際上就是對(duì)每一個(gè)樣本模型尋找適當(dāng)?shù)南嗨菩宰儞Q的過(guò)程，包括平移、縮
放和旋轉(zhuǎn)，以使得不同樣本之間盡可能的接近，另一方面，選定的變換也應(yīng)該令 各個(gè)樣本模型與由所有樣本得到的平均模型盡可能接近。為了方便表述，假設(shè)我
們的樣本只有兩個(gè)，歸一化簡(jiǎn)化為對(duì)齊這兩個(gè)樣本。每個(gè)樣本的形狀都由一個(gè)標(biāo) 記點(diǎn)坐標(biāo)組成的向量來(lái)表示
x = (x力少/， q, ..., x , _y , z )T， (1) 其中(A乂， w是某個(gè)標(biāo)記點(diǎn)在三維圖像中的空間坐標(biāo)。如果是二維醫(yī)學(xué)圖像，則標(biāo) 記點(diǎn)用(x"力)表示。
圖片上標(biāo)記點(diǎn)的獲取有多種方式，人工手動(dòng)標(biāo)記或采用一定得自動(dòng)標(biāo)記算法
都是可以的。例如，Izard等人提出的一種MRI圖像標(biāo)記方法，使用統(tǒng)一的算法 規(guī)則在不同人腦圖片中進(jìn)行組織結(jié)構(gòu)標(biāo)記。得到一組經(jīng)過(guò)標(biāo)記的訓(xùn)練樣本后，我 們需要將它們對(duì)齊到一個(gè)統(tǒng)一的坐標(biāo)系下?？梢允褂脧V義對(duì)齊算法(GPA)來(lái)將 訓(xùn)練樣本對(duì)齊，并使得所有樣本形狀到平均模型的距離平方和最小，實(shí)際上也就 是尋找變換Ti (平移，旋轉(zhuǎn)，縮放)，使下式最小化An =》m-r'(x,)卩， (2)
其中
m-丄Z ;(x,), |m|=l, (3) 且Xi是訓(xùn)練集中的一個(gè)形狀樣本，它是一個(gè)3n維的向量表示(在三維空間中有n 個(gè)標(biāo)記點(diǎn)的情況下)。
對(duì)齊后的訓(xùn)練集組成了一個(gè)三維空間上的點(diǎn)云，可以看作是一個(gè)概率密度函 數(shù)。為了降低運(yùn)算成本和內(nèi)存占用，我們采用主成分分析法(PCA)尋找點(diǎn)云中 的主元，并以前列最大的少數(shù)分量建立模型。這些被選出的分量可以描述對(duì)象主 要的形變。以心動(dòng)周期內(nèi)左心室的模型為例，我們經(jīng)過(guò)主成分分析后的用60個(gè)特 征向量已經(jīng)能足夠描述心室的形狀及其變化了。最后，PDM模型可以表達(dá)為
x =玄+ +…+ p A = i +屯c (4)
其中，i是對(duì)齊后的平均形狀向量，O是一個(gè)3nxf的矩陣，①的每一列表示主成 分方向上的單位向量，c是由形狀參數(shù)組成的"隹列向量。可見，經(jīng)過(guò)PCA分析 后形狀的維度已經(jīng)由下降到了 L
經(jīng)過(guò)以上步驟我們得到了一個(gè)類似于點(diǎn)分布模型的統(tǒng)計(jì)學(xué)模型，該模型可應(yīng) 用于柔性體的物體建模，以及在新的圖像中定位和分割新個(gè)例。在訓(xùn)練集的學(xué)習(xí) 過(guò)程中我們可以得到形狀參數(shù)b的變化范圍，只要在該范圍內(nèi)變化，任意b的值 都可以生成合理的新樣本(仿真形體)。通常bi的變化為^ (矩陣O中對(duì)應(yīng)第i大 的特征值)或者VJXj。
從樣本集中訓(xùn)練出一個(gè)統(tǒng)計(jì)學(xué)形態(tài)模型后，就可以用PDM來(lái)解釋新的圖像個(gè) 例了。為了將模型與圖像匹配，通常采用一種主動(dòng)重復(fù)迭代的方法。該方法使模 型不斷地變形來(lái)適應(yīng)新圖像中目標(biāo)對(duì)象的輪廓。模型形狀的變化受到一些統(tǒng)計(jì)數(shù) 據(jù)的約束，也即只能在訓(xùn)練樣本集中所得出的范圍內(nèi)進(jìn)行變化。對(duì)于形狀模型， 我們還要求圖像輪廓大致位于模型點(diǎn)的附近。可以表示為某點(diǎn)的梯度與輪廓上過(guò)該點(diǎn)的法線方向的變化。令統(tǒng)計(jì)模型的輪廓沿圖像輪廓不斷移動(dòng)，同時(shí)計(jì)算兩輪 廓的馬氏距離最小化就可以最終得到真實(shí)的邊界位置。因此，要在特定實(shí)例中獲 得目標(biāo)柔性對(duì)象的形狀需要反復(fù)執(zhí)行以下兩步直到收斂(l)沿模型各點(diǎn)的法線尋 找圖像輪廓上相匹配的對(duì)應(yīng)點(diǎn)(使得模型點(diǎn)與對(duì)應(yīng)圖像點(diǎn)的馬氏距離最小)；(2) 改變模型位置及形狀參數(shù)使得模型各標(biāo)記點(diǎn)更接近圖像輪廓上找到的對(duì)應(yīng)點(diǎn)。也 就是說(shuō)，模型的擬合過(guò)程實(shí)際上就是在圖像上搜尋模型點(diǎn)的最相似位置及不斷修 正整體幾何變換T及形狀參數(shù)C的過(guò)程，使之最小化。數(shù)學(xué)上可以表示為<formula>formula see original document page 9</formula> 其中X是在當(dāng)前迭代時(shí)得到的臨時(shí)模型。
這是一個(gè)尋找最小值的問(wèn)題，可以通過(guò)一些非線性優(yōu)化的迭代方法求解。最
終我們可以確定整體變換參數(shù)T，以及個(gè)體實(shí)例的形狀參數(shù)b:
<formula>formula see original document page 9</formula>
近幾年，人們已經(jīng)就基于PDM的柔性體建模方法的許多相關(guān)問(wèn)題進(jìn)行了深 入的研究，如形體對(duì)齊、自動(dòng)標(biāo)記、平均形體生成、形變建模、模態(tài)分析、圖像 分割等等問(wèn)題。此外，在該思想廣泛地應(yīng)用于非剛體對(duì)象建模的基礎(chǔ)上，并可進(jìn) 一步擴(kuò)展為諸如活動(dòng)形體模型(ASM)和活動(dòng)表現(xiàn)建模(AAM)等以應(yīng)用于更深 入的圖像分析。目前，對(duì)PDM模型最成功的表達(dá)是通過(guò)PCA方法得到的。PCA 方法非常簡(jiǎn)單易懂，首先把一個(gè)形體投影到經(jīng)過(guò)學(xué)習(xí)PCA空間(由線性無(wú)關(guān)方向 形成的正交空間)，然后從訓(xùn)練樣本中提取少量與形變模式相關(guān)的參數(shù)(稱為主 要成分)，再用這些參數(shù)來(lái)控制形體的變形。顯然，這種形變是沿著樣本集所確 定的最大形變方向進(jìn)行的。由于只研究幾個(gè)主要方向上的形變，無(wú)論是計(jì)算復(fù)雜 度還是對(duì)存儲(chǔ)空間的要求都大大降低了。從數(shù)學(xué)角度講，首先是從訓(xùn)練集中估計(jì) 得到平均形狀7和形體協(xié)方差矩陣S:<formula>formula see original document page 10</formula>這樣， 一個(gè)合理的形狀可以被表達(dá)為
y = ; +①c (8)
式中向量c由對(duì)應(yīng)于各個(gè)特征向量的形變參數(shù)組成。矩陣o則包含了一組由協(xié)方 差矩陣S分解得到的特征向量，這組向量可線性組合成任意的新形狀y。根據(jù)PCA 理論，c中的每個(gè)參數(shù)都控制形體在某一方向上的形變，并且，各個(gè)方向都彼此 獨(dú)立，參數(shù)按照形體在對(duì)應(yīng)方向上形變由大到小順序排列。本發(fā)明考慮在一個(gè) PDM模型僅包含平均形狀、特征值和特征向量而不包括原始訓(xùn)練集的數(shù)據(jù)的情況 下進(jìn)行插值運(yùn)算。艮P:
M = <;,<D,X> (9)
如上所述，單個(gè)PDM模型可以由一組柔性對(duì)象的形體獲得，例如，我們研 究的醫(yī)院患者的心臟圖片(一些被觀測(cè)的柔性體個(gè)例)。然而，有時(shí)我們要觀測(cè) 這些柔性體隨時(shí)間的變化情況，這樣可以得到一系列的形狀信息，例如心臟在整 個(gè)活動(dòng)周期內(nèi)的三維變形。而醫(yī)院現(xiàn)有設(shè)備只能獲得若干時(shí)刻的三維圖像，即對(duì) 應(yīng)于某些Ti時(shí)刻的模型Mi。
對(duì)一個(gè)特定的動(dòng)態(tài)形變物體，假設(shè)我們已經(jīng)在離散的時(shí)間坐標(biāo)(T1,T2,...) 上建立起一組PDM模型(Mh M2,...)。我們可以把問(wèn)題描述為給出已知的兩個(gè) PDM模型，TJ寸刻的模型Ma和Tb時(shí)刻的模型Mb (模型已知意味著我們知道了 平均形狀"和"、特征矩陣^a和屯b、以及特征值"和")，如何得到在t (TaS t$Tb)時(shí)刻的中間模型Mt的相關(guān)參數(shù)。這些參數(shù)有平均形狀、特征值和特征向 量，艮卩M產(chǎn)巧t,①t, Xt > (圖1)。這里我們假設(shè)在Ta到Tb的時(shí)間段中， 一個(gè) 形體從Yai連續(xù)變化到Y(jié)bl， 二者都是"維的，如Yai = (xail, xai2, xai3, ...， xain)T。這 里還定義了兩個(gè)常數(shù)
10a=i^ (10) 6 = 1-a (ii)
a和b反映了中間模型相對(duì)于Ma和Mb的距離。
(1) 中間模型的平均形狀
根據(jù)以上的定義和假設(shè)，用線性插值(可滿足大部分的實(shí)際應(yīng)用)得到Y(jié)ai
和Ybi之間的中間模型形狀為
yti = ayai + bybi (12) a和b由式(10)和(11)定義，即a+b=l。
插入的平均形狀為
S /=1 S ,=i S S 7=1
所以，我們得到以下結(jié)論
定理1線性插值時(shí)，插入的PDM中間模型平均形狀是相鄰模型平均形狀的 線性組合，其參數(shù)由式(13)確定。
(2) 中間模型的形變
相似的，中間PDM模型的形變可以這樣算出
~ = a - y, = d + ) - Oya + 6y6) _y。)+60h - yj
(14)
= ,+ ，
從上式可以得到以下結(jié)論
定理2線性插值時(shí)，插入的PDM中間模型的形狀變化速度是相鄰模型形變 速度的線性組合，其參數(shù)由式(14)確定。 (3)中間模型的協(xié)方差矩陣
中間模型的協(xié)方差矩陣可以由下式推導(dǎo)得出
=^7lX。,A; +^7Z~"A; +41>戸,；+41>盧《"
(15)令
因?yàn)?br> 一 1 ,=1
(16)
s 一丄,=i ^ 一 1 ,=i
,A乙
顯然有
S。6 = S6fl
一"盧
(19)
(17)
(18)
最終，得到以下推論。
推論l線性插值時(shí)，中間模型的協(xié)方差矩陣可由以下表達(dá)式得到 S, = a2S。 + 62S4 + WS。6 + a《 (20)
這里Sab表示模型MU和Mb的交變矩陣。我們知道協(xié)方差矩陣是通過(guò)對(duì)同一 時(shí)間段內(nèi)不同個(gè)體的采樣數(shù)據(jù)進(jìn)行統(tǒng)計(jì)分析得出的，而交變矩陣則描述了在不同 采樣時(shí)間得到的兩個(gè)對(duì)象模型的形變過(guò)程。
本發(fā)明的有益效果主要表現(xiàn)在能夠形成實(shí)時(shí)分布的中間過(guò)渡模型、
模型精度高、實(shí)用性良好。


圖1是從時(shí)間序列中的臨近模型插值到得的瞬時(shí)PDM模型的示意圖。 圖2是典型的心動(dòng)周期分為六個(gè)階段的示意圖。 圖3是心臟不同階段特征的六個(gè)PDM模型的示意圖。 圖4是插值生成新的PDM模型以進(jìn)行柔性體分析的示意圖。 圖5是通過(guò)插值生成的PDM心臟模型的時(shí)間平滑變化的示意圖。
具體實(shí)施例方式
下面結(jié)合附圖對(duì)本發(fā)明作進(jìn)一步描述。
參照?qǐng)D1 圖5， 一種基于插值算法的柔性體點(diǎn)分布模型的形狀、形變及協(xié)方差矩陣計(jì)算方法，包括以下步驟
1) 、載入兩個(gè)相鄰時(shí)刻Ta和Tb的柔性體分布模型Ma和Mb，其計(jì)算式為.-
Ma = <"，①a, "> ， Mb = <",①b, Xb〉；
2) 、通過(guò)所述柔性體分布模型Ma和Mb，分別計(jì)算平均形狀"和"、特征矩陣①a
和①b、以及特征值"和Xb，設(shè)定柔性體從Ta到Tb呈連續(xù)變化，
2.1).中間模型的平均形狀
用線性插值得到y(tǒng)^和ybi之間的中間模型形狀yti為
yti = ayai + bybi (12)
a和b由式(10)和(11)定義:j。<formula>formula see original document page 13</formula>(10)
(11)
a和b反映了中間模型相對(duì)于Ma和Mb的距離； 插入的平均形狀為
<formula>formula see original document page 13</formula> (13);
2.2)、中間模型的形變
中間模型的形變由下式算出
<formula>formula see original document page 13</formula>(14)
2.3)、中間模型的協(xié)方差矩陣 中間模型的協(xié)方差矩陣由下式推導(dǎo)得出
<formula>formula see original document page 13</formula>
(15)s^41:a戸a;, (16)
因?yàn)?br> s*-;為,-&T —tvA; (18)
s一i ,=i 丄,=i
顯然有
s。,s:。 (19);
3)、根據(jù)得到中間模型的平均形狀、形變和協(xié)方差矩陣，得到t時(shí)刻的中間模型 Mt。
在我們一個(gè)關(guān)于心臟圖像分析的研究項(xiàng)目中，可以根據(jù)心動(dòng)周期定義一個(gè)規(guī)
范化的時(shí)間坐標(biāo)軸。例如， 一個(gè)心跳周期可以分為六個(gè)生理階段心房收縮、心 室收縮、心室射血、心室舒張、心室快速填充及舒張后期休息(圖2)。也有人 進(jìn)一步將心室射血的過(guò)程分為射血前期和射血后期。左心與右心的運(yùn)動(dòng)過(guò)程相似 但不完全相同。主要的差別是，左心室與動(dòng)脈的收縮壓比右心的大三至四倍，而 且左右心泵血過(guò)程的時(shí)間方面也有顯著區(qū)別。由于整個(gè)心臟的空間形態(tài)在一個(gè)心 動(dòng)周期內(nèi)變化很大，而單一的模型很難準(zhǔn)確描述動(dòng)態(tài)的心臟運(yùn)動(dòng)過(guò)程，所以我們 需要建立一組模型來(lái)描述心動(dòng)不同階段的心臟形態(tài)。
在gated-SPECT的圖像研究中，可參考心電圖信號(hào)進(jìn)行時(shí)間規(guī)范化[19]。我們 希望通過(guò)對(duì)心臟形狀進(jìn)行統(tǒng)計(jì)后可以得到全部的六個(gè)模型，并且這些模型能用來(lái) 準(zhǔn)確地描述心臟的特定狀態(tài)，如心臟舒張末期和心臟收縮末期。圖3在標(biāo)準(zhǔn)心電 圖曲線中標(biāo)出了各個(gè)階段模型出現(xiàn)的位置。圖中我們已經(jīng)把壓力值幅度經(jīng)過(guò)了修 正，使得波峰處的壓力值為正常人體的收縮壓(120mm汞柱)，所以圖中也表示 一條典型的左心室壓力曲線。
如果在正常的心動(dòng)周期中有6個(gè)PDM模型，則式(6)中的Ti (代表某一相位 所在的時(shí)刻)可以定義為r—1000L5(f/^-L〃tfj」[ms〗 (37)
式中t是獲得gated-SPECT圖像的時(shí)間值，"是心率的的倒數(shù)，L，」表示不大于某 小數(shù)的最大整數(shù)。
然而，在實(shí)際的例子中，常會(huì)產(chǎn)生不同數(shù)目的相位(如圖4中有8個(gè))，這 是由于心臟圖像的采樣通常有固定的頻率，如每lOOins (圖4)。那么就有必要通 過(guò)插值得到采樣時(shí)刻的PDM模型，以期得到更好的三維心臟形狀的分割。
柔性體點(diǎn)分布模型插值算法，作者研究過(guò)程中分別采用計(jì)算機(jī)仿真及真實(shí)心 臟數(shù)據(jù)進(jìn)行了大量的實(shí)驗(yàn)。在仿真中，我們可以假設(shè)一組九維(或更高)的模型， yi = [Xli， x2l， ...， x9i]T。兩個(gè)PDM模型M1和M2的數(shù)據(jù)集都是隨機(jī)產(chǎn)生的。各 組數(shù)據(jù)都含有100個(gè)樣本，作為生成點(diǎn)分布模型的訓(xùn)練集，例如
一1 ioo i 100 一 —
y,區(qū)5乂' ， s'、5""')(")t
為生成一個(gè)插入模型，我們從S,計(jì)算M1的特征值及特征向量，對(duì)S2也同樣。 假設(shè)插入模型的左右距離為a=0. 3， b=0.7。為了進(jìn)行對(duì)比，我們還計(jì)算它們的交 變矩陣Sab和相應(yīng)的特征向量組。實(shí)驗(yàn)結(jié)果顯示，本中所述的方法能很好的計(jì)算出 插入模型的各種參數(shù)。由式(31)計(jì)算出的特征值誤差小于1.83%，由式(35) 計(jì)算出的特征值誤差小于0.05%。該仿真實(shí)驗(yàn)的0++程序已由所在實(shí)驗(yàn)室開發(fā)完 成。
用真實(shí)的心臟數(shù)據(jù)進(jìn)行了測(cè)試。通過(guò)灌注SPECT成像法得到的大量圖像數(shù)據(jù) 中，我們預(yù)先建立了左心室在不同時(shí)刻的五個(gè)P函模型。每個(gè)模型由246組病例 的真實(shí)圖像的訓(xùn)練得到(圖像主要由巴塞羅那醫(yī)院提供)，包括一個(gè)平均形狀的 VTK表達(dá)， 一組特征值矢量和特征向量。所有醫(yī)學(xué)圖像的處理過(guò)程和模型的插入 運(yùn)算也由都0++程序?qū)崿F(xiàn)(研究組共開發(fā)了約600MB程序代碼)。實(shí)驗(yàn)時(shí)，我們 在每?jī)蓚€(gè)相鄰的PDM之間均勻地插入四個(gè)中間模型，得到的實(shí)驗(yàn)結(jié)果令人非常滿 意。圖5顯示通過(guò)插值生成的PDM模型隨時(shí)間能夠平滑地變化左心室柔性模型。點(diǎn)分布模型作為一種描述柔性對(duì)象的統(tǒng)計(jì)學(xué)模型，已被證明在很多應(yīng)用上是 非常有效的，尤其是在醫(yī)學(xué)圖像分析領(lǐng)域。為了分析柔性體的動(dòng)態(tài)變化情況，本 發(fā)明首創(chuàng)研究了點(diǎn)分布模型的插值算法。研究結(jié)果表明，可由PDM的平均形狀、 特征值和特征向量等參數(shù)生成任意時(shí)刻的中間模型。事實(shí)上，插入模型的平均形 狀及形變就是它相鄰前后模型對(duì)應(yīng)項(xiàng)的線性組合，即卩,=^。+^,， A, ="A^+6~fa。其特征值則可以采用式"- a、a + b、b計(jì)算得到。最后，為求
出插入模型的協(xié)方差矩陣和特征向量，本發(fā)明分開了交變矩陣形式并給出了相應(yīng) 的計(jì)算方法。
權(quán)利要求
1、一種基于插值算法的柔性體點(diǎn)分布模型的形狀、形變及協(xié)方差矩陣計(jì)算方法，其特征在于所述計(jì)算方法包括以下步驟1)、載入兩個(gè)相鄰時(shí)刻Ta和Tb的柔性體分布模型Ma和Mb，其中Ma＝&lt;<overscore>y</overscore>a，Φa，λa&gt;，Mb＝&lt;<overscore>y</overscore>b，Φb，λb&gt;；2)、所述柔性體分布模型Ma和Mb，包括平均形狀<overscore>y</overscore>a和<overscore>y</overscore>b、特征矩陣Φa和Φb、以及特征值λa和λb，設(shè)定柔性體從Ta到Tb呈連續(xù)變化，則2.1)中間模型的平均形狀用線性插值得到y(tǒng)ai和ybi之間的中間模型形狀yti為yti＝ayai+bybi (12)a和b由式(10)和(11)定義<maths id="math0001" num="0001" ><math><![CDATA[ <mrow><mi>a</mi><mo>=</mo><mfrac> <mrow><mi>t</mi><mo>-</mo><msub> <mi>T</mi> <mi>a</mi></msub> </mrow> <mrow><msub> <mi>T</mi> <mi>b</mi></msub><mo>-</mo><msub> <mi>T</mi> <mi>a</mi></msub> </mrow></mfrac><mo>-</mo><mo>-</mo><mo>-</mo><mrow> <mo>(</mo> <mn>10</mn> <mo>)</mo></mrow> </mrow>]]></math></maths>b＝1-a (11)a和b反映了中間模型相對(duì)于Ma和Mb的距離；插入的平均形狀為<maths id="math0002" num="0002" ><math><![CDATA[ <mrow><msub> <mover><mi>y</mi><mo>&OverBar;</mo> </mover> <mi>t</mi></msub><mo>=</mo><mfrac> <mn>1</mn> <mi>s</mi></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>t</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><msub> <mi>y</mi> <mi>ti</mi></msub><mo>=</mo><mfrac> <mn>1</mn> <mi>s</mi></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><mrow> <mo>(</mo> <msub><mi>ay</mi><mi>ai</mi> </msub> <mo>+</mo> <msub><mi>by</mi><mi>bi</mi> </msub> <mo>)</mo></mrow><mo>=</mo><mfrac> <mi>a</mi> <mi>s</mi></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><msub> <mi>y</mi> <mi>ai</mi></msub><mo>+</mo><mfrac> <mi>b</mi> <mi>s</mi></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><msub> <mi>y</mi> <mi>bi</mi></msub><mo>=</mo><mi>a</mi><msub> <mover><mi>y</mi><mo>&OverBar;</mo> </mover> <mi>a</mi></msub><mo>+</mo><mi>b</mi><msub> <mover><mi>y</mi><mo>&OverBar;</mo> </mover> <mi>b</mi></msub><mo>-</mo><mo>-</mo><mo>-</mo><mrow> <mo>(</mo> <mn>13</mn> <mo>)</mo></mrow><mo>;</mo> </mrow>]]></math></maths>2.2)、中間模型的形變中間模型的形變由下式算出Δyn＝y(tǒng)n-<overscore>y</overscore>i＝(ayai+bybi)-(a<overscore>y</overscore>a+b<overscore>y</overscore>b)＝a(yai-<overscore>y</overscore>a)+b(ybi-<overscore>y</overscore>b) (14)＝aΔyai+bΔybi2.3)、中間模型的協(xié)方差矩陣中間模型的協(xié)方差矩陣由下式推導(dǎo)得出<maths id="math0003" num="0003" ><math><![CDATA[ <mrow><msub> <mi>S</mi> <mi>t</mi></msub><mo>=</mo><mfrac> <mn>1</mn> <mrow><mi>s</mi><mo>-</mo><mn>1</mn> </mrow></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><mrow> <mo>(</mo> <msub><mi>y</mi><mi>ti</mi> </msub> <mo>-</mo> <msub><mover> <mi>y</mi> <mo>&OverBar;</mo></mover><mi>t</mi> </msub> <mo>)</mo></mrow><msup> <mrow><mo>(</mo><msub> <mi>y</mi> <mi>ti</mi></msub><mo>-</mo><msub> <mover><mi>y</mi><mo>&OverBar;</mo> </mover> <mi>t</mi></msub><mo>)</mo> </mrow> <mi>T</mi></msup><mo>=</mo><mfrac> <mn>1</mn> <mrow><mi>s</mi><mo>-</mo><mn>1</mn> </mrow></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><msub> <mi>&Delta;</mi> <mi>yti</mi></msub><msubsup> <mi>&Delta;</mi> <mi>yti</mi> <mi>T</mi></msubsup> </mrow>]]></math></maths><maths id="math0004" num="0004" ><math><![CDATA[ <mrow><mo>=</mo><mfrac> <mn>1</mn> <mrow><mi>s</mi><mo>-</mo><mn>1</mn> </mrow></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><mrow> <mo>(</mo> <mi>a</mi> <msub><mi>&Delta;</mi><mi>yai</mi> </msub> <mo>+</mo> <mi>b</mi> <msub><mi>&Delta;</mi><mi>ybi</mi> </msub> <mo>)</mo></mrow><msup> <mrow><mo>(</mo><mi>a</mi><msub> <mi>&Delta;</mi> <mi>yai</mi></msub><mo>+</mo><mi>b</mi><msub> <mi>&Delta;</mi> <mi>ybi</mi></msub><mo>)</mo> </mrow> <mi>T</mi></msup><mo>-</mo><mo>-</mo><mo>-</mo><mrow> <mo>(</mo> <mn>15</mn> <mo>)</mo></mrow> </mrow>]]></math></maths><maths id="math0005" num="0005" ><math><![CDATA[ <mrow><mo>=</mo><mfrac> <msup><mi>a</mi><mn>2</mn> </msup> <mrow><mi>s</mi><mo>-</mo><mn>1</mn> </mrow></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><msub> <mi>&Delta;</mi> <mi>yai</mi></msub><msubsup> <mi>&Delta;</mi> <mi>yai</mi> <mi>T</mi></msubsup><mo>+</mo><mfrac> <msup><mi>b</mi><mn>2</mn> </msup> <mrow><mi>s</mi><mo>-</mo><mn>1</mn> </mrow></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><msub> <mi>&Delta;</mi> <mi>ybi</mi></msub><msubsup> <mi>&Delta;</mi> <mi>ybi</mi> <mi>T</mi></msubsup><mo>+</mo><mfrac> <mi>ab</mi> <mrow><mi>s</mi><mo>-</mo><mn>1</mn> </mrow></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><msub> <mi>&Delta;</mi> <mi>yai</mi></msub><msubsup> <mi>&Delta;</mi> <mi>ybi</mi> <mi>T</mi></msubsup><mo>+</mo><mfrac> <mi>ab</mi> <mrow><mi>s</mi><mo>-</mo><mn>1</mn> </mrow></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><msub> <mi>&Delta;</mi> <mi>ybi</mi></msub><msubsup> <mi>&Delta;</mi> <mi>yai</mi> <mi>T</mi></msubsup> </mrow>]]></math></maths>令<maths id="math0006" num="0006" ><math><![CDATA[ <mrow><msub> <mi>S</mi> <mi>ab</mi></msub><mo>=</mo><mfrac> <mi>ab</mi> <mrow><mi>s</mi><mo>-</mo><mn>1</mn> </mrow></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><msub> <mi>&Delta;</mi> <mi>yai</mi></msub><msubsup> <mi>&Delta;</mi> <mi>ybi</mi> <mi>T</mi></msubsup><mo>-</mo><mo>-</mo><mo>-</mo><mrow> <mo>(</mo> <mn>16</mn> <mo>)</mo></mrow> </mrow>]]></math></maths>因?yàn)?lt;maths id="math0007" num="0007" ><math><![CDATA[ <mrow><msub> <mi>S</mi> <mi>a</mi></msub><mo>=</mo><mfrac> <mn>1</mn> <mrow><mi>s</mi><mo>-</mo><mn>1</mn> </mrow></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><mrow> <mo>(</mo> <msub><mi>y</mi><mi>ai</mi> </msub> <mo>-</mo> <msub><mover> <mi>y</mi> <mo>&OverBar;</mo></mover><mi>a</mi> </msub> <mo>)</mo></mrow><msup> <mrow><mo>(</mo><msub> <mi>y</mi> <mi>ai</mi></msub><mo>-</mo><msub> <mover><mi>y</mi><mo>&OverBar;</mo> </mover> <mi>a</mi></msub><mo>)</mo> </mrow> <mi>T</mi></msup><mo>=</mo><mfrac> <mn>1</mn> <mrow><mi>s</mi><mo>-</mo><mn>1</mn> </mrow></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><msub> <mi>&Delta;</mi> <mi>yai</mi></msub><msubsup> <mi>&Delta;</mi> <mi>yai</mi> <mi>T</mi></msubsup><mo>-</mo><mo>-</mo><mo>-</mo><mrow> <mo>(</mo> <mn>17</mn> <mo>)</mo></mrow> </mrow>]]></math></maths><maths id="math0008" num="0008" ><math><![CDATA[ <mrow><msub> <mi>S</mi> <mi>b</mi></msub><mo>=</mo><mfrac> <mn>1</mn> <mrow><mi>s</mi><mo>-</mo><mn>1</mn> </mrow></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><mrow> <mo>(</mo> <msub><mi>y</mi><mi>bi</mi> </msub> <mo>-</mo> <msub><mover> <mi>y</mi> <mo>&OverBar;</mo></mover><mi>b</mi> </msub> <mo>)</mo></mrow><msup> <mrow><mo>(</mo><msub> <mi>y</mi> <mi>bi</mi></msub><mo>-</mo><msub> <mover><mi>y</mi><mo>&OverBar;</mo> </mover> <mi>b</mi></msub><mo>)</mo> </mrow> <mi>T</mi></msup><mo>=</mo><mfrac> <mn>1</mn> <mrow><mi>s</mi><mo>-</mo><mn>1</mn> </mrow></mfrac><munderover> <mi>&Sigma;</mi> <mrow><mi>i</mi><mo>=</mo><mn>1</mn> </mrow> <mi>s</mi></munderover><msub> <mi>&Delta;</mi> <mi>ybi</mi></msub><msubsup> <mi>&Delta;</mi> <mi>ybi</mi> <mi>T</mi></msubsup><mo>-</mo><mo>-</mo><mo>-</mo><mrow> <mo>(</mo> <mn>18</mn> <mo>)</mo></mrow> </mrow>]]></math></maths>顯然有<maths id="math0009" num="0009" ><math><![CDATA[ <mrow><msub> <mi>S</mi> <mi>ab</mi></msub><mo>=</mo><msubsup> <mi>S</mi> <mi>ba</mi> <mi>T</mi></msubsup><mo>-</mo><mo>-</mo><mo>-</mo><mrow> <mo>(</mo> <mn>19</mn> <mo>)</mo></mrow><mo>;</mo> </mrow>]]></math></maths>3)、根據(jù)得到中間模型的平均形狀、形變和協(xié)方差矩陣，得到t時(shí)刻的中間模型Mt。
全文摘要
一種基于插值算法的柔性體點(diǎn)分布模型的形狀、形變及協(xié)方差矩陣計(jì)算方法，包括以下步驟1)載入兩個(gè)相鄰時(shí)刻T_a和T_b的柔性體分布模型M_a和M_b，其中M_a＝<y_a，Φ_a，λ_a>，M_b＝<y_b，Φ_b，λ_b>；2)通過(guò)所述柔性體分布模型M_a和M_b，分別計(jì)算平均形狀y_a和y_b、特征矩陣Φ_a和Φ_b、以及特征值λ_a和λ_b，設(shè)定柔性體從T_a到T_b呈連續(xù)變化；2.1)中間模型的平均形狀；2.2)、中間模型的形變；2.3)中間模型的協(xié)方差矩陣；3)根據(jù)得到中間模型的平均形狀、形變和協(xié)方差矩陣，得到t時(shí)刻的中間模型M_t。本發(fā)明能夠形成實(shí)時(shí)分布的中間過(guò)渡模型、模型精度高、實(shí)用性良好。
文檔編號(hào)G06T17/00GK101639949SQ200910102188
公開日2010年2月3日 申請(qǐng)日期2009年8月20日 優(yōu)先權(quán)日2009年8月20日
發(fā)明者陳勝勇 申請(qǐng)人:浙江工業(yè)大學(xué)