楊麗娟, 田錚,2, 溫金環, 延偉東
(1.西北工業大學 應用數學系, 陜西 西安 710129; 2.中國科學院 遙感科學國家重點實驗室, 北京 100101)
點集匹配在圖像分析、形狀匹配等許多領域中已獲得廣泛的應用[1-2],可粗略的分為剛性和非剛性,其中非剛性匹配問題備受關注。但由于點集中通常包含大量異常值,這嚴重影響了匹配性能,因此需要研究更穩健的點集匹配算法。
近些年,基于混合高斯模型(Gaussian mixture model,GMM)的匹配方法更是受到了研究者的廣泛關注[3]。在最大似然估計框架下,利用最大化期望(expectation maximization,EM)算法求解參數的最大似然估計。Jian等[4]用混合高斯模型表示點集,通過最小化2個混合高斯之間的統計差異,進而確定點集之間的對應關系。Myronenko等[5]假設點集做一致性運動,提出了基于一致性點漂移(coherent point drift,CPD)的點集匹配方法,通過增加額外的均勻分布模型化異常值。但它本質上對應于最小二乘方法,對異常值不穩健。為了提高對異常值的穩健性,一種可選的方法是選用更合適的概率模型,如具有厚尾性質的t分布[3]。Zhou等[6]提出用混合t分布模型(student′s t mixture model,SMM)處理點集之間的匹配問題,通過最大化似然函數的對數期望,估計變換參數。但是它和CPD方法以及多數方法一樣,只允許各向同性協方差,而且需要手動確定正則化參數。
針對最大似然估計框架中存在的奇異性問題,可選的解決方式是貝葉斯分析[7-8]。Qu等[9]提出了在變分貝葉斯框架下的點集匹配方法(VBPSM),將點集匹配過程分為聚類分析和回歸分析兩部分,并分別利用2個GMM模型化,通過引入轉移變量將兩部分統一起來。但對于非剛性變換,它沒有包含對變換函數的平滑性約束。另一個問題是估計的變換函數沒有直接反映待對準點集之間的變換關系。在Simpson等[10]提出的非剛性配準的概率配準框架中,通過給變換參數指派一個合適的先驗分布可將正則化并入概率框架中,根據數據質量自適應地確定正則化參數。
貝葉斯SMM已經廣泛應用到了密度估計、聚類和模型選擇中[11-12]。Archambeau等[12]在估計變量的后驗分布過程中,考慮了變量之間的相關。Baldacchino等[13]將貝葉斯SMM應用于非線性系統識別中的穩健回歸。
綜上,針對存在異常值的非剛性點集匹配問題,本文將點集匹配問題表示為SMM對應的概率密度估計問題。通過引入真實后驗分布的變分近似,匹配問題轉化為最大化變分下界,同時允許任何形式的協方差。在參數推斷過程中,利用變分貝葉斯最大化期望(variational Bayesian EM,VBEM)算法確定變換參數。通過先驗模型,將空間正則化約束并入貝葉斯SMM模型中,可自適應地確定正則化參數。模擬點集和真實圖像上的實驗結果證明所提方法可大大提高匹配性能。
(1)
式中,Θ={A,B,{Λm},{πm},{vm}},t分布的概率密度函數為:
(2)
式中,μ和Λ分別為均值向量和精度矩陣,Γ(·)為伽瑪函數,當自由度參數v→∞時,t分布收斂到高斯分布。t分布可表示為無限個具有相同均值、不同精度的尺度高斯分布的混合,即:
(3)

于是,似然函數可寫為:
利用指數族分布的共軛先驗特性,選擇模型變量的先驗:①混合比例系數π的先驗分布服從狄利克雷分布Dir(π|κ0),其中κ0是先驗樣本大小(稱為超參數,下同);②每個混合成分的精度矩陣Λm的先驗分布服從Wishart分布W(Λm|γ0,S0),其中γ0和S0分別為自由度和尺度矩陣;③仿射變換矩陣A和非剛性形變系數矩陣B的每列元素的先驗分布分別服從零均值向量、精度參數υ和η的高斯分布;④精度參數向量υ和η均服從伽瑪分布,其中a0,b0和c0,d0分別是先驗形狀和尺度參數對;⑤正則項以先驗形式并入貝葉斯SMM模型中[10],即:
(8)
式中,λ>0是標量正則化參數;⑥正則化參數λ服從形狀參數為τ0和尺度參數為ζ0的伽瑪分布。可根據經驗確定超參數的值。通過變分推斷,可獲得變量{ΘL,ΘVB}={ΘL,{A,B,Λ,π,υ,η,λ}}的后驗分布,由于沒有自由度參數v={vm}的先驗,利用最大似然方法可確定它的值,記為ΘML={v},詳見下一節。
變量之間的聯合分布可表示為:
(9)
圖1中給出了非剛性點集匹配對應的概率圖模型,即有向無環圖[7]。

圖1 非剛性點集匹配的概率圖模型
考慮近似后驗分布q(U,Z,ΘVB),對數似然函數可表示為:
lnp(X|Y)=L(q)+
(10)
所有變量的后驗分布可分解為:
q(U,Z,π,Λ,A,υ,B,η)=q(U,Z)q(ΘVB)
(11)
記q(ΘVB)=q(π)q(Λ)q(A)q(υ)q(B)q(η)
q(λ)。關于q(U,Z)和q(ΘVB)分別最大化L(q)可獲得如下VBEM算法,在第k次迭代:
VBE-step:
q(U,Z)k∝exp(〈lnp(X,U,Z|Y,ΘVB)〉q(ΘVB)k-1)
(12)
VBM-step:
(13)
式中,〈·〉q(·)表示關于后驗分布q(·)的期望,后面下標省略。在變分貝葉斯期望步驟(VBE-step)中,固定變分變量ΘVB,計算隱變量ΘL的后驗分布;在變分貝葉斯最大化步驟(VBM-step)中,固定隱變量ΘL,重新計算ΘVB的后驗分布。下面列出部分變量的后驗分布:
1) 考慮隱變量之間的相關性,其后驗分布為:
(14)

q(znm=1)=
(15)
規范化后,示性變量元素表示為:
?n,?m
(16)
以示性變量為條件,尺度變量的后驗分布服從伽瑪分布:q(unm|znm=1)~G(unm|αnm,βnm),其中
(17)
(18)

(19)
(20)
式中,diag(〈υ〉)表示〈υ〉中元素構成的方陣,矩陣Ξ的元素為Ξnm:=〈znmunm〉,Λm(q)表示Λm的第q個對角元素。

(21)
(22)
4) 正則化參數λ的后驗分布服從伽瑪分布:q(λ)~G(λ|τ,ζ),其中
(23)
(24)
根據伽瑪分布性質,有〈λ〉=τ/ζ,可根據不同的應用自適應地確定λ的值。
其余變量的后驗分布見附錄。由于沒有自由度參數{vm}的先驗,可直接最大化下界L(q),獲得下面非線性方程:
(25)
(26)

算法1: 基于SMM的非剛性點集匹配算法
2: 初始化:
先驗:a0,b0,c0,d0,κ0,γ0,S0,τ0,ζ0
變換:A,B
3: forS0上由粗到精do
4: repeat
VBE-step:
根據(15~16)推斷示性變量Z
根據(17~18)推斷尺度變量U
VBM-step:
根據(27)更新混合比例系數π的κ
根據(19~20)推斷仿射變換矩陣A
根據(30~31)推斷精度υ的a和b
根據(21~22)推斷非剛性形變系數矩陣B
根據(32~33)推斷精度η的c和d
根據(28~29)更新精度Λ的γ和S
根據(23~24)更新正則化參數λ的τ和ζ
根據(26)更新自由度v
計算第k步的L′(q)
6: end for
算法的時間復雜度主要由Z,Λ和B的后驗更新決定。在每次迭代中,更新這些變量的時間復雜度分別近似為O(DNM2),O(D2(N+D)M)和O(DM3),一般情況下,D?N,M,故算法的總時間復雜度可近似為O(KM3)或O(KNM2),其中K為迭代次數。
為了驗證本文算法的有效性,首先測試具有形狀差異的模擬數據集,然后將本文方法應用于包含非剛性形變的圖像。并與5種算法進行了對比:GMM-L2[4],CPD[5], SMM[6]和VBPSM (默認各向異性協方差)[8]和TPS-RPM[15]。超參數的值設置如下:①超參數a0,b0,c0,d0,τ0和ζ0設定為0.01;②κ0的元素設定為M;③γ0設定為2×(D+1),若選擇各向異性協方差,則S0設定為單位矩陣I,否則設定為標量1;④高斯徑向基函數的寬度參數設定為1。針對各向異性和各向同性協方差2種選擇,本文方法分別記為VBSMM(ANI)和VBSMM(ISO)。
為了檢驗自適應確定正則化參數的有效性,選用文獻[15]中公共數據集的3組存在異常值的漢字福形狀點集,記錄算法1在固定正則化參數時獲得的匹配召回率(recall),其中正則化參數的固定取值為0.5k,k=1,2,…,40,并與算法1在自適應確定正則化參數時獲得的匹配召回率做對比。圖2給出了召回率關于正則化參數的變化曲線,其中空心對應于算法1在固定正則化參數時的匹配召回率,實心對應于算法1在自適應確定正則化參數時的匹配召回率,其中自適應確定的正則化參數分別為13.989 8,9.870 0和7.102 9。從圖中曲線可以看出,算法1通過自適應方法總是能“自動”找到相對合適的正則化參數。這在批量處理點集匹配問題時,可避免手動確定不合適的正則化參數。

圖2 召回率關于正則化參數的變化曲線(其中實心位置對應于算法1自適應確定的正則化參數)
為評價對異常值的穩健性,選用文獻[15]中的公共數據集,包含魚和漢字福兩種形狀點集,分別包含96和108個模型點。數據點樣本包含5個異常值比例水平(從0.0到2.0),對于每一個比例,均有100個數據點集樣本,采用平均的召回率評價匹配結果。圖3給出了6種算法應用于2種形狀點集的平均召回率統計結果。


圖3 6種算法應用于魚和漢字福2種形狀點集的召回率對比結果
從圖中可看到,隨著異常值比例的增加,所有的算法的匹配效果均迅速衰減,其中GMM-L2方法最嚴重,這主要是由于異常值破壞了點集的形狀結構。但相比其他5種算法,本文提出的VBSMM算法總是能給出較好的匹配結果。圖4給出了算法1在2種形狀點集上的匹配結果示例。

圖4 本文算法在不同異常值比例下的匹配結果示例
為驗證本文算法對模型點集和數據點集中均存在異常值時的穩健性,選用CMU房子序列(http:∥vasc.ri.cmu.edu/idb/html/motion/index.html),其包含不同視角的111幅房子圖像,每個圖像中標記出了30個已知對應關系的特征點。通過給每個點集增加不同數量服從均勻分布的異常值以生成待匹配點集,其中異常值比例分別為0.5和1.0。基于采樣區間10∶10∶100對圖像序列進行采樣,可獲得10組圖像對。因為30個特征點的對應關系已知,故可直接統計每個采樣區間內的平均召回率。圖5比較了6種算法的平均召回率統計結果對比。從圖中可以看出,不同方法的性能出現了不同程度的震蕩,這主要是由于異常值生成的隨機性。整體來看,相比于其他5種算法,本文提出的VBSMM算法更穩健,這與3.2部分的實驗結果一致。圖6給出了本文算法在不同異常值比例下的匹配結果示例(正確匹配對數均為30/30)。

圖5 6種算法應用于房子序列圖像的召回率對比結果

圖6 本文算法在不同異常值比例下得匹配結果示例
本文針對存在異常值的非剛性點集匹配問題提出了貝葉斯SMM模型方法。在變分貝葉斯框架中,該方法將點集匹配問題模型化為SMM對應的概率密度估計問題,利用VBEM算法迭代更新變量的后驗分布,而且允許各向異性和各向同性協方差2種選擇。通過先驗模型引入正則化項,可自適應地確定正則化參數。最后,模擬點集和真實圖像上的實驗結果驗證了本文所提方法的有效性。
其他變量的后驗分布如下:
1) 對于混合比例系數π,其后驗分布仍服從狄利克雷分布:q(π)~Dir(π|κ),其中:
(27)
2) 對于每個混合成分的精度矩陣Λm,其后驗分布仍服從Wishart分布:q(Λm)~W(Λm|γm,Sm),其中:
(28)
(29)
當選擇各向同性協方差時,混合成分的精度服從伽瑪分布。
3) 對于精度參數υ,其后驗分布仍服從伽瑪分布:q(υl)~G(υl|al,bl),其中:
(30)
(31)
4) 對于精度參數η,其后驗分布仍服從伽瑪分布,q(ηl)~G(ηl|cl,dl),其中:
(32)
(33)