By RaySaint 2011/10/12

本文的主要目的是講一下如何建立k-d tree對目標物體的特徵點集合進行資料組織和使用k-d tree最近鄰搜尋來加速特徵點匹配。上面已經講了特徵點匹配的問題其實上是一個最近鄰（K近鄰）搜尋的問題。所以為了更好的引出k-d tree，先講一講最近鄰搜尋。

最近鄰搜尋

先給出一個最近鄰的數學形式的定義。給定一個多維空間，把中的一個向量成為一個樣本點或資料點。中樣本點的有限集合稱為樣本集。給定樣本集E，和一個樣本點d，d的最近鄰就是任何樣本點d’∈E滿足None-nearer(E,d,d’)。

None-nearer如下定義：

上面的公式中距離度量是歐式距離，當然也可以是任何其他Lp-norm。

其中di是向量d的第i個分量。

現在再來說最近鄰搜尋，如何找到一個這樣的d’，它離d的距離在E中是最近的。

很容易想到的一個方法就是線性掃描，也稱為窮舉搜尋，依次計算樣本集E中每個樣本點到d的距離，然後取最小距離的那個點。這個方法又稱為樸素最近鄰搜尋。當樣本集E較大時（在物體識別的問題中，可能有數千個甚至數萬個SIFT特徵點），顯然這種策略是非常耗時的。

因為實際資料一般都會呈現簇狀的聚類形態，因此我們想到建立資料索引，然後再進行快速匹配。索引樹是一種樹結構索引方法，其基本思想是對搜尋空間進行層次劃分。k-d tree是索引樹中的一種典型的方法。

k-d tree的簡介及表示

k-d tree是英文K-dimension tree的縮寫，是對資料點在k維空間中劃分的一種資料結構。k-d tree實際上是一種二叉樹。每個結點的內容如下：

樣本集E由k-d tree的結點的集合表示，每個結點表示一個樣本點，dom_elt就是表示該樣本點的向量。該樣本點根據結點的分割超平面將樣本空間分為兩個子空間。左子空間中的樣本點集合由左子樹left表示，右子空間中的樣本點集合由右子樹right表示。分割超平面是一個通過點dom_elt並且垂直於split所指示的方向軸的平面。舉個簡單的例子，在二維的情況下，一個樣本點可以由二維向量(x,y)表示，其中令x維的序號為0，y維的序號為1。假設一個結點的dom_elt為(7,2) ，split的取值為0，那麼分割超面就是x=dom_elt(0)=7，它垂直與x軸且過點(7,2)，如下圖所示：

（紅線代表分割超平面）

於是其他資料點的x維（第split=0維）如果小於7，則被分配到左子空間；若大於7，則被分配到右子空間。例如，（5,4）被分配到左子空間，（9，6）被分配到右子空間。如下圖所示:

從上面的表也可以看出k-d tree本質上是一種二叉樹，因此k-d tree的構建是一個逐級展開的遞迴過程。

其演算法的虛擬碼如下：

1. 演算法：createKDTree 構建一棵k-d tree 
2. 輸入：exm_set 樣本集 
3. 輸出 : Kd, 型別為kd-tree 
4. 1. 如果exm_set是空的，則返回空的kd-tree 
5. 2.呼叫分裂結點選擇程式（輸入是exm_set），返回兩個值 
6.        dom_elt:= exm_set中的一個樣本點 
7.        split := 分裂維的序號 
8. 3.exm_set_left = {exm∈exm_set – dom_elt && exm[split] <= dom_elt[split]} 
9.    exm_set_right = {exm∈exm_set – dom_elt && exm[split] > dom_elt[split]} 
10. 4.left = createKDTree(exm_set_left) 
11. right = createKDTree(exm_set_right) 

現在來解釋一下分裂結點選擇程式。分裂結點的選擇通常有多種方法，最常用的是一種方法是：對於所有的樣本點，統計它們在每個維上的方差，挑選出方差中的最大值，對應的維就是split域的值。資料方差最大表明沿該維度資料點分散得比較開，這個方向上進行資料分割可以獲得最好的解析度；然後再將所有樣本點按其第split維的值進行排序，位於正中間的那個資料點選為分裂結點的dom_elt域。

下面以一個簡單的例子來解釋上述k-d tree的構建過程。假設樣本集為：{(2,3), (5,4), (9,6), (4,7), (8,1), (7,2)}。構建過程如下：

（1）確定split域，6個數據點在x,y維度上的資料方差分別為39, 28.63。在x軸上方差最大，所以split域值為0（x維的序號為0）

（2）確定分裂節點，根據x維上的值將資料排序，則6個數據點再排序後位於中間的那個資料點為(7,2)，該結點就是分割超平面就是通過(7,2)並垂直於split=0(x)軸的直線x=7

（3）左子空間和右子空間，分割超面x=7將整個空間分為兩部分，x<=7的部分為左子空間，包含3個數據點{(2,3), (5,4), (4,7)}；另一部分為右子空間，包含2個數據點{(9,6), (8,1)}。如下圖所示

（4）分別對左子空間中的資料點和右子空間中的資料點重複上面的步驟構建左子樹和右子樹直到經過劃分的子樣本集為空。下面的圖從左至右從上至下顯示了構建這棵二叉樹的所有步驟：

k-d tree的最近鄰搜尋演算法

如前所述，在k-d tree樹中進行資料的k近鄰搜尋是特徵匹配的重要環節，其目的是檢索在k-d tree中與待查詢點距離最近的k個數據點。

最近鄰搜尋是k近鄰的特例，也就是1近鄰。將1近鄰改擴充套件到k近鄰非常容易。下面介紹最簡單的k-d tree最近鄰搜尋演算法。

基本的思路很簡單：首先通過二叉樹搜尋（比較待查詢節點和分裂節點的分裂維的值，小於等於就進入左子樹分支，大於就進入右子樹分支直到葉子結點），順著“搜尋路徑”很快能找到最近鄰的近似點，也就是與待查詢點處於同一個子空間的葉子結點；然後再回溯搜尋路徑，並判斷搜尋路徑上的結點的其他子結點空間中是否可能有距離查詢點更近的資料點，如果有可能，則需要跳到其他子結點空間中去搜索（將其他子結點加入到搜尋路徑）。重複這個過程直到搜尋路徑為空。下面給出k-d tree最近鄰搜尋的虛擬碼：

1. 演算法：kdtreeFindNearest  
2. 輸入：Kd  
3. target  
4. 輸出 : nearest  
5. dist  
6. 1. 如果Kd是空的，則設dist為無窮大返回 
7. 2. 向下搜尋直到葉子結點 
8. pSearch = &Kd 
9. while(pSearch != NULL)  
10. {  
11. pSearch加入到search_path中;  
12. if(target[pSearch->split] <= pSearch->dom_elt[pSearch->split])   
13. {  
14. pSearch = pSearch->left;  
15. }  
16. else  
17. {  
18. pSearch = pSearch->right;  
19. }  
20. }  
21. 取出search_path最後一個賦給nearest 
22. dist = Distance(nearest, target);  
23. 3. 回溯搜尋路徑 
24. while(search_path不為空)  
25. {  
26. 取出search_path最後一個結點賦給pBack 
27. if(pBack->left為空 && pBack->right為空)  
28. { 
29. if( Distance(nearest, target) > Distance(pBack->dom_elt, target) )  
30. {  
31. nearest = pBack->dom_elt;  
32. dist = Distance(pBack->dom_elt, target);  
33. } 
34. } 
35. else 
36. { 
37. s = pBack->split;  
38. if( abs(pBack->dom_elt[s] - target[s]) < dist)   
39. {  
40. if( Distance(nearest, target) > Distance(pBack->dom_elt, target) )  
41. {  
42. nearest = pBack->dom_elt;  
43. dist = Distance(pBack->dom_elt, target);  
44. }  
45. if(target[s] <= pBack->dom_elt[s])   
46. pSearch = pBack->right;  
47. else  
48. pSearch = pBack->left;   
49. if(pSearch != NULL)  
50. pSearch加入到search_path中  
51. } 
52. }  
53. } 

OK，現在舉一些例子來說明上面的最近鄰搜尋演算法會比較直觀。

假設我們的k-d tree就是上面通過樣本集{(2,3), (5,4), (9,6), (4,7), (8,1), (7,2)}建立的。將上面的圖轉化為樹形圖的樣子如下：

我們來查詢點(2.1,3.1)，在(7,2)點測試到達(5,4)，在(5,4)點測試到達(2,3)，然後search_path中的結點為<(7,2), (5,4), (2,3)>，從search_path中取出(2,3)作為當前最佳結點nearest, dist為0.141；

然後回溯至(5,4)，以(2.1,3.1)為圓心，以dist=0.141為半徑畫一個圓，並不和超平面y=4相交，如下圖，所以不必跳到結點(5,4)的右子空間去搜索，因為右子空間中不可能有更近樣本點了。

於是在回溯至(7,2)，同理，以(2.1,3.1)為圓心，以dist=0.141為半徑畫一個圓並不和超平面x=7相交，所以也不用跳到結點(7,2)的右子空間去搜索。

至此，search_path為空，結束整個搜尋，返回nearest(2,3)作為(2.1,3.1)的最近鄰點，最近距離為0.141。

再舉一個稍微複雜的例子，我們來查詢點(2,4.5)，在(7,2)處測試到達(5,4)，在(5,4)處測試到達(4,7)，然後search_path中的結點為<(7,2), (5,4), (4,7)>，從search_path中取出(4,7)作為當前最佳結點nearest, dist為3.202；

然後回溯至(5,4)，以(2,4.5)為圓心，以dist=3.202為半徑畫一個圓與超平面y=4相交，如下圖，所以需要跳到(5,4)的左子空間去搜索。所以要將(2,3)加入到search_path中，現在search_path中的結點為<(7,2), (2, 3)>；另外，(5,4)與(2,4.5)的距離為3.04 < dist = 3.202，所以將(5,4)賦給nearest，並且dist=3.04。

回溯至(2,3)，(2,3)是葉子節點，直接平判斷(2,3)是否離(2,4.5)更近，計算得到距離為1.5，所以nearest更新為(2,3)，dist更新為(1.5)

回溯至(7,2)，同理，以(2,4.5)為圓心，以dist=1.5為半徑畫一個圓並不和超平面x=7相交, 所以不用跳到結點(7,2)的右子空間去搜索。

至此，search_path為空，結束整個搜尋，返回nearest(2,3)作為(2,4.5)的最近鄰點，最近距離為1.5。

兩次搜尋的返回的最近鄰點雖然是一樣的，但是搜尋(2, 4.5)的過程要複雜一些，因為(2, 4.5)更接近超平面。研究表明，當查詢點的鄰域與分割超平面兩側的空間都產生交集時，回溯的次數大大增加。最壞的情況下搜尋N個結點的k維kd-tree所花費的時間為：

後記

到此為止，k-d tree相關的基本知識就說完了。關於k-d tree還有很多擴充套件。由於大量回溯會導致kd-tree最近鄰搜尋的效能大大下降，因此研究人員也提出了改進的k-d tree近鄰搜尋，其中一個比較著名的就是 Best-Bin-First，它通過設定優先順序佇列和執行超時限定來獲取近似的最近鄰，有效地減少回溯的次數。這裡就不詳細講了，如果想知道可以查詢後面的參考資料。

參考資料

2.《影象區域性不變特性特徵與描述》王永明 王貴錦 編著 國防工業出版社

3.kdtree A simple C library for working with KD-Trees