KNN(k-nearest neighbor algorithm)–從原理到實現

KNN(k-nearest neighbor algorithm)–從原理到實現

零.廣告

        本文所有程式碼實現均可以在
DML
 找到,不介意的話請大家在github裡給我點個Star大笑

一.引入 

K近鄰演算法作為資料探勘十大經典演算法之一,其演算法思想可謂是intuitive,就是從訓練集裡找離預測點最近的K個樣本來預測分類

        因為演算法思想簡單,你可以用很多方法實現它,這時效率就是我們需要慎重考慮的事情,最簡單的自然是求出測試樣本和訓練集所有點的距離然後排序選擇前K個,這個是O(nlogn)的,而其實從N個資料找前K個資料是一個很常見的演算法題,可以用最大堆(最小堆)實現,其效率是O(nlogk)的,而最廣泛的演算法是使用kd樹來減少掃描的點,這也就是這篇文章的主要內容,本文偏實現,詳細理論教程見july的文章 ,不得不服,july這篇文章鉅細無遺!

二.前提:堆的實現

        堆是一種二叉樹,用一個陣列儲存,對於k號元素,k*2號是其左兒子,k*2 1號是其右兒子

        而大根堆就是跟比左兒子和右兒子都大,小根堆反之。

        要滿足這個條件我們需要通過up( index )操作和down( index )維護它的結構

        當然講這個的文章實在有些多了,隨便搜一篇大家看看:點選開啟連結

        大小根堆的作用是

              a) 優先佇列:因為第一個元素是最大或者最小的元素,所以可以實現優先佇列

              b) 前K個最大(最小)值:這裡限制堆的大小為k,來獲得O( n log k)的效率,但注意此時小根堆是獲得前K個最大值,大根堆是獲得前K個最小值,插入的時候先把元素和堆頂比較再決定是否插入。

        因為事先KD-tree BBF 要同時用到這兩個東西,所以把它們實現在了同一個類裡,感覺程式碼略漂亮,貼出來觀賞一下:

        此程式碼是dml / tool / heap.py

  

from __future__ import division
import numpy as np
import scipy as sp
def heap_judge(a,b):
return a>b
class Heap:
def __init__(self,K=None,compare=heap_judge):
'''
'K'                 is the parameter to restrict the length of Heap
!!! when K is confirmed,the Min heap contain Max K elements
while Max heap contain Min K elements
'compare'         is the compare function which return a BOOL when pass two variable
default is Max heap
'''
self.K=K
self.compare=compare
self.heap=['#']
self.counter=0
def insert(self,a):
#print self.heap
if self.K!=None:
print a.x,'==='
if self.K==None:
self.heap.append(a)
self.counter =1
self.up(self.counter)
else:
if self.counter<self.K:
self.heap.append(a)
self.counter =1
self.up(self.counter)
else:
if (not self.compare(a,self.heap[1])):
self.heap[1]=a
self.down(1)
return
def up(self,index):
if (index==1):
return
'''
print index
for t in range(index 1):
if t==0:
continue
print self.heap[t].x
print 
'''
if self.compare(self.heap[index],self.heap[int(index/2)]):
#fit the condition
self.heap[index],self.heap[int(index/2)]=self.heap[int(index/2)],self.heap[index]
self.up(int(index/2))
return
def down(self,index):
if 2*index>self.counter:
return
tar_index=0
if 2*index<self.counter:
if self.compare(self.heap[index*2],self.heap[index*2 1]):
tar_index=index*2
else:
tar_index=index*2 1
else:
tar_index=index*2
if not self.compare(self.heap[index],self.heap[tar_index]):
self.heap[index],self.heap[tar_index]=self.heap[tar_index],self.heap[index]
self.down(tar_index)
return
def delete(self,index):
self.heap[index],self.heap[self.counter]=self.heap[self.counter],self.heap[index]
self.heap.pop()
self.counter-=1
self.down(index)
pass
def delete_ele(self,a):
try:
t=self.heap.index(a)
except ValueError:
t=None
if t!=None:
self.delete(t)
return t

          傳入的時候不設定K就是正常的優先佇列,設定了K就是限制堆的大小了

          compare引數是比較大小的,預設是“數”的大根堆,你可以往堆裡傳任何類,只要有相適應的compare引數,比如我們KD-tree傳的就是KD-Node

        

三.KD-BFF的原理:

          首先從KD-Tree的建立說起:(直接貼《統計學習方法》的內容了)

         

          事實上從選擇哪一個feature開始切割,還可以選擇方差最大的那個引數,但是考慮到簡便,以及我們可以選擇更多的相似性度量方法,還是用《統計學習方法》裡面的選擇方式了。

          然後是KD-tree搜尋的方法:(來自《統計學習方法》,但注意這裡是最近鄰,也就是k=1的時候)

         

          那麼我們要K近鄰要怎麼做呢?就是用堆的第二個應用,用大根堆保持K個最小的距離,然後用根的距離(也就是其中最大的一個)來作為判斷的依據是否有更近的點不在結果中,這一點很重要!

          同時摘錄july部落格的一段讀者留言講得非常好的:

              在某一層,分割面是第ki維,分割值是kv,那麼 abs(q[ki]-kv) 就是沒有選擇的那個分支的優先順序,也就是計算的是那一維上的距離; 同時,從優先佇列裡面取節點只在某次搜尋到葉節點後才發生,計算過距離的節點不會出現在佇列的,比如1~10這10個節點,你第一次搜尋到葉節點的路徑是1-5-7,那麼1,5,7是不會出現在優先佇列的。換句話說,優先佇列裡面存的都是查詢路徑上節點對應的相反子節點,比如:搜尋左子樹,就把對應這一層的右節點存進佇列。

          大致這就是我們實現的基本思路了

四.KD-BFF的實現:

        知道原理了,並且有了堆這個工具之後我們就可以著手實現這個演算法了:(終於要貼程式碼了)

       程式碼~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~此程式碼是
dml / KNN / kd.py

         

from __future__ import division
import numpy as np
import scipy as sp
from operator import itemgetter
from scipy.spatial.distance import euclidean
from dml.tool import Heap
class KDNode:
def __init__(self,x,y,l):
self.x=x
self.y=y
self.l=l
self.F=None
self.Lc=None
self.Rc=None
self.distsToNode=None
class KDTree:
def __init__(self,X,y=None,dist=euclidean):
self.X=X
self.k=X.shape[0] #N
self.y=y
self.dist=dist
self.P=self.maketree(X,y,0)
self.P.F=None
def maketree(self,data,y,deep):
if data.size==0:
return None
lenght = data.shape[0]
case = data.shape[1]
p=int((case)/2)
l = (deep%self.k)
#print data
data=np.vstack((data,y))
data=np.array(sorted(data.transpose(),key=itemgetter(l))).transpose()
#print data
y=data[lenght,:]
data=data[:lenght,:]
v=data[l,p]
rP=KDNode(data[:,p],y[p],l)
#print data[:,p],y[p],l
if case>1:
ldata=data[:,data[l,:]<v]
ly=y[data[l,:]<v]
data[l,p]=v-1
rdata=data[:,data[l,:]>=v]
ry=y[data[l,:]>=v]
data[l,p]=v
rP.Lc=self.maketree(ldata,ly,deep 1)
if rP.Lc!=None:
rP.Lc.F=rP
rP.Rc=self.maketree(rdata,ry,deep 1)
if rP.Rc!=None:
rP.Rc.F=rP
return rP
def search_knn(self,P,x,k,maxiter=200):
def pf_compare(a,b):
return self.dist(x,a.x)<self.dist(x,b.x)
def ans_compare(a,b):
return self.dist(x,a.x)>self.dist(x,b.x)
pf_seq=Heap(compare=pf_compare)
pf_seq.insert(P)    #prior sequence
ans=Heap(k,compare=ans_compare)  #ans sequence
while pf_seq.counter>0:
t=pf_seq.heap[1]
pf_seq.delete(1)
flag=True
if ans.counter==k:
now=t.F
#print ans.heap[1].x,'========'
if now != None:
q=x.copy()
q[now.l]=now.x[now.l]
length=self.dist(q,x)
if length>self.dist(ans.heap[1].x,x):
flag=False
else:
flag=True
else:
flag=True
if flag:
tp,pf_seq,ans=self.to_leaf(t,x,pf_seq,ans)
#print "============="
#ans.insert(tp)
return ans
def to_leaf(self,P,x,pf_seq,ans):
tp=P
if tp!=None:
ans.insert(tp)
if tp.x[tp.l]>x[tp.l]:
if tp.Rc!=None:
pf_seq.insert(tp.Rc)
if tp.Lc==None:
return tp,pf_seq,ans
else:
return self.to_leaf(tp.Lc,x,pf_seq,ans)
if tp.Lc!=None:
pf_seq.insert(tp.Lc)
if tp.Rc==None:
return tp,pf_seq,ans
else:
return self.to_leaf(tp.Rc,x,pf_seq,ans)

       

        然後KNN就是對上面這個類的一個包裝:

        程式碼~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~此程式碼是
dml / KNN / knn.py

    

#coding:utf-8 
import numpy as np
import scipy as sp
from scipy.spatial.distance import cdist
from scipy.spatial.distance import euclidean
from dml.KNN.kd import KDTree
#import pylab as py
class KNNC:
"""docstring for KNNC"""
def __init__(self,X,K,labels=None,dist=euclidean):
'''
X is a N*M matrix where M is the case 
labels is prepare for the predict.
dist is the similarity measurement way,
The distance function can be ‘braycurtis’, ‘canberra’, 
‘chebyshev’, ‘cityblock’, ‘correlation’, ‘cosine’, 
‘dice’, ‘euclidean’, ‘hamming’, ‘jaccard’, ‘kulsinski’, 
‘mahalanobis’, 
'''
self.X = np.array(X)
if labels==None:
np.zeros((1,self.X.shape[1]))
self.labels = np.array(labels)
self.K = K
self.dist = dist
self.KDTrees=KDTree(X,labels,self.dist)
def predict(self,x,k):
ans=self.KDTrees.search_knn(self.KDTrees.P,x,k)
dc={}
maxx=0
y=0
for i in range(ans.counter 1):
if i==0:
continue
dc.setdefault(ans.heap[i].y,0)
dc[ans.heap[i].y] =1
if dc[ans.heap[i].y]>maxx:
maxx=dc[ans.heap[i].y]
y=ans.heap[i].y
return y
def pred(self,test_x,k=None):
'''
test_x is a N*TM matrix,and indicate TM test case
you can redecide the k
'''
if k==None:
k=self.K
test_case=np.array(test_x)
y=[]
for i in range(test_case.shape[1]):
y.append(self.predict(test_case[:,i].transpose(),k))
return y

       因為KNN畢竟是一個分類演算法,所以我在predict是加上了分類的程式碼,如果只想檢驗Kd-tree的話,你可以直接用for_point()找最近k個點

五.測試 後記

       測試:

       我們選取《統計學習方法》上面的例子:

     

       使用程式碼:

          

X=np.array([[2,5,9,4,8,7],[3,4,6,7,1,2]])
y=np.array([2,5,9,4,8,7])
knn=KNNC(X,1,y)
print knn.for_point([[6.5],[7]],1)

這裡y是label,是用來預測的,這個例子裡沒有實際作用,這是用來分類的

       輸出中後面帶了“===”的是掃描過的點,最後的是搜尋的結果:

       

       我們可以看到的確避免掃描了(2,3),Bingo!!

       我們再knn.for_point([[2],[2]]):可以看到避免掃了很多點!!!

      

      

       後記:

       從實現寫此文前後耗時兩天,昨天寫程式碼寫到熄燈且剛好測試通過,怎一個爽字了得!!最後,再在github上求個Star

reference:

【1】從K近鄰演算法、距離度量談到KD樹、SIFT BBF演算法 http://blog.csdn.net/v_july_v/article/details/8203674

【2】《統計學習方法》 李航

【3】最大堆的插入/刪除/調整/排序操作(圖解 程式)   http://www.java3z.com/cwbwebhome/article/article1/1362.html?id=4745