KNN演算法

用NumPy庫實現K-nearest neighbors迴歸或分類。

鄰近演算法，或者說K最近鄰(kNN，k-NearestNeighbor)分類演算法是資料探勘分類技術中最簡單的方法之一。所謂K最近鄰，就是k個最近的鄰居的意思，說的是每個樣本都可以用它最接近的k個鄰居來代表。

kNN演算法的核心思想是如果一個樣本在特徵空間中的k個最相鄰的樣本中的大多數屬於某一個類別，則該樣本也屬於這個類別，並具有這個類別上樣本的特性。該方法在確定分類決策上只依據最鄰近的一個或者幾個樣本的類別來決定待分樣本所屬的類別。 kNN方法在類別決策時，只與極少量的相鄰樣本有關。由於kNN方法主要靠周圍有限的鄰近的樣本，而不是靠判別類域的方法來確定所屬類別的，因此對於類域的交叉或重疊較多的待分樣本集來說，kNN方法較其他方法更為適合。

簡單的理解，我有一組資料，比如每個資料都是n維向量，那麼我們可以在n維空間表示這個資料，這些資料都有對應的標籤值，也就是我們感興趣的預測變數。那麼當我們接到一個新的資料的時候，我們可以計算這個新資料和我們已知的訓練資料之間的距離，找出其中最近的k個數據，對這k個數據對應的標籤值取平均值就是我們得出的預測值。簡單粗暴，誰離我近，就認為誰能代表我，我就用你們的屬性作為我的屬性。具體的簡單程式碼實現如下。

1. 資料

這裡例子來自《集體智慧程式設計》給的葡萄酒的價格模型。葡萄酒的價格假設由酒的等級與儲藏年代共同決定，給定rating與age之後，就能給出酒的價格。

def wineprice
 
(rating,age):
    """
    Input rating & age of wine and Output it's price.
    Example:
    ------
    input = [80.,20.] ===> output = 140.0
    """
    peak_age = rating - 50 # year before peak year will be more expensive
    price = rating/2.
    if age > peak_age:
        price = price*(5 -(age-peak_age))
    else
 
:
        price = price*(5*((age+1)/peak_age))

    if price < 0: price=0
    return price

a = wineprice(80.,20.)
a

140.0

根據上述的價格模型，我們產生n=500瓶酒及價格，同時為價格隨機加減了20%來體現隨機性，同時增加預測的難度。注意基本資料都是numpy裡的ndarray,為了便於向量化計算，同時又有強大的broadcast功能，計算首選。

def wineset(n=500):
    """
    Input wineset size n and return feature array and target array.
    Example:
    ------
    n = 3
    X = np.array([[80,20],[95,30],[100,15]])
    y = np.array([140.0,163.6,80.0])
    """
    X,y  = [], []
    for i in range(n):
        rating = np.random.random()*50 + 50
        age = np.random.random()*50
        # get reference price
        price = wineprice(rating,age)
        # add some noise
        price = price*(np.random.random()*0.4 + 0.8) #[0.8,1.2]
        X.append([rating,age])
        y.append(price)
    return np.array(X), np.array(y)

X,y = wineset(500)

X[:3]

array([[ 88.89511317,  11.63751282],
       [ 91.57171713,  39.76279923],
       [ 98.38870877,  14.07015414]])

2. 相似度：歐氏距離

knn的名字叫K近鄰，如何叫“近”，我們需要一個數學上的定義，最常見的是用歐式距離，二維三維的時候對應平面或空間距離。

演算法實現裡需要的是給定一個新的資料，需要計算其與訓練資料組之間所有點之間的距離，注意是不同的維度，給定的新資料只是一個sample，而訓練資料是n組，n個sample，計算的時候需要注意，不過numpy可以自動broadcat，可以很好地take care of it。

def euclidean(arr1,arr2):
    """
    Input two array and output theie distance list.
    Example:
    ------
    arr1 = np.array([[3,20],[2,30],[2,15]])
    arr2 = np.array([[2,20],[2,20],[2,20]]) # broadcasted, np.array([2,20]) and [2,20] also work.
    d    = np.array([1,20,5])
    """
    ds = np.sum((arr1 - arr2)**2,axis=1)
    return np.sqrt(ds)

arr1 = np.array([[3,20],[2,30],[2,15]])
arr2 = np.array([[2,20],[2,20],[2,20]])
euclidean(arr1,arr2)

array([  1.,  10.,   5.])

提供一個簡潔的介面，給定訓練資料X和新sample v，然後返回排序好的距離，以及對應的index(我們要以此索引近鄰們對應的標籤值)。

def getdistance(X,v):
    """
    Input train data set X and a sample, output the distance between each other with index.
    Example:
    ------
    X = np.array([[3,20],[2,30],[2,15]])
    v = np.array([2,20]) # to be broadcasted
    Output dlist = np.array([1,5,10]), index = np.array([0,2,1])
    """
    dlist = euclidean(X,np.array(v))
    index = np.argsort(dlist)
    dlist.sort()
    # dlist_with_index = np.stack((dlist,index),axis=1)
    return dlist, index  

dlist, index = getdistance(X,[80.,20.])

3. KNN演算法

knn演算法具體實現的時候很簡單，呼叫前面的函式，計算出排序好的距離列表，然後對其前k項對應的標籤值取均值即可。可以用該knn演算法與實際的價格模型對比，發現精度還不錯。

def knn(X,y,v,kn=3):
    """
    Input train data and train target, output the average price of new sample.
    X = X_train; y = y_train
    k: number of neighbors
    """
    dlist, index = getdistance(X,v)
    avg = 0.0
    for i in range(kn):
        avg = avg + y[index[i]]
    avg = avg / kn
    return avg

knn(X,y,[95.0,5.0],kn=3)

32.043042600537092

wineprice(95.0,5.0)

31.666666666666664

4. 加權KNN

以上是KNN的基本演算法，有一個問題就是該演算法給所有的近鄰分配相等的權重，這個還可以這樣改進，就是給更近的鄰居分配更大的權重(你離我更近，那我就認為你跟我更相似，就給你分配更大的權重)，而較遠的鄰居的權重相應地減少，取其加權平均。需要一個能把距離轉換為權重的函式，gaussian函式是一個比較普遍的選擇，下圖可以看到gaussian函式的衰減趨勢。從下面的單例可以看出其效果要比knn演算法的效果要好，不過只是單個例子，不具有說服力，後面給出更可信的評價。

def gaussian(dist,sigma=10.0):
    """Input a distance and return it's weight"""
    weight = np.exp(-dist**2/(2*sigma**2))
    return weight

x1 = np.arange(0,30,0.1)
y1 = gaussian(x1)
plt.title('gaussian function')
plt.plot(x1,y1);

def knn_weight(X,y,v,kn=3):
    dlist, index = getdistance(X,v)
    avg = 0.0
    total_weight = 0
    for i in range(kn):
        weight = gaussian(dlist[i])
        avg = avg + weight*y[index[i]]
        total_weight = total_weight + weight
    avg = avg/total_weight
    return avg

knn_weight(X,y,[95.0,5.0],kn=3)

32.063929602836012

交叉驗證

寫一個函式，實現交叉驗證功能，不能用sklearn庫。

交叉驗證(Cross-Validation): 有時亦稱迴圈估計， 是一種統計學上將資料樣本切割成較小子集的實用方法。於是可以先在一個子集上做分析， 而其它子集則用來做後續對此分析的確認及驗證。 一開始的子集被稱為訓練集。而其它的子集則被稱為驗證集或測試集。常見交叉驗證方法如下：

Holdout Method(保留)

• 方法：將原始資料隨機分為兩組,一組做為訓練集,一組做為驗證集,利用訓練集訓練分類器,然後利用驗證集驗證模型,記錄最後的分類準確率為此Hold-OutMethod下分類器的效能指標.。Holdout Method相對於K-fold Cross Validation 又稱Double cross-validation ，或相對K-CV稱 2-fold cross-validation(2-CV)
• 優點：好處的處理簡單,只需隨機把原始資料分為兩組即可
• 缺點：嚴格意義來說Holdout Method並不能算是CV,因為這種方法沒有達到交叉的思想,由於是隨機的將原始資料分組,所以最後驗證集分類準確率的高低與原始資料的分組有很大的關係,所以這種方法得到的結果其實並不具有說服性.(主要原因是 訓練集樣本數太少，通常不足以代表母體樣本的分佈，導致 test 階段辨識率容易出現明顯落差。此外，2-CV 中一分為二的分子集方法的變異度大，往往無法達到「實驗過程必須可以被複制」的要求。)

K-fold Cross Validation(k摺疊)

• 方法：作為Holdout Methon的演進，將原始資料分成K組(一般是均分),將每個子集資料分別做一次驗證集,其餘的K-1組子集資料作為訓練集,這樣會得到K個模型,用這K個模型最終的驗證集的分類準確率的平均數作為此K-CV下分類器的效能指標.K一般大於等於2,實際操作時一般從3開始取,只有在原始資料集合資料量小的時候才會嘗試取2. 而K-CV 的實驗共需要建立 k 個models，並計算 k 次 test sets 的平均辨識率。在實作上，k 要夠大才能使各回合中的 訓練樣本數夠多，一般而言 k=10 (作為一個經驗引數)算是相當足夠了。
• 優點：K-CV可以有效的避免過學習以及欠學習狀態的發生,最後得到的結果也比較具有說服性.
• 缺點：K值選取上

Leave-One-Out Cross Validation(留一)

• 方法：如果設原始資料有N個樣本,那麼留一就是N-CV,即每個樣本單獨作為驗證集,其餘的N-1個樣本作為訓練集,所以留一會得到N個模型,用這N個模型最終的驗證集的分類準確率的平均數作為此下LOO-CV分類器的效能指標.
• 優點：相比於前面的K-CV,留一有兩個明顯的優點：
  
  • a.每一回閤中幾乎所有的樣本皆用於訓練模型,因此最接近原始樣本的分佈,這樣評估所得的結果比較可靠。
  • b. 實驗過程中沒有隨機因素會影響實驗資料,確保實驗過程是可以被複制的.
• 缺點：計算成本高,因為需要建立的模型數量與原始資料樣本數量相同,當原始資料樣本數量相當多時,LOO-CV在實作上便有困難幾乎就是不顯示,除非每次訓練分類器得到模型的速度很快,或是可以用並行化計算減少計算所需的時間。

Holdout method方法的想法很簡單，給一個train_size,然後演算法就會在指定的比例(train_size)處把資料分為兩部分，然後返回。

# Holdout method
def my_train_test_split(X,y,train_size=0.95,shuffle=True):
    """
    Input X,y, split them and out put X_train, X_test; y_train, y_test.
    Example:
    ------
    X = np.array([[0, 1],[2, 3],[4, 5],[6, 7],[8, 9]])
    y = np.array([0, 1, 2, 3, 4])
    Then
    X_train = np.array([[4, 5],[0, 1],[6, 7]])
    X_test = np.array([[2, 3],[8, 9]])
    y_train = np.array([2, 0, 3])
    y_test = np.array([1, 4])
    """
    order = np.arange(len(y))
    if shuffle:
        order = np.random.permutation(order)
    border = int(train_size*len(y))
    X_train, X_test = X[:border], X[border:]
    y_train, y_test = y[:border], y[border:]
    return X_train, X_test, y_train, y_test

K folds演算法是把資料分成k份，進行k此迴圈，每次不同的份分別來充當測試組資料。但是注意，該演算法不直接操作資料，而是產生一個迭代器，返回訓練和測試資料的索引。看docstring裡的例子應該很清楚。

# k folds 產生一個迭代器
def my_KFold(n,n_folds=5,shuffe=False):
    """
    K-Folds cross validation iterator.
    Provides train/test indices to split data in train test sets. Split dataset 
    into k consecutive folds (without shuffling by default).
    Each fold is then used a validation set once while the k - 1 remaining fold form the training set.
    Example:
    ------
    X = np.array([[1, 2], [3, 4], [1, 2], [3, 4]])
    y = np.array([1, 2, 3, 4])
    kf = KFold(4, n_folds=2)
    for train_index, test_index in kf:
        X_train, X_test = X[train_index], X[test_index]
        y_train, y_test = y[train_index], y[test_index]
        print("TRAIN:", train_index, "TEST:", test_index)
    TRAIN: [2 3] TEST: [0 1]
    TRAIN: [0 1] TEST: [2 3]
    """
    idx = np.arange(n)
    if shuffe:
        idx = np.random.permutation(idx)
    fold_sizes = (n // n_folds) * np.ones(n_folds, dtype=np.int) # folds have size n // n_folds
    fold_sizes[:n % n_folds] += 1 # The first n % n_folds folds have size n // n_folds + 1
    current = 0
    for fold_size in fold_sizes:
        start, stop = current, current + fold_size
        train_index = list(np.concatenate((idx[:start], idx[stop:])))
        test_index = list(idx[start:stop])
        yield train_index, test_index
        current = stop # move one step forward

X1 = np.array([[1, 2], [3, 4], [1, 2], [3, 4]])
y1 = np.array([1, 2, 3, 4])
kf = my_KFold(4, n_folds=2)

for train_index, test_index in kf:
    X_train, X_test = X1[train_index], X1[test_index]
    y_train, y_test = y1[train_index], y1[test_index]
    print("TRAIN:", train_index, "TEST:", test_index)

('TRAIN:', [2, 3], 'TEST:', [0, 1])
('TRAIN:', [0, 1], 'TEST:', [2, 3])

KNN演算法交叉驗證

萬事俱備只欠東風，已經實現了KNN演算法以及交叉驗證功能，我們就可以利用交叉驗證的思想為我們的演算法選擇合適的引數，這也是交叉驗證主要目標之一。

評價演算法

首先我們用一個函式評價演算法，給定訓練資料與測試資料，計算平均的計算誤差，這是評價演算法好壞的重要指標。

def test_algo(alg,X_train,X_test,y_train,y_test,kn=3):
    error = 0.0
    for i in range(len(y_test)):
        guess = alg(X_train,y_train,X_test[i],kn=kn)
        error += (y_test[i] - guess)**2
    return error/len(y_test)

X_train,X_test,y_train,y_test = my_train_test_split(X,y,train_size=0.8)

test_algo(knn,X_train,X_test,y_train,y_test,kn=3)

783.80937472673656

交叉驗證

得到平均誤差，注意這裡KFold 生成器的使用。

def my_cross_validate(alg,X,y,n_folds=100,kn=3):
    error = 0.0
    kf = my_KFold(len(y), n_folds=n_folds)
    for train_index, test_index in kf:
        X_train, X_test = X[train_index], X[test_index]
        y_train, y_test = y[train_index], y[test_index]
        error += test_algo(alg,X_train,X_test,y_train,y_test,kn=kn)
    return error/n_folds

選最好的k值

從下圖可以看出，模型表現隨k值的變化呈現一個折現型，當為1時，表現較差；當取2,3的時候，模型表現還不錯；但當k繼續增加的時候，模型表現急劇下降。同時knn_weight演算法要略優於knn演算法，有一點點改進。

errors1, errors2 = [], []
for i in range(20):
    error1 = my_cross_validate(knn,X,y,kn=i+1)
    error2 = my_cross_validate(knn_weight,X,y,kn=i+1)
    errors1.append(error1)
    errors2.append(error2)

xs = np.arange(len(errors1)) + 1
plt.plot(xs,errors1,color='c')
plt.plot(xs,errors2,color='r')
plt.xlabel('K')
plt.ylabel('Error')
plt.title('Error vs K');