TWI378402B - Image preprocessing system for 3d image database construction - Google Patents

Description

1378402 六、發明說明： 【發明所屬之技術領域】 /本發明係、有關於影像處理，特定而言係有關於用於建 構影像貧料庫之影像資料預處理系統。 【先前技術】 磁振造影（magnetic resonance imaging)係利用核磁共 ，（nuclear magnetic resonance; NMR)現象產生影像。當物 貝例如人體組織經歷一均勻磁場（極化場（ροΙΜζ— ^eld)B0)時’組、織中之各別自旋磁矩嘗試與上述極化場對 準仁以不規則順序及其特有之拉莫爾（Larm〇r)頻率繞極 化場之方向旋進。若物質或組織經歷一於χ·γ平面且近拉 f 爾（Larmor)頻率之磁場（激發場（excitati〇n 則 淨對準矩Mz可旋轉或”傾斜，，進入χ_γ平面，以產生淨橫 ^矩激發自旋係發射—信號，且於激發信號&終 後’信號可予以接收並經過處理以形成一影像。 .當利用此類信號產生影像時，磁場梯度（Gx、Gy及Gz) 予=利用。一般而言，欲造影之區域係藉由-連串之測 ^盾％:以知描’於上述測量循環中該些梯度係根據所使 ：特疋疋位方式而改變。每一測量於此領域中係稱為― =點“，且視點之數量係決定影像之品質。接收到之核 :_MR)信號或視點或k空間樣本之最終集合係經過 ^位化及處理以利用若干已知重建技術其中之一來重建影 從磁振造影（magnetic resonance imaging)之創始以 3 1378402 來，技景;^重建方法已廣為人知。投影重建方法係以一系列 之視點取樣k空間，上述視點係從k空間之中心向外延伸 之放射狀線所取樣，而非在直線或直角以傅立葉造影形成 之掃描圖樣中取樣k空間。若獲得不足之視點數量，則在 重建影像中會產生斑驳之人造影像。

除上述感應器圍繞目標旋轉以獲得用於重建之多角度 視點之影像擷取協定之外，已有新協定被採用以透過共軛 焦顯微鏡擷取多片段影像，例如，美國申請號第ιι/ΐ69,89〇 號，標題為「Bio-EXpression System _ The 〇f I S_」。此樣方法係適合於操取高解析度之生物影像，其 係細胞等級且用於三維細胞網路之重建例如腦部，以利； 研究器官之生理學及病理學。因此，藉由系統化荒集整體 腦部中之不同神經元以及於具有座標參考之標準腦模型系 統内合理地重分配其適當位置而重建出腦部中之整體神經 網路係為一似真任務。換言之’若所有之三維神經元影 ^以演异法予以配置到具有解剖學上重要性之方位，則整 體腦部中之神經元之三維影像資料庫係為可行。當2 st 1:二29曰提出申請之標題為「Aq_S Tissue Clearing J運國專利第M72,216B1號所述之技術係予以 ^傻右/ϋ，為迴路重建所需要處理之神 ，，办像有數十萬（若為果蠅）或甚至 三維影像資料預處理之一系列搭配的演算(二= 令之神經影像之分類一致性而言係為^法對於在資料庫 人類及其他生物例如昆蟲令之神經組織包含具有延伸 4 /M02 P伤之神經元，其係用以形成傳遞電化學信號之 能:纖Γ二Π: ’例如以非常高神經密^義= 結。神^ 之轴突纖維束之神經網路所連 、工、維東因此與腦部區域之功能相當有關。 關之==之診斷、腦部手術之計畫以及其他神經學相 資^動以及腦部運作之研究活動可從詳細之解^ 汛例如軸突纖維及纖維束之追蹤而獲益。 。〗予 =，本㈣提供―料之影像預處理“， 【發::内或容7之臨床應用之三維腦部影像資料庫。 元二= 之:觀點：’影像預處理系統包含處理單 段資料，其中影像預處理單元包:f!;=處理影像片 同影像片段資料間之遺漏數以估計不 只叮 <』< 1¾ /局数值，以及影德人 片段資料接合成經接合之影像資料二及資::不 轉二處理單元以儲存經預處理之影像片段資料 影像資料產生標準座標，自以為經接合之 資料對準於鮮影像，以及神經元追經接合之影像 之影像資料t之神經元之十心線:°追_組以追縱經接合 蟲例==係為本影像預處理系統可用於建立昆 維神經影像資料庫。外/β 嶋卿r)之整體腦部之三 本發明之另一優點係為本影像預處理系統可促進依據 5 I378402 細胞神經網路之腦科學研究的進展，且其可整合入功能模 組。 、 此類及其他優蹄以下較佳㈣狀敘述並伴隨後附 圖式及申請專㈣㈣使讀者得以清楚了解本發明。 【實施方式】 本發明將以較佳實施例及觀點加以敘述，此類敘述係 解釋本發明之結構及程序，僅用以說明而非用以限制本發 明之申請專利範圍。因此，除說明書中之較佳實_以外: 本發明亦可廣泛實行於其他實施例中。 本發明係揭露-種影像預處理系統。參照第 發明之影像預處理系統10包含處理單元⑻，例如：央卢 ^(central processing unit； CPU) ^ 102 理如第一圖所示’資料*102及影像預處 要-… 至處理單元101。於此系統中亦需 記憶體、快閃記憶體、輸出入介面等。 以避免模糊本發明之主題。於本發明==敘述 像預處理單元103包含 貫丨’i中，影 明2、標準化模組二3象自=^031、影像接合模組 追縱模組1035。影像預處理單元10/:用1二：及神經元 來源輸人之様太^ 1G 3係心預處理從外部 缔箱南 MSamplebrain)之影像片段資料，以浐屮 經預處理之影像片段資斜5次… 貝討以輸出 以儲存從影像預處理單 ^ 102。資料庫102係用 像片段資料。影像預=理之樣本腦影 处主早兀1〇3之影像融合模組1〇31 6 係用以估計不同之樣本腦影 影像預處理單元103 又、/ a1逍漏數值。 之樣本腦影像片段資料：杈組1032係用以接合不同 沭X η 、 使其成為經接合之影像資粗 u 迷不同之樣本腦影像片段資：“枓，上 重疊。影像預處理單元 j 疋王度上係為 像接合模組㈣所饋入之；:準化极組1〇33係響應從影 以用於u之貝料，且係用以產生一標準座释 動針進*述經接合之影像資料。影像預處理單元103之自' :=組1034係，合至標準化模組1033,且係 广接合之影像資料對準於一標 拉型’且將經對準之影 棕準細 標準座標結合。影===:°33所產生之 咖係用以追縱上述經接二追縱模組 線。於神經元之中心線被 中: 料庫102。 办傢貝#係輸出至資 於本發明之-實施财，影像融合模 位型暨形狀型插補峰一。一〜^ terpolation)，其結合影像對位及灰階形狀型插補 (gray-levelshape_basedinterp〇iati〇n)。因此，影像融合模 :且1031包含-影像對位模組以及一灰階形狀型插補模 組。首先’同一目標係從多個方向掃描，以獲得三垂直非 均勻之體積資料集。第二’此類資料集係予以對位且施加 上述方法’每一掃描之貢獻係經過確認以適當地併入最終 、:果。實驗結果將敘述於下’其顯示上述方法能有效估計 送漏數值。本發明發展出-架構以三維對位及灰階形狀型 1378402 j $ 4 ' 插補為基礎從多個非均勻體積資料集重建三維模型。於此 架構中，將會實施一用以結合從三個不同方向獲得之三體 積資料集之方法，以產生單一高解析度三維模型。此插補 架構包含一組預處理技術例如仿射轉換（affine transform)、共同資訊量（mutuai jnf〇rmati〇n)、最小平方 (leastsquares)、B-平滑曲線（B-splines)等，以將輸入體積 201初始化成經對位體積202,接著透過萃取形狀、計算最 終邊界及強度、邊界交集（b〇undary intersecti〇n)、遺漏資 •料近似（missing data approximation)形成重建體積 2〇3，如 第二圖所示。為計算插補平面中一位置之強度值，採取下 列步驟：首先，從三個不同方向獲得三異向性體積資料集。 於此架構中，此三方向通常係彼此相互垂直。第二，對此 三體積貧料集施行空間轉換（space transf〇rmati〇n)並將其 對，。第三，萃取每-體積之形狀資訊且將其合併以找尋 出每一灰強度（gray-intensity)之適當界限。接續，高解析 •度體積遂可產生。 對位 影像對位係指用以建立二或多影像之特徵之間之空間 -對應之程序。於實驗中，影像係藉由利用相同之造影方法 ·_㈨“…以）及對位方法而獲得，包含洞悉分割及對位套件 (ITK; Insight Segmentation and Registrati〇n 丁⑽脑）中之共 同資訊量（mutual inf0rmati0n)。此方法係運用共同資訊量 (m_al infc)mati()n)，以評估兩影像中之對應立體像素之 影像強度之間之料依存性。若影像錢何上係相互對 8 1378402 準，則其統計依存性則視為最大。對位架構係如第三圖所 示。影像對位模組包含度量元件、内插元件、轉換元件以 及優化元件。其係為一種自動且基於強度之測量，無需地 標或特徵之定義。共同資訊量（mutual information/认別），參照方程式⑴，係藉由利用 Kullback-Leibler 測量方法測量共同分佈（j〇int distributi〇n) (&〆〇»)以及與完全獨立相關之分佈⑻0少户#刀之間之 距離而測量得Α及Β之依存性程度。 I{A,B) = YjPAB{a,b)\Qg PAB(a,b) —Wl.PB[b)

共同資訊量係藉由方程式（2)而關聯至熵，其中H(A)、 H(B)係各別為a及B之熵。H(A，B)係指共同熵，而h(a丨b) 或H(B|A)係代表已知b下A之條件熵或已知A下b之條 件熵。共同資訊量（MI; mutual inf〇rmati〇n)對位標準係敘述 影像係藉由空間轉換而在幾何上相對準，使其共同資訊量 鲁（mutual informationX/MA)為最大。 Ι(Α,Β) = Η(Α) + Η{Β) - Η(Α,Β) H Η^\Β) = Η{Β)-Η{Β\Α) (2) 於第二圖令，度量元件係提供固定影像被經轉換之移 =影像所匹配之程度之測量。轉換元件係代表將來自於固 :影像空間之點空間映射至移動影像空間内之點。關於此 仿射轉換係予以應用至轉換元件。内插元件係經過設 強二估：於非網格位置(_-grid P〇Siti〇nS)上之移動影像 X於此係利用三線性内插。優化元件係用以搜尋特 9 1378402

I * » 定參數以重複達到最大共同資訊量（m_ai inf0rmation)(/M，万乃。藉由利用上述方法，三個體積係對位 於同一座標系統中。 改良式灰階形狀型内插 内插技術可分為一類：場景型（scene_base(j)以及物件 .型（object-based)。二元形狀型影像内插法係首先由Raya 及 Udupa 所提出（S. Ρ· Raya and j. K Udupa，“Shape based interpolation of multidimensional objects,” IEEE Trans ⑩ Med· Imag.，v〇l. 9, no. 1，pp. 32-42, Mar. 1990)。其屬於物 件型範疇。與場景型方法直接從已知場景内插場景密度值 不同，形狀型内插包含若干步驟。第一步驟係對預處理之 二元場景實施距離轉換。於上述轉換之後，每一像素係分 配一代表至物件與背景間之界限之最短距離之數值。如此 領域所熟知，正之距離數值係分配至界限内之點，而負或 零之距離數值係分配至背景。接著，於第二步驟中，此些 φ距離圖係透過一次（first_〇rder)或更高次之内插技術予以 内插。最後，經内插之距離圖係藉由將距離圖定限 (thresholding)於零而轉換回二元場景。如此，非零之點係 代表物件内部之點。六年後，Grevera及Udupa將此方法 •延伸至多維度灰階影像，請參照G. j, Grevera，j,K. Udupa， Shape-based interpolation of multidimensional grey-level images, IEEE Transactions on Medical Imaging, 15(6):881-892, 1996。 於本發明尹，Udupa之方法係藉由蒐集從其他二垂直 1378402 方向獲得之另外之體積資訊而加以修改。換言之，首先實 施Z方向中之Udupa之灰階形狀型内插法。接著，取得從 X方向及y方向掃描而得之其他二非均勻資料集，並對其 實施三維歐基里得距離轉換（3D Euclidean distance transform)，以萃取其形狀資訊。最後，此些三非均勻資料 集係予以結合以產生高解析度結果。第四圖係顯示所有之 步驟。 本發明之灰階形狀型内插演算法係概述如下。

A\gorithml :PREPROCESSING(Vx, Vy, Vz)

Comment: Vxf Vyf Vz are non-uniform samplings of object O with dim Dx, Dyt Dz

Vy and Vz are registered to Vx Vo(xtytz) <— generate a uniform sampling volume where x=sizeof(Dx)f y=sizeof(Dy) and z=sizeof(Dz) for d — z do Vo(x,y,z) <r- Vd(xtytz) do endfor for d <r- xfy register Vd to Vo do Vo(xtyfz) Vd(x,yfz) endfor return {registered volume Vo) 演算法1 (Algorithm 1)係產生尺寸與物件之原x、y及 z尺寸相同之三維陣列，以保留三非均勻體積資料中之所 有數值。於將每一體積之強度儲存至Vo之前，Vx及Vy 係對位至Vz，且數值係配置於Vo中之對應位置。

Algombml: INTERPOLAHON(Vx, Vy, Vz)

Comment: Apply 2D gray-level shape-based interpolation to Vz and get Sz. Apply 3D Euclidean distance transform to Vx and Vy and get Sx,Sy. Then merge information from Sx,SyfSz. 11 1378402 iforzL ΦζΗ Vzi(x,y,zi)<r-Vo(x,y,zi) where Vzi(x,y,zi) means lower slice Vzfj (x,v,Zfj)<-Vo(x,y,zf{) where (x,y,z^) means upper slice for all .ice[0...Dx] >e[0...Z),v] z-zi and zfj and VzH(xJ^zH) [Τ^^Τ^^νζ^χ,γ,ζι) and VzH(x,y,zH) Vx,Vy<-Resampled to Dimension (Dx,Dy,Dz) by 'mierpolalion for /e[rmjn - rmax] ΎΖΙ^ΖΗ仆resh〇W = t ffor each x,yeVzi ,Vzff

do< [Vzi<-2DDT(Vzi) Apply 3x3 chamfer distance method \Vz^<^2DDT{Vzfj) Apply 3x3 chamfer distance method '^Vzj\f i—Average of 2DDT{Vzi) and 2DDT{Vzfj) I where zi< zj^f <zfj for each x,y,z€Vx,Vy do< • \Vx<^3DDT(Vx) Apply 3D Euclidean distance do< [Vv<-3DDT(Vy) Apply 3D Euclidean distaace transform transform

Szj^f^-foT all je[0...Dj] ^e[0...D>>] extract shape of Vzj^f where Vzj^f (χ,γ,ζ^ )>〇 Sx,^<-for all A'e[0...Z)x] ve[0...D^] ze[0...Dz] extract shape of Vx and Vy where Vx(x,y,z)>0 and Vy(x,y,z)>0 for allA'e[0...DA·] ye[0...Dy] ze[Q...Dz]

Aif (SzM true \l{x,y,zM )++

Vzj^f <^I{x,y,z^f ) Calculate gray-miensity in Vz^f Vo<^Vz^ (x,y,z^ ) Accumulate Vz^f to Vo return (merged volume Vo) 演算法2(Algorithm 2)係顯示所實施之所有步驟。其係 基於Grevera及Udupa之方法，但藉由利用其他之輸入集 而延伸至較高解析度。上述方法包含若干步驟：1)選拓（Lifting):於資料集Vz中，片段VzL& VzH係 12 1378402 任意選擇，接著計算其間之中間片段Vzm。於選拓步驟中， 於片段VzL及Vzh中之灰階場景係轉換成二元場景。首 先，找出整體場景中之最小Tmin及最大Tmax灰階值。接 著’場景係以位於Tmin及Tmax之間之每一灰階值t予以 疋限。對於具有灰階值之每一像素（x，y,ZL)或及

VzH(x，y，zH) ’ 若 VzL(x，y,ZLpt 或 VzH(x，y，ZH)^，則其分配 為1 ’否則為〇。 2)距離轉換：二維距離轉換係應用至前述步驟所產生 之二元場景。二元場景係轉換成-距離圖，其中每一傻夸 之數值係代表從該像素至界限之最小距離。二 ’、 離，3 一 dlstancemet_^ 里付距離（Euchdean distance)。方程式（3)係顯: 之斜面樣板，而方程式（4)係顯示第二回合之樣板 ° σ

Template

Template

Forward Backward 14 10 14' 10 0 0 0 0 0 '0 0 0 0 0 10 14 10 14 (3) (4) 用於Si歐基里得距離轉換:三維距離轉換演算法係使 13 1378402

Algorithm 3 : DISTANCE TRANSFORM Comment: Forward Pass forz e[〇...(Dz-\)] for y e [〇...(DyA)] do jforx ^.,.(Dx-\)\ |do{ V(x,y,z) <- infVlJ.i4l 〇 + i,y + j,z + k) + m(i.j.k))

Comment: Backward Pass for z e ]^Dz-\/..〇] for y e \(Dy-\)..〇\ do Jforjc s[(DxA)...〇] °\do{v(x,y,z) mi^ijMASiX](V(x + i,y + j,z + k) + m(i.j.k)) 1 vr 'V3 V2 vr • TemplateFoward = 1 0 0 ,TemplateFonvard = 1 V2 z=0 0 0 0 z=-\ V3 '0 0 0 " 'V3 vr TemplateBackward = 0 0 1 ,TemplateBackward = 1 z=0 β 1 V2_ z=l V3 (5) (6) 演算法3(Algorithm 3)係應用至事先藉由特定灰度值 予以定限之其他二資料集Vx及Vy。方程式（5)及方程式（6) 係為所使用之順向及反向距離樣板。其係為於每一 z方向 中之三乘三樣板。若干先前技術增加其他之斜面類型，例 如五乘五斜面。其將具有較準確之結果，但較耗時。於本 發明之實作中，鄰近之26個立體像素係納入考量，以藉由 三維歐基里得距離轉換（3D Euclidean distance transform) 及定限而萃取出Vx及Vy之表面之特徵點。 4)内插：於步驟2，計算出内含許多正及負數字之 VzL(x，y,zL)及VzH(x,y,zH)之距離圖。根據Udupa之方法， 此些點係分類成二類。接續，線性内插係應用至此二距離 14 1378402 圖’且可獲得VzM(x，y，zM)之中間距離圖。 >(广’y’Zl)or} = |P〇sitivevalue (X，y)isinteri〇r h 5 2h ) J [negative value —> (x? y) is exterior ( 7 ) + 5)轉換成二元：於上述之步驟3及4中已成功以距 離圖形式獲得VzM(x,y,zM)。距離數值係接著分類成二元數 值:以代表物件或背景。正值一般係視為物件而負值一 奴係視為背景。由於直接内插VzL(x y，ZL)及VzH(x，y，ZH)之 較低及較高之距離圖，故可能會於中間圖中獲得數值〇。 為解決此問題，利用二方法決定此些點所屬為何。於方法 1中，搜尋所關切之像素之鄰近8個像素 多數鄰近像素為正，則將其分配為正。若其Π心 素為負’則其代表背景。於方法2中，此些零之點係強迫 為正，即内部物件。 少6)相父.於實施距離轉換及内插之後，形狀資訊&M 係從VzM(X,y，ZM)獲得’而同時办係從以，办係從❼。對 鲁於任何點c^)eVzM(x，y，ZM)，其數值係稱為同一位置之及 Sy之形狀屬性。若方程式心心〜）。地叫）〇抓，叫）係 成立，則點〇,_y,zM)屬於内部物件。 7)搨印（coUaPsing):步驟1至6係依序重複直到所有 灰階強度從Tmin檢驗至Tma?^藉由計算每一強度臨限值 中之此些”真實數字”，可重建灰階中間片段VzM(x,y，ZM)。 實驗結果

Grevera及Udupa提出五個以定量方式比較内插結果 及原始片段之方法。大多數方法係比較來自正確之中間片 15 1378402 段之灰階值，即解答，與來自經内插之片段之灰階 值。二方法，包含場景型線性内插、Udupa之灰階 形狀型内插以及本發明之影像融合模組1〇31所利用之改 良式内插法，係予以實施。驗證之基準係如下所列： 1)均方差（Mean Squared DifferenceXmsi/p·)): ⑽—你，細]2 ⑻ 2)點位置不正確的數目（Number 〇f 〇f disagreement)(«5i/^): nsd{i)= (W)e/ where v(k) = ^ > θ here We Set6 = l (9) 10, otherwise 此基準係計數原始片段中之像素之灰階值與經内插之 片#又中之灰階值彼此相差大於Θ數值之像素之數量。 3)最大差異（Largest difference)(/A(%)；):

鄉)κ W)|} (10) 此基準係測量來自於原始片段及經内插之片段之灰階 值之間之最大差異。 4)總有正負號差異（Total signed differenceMbt/W)): (Π) 此基準係計算總差異。若結果為負，則經内插之片段 甲之灰度值相較於原始片段中之灰度值而言係估計過高。 若結果為正，則經估計之片段中之灰度值相較於原始片段 1378402

1 A 而&係估計過小。 5)總無正負號差異（Total unsigned difference)(如心(％)): (12) tudi{i) = Y\fc (χ, y, z) - /. (x, y, z)\ {xtytz)ei

此基準係測量原始片段與經内插之片段之間之總絕對 差異。其顯示經内插之片段之灰度值（gray_values)實際上 偏離原始片段_之灰度值之程度。於上述五個方法當中， 作為影像内插中之比較量測（C〇mparative measure)而言均 方差及點位置不正確的數目係被視為最重要。 於貫驗中係利用來自於斷層掃描之二資料集，其係為 人體下顎骨資料以及人體足部資料，且均為256χ256χ256 八位元灰階影像。原始體積資料集係各別在x、y及Z方向 予以次取樣（subsampled)並平移。接著，此三資料集係加以 對位並内插，且其結果最後係藉由上述五個基準予以評估。 1)實驗1.如第五圖所示’首先三人體下顎骨資料集（來 自Z方向之資料5(H、來自χ方向之資料5〇2及來自γ方 向之資料503)係予以對位，且選擇ζ方向中之第126片严 5〇11及第130片段5012進行實驗。四個方法即線性内又 插法、Ud零之灰階形狀型内插法以及本發明之影像融人 模組顧所_之二改良式雜型内插法係加以實施。 第七圖係顯示第1次實驗之重建結果。第一行勺八 始片段。第126片段及第130片段係用作為輸入而= 至第129片段制以評估實驗結果。藉由線性内插法 异出之重建影像係置於第二行。第三行係_示_咖pa 17 1378402 之方去所得之結果，而最後—行係顯示本發明之影像融合 模組103 1所利用之改良式内插法之結果。

第八圖係顯示實驗結果根據上述五個公式之定量測 量。每一内插法之第！27、第128以及第129中間片段係 予以比較。前面二個評估公式係視為評估程序中最具代表 性者。如第八圖所示，無論應用何種評估公式，第重 建片段相較於第127及第129片段而言具有較差之結果。 主要原因在於經内插之片段離輸入越遠，則會得到越差之 重建影像。第128片段之位置係於二輸入之中間位置。由 於與輸入之間之距離，帛128片段相較於第127或第⑶ 片段而言擁有較多不準確性。本發明之影像融合模組咖 所利用之改良式内插法因從其他二資料集合併額外之資 汛，故在實驗1中較先前之方法表現稍好。 、 2)貫驗2.如第六圖所示，首先 自元一人體足部資料集（來 自Z方向之資料6〇1、央白γ 士 貝丁寸果自X方向之資料002及來自γ方 向之資料603)係予以對位，且撰摆 且送擇z方向中之第126片段 6011及第132片段6012進行實驗，1 + 〆、 逆仃貫驗其中二輸入間之距離 係大於先前實驗之距離。第 雕弟九圖係顯不第二次實驗中之重

建影像。與第一次管，驗L ’ 同的疋本;人實驗係選擇其間JL有 較大距離之二片段。笫+同各 _ > 第十圖係顯不评估之定量結果。本發 明之影像融合模組1〇31所刹田夕里禾料 此士七 所利用之改良式内插法相較於先 則之该些方法仍具有較傑出 於本發明 1378402 相鄰調整（adjacency adjustment)、整體影像之增益補償 (gain compensation for entire images)以及多波段結合 (multi-band blending)，如第 Ί--圖所示。 錐形相關 相關計算係產生一描述二隨機變數間之匹配關係程度 之單一數字。二隨機變數X及γ，其資料對係為（Χί•，乃)，ζ·= 1，2,…，η。其平均數及變異數各別為无及、歹及外。接 著’相關r則決定為：

(13) /、、憂數（象徵六張圖）係一次列入考量，將會獲得一資 料矩陣，即相關矩陣。 首先，影像係予以減少取樣（d〇Wn_Sampled)成最小尺 度（第十二a圖），其係該錐體之頂端層級。藉由對每一像 素逐一計算與其他圖像之相關（第十二b圖）並消除不合理

之結果（一臨限值係加於變數化及Sy上）可獲得最高相關， 並可知道相對位置係位於左上角。接續，影像係予以增加 取樣（up-sampled)至下一層級，並於新位置附近之合理範圍 内加以搜尋，以精化先前所得之角落座標（第十三圖）。上 述程序係予以重複直到以最精密之層級發現重疊之位置。 另則，此裎式僅計算相關矩陣之上部三角形。於每一 相關矩陣中有：個三角形部份，—位於㈣線之下方且在 其左邊（下部三角形）以及另一位於對角線之上方且在 邊（上部三角形）。相關矩陣之二個三角形總係互為鏡像、(變 19 1378402 數x對變數y之相關總係等於變數y對變數χ之相關）。 表1相關矩陣（實驗結果之一）

影像 10] 影像 [〇] 影像⑴ 影像12】 影像[3J 影像141 影像1” \ 0.930926 0.869357 0.463536 0.456217 0.898263 影像 \ \ 0.93184 0.576419 0.429544 0.581173 影像 121 ------一 \ \ \ 0.949069 0.5361 15 0.534995 影像 [3] 一 \ \ \ \ 0.917143 0.837898 彩像 14] 影像 15] 1--- \ \ \ 0.913916 \ \ \ \ \ 於搜尋整體相關矩陣後首先發現最高相關（表1)，接著 可决定第一影像配對（image pairs)，即配對<影像[2]，影像 如表1所示。上述程序係持續直到所有編號（0〜5)圖 句出現於衫像配對清單（image pair Hst)中（表2a)。每一 〜像配對不僅代表二影像係相鄰，且代表圖像之相對位 置。於程序期間，所有料所應位在結合影像中之位置均 予以決定（表2b)。 表2a影像配對清單 20 1378402 影像[2] 影像[3] 影像[1] 影像[2] 影像[〇] 影像[1] 影像[3] 影像[4] 影像[4] 影像[5] 表2b每一圖像之x-y座標（表1之結果） 2nd 3rd 4th 5th 影像[〇] (1〇,〇) 影像[1] (4, 798) 影像[2] (0,1253) 影像η] (730,1259) 影像[4] (739,789) 影像[5] (740,26) 相鄰調整 • 於錐形相關中，六片段係以相同之序列號碼計算於六 堆疊中，且假設其均位於同一 Z座標。 為解決此問題，定義一權重C，其係影像配對清單中 所有相關之平均。C係視為可顯示此結合可能在同一平面 上之程度之參數。藉由替代相鄰之圖樣以具有新結合，可 決定何者可能會接近同一平面上六影像所欲之集合。 C = Σ

„ . . . .,. correlation of image pair all pairs in a image pair list w A (14)

number of image pairs Ist 影像m 影像[3] 2nd 影像[1] 影像[2] 3rd 影像[〇] 影像[1] 4th 影像[3] 影像[4] 5th 影像[4] 影像[5] 連續性之優點可用於開始計算何種結合係所欲之結 • 果。 經結合結果中其中之一之相對位置將會被記憶，且第 五堆疊之每一片段將會於接續之影像結合（image-blending) 程序中平移。利用堆疊當中之相似關係之優點將會省下大 21 計^㈣位之每—影像配對之相關之所需時間 (第十四a圖及第十四b圖）。 了間 整體影像之增益補償 均係：像之差異係以其平均強度加以消除。每-影像之平 m η ΣΣ^ίο) mxn 灸=1L 6 而所有輸入影像之+A( _p # 务 < 十均可糟由以下計算： 6 mk nk ΣΣΣμιλ) f/(^ Λ=1 (15) (16) k=\ 每一輸入影像之強度係調整為： !k{hj) = Ik(ij)x£- A: = 1L 6 (Π) 經過若干測試後，選擇不忽略 特別情況，乃因下列二因辛 “ 視為 於此些情況中並沒有得到：良:；；果忽;背景=重叠 景或考量重疊並沒有影像補太卜’有無忽略背 g ^ 豕補彳貝之比革太多。即使二比率不 =異太大’仍選擇方程式15纟16,乃因欲具有較佳之 多波段結合 像強係利用權重函數沿著每-波線實施影 可造;若有少量對位錯誤時，此方法 成南頻㈣之模糊化。為防止上述問題，Β叫p.and 22 1378402

Adelson, E., 1983, A multiresolution spline with application to image mosaics, ACM Transactions on Graphics, 2(4):217-236所述之多波段結合演算法（muhi_band blending algorithm)係予以利用。多波段結合之概念係將大 二間範圍之低頻與小範圍之高頻結合。 實驗令完成此方法之步驟係如下所述： 1.建立一與結合相同尺寸之大型遮罩[〇]。 2.疋義有重疊之區域為 王宜

3.於/中之像素係以遮罩[〇]中與影濛μ]之指標 Whl〜6)相同之數字加以標示。 不 4·广中之像㈣予以計算出其與㈣3設定有數字之 象素之距離’並Μ遮罩附之數字為與最靠近者相同。 _5，建Λ與步驟1 —樣大小之遮罩Π]〜遮罩[6]，若遮罩 否則之料係Ζ ’貝^夺遮罩[/](/=1〜6)中之該像素填成1。 則’將该像素設定為〇。 斯_Gausslan flltedng)以不同之變異 上也、個遮罩及六個圖像，以產生不同波段。 7·分離不同波段並乘以對應之遮罩。 8 · 一同增加所有之波段。 藉由實施此些步驟，強廑 影像之門夕拉Μ 1 強度之梯度係較為平滑，且六個 竭之接縫變成非肉眼可見。 實驗結果 汁刊用之六個輸入影像俜如坌x ^經過增益補償之六輸入影輸== 叫厂/| 經 23 1378402 過多波段結合之結果係如第十五。圖所示。 於本發明之-實施例中，標準腦模型係藉由標準化模 組1033所建立，上诚p進仆柃 、 这私卓化杈組1033包含處理系統，用 貝料；平均模型產生模組’其内嵌於計算系統内， 〃：备一維個別模型部份之輸入係饋入處理系統内時平 均換型產生模Μ係響應上述二維個別模型部份之輸入，且 能，理個別模型建構及模型平均化程序，藉此從所有之輸 入資料集產生平均模刮.择進# =处 資料庫、細胞網路次資料：：物’包含生物表現次 ^丄 貝科厍以及生物超微結構次資料庫， 中標準化資料庫係輕合至處理系統以儲存至少上述平均 杈型，以及用以產生標準腦模型及座標系統之模組。 化模組聰之詳細敘述請參照2〇〇8年8月】 出 請，標題為「Bi〇-Exp⑽ion System with M〇duie Ζ Creat.ng the Standard DrosophHa Braln Model and Its c—ate System」之美國專利申請第i2/222如號，i 揭露係於此一併列入參考。 一 /於本發明之-實施例中，自動對準模組聰係實施— 糸列之自動對準程序。於建立標準腦模型之後，部份或致 體樣本腦及透過不同波長例如波長633奈米（nm)及波長正 奈米（nm)之二通道取得之其相關單一神經元係予以利 用’如第十六a圖所示。於自動對位程序中，參 b圖二首先樣本腦輪廓之三維影像之十心位置係對準於找、 準腦模型之中心位置。接著，樣本腦輪廟之前方部份係: 準於標準腦模型之前方部份’而樣本腦輪廓之後方部份係 24 13/84U2 二準於標準腦模型之後方部份。接續，樣本腦輪廓之視覺 -之方位係對準於標準腦模型之視覺區之方位。之後，參 …、第十，、c圖，二維仿射轉換演算法（3d &出ne , matl〇n alg〇nthm)，例如辅助商用軟體如AMIRA， .1、I ^至上述&序所~之結果，以將樣本腦輪腐變形成標 >細杈型’且轉換亦實施於樣本腦内之神經元上。神經纖 隹之，腐遂因此產生。若需要校正，則上述自動對位程序 °複進行上述自動對位程序可予以自動化以批次放置 來自於所有樣本腦之神經元影像。 於本發明之一實施例中’神經元追蹤模組1035係以已 之王域最小路徑（global minimal path)為基礎利用近似最 小路徑法（approximate minimal path meth〇d)萃取由雷射掃 描共輛顯微鏡所i集而得之三維影像堆疊之神經元之中心 線。 預處理 參 由於z方向中之共軛顯微鏡影像之空間解析度係較χ 及y方向中之解析度為差，故影像係首先針對每一片段逐 一進仃預處理。接著，最小能量圖係從經處理之影像堆疊 計算而得。 對母片#又’以其強度柱狀圖（intensity histogram)為 基礎選擇一臨限值’影像片段係接著予以二元化。接著， 獲得經精化（refined)之二元影像堆疊（G)。三維26鄰近連 接元素刀析（3-£> 26-neighbor connected component analysis) 係予以實施。於大多數之情況中，最大連接元素係為所欲 25 1378402 之神經元。使κ為内含二元化神經元之體積。歐基里得距 離轉換（Euclidean distance transform)係接著應用至厂中之 每一影像片段，並建構片段々之前景中每一物件之輪廓 (skeletonsKA)。對於每一片段灸，一組候選三維端點係藉 由檢驗^之每一端點之26鄰域中9個數位平面而計算而 得。演算法之詳細敘述請參照c Pudney “Distance 〇rdered

homotopic thinning： a skeletonization algorithm for 3D digital images,55 Computer Vision and Image Understanding, = (3)，PP. 404-413, 1998。每一片段中之輪廓點之集合係扮 演設計位能函數時之重要角色。候選端點之集合係表示為 E0 〇 判定函數及最小路徑演譯程序 為使路徑位於所欲結構之中心，位能係定義為下列之 判定函數（awarding function)/。F可視為頂點為立體像素且 邊緣係由F中26鄰域所定義之網格圖。對於e厂及其 •鄰近者π有一邊緣係由連接頂點0?，…之配對所定義，/(/?,分） 係滿足下列情泥： 1) 對於某些灸’若0 E A，則/(p，g) < 〇 ,否則其等於 〇 〇 2) 使 £>isEuclid ean (/?,《）為/7及g之間之歐基里得距離。 對於某些灸，2 [/(户，《）丨<D心“心_ (尸，W£ ^。 於第二限制下，可確保於p之權重網格圖中沒有負邊緣。 藉由實施判定函數（awarding functi〇n)以從已知來源 點演譯最小路徑，s係如下所列。從已知之來源點，s，已 26 1378402 知之Dijkstra演算法係予以實施以 之單-來源最短路徑。下逑端點， 母一次重複之後，若干候選端 二重複實施。於 ☆ 〇。 攸S移除以形成五…， 1. 最長路徑係從所有最小路徑〜 2. 需注意P中可能有靠近p 、出。 點传視A多铨日 f，</矣、端點。此些候選端 占係視為夕餘且可予以移除。 夕 峋 妙予以移除，其中：夕餘鳊點係藉由 ► (Minkowski addition)以及ΤΑ揭’…閔可夫斯基加法 法遂終止。）UT為樣板結構。當心時，演算 多邊形路徑近似 定函數所得収近似中^線並不平滑 ^同時八T圖）’故每—分支係製作—近似多邊形路徑， 0 ，刀支點係予以保留。近似法係描述如下。 已知夕邊形路徑卜^…’^以及誤差界卜找尋一 ^邊幵/路徑s ’其為S之£近似。若S滿足下列基準，則》 ’···，〜>選擇性e近似夕。 1. s之頂點集合係為$之次集合。 2. 使…-心及W/+/=v" /=1，，心j，多邊形路徑<〜，，以 上之任何頂點至線段（…，Wq)間之距離係小於6。 3上之頂點數量係最小可能。 此問邊可藉由利用動態規劃（dynamic programming)技 術予以解決。 Λ (’⑺上之邊界之數量係定義成其損失（cost)。所有用於s之 27 1378402 ’ &示為C(/J)。對於 則建立最佳ε近似路 ε近似當中之最小損失係為最佳損失 ζ·=/之邊界情形，使cG j)=1。若y>z·， 徑有二種情況。 情況1 : L)係為線段（Vi，Vy)。 此情況係在頂點邮，至（⑽間之所有距離均 J、於ε時發生。㈣係ε近似<Vi•〜，.，心〉且因此 1。 情況2 : 系由二或多線段組成。

於此情況中’ t，y>可分為二個次路徑左1中 以為< Vi•，…，〜>上之頂點。需注意^)及^力兩者係ε近似多 邊形路徑。最佳ε近似之損失外，乃 係為 mini</t<y (c(») +C(D))。 基於上述敘述，最佳損失可寫成遞迴。 Φ·Λ=ι c(h j) = 1 (ν:·, vj) ε -approximates <v/,...jVy> (18) ' c(i,j) = mini<y<, (C(i,k) + c(k, j)) 且最佳解答係以由下而上（b〇ttom-up)之方式得到。 實驗結果 資料集合 所有資料係於新竹國立清華大學腦研究中心中以 Zeiss LSM 510共軛顯微鏡獲得。果蠅之投影神經元影像堆 疊之平均立體像素尺寸係0.17x0.17x1.0立方微米。 提出方法之評估 28 1378402 於實驗中’設定；1=1.0、W=1.0，且樣板τ為具有5 χ5χ3之尺寸之盒體’判定函數（awarding function)係僅如 下： f (P,q) = -0.5 if Sk for some k (19) 雖此判定函數看似單純，但其結果良好。於第十七a 圖至第十七d圖中係呈現原始資料之最大強度投影 (maximum intensity projection; MIP)以及當實施判定函數 時之經追縱結果。第十七&圖至第十七d圖係顯示果蠅之 嗅覺神經（calyx)中之投影神經元。第十七a圖係顯示原始 景〉像堆豐之最大強度投影（maximuin intensity projecti〇n; MIP)。第十七b圖係顯示無使用判定函數之神經元追蹤結 果。於第十七b圖中，追蹤結果係以黑色線繪製，且與繪 製成灰色之原始資料之最大強度投影（MIP)相重疊。第十七 c圖中之黑色線係顯示有實施判定函數之追蹤結果。第十 七d圖係顯示第十七c圖中之追蹤結果之ε近似，其中 ί = 。第十七C圖及第十七d圖兩者均係藉由直接體積重 構（directed volume rendering)而繪製。立體像素尺寸係〇 18 χΟ· 1 8x 1.0立方微米（e m3)。於第十八圖中係呈現一較複雜 之結構，即果繩腦側角區中之投影神經元之ε近似中心 線。第十八圖係顯示果蠅之腦側角區中之投影神經元之ε 近似，^ = W。黑色線係為追蹤結果。第十八圖係藉由直接 體積重構（directed volume rendering)而繪製。立體像素尺 寸係 0.16x0.16x1.0 立方微米。 29 以，此’本發明係提供上述之影像預處理系統，其可用 資料庫it昆蟲腦影像資料庫系統’例如三維果繩腦影像 ’以探索昆蟲腦部例如果罐腦部之内部結構。 上述敘述係為本發明之較佳實施例。此領域 領會其係用以說明本發明而非用以限 ； =1=:，利保護範圍當視後附之申請專利 離本糞靜5…而疋。凡熟悉此領域之技藝者，在不脫 本專利精神或範_，所作之更 神下所完成之等效改變或設計，且 述之申凊專利範圍内。 【圖式簡單說明】 *後：==明書中之若干較佳實施例及詳細敘述 明中之Γ:Γ 。圖式中相同之元件符號係指本發 實而’應理解者為，本發明之所有較佳 貫施例係僅用以說明而非用以限制申請專利範圍，並中. 第一圖係根據本發明之影像預處理系統之方塊圖。 第-圖係根據本發明說明對位型形狀型内插之系 要圖。 第三圖係根據本發明說明影像對位模組之基本元件。 t圖係根據本發明說明2方向中之灰階形狀型内插 與伙X及Y方向萃取而得之二額外形狀f訊之結合。 第五圖係根據本發明說明關於影像融合模 之結果。 第六圖係根據本發明說明關於影像融合模組之實驗2 1378402 之結果。 第七圖係根據本發明說明實驗1之重建結果。 第八圖係根據本發明說明實驗結果根據五個檢驗公式 之定量測量。 第九圖係根據本發明說明實驗2之重建影像。 第十圖係根據本發明說明評估之定量結果。 第十一圖係根據本發明說明影像接合模組之整體程 序。 第十二a圖係根據本發明說明依序排列之經減少取樣 之影像。 第十二b圖係根據本發明說明逐一對每一像素之相關 s十算’其中虛線係標示b之搜尋範圍。 第十一 c圖係根據本發明說明相關中之未匹配結果 (左邊）以及理想匹配結果（右邊）。 第十三圖係根據本發明說明下一層級（虛線）之搜尋範 圍。 第十四a圖係根據本發明說明堆疊間之理想關係。 第十四b圖係根據本發明說明實驗中之堆疊間之關 係。 第十五a圖係根據本發明說明六個輸入影像。 第十五b圖係根據本發明說明經過增益補償後之六輸 入影像。 第十五c圖係根據本發明說明經過多波段結合後之結 果0 31 1J/MU2 第十六a圖係根據本發明說明自動對準 動對準程序令所利用之三維腦部影像。、，·且e 自動圖係根據本發明說明自動對準模組所實施之 之 自動::::圖係根據本發明說明自動對準模組所實施 第十七a圖係根據本發明說明、^ ^ ^ 度投影（MIP)。 以彳豕堆f之取大強 函數之神經 第十七b圖係根據本發明說明無使用 元追蹤結果。 結果 第十七e圖#、根據本發明說明有實施判定函數之追縱 第:七d圖係根據本發明說明第十七c圖中之 果之ε近似，其中f = VJ。 视、、，。 第十八圖係根據本發明說明果蠅之腦側角區 神經兀之ε近似，£ =万。 又如 【主要元件符號說明】 W影像預處理系統 1 〇 1處理單元 102資料庠 103影像預處理單元 201輸入體積 202經對位體積 203重建體積 32 1378402 »· i 501來自Z方向之資料 502來自X方向之資料 503來自Y方向之資料 601來自Z方向之資料 602來自X方向之資料 603來自Y方向之資料 1031影像融合模組 1032影像接合模組 • 1033標準化模組 1034自動對準模組 1035神經元追蹤模組 5011第126片段 5012第130片段 6011第126片段 6012第132片段

33

Claims (1)

1. 丄 3/MUZ* 七、申請專利範圍： 2012/9/19無劃線替換頁 叫年日修正J 1. 一種影像預處理系統，包含 一處理單元； 一衫像預處理單元，苴▲ u # ^ 一耦合至該處理單元以預處理影像 二 =’其中該影像預處理單元包含一影 : ==影像片段資料間之遺漏數值，-影像二 料，以及一 ^ 」二片段資料接合成經接合之影像資 中之神經元之，心線；以及㈣敲接合之景，像資料 為料庫，其耗合至該處理罝一、 像片段資料。 70以儲存該經預處理之影 2. 如請求項1所述之影像預處m 組係利用對位型形狀型内插。 其中§玄影像融合模 月求項1所述之影像預處理系統 組包含一影像對位模組。 其中該影像融合模 4.如凊求項3所述之影像預處理系統 組包含一度量元件。 其中該影像對位模 如3月求項3所述之影像預處理系統 組包含一内插元件。 其中該影像對位模 34 5. 1378402 6.如請求 組包含 影像預處理系統’其中模 轉換元件 所述，影像預處理系統’其中該影像對位 優化元件 模 8:=::=:系統’其―合模 9. 項1所述之影像預處理系統，其中該影像接合模 -系列演算法’其中該演算法包含錐形相關、相 鄰調正、整體影像之增益補償以及多波段結合。 10.=請求項1所述之影像預處理系統，其中該影像預處理 早7G更包含一標準化模組，以為該經接合之影像資料產 生一標準座標。 11.如請求項10所述之影像預處理系統，其中該標準化模 組包含： ' 一處理系統，用以處理資料； 一平均模型產生模組’其内嵌於一計算系統内，其中當 二維個別模型部份之一輸入係饋入該處理系統内時該 平均模型產生模組係響應該二維個別模型部份之該輸 入’且能處理個別模型建構及模型平均化程序，藉此從 35 U78402. 2012/9/19無劃線替換頁 所有之輸入資料集產生一平均模型； /示準化貢料庫’其包含生物表現次資料庫、細胞網路 2料庫以及生物超微結構次資料庫，其中該標準化資 一庫係輕合至該處理系統以儲存至少該平均模型；以及 -用以產生標準腦模型及座標系統之模組。 l2m所述之影像預處理系统，其中該影像預處 料對準於-標準影像。 將缝接合之影像資 13·如請求項12所述之料預處理 模組係結合該經對準之影像資料盘_ =自動對準 生之該標準座標。 胃i礼準化模組所產 14.如請求項1所述之影像預處理“， 模组俜刹田也丨宁τ缸 ’、v 5亥神經元追蹤 #j用-判疋函數以及-最小路徑演譯程序。 ΐ5·如請求項丨所述之影像預處理_， 模組係利用-多邊形路徑近似計算法、。°Λ神經70追縱 36