簡介這份資源是面向本科及以上學習者、科研人員與工程實踐者的MATLAB聚類分析工具包圍繞k-means算法解決數(shù)據(jù)分組與模式識別問題適合課程作業(yè)、論文實驗及項目原型開發(fā)等場景。壓縮包共10個文件約378KB包含2個.m主程序腳本、1個.xls與1個.xlsx數(shù)據(jù)表格以及6張jpg結果示意圖代碼完整且附有注釋數(shù)據(jù)齊全便于直接運行與后續(xù)擴展。資源已有232人學習下載說明其在教學與實踐中具有一定參考價值。讀者可獲得一套可復用的k-means實現(xiàn)流程包括數(shù)據(jù)讀取、聚類計算與結果可視化配合示例數(shù)據(jù)與運行截圖能快速理解算法參數(shù)設置與效果評估并在此基礎上修改距離度量、聚類數(shù)目或替換數(shù)據(jù)集完成創(chuàng)新性擴展。若運行中遇到疑問也可聯(lián)系作者獲取支持。1. 從一份能跑的 kmeans 聚類 MATLAB 代碼說起數(shù)據(jù)齊全到底意味著什么很多人第一次接觸聚類是在一堆沒有標簽的樣本面前發(fā)懵分類模型用不了因為沒人告訴你哪條數(shù)據(jù)屬于哪一類。kmeans 聚類分析就是干這個的——它不需要標簽只靠樣本之間的距離把相似的東西自動歸堆。MATLAB 做這件事有天然優(yōu)勢矩陣運算快、可視化順手、統(tǒng)計與機器學習工具箱里現(xiàn)成的函數(shù)拿來就能用。但真正卡住新手的往往不是算法本身而是「代碼完整、數(shù)據(jù)齊全」這六個字背后的東西數(shù)據(jù)長什么樣、維度怎么統(tǒng)一、聚類數(shù) K 怎么定、結果怎么驗證。這篇筆記就圍繞一份可直接復現(xiàn)的 kmeans 聚類 MATLAB 代碼把數(shù)據(jù)準備、參數(shù)設置、結果評估和踩坑記錄一次講透適合剛上手聚類、手里有數(shù)據(jù)但不知道怎么落地的人。2. kmeans 在 MATLAB 里到底怎么算原理、選型與最小可跑代碼2.1 算法核心與 MATLAB 的實現(xiàn)路徑kmeans 的目標很樸素把 N 個樣本分成 K 個簇讓每個樣本到它所屬簇中心的距離平方和最小。數(shù)學上就是最小化簇內(nèi)平方誤差和這個目標函數(shù)沒有解析解所以用迭代逼近。標準流程是四步初始化 K 個中心、把每個樣本分配給最近的中心、重新計算每個簇的中心、重復分配和更新直到中心不再明顯移動或達到迭代上限。MATLAB 里實現(xiàn) kmeans 有兩條路。一條是統(tǒng)計與機器學習工具箱的kmeans函數(shù)一行調(diào)用就能出結果支持距離度量、重復次數(shù)、并行等參數(shù)另一條是自己寫循環(huán)適合理解算法細節(jié)或做教學演示。實際項目里我一般先用內(nèi)置函數(shù)跑通基線確認數(shù)據(jù)沒問題、K 值合理再考慮要不要手寫改造。內(nèi)置函數(shù)底層用的是 Lloyd 算法配合 kmeans 初始化能顯著降低陷入局部最優(yōu)的概率。選型上要注意如果你的數(shù)據(jù)維度很高比如上百維歐氏距離會失效這時候要么先降維PCA、t-SNE要么換余弦距離。如果簇的形狀不是球形kmeans 本身就不合適得考慮 DBSCAN 或譜聚類。這些邊界在動手前就要想清楚否則跑出來的結果看著有模有樣實際沒法用。2.2 數(shù)據(jù)準備從原始表格到聚類矩陣「數(shù)據(jù)齊全」不是指文件多而是指數(shù)據(jù)能直接喂進算法。kmeans 要求輸入是一個 N×D 的數(shù)值矩陣每行一個樣本每列一個特征。常見的數(shù)據(jù)問題有三類缺失值、量綱不統(tǒng)一、類別型字段沒編碼。缺失值處理上我一般先看缺失比例。低于 5% 的用列均值或中位數(shù)填補高于 20% 的考慮直接刪列或換特征。量綱問題更隱蔽比如一個特征是年齡0-100另一個是年收入0-1000000不標準化的話收入會完全主導距離計算年齡等于白給。標準做法是 z-score 標準化讓每個特征均值為 0、標準差為 1。下面是一段數(shù)據(jù)準備代碼假設原始數(shù)據(jù)存在 Excel 里包含數(shù)值列和幾個類別列% 讀取原始數(shù)據(jù)第一行為表頭 rawData readtable(raw_data.xlsx); % 查看前幾行確認列名和數(shù)據(jù)類型 head(rawData); % 分離數(shù)值特征和類別特征 numFeatures rawData(:, {Age, Income, Score, Spend}); catFeatures rawData(:, {Gender, City}); % 類別特征做獨熱編碼轉成數(shù)值 catEncoded onehotencode(catFeatures, 1:width(catFeatures)); % 合并成完整特征矩陣 featureMatrix [table2array(numFeatures), catEncoded]; % 缺失值用列中位數(shù)填補 featureMatrix fillmissing(featureMatrix, constant, ... median(featureMatrix, 1, omitnan)); % z-score 標準化每列減均值除標準差 featureMatrix zscore(featureMatrix); % 確認最終矩陣尺寸 fprintf(樣本數(shù): %d, 特征數(shù): %d\n, size(featureMatrix, 1), size(featureMatrix, 2));這段代碼的邏輯是先讀表把數(shù)值列和類別列分開處理類別列用獨熱編碼變成 0/1 向量再拼回一個大矩陣。fillmissing用列中位數(shù)填補比均值更抗異常值。zscore是標準化關鍵少了這一步后面聚類結果基本不可信。參數(shù)上onehotencode的第二個參數(shù)指定對哪些列編碼fillmissing的constant配合中位數(shù)是常見組合。跑完看輸出尺寸如果特征數(shù)和你預期對不上多半是獨熱編碼把某一列拆成了多列。2.3 最小可跑的 kmeans 調(diào)用與參數(shù)含義數(shù)據(jù)準備好之后核心調(diào)用就一行。但這一行里的參數(shù)決定了結果好壞不能隨便填% 設定聚類數(shù) K先用肘部法粗定一個范圍 K 4; % 調(diào)用 kmeans關鍵參數(shù)逐個說明 [idx, C, sumd, D] kmeans(featureMatrix, K, ... Distance, sqeuclidean, ... % 距離度量默認平方歐氏 Replicates, 10, ... % 重復 10 次取最優(yōu)降低局部最優(yōu)風險 Start, plus, ... % kmeans 初始化 MaxIter, 500, ... % 單次迭代上限 Display, final); % 只輸出最終結果避免刷屏 % idx 是每個樣本的簇編號C 是 K 個簇中心sumd 是簇內(nèi)距離和 fprintf(各簇樣本數(shù): ); disp(histcounts(idx, 1:K1));idx是 N×1 的簇標簽C是 K×D 的中心矩陣sumd是每個簇內(nèi)樣本到中心的距離平方和D是每個樣本到所有中心的距離。Replicates設 10 是經(jīng)驗值數(shù)據(jù)量大或 K 大時可以加到 20代價是時間線性增長。Start用plus就是 kmeans比默認的均勻采樣穩(wěn)。MaxIter一般 300 到 500 夠用設太小可能沒收斂就停了。跑完用histcounts看各簇樣本數(shù)如果某一簇只有個位數(shù)樣本要么是 K 設大了要么是數(shù)據(jù)里有離群點。3. 聚類數(shù) K 怎么定肘部法、輪廓系數(shù)與業(yè)務約束的三方博弈3.1 肘部法的計算與讀圖K 是 kmeans 唯一需要人為指定的關鍵參數(shù)也是最容易拍腦袋的地方。肘部法的思路是隨著 K 增大簇內(nèi)距離和必然下降但下降速度會在某個點明顯變緩那個拐點就是候選 K。實現(xiàn)上就是循環(huán)跑不同 K記錄sumd總和% 測試 K 從 1 到 10 的簇內(nèi)距離和 K_range 1:10; wss zeros(length(K_range), 1); for i 1:length(K_range) [~, ~, sumd] kmeans(featureMatrix, K_range(i), ... Replicates, 5, Start, plus, Display, off); wss(i) sum(sumd); end % 畫肘部圖 figure; plot(K_range, wss, -o, LineWidth, 1.5); xlabel(聚類數(shù) K); ylabel(簇內(nèi)距離和); title(肘部法確定 K); grid on; % 計算相鄰點的下降率輔助判斷拐點 dropRate -diff(wss) ./ wss(1:end-1); disp(table(K_range(2:end), dropRate, VariableNames, {K, DropRate}));wss是 within-cluster sum of squares隨 K 單調(diào)下降?？磮D時找下降率突然變小的位置比如從 K3 到 4 降了 30%從 4 到 5 只降了 8%那 4 就是候選。代碼里額外算了dropRate比肉眼讀圖更客觀。注意Replicates這里設 5 就夠因為只是比較趨勢不需要每個 K 都跑到最優(yōu)。3.2 輪廓系數(shù)比肘部法更硬的指標肘部法主觀性強輪廓系數(shù)silhouette能給出每個樣本的聚類質(zhì)量分數(shù)范圍 -1 到 1越接近 1 說明樣本離本簇近、離其他簇遠。MATLAB 里silhouette函數(shù)直接算% 對候選 K 計算平均輪廓系數(shù) K_candidates 2:8; silScores zeros(length(K_candidates), 1); for i 1:length(K_candidates) idx kmeans(featureMatrix, K_candidates(i), ... Replicates, 10, Start, plus, Display, off); silScores(i) mean(silhouette(featureMatrix, idx)); end % 輸出對比表 disp(table(K_candidates, silScores, VariableNames, {K, Silhouette})); % 找最高分對應的 K [bestScore, bestIdx] max(silScores); fprintf(最佳 K %d, 輪廓系數(shù) %.4f\n, K_candidates(bestIdx), bestScore);輪廓系數(shù)對距離度量敏感如果前面沒做標準化這里分數(shù)會普遍偏低且不可比。一般平均輪廓系數(shù)高于 0.5 算結構清晰0.3 到 0.5 算可接受低于 0.25 就要懷疑數(shù)據(jù)本身沒有明顯簇結構。注意輪廓系數(shù)在 K2 時往往偏高這是它的已知偏向所以不能只看分數(shù)要結合肘部法和業(yè)務含義。3.3 業(yè)務約束下的 K 選擇純數(shù)學指標給的是候選最終定 K 還要看業(yè)務能不能用。比如做用戶分群分成 3 群和 5 群對應的運營策略完全不同5 群可能細到?jīng)]法針對性投放。我一般會做一張對照表把不同 K 下的簇大小、中心特征、輪廓系數(shù)列出來和業(yè)務方一起過一遍。如果某個 K 下出現(xiàn)一個超大簇加幾個極小簇通常說明 K 偏大或者數(shù)據(jù)里有離群點沒處理。這一步?jīng)]有代碼能替代但前面算出的C和idx就是討論的素材。4. 結果可視化與簇特征解讀讓聚類結果能講出人話4.1 二維和三維散點圖的畫法聚類結果如果只給一堆標簽沒人看得懂?？梢暬亲尳Y果落地的關鍵一步。高維數(shù)據(jù)沒法直接畫常規(guī)做法是先用 PCA 降到 2 維或 3 維再按簇標簽上色% PCA 降到二維用于可視化 [coeff, score, ~, ~, explained] pca(featureMatrix); score2d score(:, 1:2); % 按簇標簽畫散點圖 figure; gscatter(score2d(:,1), score2d(:,2), idx, lines(K), ., 12); xlabel(sprintf(PC1 (%.1f%%), explained(1))); ylabel(sprintf(PC2 (%.1f%%), explained(2))); title(kmeans 聚類結果PCA 二維投影); grid on; % 疊加簇中心在 PCA 空間的投影 hold on; center2d (C - mean(featureMatrix)) * coeff(:, 1:2); plot(center2d(:,1), center2d(:,2), kx, MarkerSize, 14, LineWidth, 2); hold off;gscatter按idx分組上色比手動循環(huán)scatter省事。explained告訴你前兩個主成分解釋了多少方差如果加起來不到 50%說明二維投影丟失信息太多圖只能當參考不能下結論。中心點投影那一步是把原始空間的C通過同樣的均值和coeff變換到 PCA 空間這樣中心點和樣本點在同一坐標系里方便看簇的緊致程度。4.2 簇中心反標準化與特征畫像PCA 圖看的是整體分布要解釋每個簇是什么得回到原始特征空間看中心。但前面做了 z-score中心值是標準化的得反變換回去% 保存標準化參數(shù)在 zscore 那一步之后 mu mean(featureMatrix_raw, 1); sigma std(featureMatrix_raw, 0, 1); % 反標準化簇中心 C_original C .* sigma mu; % 把中心轉成表格方便對照列名 centerTable array2table(C_original, ... VariableNames, featureMatrix_colnames); disp(centerTable); % 對每個簇找出中心值最高的三個特征 for k 1:K [~, topIdx] maxk(C_original(k,:), 3); fprintf(簇 %d 主導特征: %s\n, k, strjoin(featureMatrix_colnames(topIdx), , )); end這里的關鍵是mu和sigma必須在標準化之前從原始矩陣算出來并保存否則反變換對不上。C_original的每一行是一個簇在原始量綱下的中心比如年齡 35、收入 8000 這樣業(yè)務方一看就懂。maxk找每個簇最突出的特征快速生成畫像描述。如果某個簇在多個特征上都偏高說明這個簇的特征不單一可能需要拆或者合并。4.3 用輪廓圖定位問題樣本整體輪廓系數(shù)是平均值掩蓋了個體差異。輪廓圖能把每個樣本的分數(shù)畫出來一眼看出哪些樣本分錯了figure; [silVals, ~] silhouette(featureMatrix, idx, sqeuclidean); title(各樣本輪廓系數(shù)); % 找出輪廓系數(shù)為負的樣本這些是可能分錯的 negIdx find(silVals 0); fprintf(輪廓系數(shù)為負的樣本數(shù): %d (占比 %.1f%%)\n, ... length(negIdx), 100*length(negIdx)/length(silVals)); % 輸出這些樣本的原始索引和當前簇標簽 if ~isempty(negIdx) disp(table(negIdx, idx(negIdx), silVals(negIdx), ... VariableNames, {SampleIndex, Cluster, Silhouette})); end輪廓系數(shù)為負意味著樣本到其他簇的平均距離比到本簇還近基本可以判定分錯了。占比低于 5% 可以接受高于 10% 就要回頭檢查 K 是否合理、特征是否夠區(qū)分。這些負分樣本往往是邊界情況業(yè)務上可能正好是需要單獨關注的那批人。5. 避坑與排查kmeans 聚類 MATLAB 實現(xiàn)里最容易翻車的五件事5.1 沒標準化導致某列特征獨大現(xiàn)象聚類結果里某一簇的樣本在某個特征上高度一致其他特征完全隨機看起來像按單一維度分的。原因不同特征量綱差異大距離計算被大量綱特征主導。解決聚類前對所有數(shù)值特征做 z-score 或 min-max 標準化并在反標準化解讀中心時用對應的均值和標準差還原。判斷方法很簡單看簇中心表里各特征的數(shù)值范圍如果某一列數(shù)值比其他列大幾個數(shù)量級基本就是這個問題。5.2 K 值拍腦袋定結果沒法解釋現(xiàn)象跑出來的簇大小嚴重不均或者業(yè)務方問「為什么是 4 類不是 3 類」時答不上來。原因只跑了一次 kmeans沒做 K 的掃描和對比。解決至少用肘部法和輪廓系數(shù)各掃一遍 K 的范圍把不同 K 下的簇大小、輪廓系數(shù)、中心特征列成表選數(shù)學指標和業(yè)務含義都說得通的那個。我一般會把 2 到 8 的結果都留著業(yè)務討論時隨時調(diào)出來看。5.3 忽略 Replicates 導致結果每次不一樣現(xiàn)象同樣的數(shù)據(jù)和 K兩次運行得到的簇標簽和中心不同。原因kmeans 對初始中心敏感單次運行容易陷入局部最優(yōu)。解決Replicates設 10 以上Start用plus。如果數(shù)據(jù)量特別大導致重復太慢可以先用sample初始化跑一次看大概再用plus加Replicates精跑。另外注意即使這樣簇的編號也可能不同比較兩次結果時要看中心而不是看標簽數(shù)字。5.4 缺失值沒處理直接進 kmeans現(xiàn)象代碼報錯NaN相關或者結果里某些樣本的簇標簽異常。原因kmeans不接受含 NaN 的輸入矩陣。解決進 kmeans 之前必須fillmissing或rmmissing。填補方法上數(shù)值列用中位數(shù)比均值穩(wěn)類別列用眾數(shù)。如果某列缺失超過 30%我傾向于直接刪掉這列因為填補引入的偏差可能比丟掉這列更大。5.5 把聚類結果當分類標簽用現(xiàn)象拿 kmeans 的idx去訓練一個分類器然后在新數(shù)據(jù)上預測發(fā)現(xiàn)效果很差。原因kmeans 給出的簇編號沒有跨數(shù)據(jù)集的一致性新數(shù)據(jù)跑一遍 kmeans 得到的編號和舊數(shù)據(jù)對不上。解決如果要做預測應該用聚類中心訓練一個分類器比如最近鄰或 SVM把中心作為「偽標簽」的來源而不是直接用編號。或者用knnsearch把新樣本分配到最近的已有中心。這個坑很隱蔽因為在自己數(shù)據(jù)集上驗證時看著沒問題一上生產(chǎn)就露餡。6. 從能跑到好用kmeans 結果穩(wěn)定性驗證與增量分配的一個實用技巧代碼能跑通只是起點真正投入使用前我會做一件事驗證聚類結果的穩(wěn)定性。方法不復雜把數(shù)據(jù)隨機分成兩半各自跑 kmeans然后比較兩半得到的簇中心是否接近。如果中心差異很大說明數(shù)據(jù)本身沒有穩(wěn)定結構或者 K 選得不對。MATLAB 里可以用pdist2算兩組中心的距離矩陣看最小距離是否在可接受范圍內(nèi)。% 隨機對半切分 n size(featureMatrix, 1); halfIdx randperm(n, floor(n/2)); dataA featureMatrix(halfIdx, :); dataB featureMatrix(setdiff(1:n, halfIdx), :); % 各自跑 kmeans K 4; [~, CA] kmeans(dataA, K, Replicates, 10, Start, plus, Display, off); [~, CB] kmeans(dataB, K, Replicates, 10, Start, plus, Display, off); % 計算兩組中心的兩兩距離 centerDist pdist2(CA, CB); minDist min(centerDist, [], 2); fprintf(各中心到另一組最近中心的距離: ); disp(minDist); % 如果每個中心都能在另一組找到距離小于閾值的對應中心認為穩(wěn)定 threshold 0.5; % 標準化空間下的經(jīng)驗閾值 stable all(minDist threshold); fprintf(聚類結果穩(wěn)定: %s\n, string(stable));這段代碼的核心是pdist2算兩組中心的距離矩陣minDist是每個 A 組中心到 B 組最近中心的距離。閾值 0.5 是在標準化空間下的經(jīng)驗值因為標準化后特征標準差為 1中心距離小于半個標準差算接近。如果某個中心的最小距離超過 1說明這一簇在兩半數(shù)據(jù)里位置差異大要么是樣本太少不穩(wěn)定要么是 K 偏大。另一個實用技巧是增量分配當有新樣本進來時不需要重新跑整個 kmeans直接用已有的中心做最近鄰分配。knnsearch或者手動算距離都行% 新樣本標準化用訓練時的 mu 和 sigma newSample (newRaw - mu) ./ sigma; % 分配到最近的中心 [d, assignedCluster] min(pdist2(newSample, C), [], 2); fprintf(新樣本分配到簇 %d距離 %.4f\n, assignedCluster, d);這樣做的前提是新樣本的分布和訓練數(shù)據(jù)一致如果業(yè)務發(fā)生突變比如用戶行為模式整體偏移增量分配會失效這時候需要重新聚類。我一般會監(jiān)控新樣本到最近中心的平均距離如果持續(xù)上升就是重新訓練的觸發(fā)信號。這套流程跑下來從數(shù)據(jù)準備到結果驗證大概兩三百行代碼覆蓋了 kmeans 聚類分析在 MATLAB 里落地的完整鏈路。我自己踩過最深的坑是早期不做標準化直接跑結果對著簇中心表看了半天沒看出規(guī)律后來才發(fā)現(xiàn)是收入那一列把距離全吃掉了。希望幫到你。本文還有配套的精品資源點擊獲取