Python实现Mean Shift算法

刺骨的言语ヽ痛彻心扉 2021-09-29 12:58 328阅读 0赞

声明：代码的运行环境为Python3。Python3与Python2在一些细节上会有所不同，希望广大读者注意。本博客以代码为主，代码中会有详细的注释。相关文章将会发布在我的个人博客专栏[《Python从入门到深度学习》][Python]，欢迎大家关注~

--------------------

在K-Means算法中，聚类的类别个数需要提前指定，对于类别个数未知的数据集，K-Means算法和K-Means++算法将很难对其进行求解，所以需要一些能够处理未知类别个数的算法来处理此类问题。Mean Shift算法，又称作均值漂移算法，它跟K-Means算法一样，都是基于聚类中心的聚类算法，不同的是，它不需要提前指定聚类中心的个数，聚类中心是通过在给定区域中样本的均值来确定的，通过不断更新聚类中心，直至聚类中心不再改变为止。

**一、Mean Shift向量与核函数**

1、Mean Shift向量

对于给定的n维空间![20190713121057956.png][]中的m个样本点![20190713121156601.png][]，对于其中的一个样本X，其Mean Shift的向量为：

![20190713121445955.png][]

其中，![20190713121635531.png][]指的是一个半径为h的高维球区域，![20190713121709361.png][]定义为：

![20190713121832294.png][]

2、核函数

通过上述方式求出的Mean Shift向量时存在问题的，即在![20190713121635531.png][]区域内每一个![20190713122316225.png][]对样本X的贡献是一样的，然而实际上，每一个样本![20190713122316225.png][]对样本X的贡献是不一样的，我们可以通过核函数对每一个样本的贡献进行度量。

核函数的定义如下：

设Z是输入空间，H是特征空间，如果存在一个Z到H的映射：![20190713124429601.png][]使得所有![20190713124541100.png][]，函数![20190713124616313.png][]满足条件：![20190713124712403.png][]，则称![20190713124616313.png][]为核函数，![20190713124754156.png][]为映射函数。

我们在Mean Shift算法中使用的是高斯核函数，这也是最常用的核函数之一，高斯核函数的表达式为：

![20190713125150741.png][]

其中，h为带宽，当带宽一定时，样本点之间的距离越近，核函数的值越大；当样本点距离一定时，带宽越大，核函数的值越小。

下面我们使用Python代码实现高斯核函数：

import numpy as np
    import math
    
    def gs_kernel(dist, h):
        '''
        高斯核函数
        :param dist: 欧氏距离
        :param h: 带宽
        :return: 返回高斯核函数的值
        '''
        m = np.shape(dist)[0]  # 样本个数
        one = 1 / (h * math.sqrt(2 * math.pi))
        two = np.mat(np.zeros((m, 1)))
        for i in range(m):
            two[i, 0] = (-0.5 * dist[i] * dist[i].T) / (h * h)
            two[i, 0] = np.exp(two[i, 0])
        
        gs_val = one * two
        return gs_val

**二、Mean Shift原理**

在Mean Shift中通过迭代的方式找到最终的聚类中心，即对每一个样本点计算其漂移均值，以计算出来的漂移均值点作为新的起始点重复上述步骤，直到满足终止条件，得到的最终的漂移均值点即为最终的聚类中心。

Mean Shift算法实现过程如下：

def mean_shift(points, h=2, MIN_DISTANCE=0.000001):
        '''
        训练Mean Shift模型
        :param points: 特征点
        :param h: 带宽
        :param MIN_DISTANCE: 最小误差
        :return: 返回特征点、均值漂移点、类别
        '''
        mean_shift_points = np.mat(points)
        max_min_dist = 1
        iteration = 0  # 迭代的次数
        m = np.shape(mean_shift_points)[0]  # 样本的个数
        need_shift = [True] * m  # 标记是否需要漂移
    
        # 计算均值漂移向量
        while max_min_dist > MIN_DISTANCE:
            max_min_dist = 0
            iteration += 1
            print("iteration : " + str(iteration))
            for i in range(0, m):
                if not need_shift[i]:  # 判断每一个样本点是否需要计算偏移均值
                    continue
                point_new = mean_shift_points[i]
                point_new_start = point_new
                point_new  = shift_point(point_new, points, h)  # 对样本点进行漂移计算
                dist = distince(point_new, point_new_start)  # 计算该点与漂移后的点之间的距离
                
                if dist > max_min_dist:
                    max_min_dist = dist
                if dist < MIN_DISTANCE:
                    need_shift[i] = False
                
                mean_shift_points[i] = point_new
        # 计算最终的类别
        lb = lb_points(mean_shift_points)  # 计算所属的类别
        return np.mat(points), mean_shift_points, lb

其中，shift\_point()方法目的在于计算漂移量，lb\_points()方法的目的在于计算最终所属分类，distance()方法用于计算欧氏距离，三个方法的实现过程分别如下：

（1）shift\_point()方法

def shift_point(point, points, h):
        '''
        计算漂移向量
        :param point: 需要计算的点
        :param points: 所有的样本点
        :param h: 带宽
        :return: 返回漂移后的点
        '''
        points = np.mat(points)
        m = np.shape(points)[0]  # 样本的个数
        # 计算距离
        point_dist = np.mat(np.zeros((m, 1)))
        for i in range(m):
            point_dist[i, 0] = distince(point, points[i])
    
        # 计算高斯核函数
        point_weights = gs_kernel(point_dist, h)
    
        # 计算分母
        all_sum = 0.0
        for i in range(m):
            all_sum += point_weights[i, 0]
    
        # 计算均值偏移
        point_shifted = point_weights.T * points / all_sum
        return point_shifted

（2）lb\_points()

def lb_points(mean_shift_points):
        '''
        计算所属类别
        :param mean_shift_points: 漂移向量
        :return: 返回所属的类别
        '''
        lb_list = []
        m, n = np.shape(mean_shift_points)
        index = 0
        index_dict = {}
        for i in range(m):
            item = []
            for j in range(n):
                item.append(str(("%5.2f" % mean_shift_points[i, j])))
    
            item_1 = "_".join(item)
            if item_1 not in index_dict:
                index_dict[item_1] = index
                index += 1
    
        for i in range(m):
            item = []
            for j in range(n):
                item.append(str(("%5.2f" % mean_shift_points[i, j])))
    
            item_1 = "_".join(item)
            lb_list.append(index_dict[item_1])
        return lb_list

（3）distince()方法

def distince(pointA, pointB):
        '''
        计算欧氏距离
        :param pointA: A点坐标
        :param pointB: B点坐标
        :return: 返回得到的欧氏距离
        '''
        return math.sqrt((pointA - pointB) * (pointA - pointB).T)

**三、Mean Shift算法举例**

1、数据集：数据集含有两个特征，如下图所示：

![20190713145124459.png][]

2、加载数据集

我们此处使用如下方法加载数据集，也可使用其他的方式进行加载，此处可以参考我的另外一篇文章[《Python两种方式加载文件内容》][Python 1]。加载文件内容代码如下：

def load_data(path, feature_num=2):
        '''
        导入数据
        :param path: 路径
        :param feature_num: 特行总数
        :return: 返回数据特征
        '''
        f = open(path)
        data = []
        for line in f.readlines():
            lines = line.strip().split("\t")
            data_tmp = []
            if len(lines) != feature_num:  # 判断特征的个数是否正确，把不符合特征数的数据去除
                continue
            for i in range(feature_num):
                data_tmp.append(float(lines[i]))
            data.append(data_tmp)
        f.close()
        return data

3、保存聚类结果

通过Mean Shift聚类后，我们使用一下方法进行聚类结果的保存

def save_result(file_name, data):
        '''
        保存聚类结果
        :param file_name: 保存的文件名
        :param data: 需要保存的文件
        :return:
        '''
        f = open(file_name, "w")
        m, n = np.shape(data)
        for i in range(m):
            tmp = []
            for j in range(n):
                tmp.append(str(data[i, j]))
            f.write("\t".join(tmp) + "\n")
        f.close()

4、调用Mean Shift算法

if __name__ == "__main__":
        data = load_data("F://data", 2)
        points, shift_points, cluster = mean_shift(data, 2)
        save_result("sub", np.mat(cluster))
        save_result("center", shift_points)

5、结果展示

得到的聚类结果如下所示：

![watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2dka3l4eTIwMTM_size_16_color_FFFFFF_t_70][]

你们在此过程中遇到了什么问题，欢迎留言，让我看看你们都遇到了哪些问题。

[Python]: https://blog.csdn.net/gdkyxy2013/column/info/22367
[20190713121057956.png]: https://img-blog.csdnimg.cn/20190713121057956.png
[20190713121156601.png]: /images/20210726/62aff76a14ea4942856313e922c0afb5.png
[20190713121445955.png]: /images/20210726/b3047420c00741a1ae3ca946af37723d.png
[20190713121635531.png]: https://img-blog.csdnimg.cn/20190713121635531.png
[20190713121709361.png]: /images/20210726/836451c6d4ca47e3b749413d6d29f35f.png
[20190713121832294.png]: /images/20210726/c57dca5c3ad6408b9fcc838a8fcbd897.png
[20190713122316225.png]: /images/20210726/5663ccd091c54b34bfff18ddbc1b89de.png
[20190713124429601.png]: https://img-blog.csdnimg.cn/20190713124429601.png
[20190713124541100.png]: https://img-blog.csdnimg.cn/20190713124541100.png
[20190713124616313.png]: https://img-blog.csdnimg.cn/20190713124616313.png
[20190713124712403.png]: https://img-blog.csdnimg.cn/20190713124712403.png
[20190713124754156.png]: /images/20210726/d08c051c037f45b0870af765d8c887aa.png
[20190713125150741.png]: /images/20210726/ecfbe362274f4f0eb3e2331a9c9bcc55.png
[20190713145124459.png]: /images/20210726/ce81f2f06aa3450daff6a51bd4a10a83.png
[Python 1]: https://blog.csdn.net/gdkyxy2013/article/details/81979214
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2dka3l4eTIwMTM_size_16_color_FFFFFF_t_70]: /images/20210726/a036e7052150480ab3fb246418ceecd5.png