matlab練習程序（馬爾可夫聚類MCL）

本文主要參考：

https://wenku.baidu.com/view/b7907665caaedd3383c4d31b.html

https://blog.csdn.net/u010376788/article/details/50187321

雖然參考上面兩篇文章，不過這里和他給出的算法步驟不完全一致。

因為上面文章是針對Graphs的，矩陣中有邊為0的數據，而我的數據是Vector的，邊是根據距離計算出來的，應該還是有點區別的，而且我還用了kmeans。

首先給出他的算法步驟：

　　1.輸入一個無向圖，Expansion的冪e和Inflation的參數r。

　　2.創建鄰接矩陣。

　　3.對每個結點添加自循環（可選的）。

　　4.標准化矩陣（每個元素除以所在列的所有元素之和）。

　　5.計算矩陣的第e次冪。

　　6.用參數r對求得的矩陣進行Inflation處理。

　　7.重復第5步和第6步，直到狀態穩定不變（收斂）。

　　8.把最終結果矩陣轉換成聚簇。

然后是這里的算法步驟：

　　1.輸入數據，創建鄰接矩陣。

　　2.標准化矩陣（每個元素除以所在列的所有元素之和）。

　　3.對標准化后的矩陣進行馬爾可夫狀態轉移。

　　4.重復第2步和第3步，直到狀態穩定不變（收斂）。

　　5.使用kmeans把最終結果矩陣轉換成聚簇。（不明白kmeans的可以看看這篇文章）

代碼如下：

clear all;
close all;
clc;

K=3;
theta=0:0.01:2*pi;
p1=[3*cos(theta) + rand(1,length(theta))/2;3*sin(theta)+ rand(1,length(theta))/2];
p2=[2*cos(theta) + rand(1,length(theta))/2;2*sin(theta)+ rand(1,length(theta))/2];
p3=[cos(theta) + rand(1,length(theta))/2;sin(theta)+ rand(1,length(theta))/2];
p=[p1 p2 p3]';

randIndex = randperm(length(p))';       %打亂數據順序
p=p(randIndex,:);
plot(p(:,1),p(:,2),'.')

for i = 1:length(p)
    for j =1:length(p)
        W(i,j) = sqrt(sum((p(i,:)-p(j,:)).^2));　　%根據距離初始化無向圖的邊
    end
end

preW=W;
while 1
    x=repmat(sum(W),length(p),1);
    W=W./x;
    W=W*W;                      %馬爾科夫狀態轉移
    
    if sum(sum(preW-W))<1e-15
        break;
    end
    
    preW=W;
end

[idx,ctrs] = kmeans(W(:,1),K);          %用kmeans將收斂矩陣轉換為聚簇
figure;
plot(p(idx==1,1),p(idx==1,2),'r.') 
hold on;
plot(p(idx==2,1),p(idx==2,2),'g.')    
plot(p(idx==3,1),p(idx==3,2),'b.') 

原始數據：

聚類后：