基于MATLAB的K-means算法实现及评价指标计算代码

基于MATLAB的K-means算法实现及评价指标计算代码

---

一、完整代码实现

function kmeans_with_metrics()
    %% 数据生成与预处理
    rng(42); % 固定随机种子
    data = [mvnrnd([1,2], eye(2), 100); 
            mvnrnd([5,6], eye(2), 100); 
            mvnrnd([9,4], eye(2), 100)]; % 生成3类数据
    
    % 数据标准化
    data = zscore(data);
    
    %% 参数设置
    maxK = 5; % 最大簇数
    numRuns = 10; % 每个K值运行次数
    distance = 'sqeuclidean'; % 距离度量
    
    %% 运行K-means并计算指标
    silhouetteScores = zeros(1, maxK-1);
    chScores = zeros(1, maxK-1);
    dbScores = zeros(1, maxK-1);
    sse = zeros(1, maxK);
    
    for k = 2:maxK
        [bestIdx, bestC, bestSSE] = run_kmeans(data, k, numRuns, distance);
        silhouetteScores(k-1) = mean(silhouette(data, bestIdx));
        chScores(k-1) = calinski_harabasz(data, bestIdx);
        dbScores(k-1) = davies_bouldin(data, bestIdx);
        sse(k) = bestSSE;
    end
    
    %% 可视化结果
    figure;
    subplot(2,2,1);
    gscatter(data(:,1), data(:,2), bestIdx);
    hold on; plot(bestC(:,1), bestC(:,2),'kx','MarkerSize',15); 
    title('聚类结果可视化');
    
    subplot(2,2,2);
    plot(2:maxK, silhouetteScores,'bo-','LineWidth',2);
    hold on; plot(2:maxK, sse(2:end)/max(sse),'r--');
    title('轮廓系数与SSE对比');
    legend('轮廓系数','SSE');
    
    subplot(2,2,3);
    bar([mean(silhouetteScores), mean(chScores), mean(dbScores)]);
    set(gca,'XTickLabel',{'轮廓系数','CH指数','DB指数'});
    ylabel('平均得分');
    
    %% 输出最佳K值
    [~, bestK] = max(silhouetteScores);
    fprintf('推荐最佳簇数: %d
', bestK+1);
end

%% K-means核心算法
function [idx, centroids, sse] = run_kmeans(data, k, numRuns, distance)
    [n, ~] = size(data);
    bestSSE = inf;
    bestIdx = [];
    bestC = [];
    
    for run = 1:numRuns
        % 初始化质心
        centroids = data(randperm(n,k), :);
        
        % 迭代优化
        prevC = centroids;
        for iter = 1:100
            % 分配样本
            distances = pdist2(data, centroids, distance);
            [~, idx] = min(distances, [], 2);
            
            % 更新质心
            for i = 1:k
                centroids(i,:) = mean(data(idx==i,:), 1);
            end
            
            % 收敛判断
            if norm(centroids - prevC) < 1e-5
                break;
            end
            prevC = centroids;
        end
        
        % 计算SSE
        sse_current = sum(sum((data - centroids(idx,:)).^2));
        if sse_current < bestSSE
            bestSSE = sse_current;
            bestIdx = idx;
            bestC = centroids;
        end
    end
end

%% 评价指标计算
function score = silhouette(data, idx)
    n = size(data,1);
    a = zeros(n,1);
    b = zeros(n,1);
    
    for i = 1:n
        cluster = idx(i);
        sameCluster = data(idx==cluster,:);
        diffCluster = data(idx~=cluster,:);
        
        a(i) = mean(pdist2(data(i,:), sameCluster, 'euclidean'));
        if ~isempty(diffCluster)
            b(i) = min(mean(pdist2(data(i,:), diffCluster, 'euclidean')));
        else
            b(i) = Inf;
        end
    end
    
    score = mean((b - a)./max(a, b));
end

function score = calinski_harabasz(data, idx)
    k = max(idx);
    n = size(data,1);
    ssb = 0;
    ssw = 0;
    
    for i = 1:k
        cluster = data(idx==i,:);
        mu = mean(cluster);
        ssb = ssb + size(cluster,1)*(sum((mu - mean(data)).^2));
        ssw = ssw + sum(sum((cluster - mu).^2));
    end
    
    score = (ssb/(k-1))/(ssw/(n-k));
end

function score = davies_bouldin(data, idx)
    k = max(idx);
    n = size(data,1);
    distances = pdist2(data, data, 'euclidean');
    
    maxRatio = zeros(k,1);
    
    for i = 1:k
        cluster = data(idx==i,:);
        mu = mean(cluster);
        for j = 1:k
            if i ~= j
                otherCluster = data(idx==j,:);
                ratio = (mean(pdist2(cluster, mu, 'euclidean')) + ...
                        mean(pdist2(otherCluster, mu, 'euclidean'))) / ...
                        mean(pdist2(cluster, otherCluster, 'euclidean'));
                maxRatio(i) = max(maxRatio(i), ratio);
            end
        end
    end
    
    score = mean(maxRatio);
end

---

二、关键功能说明

1. 算法实现

• K-means核心：包含初始化、分配、更新三阶段，支持多次运行避免局部最优
• 距离度量：支持欧氏距离、曼哈顿距离等多种选项
• 收敛条件：质心移动小于1e-5或达到100次迭代

2. 评价指标

指标名称	计算公式	代码函数
轮廓系数	$s(i)=max(a(i),b(i))b(i)−a(i)$	silhouette
Calinski-Harabasz指数	$CH=SSW/(n−k)SSB/(k−1)$	calinski_harabasz
Davies-Bouldin指数	$DB=k1∑i=1kmaxj\\=id(i,j)Si+Sj$	davies_bouldin

3. 可视化模块

• 聚类结果：显示数据分布与质心位置
• 指标对比：绘制轮廓系数与SSE曲线
• 指标评分：柱状图展示各指标平均值

---

三、使用示例

% 生成测试数据
data = [mvnrnd([1,2], eye(2), 100); 
        mvnrnd([5,6], eye(2), 100); 
        mvnrnd([9,4], eye(2), 100)];

% 运行完整流程
kmeans_with_metrics();

% 自定义参数运行
kmeans_with_metrics('maxK', 6, 'numRuns', 20);

参考代码 用matlab 实现了kmeans算法还附有评价指标计算 www.youwenfan.com/contentcsp/98288.html

四、应用场景

1. 客户分群：分析消费行为特征
2. 图像分割：处理RGB图像颜色聚类
3. 异常检测：识别离群点分布特征
4. 生物信息学：基因表达谱聚类分析

---

五、扩展功能

1. 动态聚类：添加实时数据更新功能
2. 三维可视化：支持三维数据降维展示
3. 交互式分析：通过GUI调整参数实时查看结果
4. 多指标融合：结合轮廓系数与CH指数进行综合评估