基于因子分析法的超聲肝灌注定量分析初步研究

張 冀 丁明躍* 孟 燔，2 張旭明 尉遲明

1(華中科技大學生命科學與技術學院 圖像處理和智能控制國家教育部重點實驗室，武漢 430074)

2(武漢市第一醫(yī)院設備科，武漢 430070)

引言

超聲造影成像是一種通過注入造影劑來動態(tài)、清晰地顯示微細血管的新型成像技術。利用該技術可以獲得組織豐富的供血及血流灌注信息，這對臨床診斷組織病變，特別是惡性腫瘤具有重要的臨床意義。與增強CT和增強MRI相比，它具有實時性。目前該技術已被廣泛應用于局灶性肝結節(jié)(focal liver lesions，F(xiàn)LLs)的鑒別診斷［1］。

由于不同類型的FLLs具有各自獨特的增強模式，F(xiàn)LLs通常能被從灌注序列圖像中提取的血流時間-強度曲線(time-intensity curves，TICs)量化評價［2］。因此，分析 TICs對 FLLs的輔助診斷具有良好應用的前景，尤其對肝癌的早期探查和治療。提取TICs的最常見的方法是在圖像序列中計算感興趣區(qū)域(region of interest，ROI)中的平均像素強度。目前，醫(yī)院大多數(shù)超聲機配置了基于該方法的灌注量化商業(yè)軟件。然而，由于感興趣區(qū)域需要人工選定，因此，該方法的可靠性和準確性受到質疑。最近有實驗表明，ROI的在肝臟中的位置深度可以引起TICs參數(shù)的不穩(wěn)定性［3］。同時，由于肝實質位置的差異，所提取到的TICs的參數(shù)也會發(fā)生變化［4］。

為了克服人工操作的主觀性，近年動態(tài)結構因子分析法(factor analysis of dynamic structures，F(xiàn)ADS)被嘗試用于自動地提取肝灌注序列圖像中的特征性 TICs［5－6］。FADS是一種基于圖像序列的線性模型，一般用于提取灌注序列圖像中的幾個生理曲線和相應的因子圖［7－8］。每條因子曲線對應一張因子圖并具有因子圖中生理結構的時間變化特點。因此，該因子曲線可以被認為是相應因子結構的TIC。

確定因子的個數(shù)對于FADS是關鍵，一般可根據(jù)動態(tài)灌注過程中的生理成分數(shù)目進行估計。肝臟具有肝動脈和門靜脈雙重供血，則一般肝灌注模型包含3個在時間區(qū)域重疊的增強相，即動脈相(5～30 s)、門靜脈相(31～120 s)和延遲相(2 min后)［5］。大多數(shù) FLLs在動脈相和門靜脈相中具有獨特的增強模式［1］。有的腫瘤可以在早期動脈相出現(xiàn)離心或向心式充填的增強變化。有的惡性腫瘤部分結構也可以在早期動脈相出現(xiàn)提早增強，然后，增強結構的回聲在門靜脈相又快速減弱。因此，動脈相可以進一步劃分為早期動脈相和晚期動脈相［5］。文獻［5－6］認為 3個因子的分析可以較全面地反映大多數(shù)肝結節(jié)的灌注過程，則文中將重點探討FADS方法在FLLs灌注中進行3個因子分析的應用。

解的不惟一性是 FADS模型主要問題［7］，傳統(tǒng)FADS方法往往通過因子旋轉［8-9］和運用與空間分布或時間曲線基本函數(shù)相關的附加先驗知識來獲取唯一解［10］。Sitek等提出一種運用懲罰性最小方差目標函數(shù)法來校正對動態(tài)SPECT結構進行因子分析中的不確定解［11］。雖然這種方法不需要因子旋轉和先驗知識，但是目標函數(shù)中關鍵參數(shù)的估計限制了該方法應用于其他成像模式。針對目前FADS方法的不足，本研究提出一種新型的基于頂點搜尋的替代-近似(replace-approximation，RA)法。該算法不需要任何特定參數(shù)的估計和傳統(tǒng)的因子旋轉。RA算法已在體模實驗中進行初步驗證［12］，為了進一步驗證該法的有效性，采用該法對6個臨床病例進行因子分析，并初步探索了其臨床應用性。

1 FADS模型

N張序列圖像可組成一個N×M 2-D矩陣X。其中，N個時間點為列，M個像素排成行(M是圖片大小)。研究中某像素的強度隨時間變化構成該像素的 TIC，被稱作動態(tài)研究像素(dynamic study pixel，dixel)［8］。向量 xi(i=1，2，…，M)表示 i的動態(tài)研究像素。xi(t)(t=1，2，…，N)是像素i在時間 t時的信號強度。FADS 模型如下［8－9］:

或

其中

式(1)中，N×q矩陣 F由 q條因子曲線組成，而A為一個q×M的因子系數(shù)矩陣。ε是誤差矩陣。式(2)中，fk是第k條因子曲線，ak表示為第k張因子圖。每張因子圖代表了具有相應因子曲線的一個生理結構。由于殘差太小可以被忽略，式(1)可近似為

xi(t)可由平均曲線(t)和(q－1)個較大的主成分(principle components，PCs)構成，表示為

式中，uk是第 k條主成分曲線，vk是相應的系數(shù)向量，λk是對應第k個特征向量的特征值uk［8］。

所有的dixel通常可被不同方法進行預處理，X矩陣經過預處理后為矩陣Xn。當F確定后，空間系數(shù)矩陣A可以由F的偽逆算子計算出

A 可被顯示成圖像［5］。

從幾何角度看，當q=3時，所有的dixel必須位于一個三角形內，而此三角形的3個頂點為處于二維空間的正值域中的因子曲線。當q=4時，4條因子曲線對應了一個三維空間中四面體的四個頂點，在高維正值域中以此類推，可以找到多個頂點。然而，在正值域中，不止一個單純形(例如:三角形)包圍這些dixel，這導致了解的不惟一性［7］。

2 替代-近似算法

FADS法可用于分離具有不同時間曲線特性且部分重疊的生理結構［13］。重疊結構主要包含混合信號而未重疊結構含有單一信號，因此，未重疊結構對應的TIC吻合于真實的因子曲線。RA算法的主要思想為在這些未重疊區(qū)域中，找到最佳的一個dixel以近似于真實的因子曲線。然而，如何找到最佳曲線是關鍵。

該算法是基于頂點搜尋法，因此，其目標是獲取單純形的頂點。通常，因子個數(shù)不會超過3。雖然在二維空間中不止一個三角形圍繞點云，但希望找在正值域中找到一個包圍最緊的三角形。然而，實際中，由于噪聲影響，有些點會處于找到的三角形的外部。因此，設定所找的三角形應該包含大多數(shù)點。在本實驗中，所找的三角形包含了超過80%的點。

雖然從視覺上找頂點不難，但是提供一個自動頂點搜尋程序并不容易［8］。本方法總是從一維空間開始搜尋頂點，而不像其他方法從q-1維空間開始搜尋頂點［5－6，8－10］。RA 法的流程圖見文獻［12］，具體步驟如下:

步驟1:預處理:每個像素xi(t)先按以下公式標準化處理

然后，每個dixel xni中心化，生成一個新的矩陣Xn。

步驟2:第一個頂點搜尋

1)運用主成分分析法(principal components analysis，PCA)將Xn映射至一維空間;

2)在一維空間中，只有兩個方向:正和負方向。兩個最優(yōu)解應該分別位于這兩個方向上。假設AP是在正方向上Xn中k個dixel pi(i=1，2，…，k)的集合，而 AN是在負方向上Xn中 h個dixel gi(i=1，2，…，h)的集合。首先，參照 Bazin選擇法［9］找到頂點fe，fe∈AN。該點是在負方向上離原點最遠點，并且將fe作為第一個頂點估計初值。通過式(4)，利用fe的坐標合成曲線 fr。其次，建立集合 C含有 u個dixel ci，ci∈C(i=1，2，…，u)，C ? AP，‖ci‖2≥向量 c和 f合ir并組成矩陣F。然后，當A中負值的系數(shù)個數(shù)達到最小時，由式(5)得到C中最佳頂點。該點為第一個頂點 f1。

步驟3:第二個頂點搜尋:F由f1和AN中任一dixel組成。然后，第二個頂點 f2可以由相同方法得到。

步驟4:第3個頂點搜尋:運用PCA將Xn再次映射到二維空間。同理，在獲得兩個頂點后，運用步驟2中相同準則確定第3個頂點f3。假設Z是該空間中所有dixel的集合，d是原點到直線的距離。建立集合 W，wi∈W(i=1，2，…，v)，W ? Z，‖wi‖2≥d。此時，F(xiàn)由f1，f2和 wi組成。因此，當 A中負值的系數(shù)個數(shù)達到最小時，由式(5)得到W中最佳頂點。該點為第3個頂點f3。

3 材料和方法

3.1 臨床病例分析

實驗中所采用的6個肝細胞癌(hepatocellular carcinomas，HCC)病例由華中科技大學附屬醫(yī)院協(xié)和醫(yī)院提供，均通過穿刺或手術標本進行病理診斷。HCC是一種最常見的原發(fā)性肝細胞癌。其典型增強表現(xiàn)為造影劑在腫瘤內“快進快出”。即腫瘤在早期動脈相期間明顯強化，隨后腫瘤區(qū)域內的造影劑在門靜脈相或延遲相時迅速流出，病灶內回聲明顯低于周圍正常的肝組織［1］。腫瘤大部分區(qū)域的平均通過時間(mean transit time，MTT)短于周圍正常組織。病例由Acuson Sequoia 512超聲探頭及Sequoia超聲系統(tǒng)(Siemens ultrasound，Mountain View，CA)采集。選擇最佳切面后進入實時超聲造影模式，調節(jié)聲輸出功率達低機械指數(shù)狀態(tài)(MI＜0.20)。經患者肘部團注造影劑約2.4 mL SonoVue(Brraco，Milan，Italy)。注射后，動態(tài)灌注過程錄制40～60 s，幀率7～10 fps。在檢查中，引導病人盡量保持規(guī)律地淺呼吸狀態(tài)。實驗采用文獻［14］中的呼吸運動校正策略對每個臨床數(shù)據(jù)進行預處理，減少圖像序列中的呼吸運動影響，從而可增加FADS分析的準確性。該校正策略通過模板匹配后，采用一個雙選擇法來選取待分析的圖片。此方法屬于一種幀選擇法，它從大量的圖像序列中去除許多與模板偏離較遠的圖像，因此，每個圖像序列在選擇后的最終圖像數(shù)量遠遠少于最初數(shù)量。實驗中，校正后的圖像數(shù)量不超過最初的15%，大約為10～20幀。所有程序均采用 Matlab(The Mathworks Inc.，Natick，MA)編寫。

3.2 ROI測定標準

當要分離的生理結構區(qū)域可通過先驗知識估計時，ROI通常置于該結構區(qū)域內，此時提取的ROI測定曲線往往被視作因子曲線的參考標準［11，13］，。在參照文獻［5］中，大多數(shù)的FLLs增強變化依次發(fā)生在早期動脈相、晚期動脈相和門靜脈相，因此，本實驗中所進行分析的3個因子也確定為上述三相。它們所對應的主要結構分別為腫瘤成分，肝臟的動脈成分和肝實質成分。本實驗對這3個生理結構采用先驗知識進行估計，然后將ROI置于這3種成分中，則可提取到相對應的標準曲線。將ROI置于腫瘤中心則提取到早期動脈相的標準曲線。由于肝動脈血流先流入肝臟，將ROI置于動脈血管成分中或最先增強的腫瘤以外的區(qū)域則提取到肝晚期動脈相的標準曲線。肝臟中四分之三的血流是由肝門靜脈供應的［15］，因此，將 ROI置于肝實質中，并且避開明顯的血管，提取到肝門靜脈相的標準曲線。位于腫瘤中心的ROI稍大，取為10像素×10像素，使得提取到的曲線能充分反應腫瘤區(qū)域增強變化。其他兩處的ROI取為4像素×4像素。

3.3 因子曲線的量化評價

為了評價因子曲線的準確性，本實驗測定了因子曲線和ROI標準曲線間的相關系數(shù)。實際病例中肝灌注的曲線滿足數(shù)學灌注模型，因此，采用一個修正的對數(shù)正態(tài)分布模型［16］擬合病例中提取的因子曲線和ROI曲線，為

O，A，m和 s是擬合參數(shù)。O是偏移值，A是幅度參數(shù)。m和s分別是均值和正態(tài)分布自然對數(shù)t的標準差。當曲線擬合后，實驗計算了因子曲線與ROI曲線之間3個參數(shù)差值。這3個參數(shù)分別為峰值時間(peak time，PT)，注入時間(wash-in time，WIT)和MTT。其對應的參數(shù)差值分別為 PT-D，WIT-D，MTT-D。

3.4 臨床病例中因子分析的區(qū)域選擇和曲線標準化比較

肝臟的超聲圖像一般含有不均一結構，如肝實質，腎臟，膽囊或肋骨陰影等。不同類型的組織結構具有獨特的增強模式。過多的區(qū)域會增加能夠充分描述分析數(shù)據(jù)的因子個數(shù)，這會降低待分析因子個數(shù)的準確性［13］。因此，在 FADS分析前，本實驗只選擇圖像局部區(qū)域進行分析，而非肝臟組織被遮罩。

為了進行客觀評價，將實驗中提取的因子曲線和ROI曲線縮放至［0，1］，并且每條曲線的基線值被置為0。

4 結果

在運用RA法對HCC病例進行3個因子的FADS分析中，均可提取到3條不同的因子曲線。這些曲線均由一個較急速的上升部分，一個明顯的峰和一個較緩慢的下降部分構成。這3條曲線峰值時間均有一定延遲。根據(jù)曲線上峰值出現(xiàn)時間的順序，可區(qū)分為與早期動脈相，晚期動脈相和門靜脈相因子結構相對應的因子曲線。圖1顯示了一例HCC的灰度圖。該病例的三相不同的因子曲線如圖2所示。圖3顯示了與因子曲線相對應的因子圖。早期動脈相因子圖中腫瘤的部分區(qū)域開始增強，其對應的因子曲線具有較短的 MTT，表明腫瘤中造影劑快速注入并且快速流出。在晚期動脈相因子圖中，腫瘤完全增強，并且周圍組織也增強，表明有肝動脈血流開始進入周圍組織區(qū)域。其對應的因子曲線峰值時間延遲于早期動脈相因子曲線，而且因子曲線的下降部分也相對變緩。最后，在門靜脈相的因子圖中腫瘤回聲明顯下降，而周圍組織回聲并未完全減弱，有部分組織開始增強，表明有門靜脈血流進入周圍組織，而腫瘤中的造影劑已經快速流出。上述病灶的增強變化呈典型的“快進快出”的HCC增強模式。該因子圖對應的因子曲線的峰值時間延遲于晚期動脈相因子曲線，而且因子曲線的下降部分變將更加緩慢。3條曲線在因子圖中對應的權重分別為28%、38%和34%。

本實驗測定了提取的因子曲線與ROI曲線的相關系數(shù)和幾個參數(shù)差值。如表1所示，RA法提取的因子曲線與對應的ROI測定曲線的相關系數(shù)均值達0.91左右，顯示出較高的相關性。在表2中，PT和WIT的參數(shù)差值在5 s以內，MTT參數(shù)差值偏大，尤其是門靜脈相。總體來說，RA法提取的因子曲線和其對應的因子圖具有一定的生理意義。

圖1 HCC樣例的一張超聲造影圖和對應的灰度圖(在右邊灰度圖中，腫瘤由一個ROI勾畫)Fig.1 One of CEUS images and the corresponding greylevel image for one sample from HCC(in the right greylevel image，the tumor is delineated with a ROI)

圖2 RA法提取的因子曲線(HCC樣例)。(a)早期動脈相;(b)晚期動脈相;(c)門靜脈相。所有的曲線經過擬合Fig.2 Factor curves extracted by the RA method(one sample from HCC).(a)early arterial phase;(b)delayed arterial phase;(c)portal vein phase.All the curves are fitted

圖3 RA法提取的因子圖(HCC樣例)。(a)早期動脈相;(b)晚期動脈相;(c)門靜脈相腫瘤位于因子圖的中心Fig.3 Factor images obtained by the RA method(the sample from HCC).(a)early arterial phase;(b)delayed arterial phase;(c)portal vein phase.The tumors are in the center of the factor images

表1 采用RA法從所有病例中提取的因子曲線的相關系數(shù)值(平均值±標準差)Tab.1 Correlation coefficients of the factor curves extracted from all the cases by the RA method(mean±SD)

表2 采用RA法從所有病例中提取的因子曲線與對應ROI曲線的參數(shù)差值(平均值±標準差)Tab.2 Parameter differences between the factor curves extracted from all the cases by the RA method and the corresponding ROI curves(mean±SD)

5 討論和結論

FADS方法最常見的問題是難以找到一個金標準來評價FADS的結果。該評價通常依賴于結果的合理性和其具有的生理意義。因此，F(xiàn)ADS的目標是獲取容易解釋的結果。實驗結果表明新型的RA算法從HCC超聲灌注病例中提取到了3條具有生理意義的因子曲線和因子圖，初步證明了RA算法在FLLs灌注的量化分析具有可行性，真實的因子曲線難以提出，因此，為了評價RA法提取因子曲線的準確性，實驗采用ROI曲線作為參考曲線。由于提取ROI曲線的生理結構需通過先驗知識估計，而且人工選擇具有一定主觀性，因此，表2結果顯示提取的因子曲線與ROI曲線的主要參數(shù)還存在一定差異。這種差異是否影響腫瘤的量化分析有待進一步研究。

實際中，肝臟結構比較復雜，雖然除了主要的幾個因子結構外還存有其他一些非均一性結構，使得臨床數(shù)據(jù)中存在一定量的噪聲成分。噪聲成分的多少與所選取的區(qū)域結構有關。此外，由于實驗采用的數(shù)據(jù)是在病人自由呼吸運動下采集的，雖然通過呼吸運動校正方法可盡量減少數(shù)據(jù)中的呼吸運動影響，但是這種影響不可能被徹底消除。因此，呼吸運動的影響也增加了臨床數(shù)據(jù)中的噪聲成分。上述因素均導致因子圖中的獨立結構不可避免地被其他幾個結構所影響。然而，實驗結果顯示，三相因子圖中腫瘤的增強模式變化比較明顯，符合典型的惡性腫瘤造影表現(xiàn)特點，這表明主要的因子結構仍然可被分離。而且，RA法提取的因子圖能反應腫瘤的灌注特點，有助于醫(yī)生對FLLs的診斷。Renault等利用HCC的3張因子圖合成一張參數(shù)圖［6］，從而形象地描述了該腫瘤動態(tài)灌注的過程。由此，值得進一步探索因子圖的臨床應用價值。

考慮檢驗RA法的臨床實用性，本實驗采用自由呼吸的病例數(shù)據(jù)。實驗結果發(fā)現(xiàn)RA方法能抵抗一定的呼吸運動影響，提取到了具有生理意義的因子曲線和因子圖，體現(xiàn)出較強的魯棒性。本方法不需要任何特定參數(shù)的估計或因子旋轉，因此，這種簡單實用的方法對于FLLs灌注的量化分析具有一定潛力，但是今后還需在更多數(shù)量以及其他類型的FLLs中進一步驗證RA法的有效性。

［1］Bartolotta TV，Taibbi A，Midiri M，et al.Focal liver lesions:contrast-enhanced ultrasound［J］.Abdom Imaging，2009，34:193－209.

［2］Huang-Wei C，Bleuzen A，Bourlier P，et al.Differential diagnosis of focal nodular hyperplasia with quantitative parametric analysis in contrast-enhanced sonography［J］.Invest Radiol，2006，41:363 －368.

［3］Ignee A，Jedrejczyk M，Schuessler G，et al.Quantitative contrast enhanced ultrasound of the liver for time intensity curves-Reliability and potential sources of errors［J］.Eur J Radiol，2010，73(1)，:153 － 158.

［4］李杰，董寶瑋，于曉玲，等.肝不同部位灰階超聲造影定量分析的實驗研究［J］.中華超聲影像學雜志，2004，13(6):460－462.

［5］Lueck GJ，Kim TK， Burns PN，et al. Hepatic Perfusion imaging using factor analysis of contrast enhanced ultrasound［J］.IEEE Trans Med Imaging，2008，27(10):1449 －1457.

［6］Renault G， Tranquart F，Perlbarg V， et al. A posteriori respiratory gating in contrast ultrasound for assessment of hepatic perfusion［J］.Phys Med Biol，2005，50(19):4465 －4480.

［7］Houston AS，Nijran KS.Constraint problems in factor analysis of dynamic structures in nuclear medicine［C］//Proceedings of IEE 3rd International Conference on Image Processing and its Applications.London:IEE，1989:333－337.

［8］Barber DC.The use of principal components in the quantitative analysis of gamma camera dynamic studies［J］.Phys Med Biol，1980，25:283－292.

［9］Di Paola R，Bazin JP，Aubry F，et al.Handling of dynamic sequences in nuclear medicine ［J］.IEEE Trans Nucl Sci，1982，29:1310－1321.

［10］Buvat I，Benali H，Rouin F，et al.Target apex-seeking in factor analysis of medical image sequences ［J］.Phys Med Biol，1993，38:123－138.

［11］Sitek A，Gullberg GT，Huesman RH.Correction for ambiguous solutions in factor analysis using a penalized least squares objective［J］.IEEE Trans Med Imaging，2002，21(3):216 －225，.

［12］Zhang Ji， Ding Mingyue， YuchiMing， etal. Replaceapproximation method for ambiguous solutions in factor analysis of ultrasonic hepatic perfusion［C］//D'hooge J，McAleavey S A.，eds.Proceedings of SPIE Medical Imaging.San Diego:SPIE，2010，7629:1－11.

［13］Sitek A，Di Bella EVR，Gullberg GT.Factor analysis with a priori knowlege—application in dynamic cardiac SPECT ［J］.Phys Med Biol，2000，45:2619 －2638.

［14］Zhang Ji，Ding Mingyue，Meng Fan，et al.Fast algorithm for respiratory motion correction in free-breathing contrast-enhanced ultrasound imaging［C］.//D'hooge J，Doyley M M，eds.Proceedings of SPIE Medical Imaging.Orlando:SPIE，2011.7968－40.

［15］Schenk WG Jr，McDonald DJ，McDonald DK，et al.Direct measurement of hepatic blood flow in surgical patients:with related observations on hepatic flow dynamics in experimental animals［J］.Annals of Surgery，1962，156:463 －471.

［16］Rognin N，F(xiàn)rinking P，Messager T，et al.A new method for enhancing dynamic vascular patterns of focal liver lesions in contrast ultrasound ［C］//Marjorie Passini Y，eds.IEEE Ultrason Symp，New York:IEEE，2007:546－549.