流形學習算法中鄰域大小參數的遞增式選取

邵 超，萬春紅，趙靜玉

(河南財經政法大學計算機與信息工程學院，鄭州 450002)

1 概述

近年來，人們陸續發現實際中的很多高維數據都分布在某個低維非線性流形上［1-2］，為此提出了多種能有效學習這種結構的流形學習算法如等距映射［3-4］(Isometric Mapping，ISOMAP)和局部線性嵌入［5-6］(Locally Linear Embedding，LLE)。現有的大多數流形學習算法能否成功應用還依賴于其鄰域大小參數的選取是否合適，然而，目前該參數在實際中通常還難以高效選取，另外，數據中的噪音對鄰域大小參數的合適性也會產生一定的影響［7］，從而限制了它們的實際應用［7-10］。

針對流形學習算法敏感于鄰域大小參數的這一問題，目前大多數方法都基于殘差來選取合適的鄰域大小參數［3，11-12］，但這需要多次運行整個流形學習算法以計算相應的殘差(盡管文獻［11 -12］都只對其中具有極小重建誤差或極小成本函數的少數幾個鄰域大小參數分別計算相應的殘差);另外，殘差只是維數約簡過程中產生的誤差，難以正確指導鄰域大小參數的合適選取［13］。

通過對采用不同鄰域大小參數得到的多個運行結果進行集成可以在一定程度上減輕流形學習算法對鄰域大小參數的敏感程度［9］，但算法的時間復雜度會顯著增加。通過識別和刪除鄰域圖中可能出現的“短路”邊也能在一定程度上避免鄰域大小參數難以高效選取的問題，例如，文獻［14］刪除每一個數據點中具有最小重建權重的幾個近鄰點，文獻［15］刪除每一個數據點中具有最大圖距離的幾個近鄰點，但它們又引入了一個同樣難以選取的額外參數。文獻［16 -17］通過組合若干個最小生成樹來創建鄰域圖，沒有了鄰域大小參數，但最小生成樹的數量同樣難以有效確定。文獻［18 -19］可以自適應地為每一個數據點選取不同的鄰域大小參數，但計算比較復雜。

根據流形的局部歐氏性，鄰域圖上的所有鄰域都呈線性或近似線性，是鄰域大小參數被認為合適的基礎，此時所有鄰域的線性度量可聚成一類;而鄰域大小參數一旦變得不合適，鄰域圖上的某些鄰域就不再呈線性或近似線性［20］，其線性度量也不再能聚成一類。本文用加權主成分分析［21］(Principal Component Analysis，PCA)得到的重建誤差對鄰域圖上所有鄰域的線性程度進行度量，然后用貝葉斯信息準則［2，22］(Bayesian Information Criterion，BIC)對其聚類個數進行探測，從而能遞增式地選取合適的鄰域大小參數。

2 鄰域大小的遞增式選取

具有合適鄰域大小參數的鄰域圖應能正確表達數據的鄰域結構，按照流形的局部歐氏性，所有鄰域都應呈線性或近似線性;而鄰域大小參數一旦變得不合適，鄰域圖中開始出現“短路”邊，從而使某些鄰域不再呈線性或近似線性。線性度量有很多，考慮到魯棒性，本文采用的是加權PCA［21］重建誤差。當鄰域大小k 合適時，鄰域圖上所有鄰域的加權PCA重建誤差都相對較小，可聚成一類，如圖1 所示;而鄰域大小k 一旦變得不合適，某些不再線性的鄰域其加權PCA 重建誤差就會變得很大，這樣就使鄰域圖上所有鄰域的加權PCA 重建誤差不再能聚成一類，如圖2 所示。因此，根據鄰域圖上所有鄰域的加權PCA 重建誤差所聚成的類別個數即可遞增式地選取合適的鄰域大小k。本文采用BIG［2，22］來探測其聚類個數。

圖1 當鄰域大小k 合適時的加權PCA 重建誤差

圖2 當鄰域大小k 不合適時的加權PCA 重建誤差

2.1 加權PCA 重建誤差的計算

PCA 重建誤差可用來度量數據的線性程度，一般而言，線性程度越大，PCA 重建誤差就越小。然而，PCA 對數據中的異常點(outlier)比較敏感［21］，異常點的加入會極大地改變主成分，從而使PCA 重建誤差難以準確度量數據的線性程度。如圖3 所示，異常點(圖中標注為“* ”)的加入使原本以實線標注的主成分變成了以虛線標注的主成分，從而降低了該數據集合的最大PCA 重建誤差。

圖3 異常點對PCA 的影響

鄰域大小參數一旦變得不合適，鄰域圖上就會有某些鄰域由線性變成非線性，和之前線性鄰域中的那些數據點相比，新加入的導致鄰域非線性的數據點就類似于異常點，其PCA 重建誤差(通常也是該鄰域中的最大PCA 重建誤差)用來度量該鄰域線性程度隨鄰域大小的變化趨勢是比較合適的。然而，為準確計算這些“異常點”的PCA 重建誤差，需要盡量降低異常點對PCA 的影響，因此，本文采用加權PCA 算法［21］，通過迭代給異常點賦以較小的權重，從而得到盡可能接近于原始主成分(圖3 中標注為實線)的主成分(圖3 中標注為點劃線)。

加權主成分分析算法可簡要描述如下［21］:

(1)采用標準PCA 算法獲得數據集合X 的初始數據中心μ(0)和初始主成分p(0);

(2)t=0;

(3)DO

1)t=t+1;

2)按照式(1)計算每一個數據點Xi在數據中心μ(t-1)和主成分p(t-1)下的PCA 重建誤差

3)按照式(2)計算每一個數據點Xi的權重wi:

4)重新計算數據中心μ(t)=

Until μ(t)和p(t)都相對穩定下來。

如上所述，本文用每一個鄰域Nk(Xi)(即數據點Xi及其k 個最近鄰數據點的集合，也就是加權PCA 算法中的數據集合X，因此有s=k +1)的最大加權PCA 重建誤差(一個標量值，維數為1)來度量其線性程度，進而探測該鄰域中是否出現了“短路”邊。

2.2 貝葉斯信息準則的計算

由于具有不同聚類個數的聚類模型可以看作是不同的模型，因此BIC 作為模型選擇準則就可以用來探測數據的聚類個數［2，22］。BIC 使用一個帶修正項的對數似然統計量對不同的模型進行打分，BIC分值越大的模型就越適合于對應的數據集合。

模型Mj的BIC 分值BIC(Mj)如式(3)所示［2，22］:

本文采用K-均值聚類算法對鄰域圖上所有鄰域的線性度量(即最大加權PCA 重建誤差max(norm進行聚類，因此，數據集合D 就由鄰域圖上所有鄰域的最大加權PCA 重建誤差組成，即D=，則有m=n，pj=K(R +1)，其中，K 為聚類個數;R=1 為D 的維數。

不同的聚類個數K 就對應了不同的K-均值聚類模型，因此，分別用BIC(K=1)和BIC(K=2)表示聚類個數為1 和2 時對應K-均值聚類模型的BIC 分值。例如，圖4 所示的2 個數據集合，其聚類個數分別為1和2，與之相適應的K-均值聚類模型的BIC 分值都相對較大(圖4(a)中BIC(K=1)=-354.373 6，BIC(K=2)=-458.869 9，圖4(b)中BIC(K=1)=-558.274 3，BIC(K=2)=-507.237 9)。

圖4 BIC 用于探測聚類個數的數據集合

2.3 算法流程

該方法可簡要描述如下:

(1)使用廣度優先搜索算法獲得能使鄰域圖連通的最小鄰域大小參數kmin，作為遞增式選取鄰域大小參數的起點;

(2)k=kmin-1;

(3)DO

1)k=k+1;

2)對鄰域圖上的每一個鄰域Nk(Xi)執行加權PCA 算法，獲得鄰域圖上每一個鄰域的最大加權PCA 重建誤差，即式(4)中的數據集合D;

3)計算BIC(K=1)(由于只有一類，所以數據中心即為聚類中心);

4)計算BIC(K=2)，這需要對數據集合D 執行K-均值聚類算法將其聚成2 類;

Until BIC(K=1)＜BIC(K=2)

(4)kselected=k-1。

通過該方法的遞增式選取，最終的kselected即為合適的鄰域大小參數。和傳統的基于殘差的參數選取方法相比，該方法只需遞增式地執行加權PCA 算法和K-均值聚類算法，并計算相應的BIC(K=1)和BIC(K=2)，而無需運行整個流形學習算法和計算相應的殘差。

2.4 時間復雜度分析

該方法主要涉及以下3 個部分的計算:

(1)為確定kmin，需要從k=1 開始反復執行廣度優先搜索算法(共執行kmin次)。廣度優先搜索算法的時間復雜度為O(n +e)=O(kn)(e 為鄰域圖的邊數)，因此，這一部分的時間復雜度為通常情況下，kmin都比較小。

(2)對于鄰域大小k 的每一次迭代(共迭代(kselected-kmin+2)次)，在鄰域圖的每一個鄰域(包括k+1 個d 維數據點，d 即為數據的維數)上分別執行加權主成分分析算法(共執行n 次)。加權主成分分析算法的時間復雜度為O(((k +1)d2+d3)c)(c為加權主成分分析算法的迭代次數，通常都比較小)，因此，這一部分的時間復雜度為O(((k +1)d2+d3)cn(kselected-kmin+2))。

(3)對于鄰域大小k 的每一次迭代(共迭代(kselected-kmin+2)次)，對數據集合D(包括n 個1 維數據點)執行K-均值聚類算法將其聚成兩類。K-均值聚類算法的時間復雜度為O(nKpl)=O(2nl)(n和p=1 分別為D 中數據點的個數及其維數，K=2為聚類個數，l 為K-均值聚類算法的迭代次數，通常都比較小)，因此，這一部分的時間復雜度為O(2nl(kselected-kmin+2))。

3 實驗結果與分析

為驗證該方法的有效性，本文采用具有2 000 個隨機樣本點的Swiss roll［3］和S-curve［5］數據集合(為了降低計算量，本文采用K-均值算法分別從中選取了n=500 個代表點，如圖5 所示)，該方法的運行結果如圖6 所示。

從圖6 可以看出，對于Swiss roll 數據集合，該方法選取的鄰域大小參數為k=8;對于S-curve 數據集合，該方法選取的鄰域大小參數為k=16。它們的合適性可通過ISOMAP 算法的運行結果得以證實。

眾所周知，數據中的噪音對鄰域大小參數的合適性也會產生一定的影響［7］。為了驗證該方法的魯棒性，本文在以上2 個數據集合上分別加入一定的噪音［7］(如圖7 所示)，該方法的運行結果如圖8 所示。

圖5 實驗中的2 個數據集合

圖6 在不同數據集合上的運行結果

圖7 加入了噪音的2 個數據集合

圖8 本文方法在不同數據集合上的運行結果

從圖8 可以看出，對于加入了噪音的Swiss roll數據集合，該方法選取的鄰域大小參數為k=6;對于加入了噪音的S-curve 數據集合，該方法選取的鄰域大小參數為k=11。它們的合適性同樣可通過ISOMAP 算法的運行結果得以證實。

此外，該方法基于的是流形的局部歐氏性，那么在由于采樣不均勻而出現“空洞”的數據集合上是否依然有效?為了驗證這一點，本文采用帶有“空洞”的Swiss roll 和S-curve 數據集合(如圖9 所示)，該方法的運行結果如圖10 所示。

圖9 帶有“空洞”的2 個數據集合

從圖10 可以看出，對于帶有“空洞”的Swiss roll數據集合，該方法選取的鄰域大小參數為k=8;對于帶有“空洞”的S-curve 數據集合，該方法選取的鄰域大小參數為k=11。它們的合適性同樣可通過ISOMAP 算法的運行結果得以證實。

與傳統的基于殘差的參數選取方法(括號內指定的是鄰域大小的選取范圍)相比，本文方法不但能選取合適的鄰域大小(分別如圖6、圖8 和圖10 所示)，而且還具有較高的運行效率，如表1 所示(計算機配置為Intel i3-2120 CPU，主頻為3.30 GHz，內存容量為4 GB)。

圖10 本文方法在不同數據集合上的運行結果

表1 2 種方法在運行時間上的對比 s

4 結束語

按照流形的局部歐氏性，本文方法通過對鄰域圖上每一個鄰域的線性程度進行度量，并對其聚類個數進行探測，能夠高效地選取合適的鄰域大小參數，且無需任何額外參數(在數據維數不是特別高的情況下)。下一步的研究方向是當數據的維數非常高時，如何高效地選取合適的鄰域大小。

［1］Seung H S，Lee D D.The Manifold Ways of Perception［J］.Science，2000，290(5500):2268-2269.

［2］楊 劍，李伏欣，王 玨.一種改進的局部切空間排列算法［J］.軟件學報，2005，16(9):1584-1590.

［3］Tenenbaum J B，de Silva V，Langford J C.A Global Geometric Framework for Nonlinear Dimensionality Reduction［J］.Science，2000，290(5500):2319-2323.

［4］王耀南，張 瑩，李春生.基于核矩陣的Isomap 增量學習算法研究［J］.計算機研究與發展，2009，46(9):1515-1522.

［5］Roweis S T，Saul L K.Nonlinear Dimensionality Reduction by Locally Linear Embedding［J］.Science，2000，290(5500):2323-2326.

［6］Zhang S.Enhanced Supervised Locally Linear Embedding［J］.Pattern Recognition Letters，2009，30(13):1208-1218.

［7］Balasubramanian M，Shwartz E L，Tenenbaum J B，et al.The ISOMAP Algorithm and Topological Stability［J］.Science，2002，295(5552):7-17.

［8］Saul L K，Roweis S T.Think Globally，Fit Locally:Unsupervised Learning of Low Dimensional Manifolds［J］.Journal of Machine Learning Research，2003，4(1):119-155.

［9］詹德川，周志華.基于集成的流形學習可視化［J］.計算機研究與發展，2005，42(9):1533-1537.

［10］邵 超，黃厚寬，趙連偉.一種更具拓撲穩定性的ISOMAP 算法［J］.軟件學報，2007，18(4):869-877.

［11］Kouropteva O，Okun O，Pietikainen M.Selection of the Optimal Parameter Value for the Locally Linear Embedding Algorithm［C］//Proc.of the 1st International Conference on Fuzzy Systems and Knowledge Discovery，Orchid Country Club.Singapore:IEEE Press，2002:359-363.

［12］Samko O，Marshall A D，Rosin P L.Selection of the Optimal Parameter Value for the Isomap Algorithm［J］.Pattern Recognition Letters，2006，27(1):968-979.

［13］黃啟宏，劉 釗.流形學習中非線性維數約簡方法概述［J］.計算機應用研究，2007，24(11):19-25.

［14］Saxena A，Gupta A，Mukerjee A.Non-linear Dimensionality Reduction by Locally Linear Isomaps［C］//Proc.of the 11th International Conference on Neural Information Processing.Calcutta，India:Springer，2004:1038-1043.

［15］Wen G，Jiang L，Shadbolt N R.Using Graph Algebra to Optimize Neighborhood for Isometric Mapping［C］//Proc.of the 20th International Joint Conference on Artificial Intelligence.Hyderabad，India:AAAI Press，2007，2398-2403.

［16］Carreira-Perpinan M A，Zemel R S.Proximity Graphs for Clustering and Manifold Learning［C］//Proc.of the 18th Annual Conference on Neural Information Processing Systems.Vancouver，Canada:MIT Press，2004:225-232.

［17］Yang L.Building k Edge-disjoint Spanning Trees of Minimum Total Length for Isometric Data Embedding［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，2005，27(10):1680-1683.

［18］Zhang Z，Wang J，Zha H.Adaptive Manifold Learning［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，2012，34(2):1473-1480.

［19］文貴華，江麗君，文 軍.鄰域參數動態變化的局部線性嵌入［J］.軟件學報，2008，19(7):1666-1673.

［20］邵 超，張 斌，萬春紅.流形學習中鄰域大小參數的合適性判定［J］.計算機工程與應用，2010，46(20):172-175.

［21］Chang H，Yeung D.Robust Locally Linear Embedding［J］.Pattern Recognition，2006，39(6):1053-1065.

［22］Pelleg D，Moore A.X-means:Extending K-means With Efficient Estimation of the Number of Clusters［C］//Proc.of the 17th International Conference on Machine Learning.San Francisco，USA:Morgan Kaufmann Publishers，2000:727-734.