%% Unsupervised Learning
%
%% k-Means Clustering
%
% Read IRIS Data

[X, y] = iris_datamatrix();

%
%% Create 3 clusters
%

% Form three clusters
[idx3, centroids3] = kmeans (X, 3);

%
% A silhouette plot displays a measure of how close each point in one cluster is to points in the neighboring clusters. This measure ranges from 1 (indicating points that are very distant from neighboring clusters) through 0 (points that are not distinctly in one cluster or another) to –1 (points that are probably assigned to the wrong cluster)
%

[silh3, h] = silhouette(X, idx3);
xlabel('Silhouette Value')
ylabel('Cluster')
fprintf("silhoutte mean score with 3 clusters is %4.3f%", mean(silh3));

fprintf("Purity of k-Means wwith 3 clusters is %5.3f%%", ...
        compute_purity(idx3, y)*100);

%
%% Create 2 clusters
%

[idx2, centroids2] = kmeans (X, 2);
[silh2, h] = silhouette(X, idx2);
xlabel('Silhouette Value')
ylabel('Cluster')
fprintf("silhoutte mean score with 2 clusters is %4.3f%", mean(silh2))

fprintf("Purity of k-Means with 2 clusters is %5.3f%%", ...
        compute_purity(idx2, y)*100);

%

hold off;

%

scatter(X(:, 1), X(:,3), 50, idx3, "filled")
hold on;

%

plot(centroids3(:,1), centroids3(:,3), 'ro', ...
            'MarkerSize', 10, 'LineWidth', 4 );

%
%% Elbow Analysis
%

hold off;

%

rng(0);
SSE = [];
for k = 1:10
    [~, ~, sumd] = kmeans (X, k);
    SSE = [SSE; sum(sumd)];
end
plot (1:10, SSE, '-o');
grid;
xlabel('Number of Clusters (K)');
ylabel('SSE');
title ('Elbow Plot');

%
%% Agglomerative Clustering
%

%dist_mat = pdist(X);       %default
dist_mat = pdist(X, "cityblock");
Z = linkage(dist_mat, 'average');
dendrogram(Z)

%
%% Cophenet Similarity
% You can use the cophenetic correlation coefficient to compare the results of clustering the same data set using different distance calculation methods or clustering algorithms. For example, you can use the cophenet function to evaluate the clusters created for the sample data set.
%

c = cophenet(Z, dist_mat);
fprintf("Cophenet Similarity is %4.3f", c);

%
%% Inconsistency Matrix
% Column
% Description
% 1
% Mean of the heights of all the links included in the calculation
% 2
% Standard deviation of all the links included in the calculation
% 3
% Number of links included in the calculation
% 4
% Inconsistency coefficient
%

I = inconsistent(Z);
fprintf("Average Inconsistency Value  %4.3f", mean(I(:,4)));

%
%% 3 Clusters
%

idx3 = cluster(Z, 'maxclust', 3);
fprintf("Purity of Heirarchical Method with 3 Clusters is %5.3f%%", ...
            compute_purity(idx3, y)*100);

%
%% 2 Clusters
%

idx2 = cluster(Z, 'maxclust', 2);
fprintf("Purity of Heirarchical Method with 2 Clusters is %5.3f%%", ...
            compute_purity(idx2, y)*100);

%
%% Elbow Analysis
%

SSE = [];
for k = 1:10
    %idx = cluster(Z, "cutoff", cut_off);
     idx = cluster(Z, "maxclust", k);
    centroids = zeros(k, size(X,2));
    for i = 1:k
        centroids(i, :) = mean(X(idx==i, :));
    end

    SSE_k = 0;
    for i =1:k
        pts = X(idx == i, :);
        SSE = SSE + sum(sum((pts - centroids(i, :)).^2));
    end

    SSE = [SSE; SSE_k];
end

plot (1:10, SSE, '-o');
grid;
xlabel('Number of Clusters');
ylabel('SSE');
title ('Elbow Plot');

%

idx = cluster(Z, "cutoff", 1);
fprintf("Purity of Heirarchical Method with 1.5 Cut-Off is %5.3f%%", ...
            compute_purity(idx, y)*100);

%
%% Functions
%

function purity = compute_purity(idx, y)
    unique_clusters = unique(idx);
    purity = 0;
    for i=1:length(unique_clusters);
        cluster_points = y(idx == unique_clusters(i));
        mode_label = mode(cluster_points);
        num_correct = sum(cluster_points == mode_label);
        purity = purity + num_correct;
    end
    purity = purity / length(y);
end
