寫在前面

前一段時間看到一篇文章，建議學生時代寫程式碼不要光呼叫庫和複製貼上，而是要儘量每一行程式碼都自己寫。因為以後工作的時候都主要是用別人寫好的東西，就沒有這樣鍛鍊基本功的機會了。

筆者最近入門python，希望能夠通過這些重複造輪子的簡單工作來加強基本功，讀者如果有什麼建議可以在評論區指出。

筆記

matrix, array 還是 list?

numpy包中含有matrix和array兩種常見的資料型別，python內建的list也和他們非常相像。筆者一開始將這三種資料型別混合使用，結果發現要不停地轉換，非常地麻煩，那麼如何解決呢？

首先來比較Matrix和Array兩種資料型別，直接看numpy官網給的結論。

The array is thus much more advisable to use. Indeed, we intend to deprecate matrix eventually. 我們更建議使用array。事實上，我們最終準備棄用matrix這個型別。

那麼在array和list之間又應該如何選擇呢？ 事實上numpy針對資料的運算做了很多優化，並且大多數的第三方包函式返回的也都是array。因此至少當你在使用Numpy包處理資料的時候，array是一個更好的選擇。

可以看到array是一種更好的資料型別，筆者因此將程式碼中的list和marix全改成了array，程式碼確實順暢了很多。

陣列排序：argsort

這個排序函式筆者一開始老搞錯，在這裡做一筆記：

argsort:

文件如下：

argsort(a, axis=-1, kind=‘quicksort’, order=None) Perform an indirect sort along the given axis using the algorithm specified by the ‘kind’ keyword. It returns an array of indices of the same shape as ‘a’ that index data along the given axis in sorted order. 沿給定的軸執行一種非直接的排序，排序方法可以由kind指定。返回一個和待排陣列’a’相同大小的陣列，該陣列給出了演給定座標軸和排序順序排序後的下標。

axis = 0：按行排列，axis=1：按列排列

例子：

>> x = np.array([6, 4, 2, 5, 5])
>> ind = np.argsort(x) 
>> print(ind)
>> print(x[ind])
[2 1 3 4 0]
[2 4 5 5 6]

對於一維陣列，可以看到argsort返回的是原陣列從小到大排列的下標，ind[0]是最小數字的下標，ind[-1]是最大數字的下標，可以通過x[ind]取得陣列x從小到大的排序序列。

對於多維陣列有：

>> x = np.array([[0, 3, 5], [2, 2, 1], [5, 0, 3]])
>> ind = np.argsort(x)
>> print(x)
>> print(ind)
[[0 3 5]
 [2 2 1]
 [5 0 3]]
[[0 1 2]
 [2 0 1]
 [1 2 0]]

可以看到對於二維陣列不能再用x[ind]，但可以用x[ind[i, :]]取得每個一維陣列第i大的排序值。

程式碼

優化前

import csv
import numpy as np
import cmath
import matplotlib.pyplot as plt
from sklearn import preprocessing as prp
import scipy


def k_means(k, x):
    init_center = np.random.randint(0, len(x) - 1, k)
    center = [x[i, :] for i in init_center]
    cost_new = 1
    cost_old = 0
    while abs(cost_old - cost_new) >= 0.0001:
        # 將各樣本分配到離他們最近的center
        d = []
        for j in range(k):
            d_sam_cen1 = [np.linalg.norm(center[j] - x[i, :]) for i in range(len(x))]
            d.append(d_sam_cen1)
        h = np.mat(d).T.argsort()
        flag = np.argwhere(h == 0)[:, 1]
        index = [np.argwhere(flag == i) for i in range(k)]
        # 將各樣本的平均值作為新的center
        group = []
        for i in range(k):
            group_i = [x[index[i][j]] for j in range(len(index[i]))]
            group.append(group_i)
        center = [np.mean(group[i], axis=0) for i in range(k)]
        d2 = []
        for i in range(k):
            for j in range(len(group[i])):
                d2.append(np.linalg.norm(center[i] - group[i][j]))
        cost_old = cost_new
        cost_new = np.array(d2).var()
    print(cost_new)
    return flag, group

# 資料匯入
filename = 'E:\Administrator\OneDrive\文件\學校課程\模式識別\pr_homework\homework_03_kmeans\dataset_circles.csv'
with open(filename) as f:
    reader = csv.reader(f)
    data = list(reader)
# 資料預處理
data = np.mat(data).astype(np.float)
dataRec = data[:, 0:2]
dataFlag = data[:, 2]
dataComplex = [complex(dataRec[i, 0], dataRec[i, 1])for i in range(len(dataRec))]
dataPolar = np.mat([cmath.polar(dataComplex[i]) for i in range(len(dataComplex))])
ss = prp.StandardScaler().fit(dataPolar)
staPolar = ss.transform(dataPolar)
# 資料計算
k = 2
staFlag, staGro = k_means(k, staPolar)
unStaFlag, unStaGro = k_means(k, dataPolar)
# 資料恢復
staGro = [np.matrix(np.array(staGro[i])) for i in range(k)]
unStaGro = [np.matrix(np.array(unStaGro[i])) for i in range(k)]
unStaRec = []
staGroRec = []
for i in range(k):
    unStaRec.append([cmath.rect(unStaGro[i][j, 0], unStaGro[i][j, 1])
                     for j in range(len(unStaGro[i]))])
    unStaRec[i] = np.matrix([[unStaRec[i][j].real, unStaRec[i][j].imag]
                             for j in range(len(unStaRec[i]))])
    temp_i = ss.inverse_transform(staGro[i])
    staGroRec.append([cmath.rect(temp_i[j, 0], temp_i[j, 1])
                     for j in range(len(temp_i))])
    staGroRec[i] = np.matrix([[staGroRec[i][j].real, staGroRec[i][j].imag]
                             for j in range(len(staGroRec[i]))])
# print(cmath.rect(unStaGro[0][1, 0], unStaGro[0][1, 0]))
# print(un)

# np.mat(unstaflag)
plt.figure()
plt.title('without standardizing')
for i in range(k):
    style = ['ro', 'bo', 'yo', 'go', 'co']
    plt.plot(unStaRec[i][:, 0], unStaRec[i][:, 1], style[i])
plt.show()
plt.figure()
plt.title('with standardizing')
for i in range(k):
    style = ['ro', 'bo', 'yo', 'go', 'co']
    plt.plot(staGroRec[i][:, 0], staGroRec[i][:, 1], style[i])
plt.show()

優化後

import numpy as np
import cmath
import matplotlib.pyplot as plt
from sklearn import preprocessing as prp


def k_means(k, x):
    # x是一個(m,2)維的陣列
    init_center = np.random.randint(0, len(x) - 1, k)
    center = [x[i, :] for i in init_center]
    cost_new = 1
    cost_old = 0
    flag = 0
    while cost_old != cost_new:
        # 將各樣本分配到離他們最近的center
        d = []
        for j in range(k):
            d_sam_cen1 = [np.linalg.norm(center[j] - x[i, :]) for i in range(len(x))]
            d.append(d_sam_cen1)
        flag = np.array(d).T.argsort()[:, 0]
        group = [x[np.argwhere(flag == i)[:, 0]] for i in range(k)]
        # 將各樣本的平均值作為新的center
        center = [np.mean(group[i], axis=0) for i in range(k)]
        # 計算新center的代價函式
        d2 = [np.linalg.norm(center[i] - group[i], axis=1) for i in range(k)]  # [array*3]
        cost_old = cost_new
        cost_new = np.mean([np.mean(np.array(d2[i])) for i in range(k)])
    return flag


def rectangular_polar(x):
    # x 位(m, 2)維的陣列，第一列維x座標，第二列為y座標
    x_complex = [complex(x[i, 0], x[i, 1]) for i in range(len(x))]
    x_polar = np.array([cmath.polar(x_complex[i]) for i in range(len(x_complex))])
    return x_polar


if __name__ == '__main__':
    # 資料匯入和預處理
    data = np.loadtxt('dataset_circles.csv', delimiter=',')
    dataRect = data[:, 0:2]
    dataPolar = rectangular_polar(dataRect)
    staPolar = prp.StandardScaler().fit(dataPolar).transform(dataPolar)
    # 資料計算
    k1 = 2
    staFlag = k_means(k1, staPolar)
    unStaFlag = k_means(k1, dataPolar)
    # 資料plot
    plt.scatter(data[:, 0], data[:, 1], c=staFlag)
    plt.show()
    plt.scatter(data[:, 0], data[:, 1], c=unStaFlag)
    plt.show()

除了k_means這個函式本身，我對載入資料和畫圖部分也都做了改進，感覺還是蠻有收穫的。