python学习之旅 | k_means算法python实现

写在前面

前一段时间看到一篇文章，建议学生时代写代码不要光调用库和复制粘贴，而是要尽量每一行代码都自己写。因为以后工作的时候都主要是用别人写好的东西，就没有这样锻炼基本功的机会了。

笔者最近入门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这个类型。

numpy官网上对这两种数据类型的特性有详细的对比，可以参考网址：
https://docs.scipy.org/doc/numpy/user/numpy-for-matlab-users.html#array-or-matrix-which-should-i-use

那么在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这个函数本身，我对载入数据和画图部分也都做了改进，感觉还是蛮有收获的。