带约束的K-means聚类算法

带约束的K-means 聚类算法

1. 前言

上一期学习了K-means聚类算法，聚类是不受到限制的，单纯的无监督学习，但是当存在一些约束时，比如对每一簇的聚类样本点数量有限制，或者每个样本点带需求，每一簇存在满足约束的最大能力，这时候应该怎么添加约束，并且不破坏聚类的效果。

2. 场景构建

假设场景：有一定数量的区域需要分配邮递员配送信件，怎么将区域划分给邮递员？每个区域有一定的需求量，一个邮递员一天的工作量只能满足一定的需求，那么应该怎么将区域划分给邮递员？考虑到成本问题，这里设置一个目标——在满足需求的情况下最小化邮递员数量。

这样的场景在快递行业应该很常见，不太清楚实际业务中是怎么分配的，但我觉得带约束的聚类算法应该能求解一个可行解，为了给K-means加约束构想的场景，场景简化了很多，可能不太符合实际，欢迎提建议。

3. 带约束的K-means 聚类设计

K-means聚类算法流程

①算法说明：由于存在需求约束，那么就存在最小的K（总需求/簇的能力），可以从最小值K开始遍历。而约束就加在（3）将对象分配到最相似的簇中，具体修改为：在满足约束的情况下，将对象分配到最相似的簇中，如果最相似的簇无法容纳，就找第二相似的簇，依次类推。

②遍历改进：对遍历对象进行一个小的修改，通常K-means是随机/按列表顺序遍历对象，将对象分配到最相似的簇，这里改进为：取所有对象中距离某个簇最近的点，优先分配（这算是随机/按列表顺序遍历对象的一个特殊情况）。

③初始解敏感：为了减小K-means对初始解敏感的问题，这里采用多次求解的方式保留最优的方式降低对初始解的依赖。

④效果评价：聚类效果使用轮廓系数进行评价。

3. python实现

# -*- coding: utf-8 -*-
import pandas as pd
import numpy as np
import random
import math
from sklearn.metrics import silhouette_score
import matplotlib.pyplot as plt
from matplotlib.pylab import mpl
mpl.rcParams['font.sans-serif'] = ['SimHei']  # 添加这条可以让图形显示中文


def calDistance(centers,nodes):
    ''' 
    计算节点间距离
    输入：centers-中心，nodes-节点；
    输出：距离矩阵-dis_matrix
    '''
    dis_matrix = pd.DataFrame(data=None,columns=range(len(centers)),index=range(len(nodes)))
    for i in range(len(nodes)):
        xi,yi = nodes[i][0],nodes[i][1]
        for j in range(len(centers)):
            xj,yj = centers[j][0],centers[j][1]
            dis_matrix.iloc[i,j] = round(math.sqrt((xi-xj)**2+(yi-yj)**2),2)
    
    return dis_matrix


def scatter_diagram(clusters,nodes):
    '''
    #画路径图
    输入：nodes-节点坐标；
    输出：散点图
    '''
    for cluster in clusters:
        x,y= [],[]
        for Coordinate in cluster:
            x.append(Coordinate[0])
            y.append(Coordinate[1])
        plt.scatter(x, y,alpha=0.8,)
    plt.xlabel('x')
    plt.ylabel('y')
    plt.show()


def distribute(center,nodes,K,demand,d_limit):
    '''
    将节节点分给最近的中心
    输入：center-中心,nodes-节点,K-类数量,demand-需求,d_limit-一个簇的满足需求的能力
    输出：新的簇-clusters，簇的需求-clusters_d
    '''
    clusters = [[] for i in range(K)]#簇
    label = [None for i in range(len(nodes))]
    clusters_d = [0 for i in range(K)]#簇需求
    
    dis_matrix = calDistance(center,nodes).astype('float64')
    
    for i in range(len(dis_matrix)):
        row, col = dis_matrix.stack().idxmin()
        node = nodes[row]
        j = 1
        while clusters_d[col] + demand[row] > d_limit:#检验约束
            if j < K: 
                j += 1
                dis_matrix.loc[row,col] = math.pow(10,10)
                col = dis_matrix.loc[row,:].idxmin()
            else:#所有类都不满足需求量约束
                print("K较小")
                scatter_diagram(clusters,nodes)
                return None
            
        #满足约束正常分配
        clusters[col].append(node)
        label[row] = col
        clusters_d[col] += demand[row]
        dis_matrix.loc[row,:] = math.pow(10,10)
        
    return clusters,clusters_d,label


def cal_center(clusters):
    '''
    计算簇的中心
    输入：clusters-类坐标；
    输出：簇中心-new_center
    '''
    new_center = []
    for cluster in clusters:
        x,y= [],[]
        for Coordinate in cluster:
            x.append(Coordinate[0])
            y.append(Coordinate[1])
        x_mean = np.mean(x)
        y_mean = np.mean(y)
        new_center.append((round(x_mean,2),round(y_mean,2)))
    
    return new_center
    

if __name__ == "__main__":
    #输入数据
    nodes = [(11, 36),(13, 35),(14, 12),(25, 22),(15, 11),(16, 24),(14, 30),(16, 32),(26, 10),(27, 31),(22, 29),(1, 29),(20, 29),(17, 15),(5, 38),
             (3, 15),(22, 20),(4, 13),(4, 22),(12, 18),(23, 12),(18, 15),(5, 13),(9, 27),(17, 36),(2, 14),(1, 16),(7, 15),(25, 15),(19, 26),(25, 40),
             (26, 34),(25, 35),(2, 36),(29, 24),(17, 17),(8, 26),(4, 14),(5, 25),(6, 37),(1, 14),(6, 39),(11, 13),(10, 20),(21, 11),(5, 19),(5, 35),(1, 34),
             (16, 39),(19, 24),(39, 31),(49, 31),(41, 50),(31, 33),(32, 40),(35, 30),(31, 39),(34, 48),(42, 32),(32, 35),(35, 33),(35, 34),(43, 41),(35, 47),
             (49, 36),(37, 41),(43, 46),(41, 41),(45, 50),(41, 35),(45, 44),(41, 30),(43, 33),(31, 45),(48, 32),(39, 49),(38, 42),(33, 39),(49, 33),(43, 44),
             (32, 30),(40, 47),(36, 46),(47, 47),(37, 33),(35, 31),(42, 38),(43, 47),(30, 47),(30, 30),(37, 34),(41, 45),(27, 33),(42, 39),(43, 43),(50, 43),
             (28, 40),(35, 41),(32, 41),(31, 30)]
    demand = [3,3,5,9,6,1,4,8,8,3,6,5,6,4,1,4,7,1,3,6,2,5,7,2,2,1,9,3,7,8,4,3,1,2,5,1,3,4,8,5,2,9,3,10,6,1,9,3,5,3,3,3,9,9,6,2,5,2,4,
              10,4,10,10,6,3,7,9,4,2,9,7,4,5,3,3,8,2,3,2,9,1,3,3,3,9,5,7,8,8,5,2,8,3,4,2,4,6,3,7,4]
    
    scatter_diagram([list(nodes)],nodes)#初始位置分布图
    
    #参数
    d_limit = 150 # 一个区的约束
    K = math.ceil(sum(demand)/d_limit)  #簇的最小边界
    
    #历史最优参数
    best_s = -1
    best_center = 0
    best_clusters = 0
    best_labels = 0
    best_clusters_d = 0
    
    for n in range(20):#多次遍历，kmeans对初始解敏感
        #初始化随机生成簇中心
        index = random.sample(list(range(len(nodes))),K)
        new_center = [nodes[i] for i in index]
        i = 1
        center_list = []
        while True:
        	#节点——> 簇
            center = new_center.copy()
            center_list.append(center)#保留中心，避免出现A生成B，B生成C，C有生成A的情况
            clusters,cluster_d,label = distribute(center,nodes,K,demand,d_limit)
            new_center = cal_center(clusters)
            if (center == new_center) | (new_center in center_list):
                break
            i += 1
            
        s = silhouette_score(nodes,label, metric='euclidean') # 计算轮廓系数，聚类的评价指标
        # scatter_diagram(clusters,nodes)#当前最优图
        
        if best_s < s:
            best_s = s
            best_clusters_d = cluster_d
            best_center = center.copy()
            best_clusters = clusters.copy()
            best_labels = label.copy()
    
    scatter_diagram(best_clusters,nodes)#历史最优图
    # print("各类需求量：",best_clusters_d)
    

【结果展示】

图1 需求点位置

图2 聚类效果

4 谈论

从结果上看，带约束的聚类效果还是可以的，但是求解过程中会出现图3的聚类效果，蓝色的点距离黄色的点反而更近一点。而且还存在另外一个问题，假设需求量是500，而每簇满足需求的能力是100，理论上最少需要5个簇就可以，但是由于约束的存在和需求不可拆分，可能会导致无法刚好满足五个簇的容量，遇到这种情况我写的算法中是直接报错，虽然这种情况下大概率是出现如图3蓝色点的情况存在，但是也会减小求解最优解的可能性，程序的健壮性较差。

图3 聚类效果