DSMをJuliaでクラスタリングし、自動で整列させる

前回の記事の続きです。

otepipi.hatenablog.com

前回はJuliaでDSMのクラスタリングのみしましたが、今回はそのクラスタリング結果に基づいてDSMを並べ替えるところまでやります。

実装

実装したコードは以下のようになりました。

using LinearAlgebra
using Seaborn

const powcc = 1;
const powbid = 0;
const powdep = 1;
const max_cluster_size = 61;
const rand_accept = 30;
const rand_bid = 30;
const times = 5;
const stable_limit = 2;

original_DSM=[ 0 1 0 0 1 1 0; 0 0 0 1 0 0 1; 0 1 0 1 0 0 1; 0 1 1 0 1 0 1; 0 0 0 1 0 1 0; 1 0 0 0 1 0 0; 0 1 1 1 0 0 0];

original_label = ["A","B","C","D","E","F","G"];
(DSMsize,) = size(original_DSM);


Clustermatrix = zeros(DSMsize,DSMsize) + Diagonal(ones(DSMsize,DSMsize)) ;

DSM = original_DSM;


function ClusterBid(Clustermatrix,DSM,DSMsize,powdep,powbid,element)

    bestCluster = 0;
    bestBid = -1;

    for i = 1 : DSMsize

        Clusterlist= findall(x-> x == 1 , Clustermatrix[i,:]);
        (Clustersize,) = size(Clusterlist);
        bid = 0;

        for j = 1 : Clustersize

            bid += DSM[element,Clusterlist[j]];

        end

        bid = bid ^ powdep / Clustersize ^ powbid;

        if bid > bestBid

             bestBid = bid;
            bestCluster = i;
        end

    end


    return bestCluster

end

function TCC(DSM,DSMsize,Clustermatrix,powcc)

    intraCost=0;
    extraCost=0;


    for i = 1: DSMsize
    
        for j = 1 : DSMsize

            Clusterofi = findall(x-> x == 1 , Clustermatrix[:,i]);
            Clusterofj = findall(x-> x == 1 , Clustermatrix[:,j]);
            (Clustersizeofi,) = size(Clustermatrix[Clusterofi,:]);

            cost = DSM[i,j] + DSM[j,i];

            if Clusterofi == Clusterofj

                intraCost += cost*Clustersizeofi^powcc;

            else

                extraCost += cost*DSMsize^powcc;
            end


        end
    
    end

    totalCost = intraCost + extraCost;

    return totalCost
end


TCCost = TCC(DSM,DSMsize,Clustermatrix,powcc)

costlist = zeros(DSMsize*times,1);

for i = 1 : DSMsize * times;

    preClustermatrix = Clustermatrix;
    element = rand(1:DSMsize);
    bestCluster = ClusterBid(Clustermatrix,DSM,DSMsize,powdep,powbid,element);
    Clustermatrix[:,element] = zeros(DSMsize,1);
    Clustermatrix[bestCluster,element] = 1;
    newCost = TCC(DSM,DSMsize,Clustermatrix,powcc);
    if newCost <= maximum([TCCost,rand_accept])
        global TCCost = newCost;
    else
       global Clustermatrix = preClustermatrix;
    end

    global costlist[i] = newCost
end

Clustermatrix
costlist'

orderCluster = [sum(Clustermatrix,dims=2) Clustermatrix]
orderedCluster = sortslices(orderCluster,dims = 1, rev=true)
reorderedCluster = orderedCluster[:,2:DSMsize+1]

Order = [];

for i = 1 : DSMsize
    global Order = [Order ; findall(x-> x == 1 , reorderedCluster[i,:])]
end

I =zeros(DSMsize,DSMsize) + Diagonal(ones(DSMsize,DSMsize)) #単位行列の作成
P =zeros(DSMsize,DSMsize) + Diagonal(ones(DSMsize,DSMsize)) #置換行列の作成 
Clusteredlabel = copy(original_label) #ラベルも並べ替える
copylabel= copy(original_label)

for i = 1 : DSMsize
    global P[i,:] = I[Order[i],:]
    global Clusteredlabel[i] = copylabel[Order[i]]
end

ClusteredDSM = P*original_DSM*P'    

heatmap(original_DSM,annot=true, cbar=false,square=true,cmap="Blues",linewidths="0.5",linecolor="grey",xticklabels = original_label, yticklabels = original_label )
plt.savefig("original_DSM.png",bbox_inches="tight");

heatmap(ClusteredDSM,annot=true, cbar=false,square=true,cmap="Blues",linewidths="0.5",linecolor="grey",xticklabels = Clusteredlabel, yticklabels = Clusteredlabel )
plt.savefig("ClusteredDSM.png",bbox_inches="tight");