【机器学习】快速有效理解 K-Means 算法

2019-01-14 18:06:21 浏览数 (1)

什么是 K-Means ?

学习 K-Means 之前,大家首先需要对聚类有一个概念.

我们都知道,机器学习可以划分为 3 类:监督学习、无监督学习、强化学习.

无监督学习指的是数据没有标签,也就是说我们只有数据的特征,但并不知道这些数据都是什么,无监督学习算法或者是模型需要从这样的数据中学习给数据按照某种规律进行分类的能力,或者是找出不同特征之间的关联性等等.

聚类(Clustering)也归属于无监督学习,通过对所有数据进行特征分析,然后把相似的对象划分到同一个集合中,这个集合我们称之为簇.

它需要给一堆没有数据标签的数据找到给定数量的簇,每个簇的周围聚拢了许多数据,所以就可以近似认为聚焦在同一个簇周围的数据就是一个类别,这样达到了分类的效果.

比如,给一张电子表格给你,里面有 10000 个成员的信息.这些信息包含了 5 个特征,分别是性别年龄身高体重颜值.

姓名

性别

年龄

身高

体重

颜值

person0

0

24

177

60

100

person1

0

20

177

70

98

person2

1

29

177

70

89

person3

0

24

177

70

60

person4

1

20

170

70

78

现在,给你的任务就是强行给这 10000 人划分为 3 组,划分的过程就叫做聚类.

而 K-means 就是一种可以完成这样任务的算法.

如何理解 K-Means?

K 是指目标的数量,也是刚刚提到簇的数量,比如你要分 3 个类别的话,K = 3. K 的大小不固定,按照实际需求而定.

Means 是均值的意思,这说明 K-Means 和数据的均值息息相关.

在 K-means 中,每个簇由它对应的质心(centroid)表示,质心是这个簇中所有点的中心,所谓均值也就是这个意思.

值得注意的是我们这里讲的点是指代数据集中的每一条记录的,每条记录其实就是一个向量,所以质心也不是真正的一个点,它也是一个向量.

算法思想

K-means 具体算法是什么呢?

其实非常的简单和容易被人理解.

  1. 随机生成 k 个质心.
  2. 将数据集中每条记录,在这里可以当成一个点,分配到一个簇当中.分配的原则是将点对照 k 个质心,分配到离它最近的质心所在的那个簇.
  3. 因为初始的质心取值是随机的,所以结果可能不正确.所以,算法需要不停地更新来修正.修正的依据是步骤 2 分配完所有的数据后,将每个簇的质心的取值更新为整个簇的平均值.
  4. 2 ~ 3 是一个周期,整个算法需要多个周期才能完成.完成的条件是步骤 2 中,数据集中所有的点都依附在当前质心上,也可以理解为,每次对整个簇进行平均值时,得到的结果都是一样的.

代码演练

依照前面的设想,假如我们手里有一份清单,清单里面是每个成员的信息,但为了演示的方便,我只保留了 2 个特征.

姓名

身高

颜值

frank

177

100

person1

177

98

person2

177

89

现在有 10 个成员,目标把他们分成 3 类,那么怎么编写代码?

我们用 Python 实现它.

生成 K 个随机质心

虽然说是随机,但我们应该都能理解这个随机也应该是有约束条件的.

比如下面的图,方块代表数据集合中的数据,圆形代表生成的质点,显然生成的质点就不合理.

这里抛出一个问题,大家思考一下为什么这样不合理?

正确的质点应该在所有的数据集合中随机生成.

我们可以让它在数据值中最小值和最大值中生成随机数.

代码语言:javascript复制
def rand_centroid(dataset,k):
    n = dataset.shape[1]
    centroid = np.zeros((k,n))

    for i in range(n):
        min = np.min(dataset[:,i])
        max = np.max(dataset[:,i])
        stride = max - min

        print(" min ",min," max ",max)
        centroid[:,i] = min   np.random.rand(k).T*stride

    return centroid

我们可以测试一下这个方法.

代码语言:javascript复制
def test_centroid:

    dataset = np.array([[177,99],
            [169,80],
            [170,88],
            [190,86]],dtype=np.float)

    print(" centroid ",rand_centroid(dataset,2))

结果如下:

代码语言:javascript复制
(' min ', 169.0, ' max ', 190.0)
(' min ', 80.0, ' max ', 99.0)
(' centroid ', array([[ 189.43140332,   87.76774545],
       [ 175.81451753,   89.8356134 ]]))

为数据集中每一个点分配簇

依照前面的流程图,这一步中有几个关键信息.

将数据集中的点分配到簇的过程当中,要计算它们到质心的距离.

在机器学习领域,距离的表示有许多种,这里采用欧氏距离.

把每个点分配到最近的质心,这就要求再建立一张专门的表,用来保存对应的信息.这张表分为 2 列,第一列是对应到数据集中每一行数据所分配到的簇的序号,第二列对应每个数据到质心的距离.

迭代与更新

每一轮簇的更新之后,需要将质心重新调整为整个簇的数据平均值,直到所有的点都不再需要更新.

前面有提到过,更新的终止条件是数据集中所有的点都依附在当前质心上,也可以理解为,每次对整个簇进行平均值时,得到的结果都是一样的.

代码语言:javascript复制
def get_dist(A,B):
    return np.sqrt(np.sum(np.square(A - B)))
  
def inference(dataset,k):
    record_n = dataset.shape[0]
    cluster_table = np.zeros((record_n,2))

    centroid = rand_centroid(dataset,k)

    print("centroid ",centroid)

    need_update = True

    while need_update:
        need_update = False

        for i in range(record_n):
            min_index = 0
            min_dist = np.inf;
            for j in range(k):
                dist = get_dist(dataset[i,:],centroid[j,:])
                print("dist ",dist)
                if dist < min_dist : 
                    min_dist = dist
                    min_index = j
                   
            print('---------------------- index ',min_index)
            if cluster_table[i,0] != min_index:
                cluster_table[i,0] = min_index
                need_update = True
            cluster_table[i,0] = min_index
            cluster_table[i,1] = min_dist

        print("cluster_table ",cluster_table)
        
        for j in range(k):
            tmp = cluster_table[:,0] == j
            tmp = dataset[np.nonzero(tmp)[0]]
            # print ' cluster tmp ',tmp
            centroid[j,:] = np.mean(tmp,axis=0)

        # need_update = False

    return cluster_table,centroid 

代码也许很抽象,希望大家结合前面的流程图来加深理解。

一定要搞懂 cluster_tablecentroid 这两张表的意义。

  1. cluster 每一行与 dataset 每一行对应
  2. centroid 存放了 k 个质心的特征。

最终结果如何,我们可以编写测试代码验证.

代码语言:javascript复制
def test_kmeans():
    dataset = np.array([[177,99],
                [177,70],
                [180,72],
                [169,80],
                [170,88],
                [190,86]],dtype=np.float)

    # tmp = dataset[:,0] == 177
    # tmp = dataset > 170

    # print(" tmp ",np.nonzero(tmp)," ",np.nonzero(tmp)[0]," d ",dataset[np.nonzero(tmp)[0]])

    # print(" mean ",np.mean(dataset[np.nonzero(tmp)[0]],axis=0))

    print(inference(dataset,3))


test_kmeans()

输出的结果如下.

代码语言:javascript复制
(array([[ 2.        ,  9.19238816],
       [ 1.        ,  1.80277564],
       [ 1.        ,  1.80277564],
       [ 0.        ,  4.03112887],
       [ 0.        ,  4.03112887],
       [ 2.        ,  9.19238816]]), array([[ 169.5,   84. ],
       [ 178.5,   71. ],
       [ 183.5,   92.5]]))

当然,你可以用可视化的手段将算法过程直观表现出来,如下图:

上面的动画中,红色的就是质心。可以观察到,它是不停变动的。

其他 3 中颜色,代表三种类别的数据,数据会有少许变动,代表 K-Means 算法更新过程中,数据的分类发生了变化。

修改后的参考代码如下:

代码语言:javascript复制
import matplotlib.pyplot as plt
import numpy as np
from moviepy.editor import VideoClip
from moviepy.video.io.bindings import mplfig_to_npimage


global duration
global size
global fps
global img

def rand_centroid(dataset,k):
    n = dataset.shape[1]
    centroid = np.zeros((k,n))

    for i in range(n):
        min = np.min(dataset[:,i])
        max = np.max(dataset[:,i])
        stride = max - min

        print(" min ",min," max ",max)
        centroid[:,i] = min   np.random.rand(k).T*stride

    return centroid



def test_centroid():

    dataset = np.array([[177,99],
            [169,80],
            [170,88],
            [190,86]],dtype=np.float)

    print(" centroid ",rand_centroid(dataset,2))


def get_dist(A,B):
    return np.sqrt(np.sum(np.square(A - B)))
    # return np.sqrt(np.sum(np.square(A))) - np.sqrt(np.sum(np.square(B)))






def inference(dataset,k):
    record_n = dataset.shape[0]
    cluster_table = np.zeros((record_n,2))

    centroid = rand_centroid(dataset,k)

    print("centroid ",centroid)

    need_update = True

    while need_update:
        need_update = False

        for i in range(record_n):
            min_index = 0
            min_dist = np.inf;
            for j in range(k):
                dist = get_dist(dataset[i,:],centroid[j,:])
                print("dist ",dist)
                if dist < min_dist : 
                    min_dist = dist
                    min_index = j
                   
            print('---------------------- index ',min_index)
            if cluster_table[i,0] != min_index:
                cluster_table[i,0] = min_index
                need_update = True
            cluster_table[i,0] = min_index
            cluster_table[i,1] = min_dist

        print("cluster_table ",cluster_table)
        
        for j in range(k):
            tmp = cluster_table[:,0] == j
            tmp = dataset[np.nonzero(tmp)[0]]
            # print ' cluster tmp ',tmp
            centroid[j,:] = np.mean(tmp,axis=0)

        # need_update = False

    return cluster_table,centroid 

def make_frame(t):
    global duration
    global size
    global fps
    global img

    index = int(t * fps) - 1
    if index == -1 : index = 0
    print(" index ",index," t ",t)

    return img[index]

def inference_log(dataset,k):
    global img

    img = []

    record_n = dataset.shape[0]
    cluster_table = np.zeros((record_n,2))

    centroid = rand_centroid(dataset,k)

    print("centroid ",centroid)

    need_update = True

    duration = 2

    fig, ax = plt.subplots()

    color_list = ['m','y','c']
    marker_list = ['h','v','^']

    while need_update:
        need_update = False

        for i in range(record_n):
            min_index = 0
            min_dist = np.inf;
            for j in range(k):
                dist = get_dist(dataset[i,:],centroid[j,:])
                print("dist ",dist)
                if dist < min_dist : 
                    min_dist = dist
                    min_index = j
                   
            print('---------------------- index ',min_index)
            if cluster_table[i,0] != min_index:
                cluster_table[i,0] = min_index
                need_update = True
            cluster_table[i,0] = min_index
            cluster_table[i,1] = min_dist

        print("cluster_table ",cluster_table)
        
        ax.clear()
        for j in range(k):
            tmp = cluster_table[:,0] == j
            tmp = dataset[np.nonzero(tmp)[0]]
            ax.scatter(tmp[:,0],tmp[:,1],c=color_list[j],marker=marker_list[j])
            ax.set_ylim(60.0,100.0)
            # print ' cluster tmp ',tmp
            centroid[j,:] = np.mean(tmp,axis=0)

        print(" dataset ",dataset[:,0].T)

        
        
        ax.scatter(centroid[:,0],centroid[:,1],c='r',lineWidth=3)
        
        img_figure = mplfig_to_npimage(fig)
        img.append(img_figure)
        # need_update = False
    plt.close()

    return cluster_table,centroid 




# a = np.array([1,2,3],dtype=np.float)
# b = np.array([4,7,9],dtype=np.float)

# print(" a ",np.square(a),"sum ",np.sqrt(np.sum(np.square(a))))
# print ("dist ",get_dist(a,b))

def test_kmeans():
    global duration
    global size
    global fps
    global img

    dataset = np.zeros((100,2),dtype=np.float)

    dataset[:,0] = np.random.rand(100)*40.0 150.0
    dataset[:,1] = np.random.randint(0,40,size=100) 60.0

    # dataset = np.array([[177,99],
    #             [177,70],
    #             [172,70],
    #             [173,79],
    #             [174,70],
    #             [175,70],
    #             [180,72],
    #             [169,80],
    #             [170,88],
    #             [190,86]],dtype=np.float)

    # tmp = dataset[:,0] == 177
    # tmp = dataset > 170

    # print(" tmp ",np.nonzero(tmp)," ",np.nonzero(tmp)[0]," d ",dataset[np.nonzero(tmp)[0]])

    # print(" mean ",np.mean(dataset[np.nonzero(tmp)[0]],axis=0))

    # print(inference(dataset,3))

    inference_log(dataset,3)

    

    fps = 2

    time_span = 1.0 / fps



    size = img.__len__()

    duration = size * time_span

    print(" duration ",duration," time_span ",time_span," size ",size)

    animation = VideoClip(make_frame, duration=duration)
    animation.write_gif('kmeans.gif', fps=fps)

    img = None





test_kmeans()

思考

K-means 算法最关键的地方是什么?

我认为是距离,因为在机器学习领域内距离的表征有许多种,这里采用的是欧式距离. 但比如曼哈顿距离、余弦相似、马氏距离等等.

采用的距离表征不同,结果也不一样.

如何验证 K-means 结果?

K-means 运行结束后,但我们还得给这个算法进行评估,这样才能确保这个算法模型结果更精确.

因此,我们需要一些判断指标来验证 K-means 划分的每个簇是否合理.

SSE(误差的平方和)就是一种可用的判断表.

我们回到上面例子的输出结果,观察 cluster_table

代码语言:javascript复制
array([[ 2.        ,  9.19238816],
       [ 1.        ,  1.80277564],
       [ 1.        ,  1.80277564],
       [ 0.        ,  4.03112887],
       [ 0.        ,  4.03112887]

数组的第 0 列是簇的序号,而第 1 列存储的是数据集中每个点到质心的距离,我们称之为误差,误差越小,结果就越准确.

把每个簇中所有的误差取平方,然后相加的结果就是 SSE.如果一个簇的 SSE 太大了,就说明这个簇的数据太散了,不够紧凑.

如果要加以改进的话,我们可以找出 SSE 最大的簇,然后运用 K-means 对它进行再次划分.这个时候 k 一般取值 2.

但这种操作之后,簇的个数会增加 1 个,怎么办呢?

我们可以将距离很近的两个质心所代表的簇合并.

K-means 的变种二分 K-means

由于原始的 K-means 不够理想,有人就提出了另外一种叫做二分 k-means 的算法.

思路如下:

  1. 将所有数据在开始时归为一个簇.
  2. 将这个簇一分为二.
  3. 找出 SSE 最大或者最小的簇,继续划分.直到簇的个数等于 K.

如何挑选一个簇进行划分呢?依据是划分之后总的 SSE 能够最大程度降低.

当然,也可以暴力一点,每次划分找 SSE 最大的那个簇.

0 人点赞