Shandong University machine learning experiment VI k-means

Machine learning experiment of Shandong University 6 The report

Experimental hours ： 4 Date of experiment ：2021.11.29

Experimental topics ：Experiment 6 : K-Means

The experiment purpose

In this exercise, you will use K-means to compress an image by reducing the
number of colors it contains

Use K-means Compress the image by reducing the number of colors

Experimental environment

Software environment

​ Windows 10 + Matlab R2020a

​ CPU: Intel Core i5-8300H CPU @ 2.30GHz 2.30 GHz

Experimental steps and contents

understand K-Means

K-Means (Lloyd 1957)

1. Update for each point Class K: By calculating each data point $X_i$ To Center point of each category ${\mu }_k$ Distance of , To choose which category a point belongs to

   

2. Update each cluster means : That is, according to the current point of this category , Carry out the search mean operation , Come to the center

   

3. When ${\mu }_k$ perhaps Loss Constant time , Think the iteration is complete .

   Loss:

   

Loss Function

${\mu }_k$ Think it is the K The centroid of a cluster ,$z_{i,k}$ Think it's $x_i$ Whether to belong to $C_k$ An indication of , Then for each piece of data, there will be a $z_i$, So define a data $x_i$ Of $loss$ by ：

The total $LossFunction:$ among $XisN\times D$ then $ZisN\times KAnd\mu isK\times D$

It means X Express all D D data ,Z Express all N The cluster category of data indicates ,$\mu$ Express K A cluster of D Dimensional center point .K-Means Is to minimize this Loss Function

K-Means Objective

​ K-Means Is a heuristic method , To solve this NP-hard problem , And he is a Non-Convex problem, There is a lot of local minima

The algorithm is described as ：

Choosing K

​ Choosing the number of clusters is also a problem .

• ​ You can try different K value , Painting for different K Under the Loss chart , choice elbow point that will do .

• ​ You can also use AIC To solve the calculation

---

Task Description:

​ Your task is to compute 16 cluster centroids from this image, with each centroid being a vector of length three that holds a set of RGB values.

​ Calculation 16 Clusters , Every central point is a $size=3\times 1$ Of RGB $pixel$ Pixel point $value$ Express RGB Value

​ Whereas calculation $538\times 538\times 3$ The center point of the picture will take a lot of time , So first Little picture $(128\times 128\times 3)$ After training , Apply to large pictures , Get one that only uses 16 individual color It means the new picture of

Experimental steps :

K-Means Algorithm

1. Random initialization ： From the whole picture Randomly choose 16 individual pixel As the initial iteration center point . The random number scheme I adopted can realize no repeated sampling .

   function code:

   %%  Select the initial point randomly  return K * 3 
   function sample_Datas = sample_random(num,datas,dimx,dimy)
       % datas  For raw data  num  Is the target number  
       sample_Datas = zeros(num,3);
       a = rand(dimx,1);
       b = rand(dimy,1);
       [vala,posa] = sort(a);
       [valb,posb] = sort(b);
       chose_x = posa(1:num,1);
       chose_y = posb(1:num,1);
       for i=1:num
           sample_Datas(i,:) = datas(chose_x(i),chose_y(i),:);
       end
   end
   

2. Calculate each Pixel The nearest neighbor of ： Traverse every pixel , Every pixel There is one. $RGBvector$ also $size=3\times 1$ And then Find the distance of all center points , Select the center point corresponding to the minimum distance K As this pixel Of class

   

   function code:

   %%  Calculate each pixel  Of   Category  return Clusters: N * N val  by   Category  
   function Clusters = calculate_Ci(centroids,dimx,dimy,datas,K)
       Clusters = [];
       %  Through each pixel  Calculate a z(i,j)
       for i= 1:dimy
           for j = 1:dimx
               %  obtain xi
               pixel_rgb = [datas(i,j,1),datas(i,j,2),datas(i,j,3)];
               diff = ones(K,1)*pixel_rgb - centroids;
               distance = sum(diff.^2,2);
               [~,class] = min(distance); %  Get the smallest corresponding category index
               Clusters(i,j) = class;
           end
       end
   end
   

3. Update the center point $\mu$ adopt

​ Function code:

​

%%  Update the center point  return 16 * 3
function newcenters = updatemeans(Clusters,dimx,dimy,dimz,datas,K)
    newcenters = zeros(K,dimz);
    nums =  zeros(K);
    sumred = zeros(K,1);
    sumgreen = zeros(K,1);
    sumblue = zeros(K,1);
    for i=1:dimx
        for j=1:dimy
            class = Clusters(i,j);
            nums(class) = nums(class) + 1;
            sumred(class) = sumred(class) + datas(i,j,1);
            sumgreen(class) = sumgreen(class) + datas(i,j,2);
            sumblue(class) = sumblue(class) + datas(i,j,3);
        end
    end
    for i=1:K
        if nums(i) ~= 0
            sumred(i) = sumred(i) / nums(i);
            sumgreen(i) = sumgreen(i) / nums(i);
            sumblue(i) = sumblue(i) / nums(i);
        end
    end
    newcenters = [sumred,sumgreen,sumblue];
    newcenters = round(newcenters);
end

4. Judge convergence ：

   I set the maximum number of convergence to 150 Time , Or the sum of the squares before and after updating according to all the center points if <$10^{-5}$ I think it has not changed , Convergence .

Programming to achieve the above process , Finally get Corresponding to each pixel **$Clusters,size=N\times N$** as well as Of each central point RGB value ：$Centroids,size=K\ast 3$

Reassigning Colors to The Large Image

​ ad locum , I have processed both small and large graphs respectively .

​ Little picture ： For each pixel Turn into Corresponding class Of RGB value

​ Big picture ： First calculate each pixel and The nearest point to the center point after training , Then replace RGB value , Get new pictures .

​ Of a random initial value effect ：

​ Little picture ：

Big picture ： Higher resolution

Choosing A Good Initial Point

Choose a good starting point , To this end, I cycle through the above process , Set times . It is intended to choose a better starting point .

After each calculation of the last picture , Compare with the original drawing one by one pixel Of rgb Difference solution , Then on RGB Perform the square mean , To calculate the Loss, Get the lowest loss Pictures of the

As Best New Small / Large Image

The effect is as follows ：

Little picture ：

Big picture ：

Compare with the original picture , The result is right ：

At the same time, record the following Loss and Different times of random initial relations , Draw a chart ：

It can be seen that ： The overall effect is related to the initial of the initial point , Possible Local Mininum

Conclusion analysis and experience

K-Means Overall algorithm architecture

​ namely ：

​ Fix the center point each time , Then update the category of each point

​ Or fix the category of each point , According to this, the position information of the center point is updated

Some Limitations Of K-Means

K-Means With applicable scope ：

• Apply to every cluster It's about the same size , At this time, the effect will be good
• stay round-shaped The effect is good

At the same time non-convex shaped Bad performance

Kernel K-Means

That is, high-dimensional mapping , To solve that non-convex problem . While still not explicitly projecting , It's about using $K$

Hierarchical Clustering

Hierarchical clustering ,bottom-up、top-down

Experiment source code

Lab61.m

clear,clc
%% K-means
% small_img : 128 * 128 * 3
% large_img : 538 * 538 * 3
small_img = imread('bird_small.tiff');
large_img = imread('bird_large.tiff');
%  Turn into double  matrix  A(50,33,3)  Express  y = 50 x =33  Of Blue  Of val
small_img_matrix = double(small_img);
large_img_matrix = double(large_img);
%% show image 
% imshow(large_img)
%%
K = 16; % 16 individual class
[small_dimx,small_dimy,small_dimz] = size(small_img_matrix);
[large_dimx,large_dimy,large_dimz] = size(large_img_matrix);
max_iter_times = 200;
convergence_condition = 10^-5;
[centroids,Clusters] = kmeans(K,small_img_matrix,max_iter_times,convergence_condition);
%%  Redisplay the small picture 
new_small_pic=[];
new_large_pic=[];
%  Small picture processing 
for i=1:small_dimx
    for j=1:small_dimy
        new_small_pic(i,j,:) =  centroids(Clusters(i,j),:);
    end
end

%  Big picture processing 
for i=1:large_dimy
    for j=1:large_dimy
        pixel_rgb = [large_img_matrix(i,j,1),large_img_matrix(i,j,2),large_img_matrix(i,j,3)];
        diff = ones(K,1)*pixel_rgb - centroids;
        distance = sum(diff.^2,2);
        [~,class] = min(distance); %  Get categories 
        new_large_pic(i,j,:) =  centroids(class,:);
    end
end
%%  Show small pictures   and   Big picture 
figure
subplot(2,1,1);
imshow(uint8(small_img_matrix)),title('Origin Small Image');
subplot(2,1,2);
imshow(uint8(new_small_pic)),title('New Small Image')
%% 
figure
subplot(2,1,1);
imshow(uint8(large_img_matrix)),title('Origin Large Image');
subplot(2,1,2);
imshow(uint8(new_large_pic)),title('New Large Image')
%% K-Means Fuction Return ?
function [centroids,Clusters] = kmeans(K,datas,max_iter_times,convergence_condition)
    %  Return to the center point  K * dimz
    centroids = [];
    [dimy,dimx,dimz] = size(datas);
    %  Initialize the center point randomly  K*3
    centroids = sample_random(K,datas,dimx,dimy);
    for it = 1 : max_iter_times
        %  Calculation its nearest mean
        Clusters = calculate_Ci(centroids,dimx,dimy,datas,K);
        %  Update the center point 
        new_centroids = updatemeans(Clusters,dimx,dimy,dimz,datas,K);
        %  See if it converges 
        convergence = judge(new_centroids,centroids,convergence_condition);
        centroids = new_centroids;
        if convergence
            break
        end
    end
end

%%  Select the initial point randomly  return K * 3 
function sample_Datas = sample_random(num,datas,dimx,dimy)
    % datas  For raw data  num  Is the target number  
    sample_Datas = zeros(num,3);
    a = rand(dimx,1);
    b = rand(dimy,1);
    [vala,posa] = sort(a);
    [valb,posb] = sort(b);
    chose_x = posa(1:num,1);
    chose_y = posb(1:num,1);
    for i=1:num
        sample_Datas(i,:) = datas(chose_x(i),chose_y(i),:);
    end
    
end

%%  Calculate each pixel  Of   Category  return Clusters: N * N val  by   Category  
function Clusters = calculate_Ci(centroids,dimx,dimy,datas,K)
    Clusters = [];
    %  Through each pixel  Calculate a z(i,j)
    for i= 1:dimy
        for j = 1:dimx
            %  obtain xi
            pixel_rgb = [datas(i,j,1),datas(i,j,2),datas(i,j,3)];
            diff = ones(K,1)*pixel_rgb - centroids;
            distance = sum(diff.^2,2);
            [~,class] = min(distance); %  Get the smallest corresponding category index
            Clusters(i,j) = class;
        end
    end
end

%%  Update the center point  return 16 * 3
function newcenters = updatemeans(Clusters,dimx,dimy,dimz,datas,K)
    newcenters = zeros(K,dimz);
    nums =  zeros(K);
    sumred = zeros(K,1);
    sumgreen = zeros(K,1);
    sumblue = zeros(K,1);
    for i=1:dimx
        for j=1:dimy
            class = Clusters(i,j);
            nums(class) = nums(class) + 1;
            sumred(class) = sumred(class) + datas(i,j,1);
            sumgreen(class) = sumgreen(class) + datas(i,j,2);
            sumblue(class) = sumblue(class) + datas(i,j,3);
        end
    end
    for i=1:K
        if nums(i) ~= 0
            sumred(i) = sumred(i) / nums(i);
            sumgreen(i) = sumgreen(i) / nums(i);
            sumblue(i) = sumblue(i) / nums(i);
        end
    end
    newcenters = [sumred,sumgreen,sumblue];
    newcenters = round(newcenters);
end

%%  Judge whether convergence or not 
function convergence = judge(newcenter,oldcenter,condition)
    convergence = 0;
     d = sum(sqrt(sum((newcenter - oldcenter).^2, 2)));
    if d < condition
        convergence = 1;
    end
end

%%  Return row and column coordinates 
function [row,col] = findrc(Clusters,val)
    [dimx,dimy] = size(Clusters);
    row = [];
    col = [];
    for i=1:dimx
        for j=1:dimy
            if Clusters(i,j) == val
                row = [row;i];
                col = [col;j];
            end
        end
    end
end

Lab62.m ( Select the best effect through multiple iterations )

clear,clc
%%  Do it a few more times , The aim is to have a different starting point , Find the best .
%% K-means
% small_img : 128 * 128 * 3
% large_img : 538 * 538 * 3
small_img = imread('bird_small.tiff');
large_img = imread('bird_large.tiff');
%  Turn into double  matrix  A(50,33,3)  Express  y = 50 x =33  Of Blue  Of val
small_img_matrix = double(small_img);
large_img_matrix = double(large_img);
Best_new_small_pic=[];
Best_new_large_pic=[];
%% show image 
% imshow(large_img)
%%
inital_times = 10;
for inital_time = 1 : inital_times
    K = 16; % 16 individual class
    [small_dimx,small_dimy,small_dimz] = size(small_img_matrix);
    [large_dimx,large_dimy,large_dimz] = size(large_img_matrix);
    max_iter_times = 200;
    convergence_condition = 10^-5;
    [centroids,Clusters] = kmeans(K,small_img_matrix,max_iter_times,convergence_condition);
    %%  Redisplay the small picture 
    loss_small = 1e9;
    loss_large = 1e9;
    new_small_pic=[];
    new_large_pic=[];
    %  Small picture processing 
    for i=1:small_dimx
        for j=1:small_dimy
            new_small_pic(i,j,:) =  centroids(Clusters(i,j),:);
        end
    end
    new_loss_small = calculate_Kmeans_Loss(small_img_matrix,new_small_pic,small_dimx,small_dimy);
    if new_loss_small < loss_small
        loss_small = new_loss_small;
        Best_new_small_pic = new_small_pic;
    end

    %  Big picture processing 
    for i=1:large_dimy
        for j=1:large_dimx
            pixel_rgb = [large_img_matrix(i,j,1),large_img_matrix(i,j,2),large_img_matrix(i,j,3)];
            diff = ones(K,1)*pixel_rgb - centroids;
            distance = sum(diff.^2,2);
            [~,class] = min(distance); %  Get categories 
            new_large_pic(i,j,:) =  centroids(class,:);
        end
    end
    
    new_loss_large = calculate_Kmeans_Loss(large_img_matrix,new_large_pic,large_dimy,large_dimx);
    if new_loss_large < loss_large
        loss_large = new_loss_large;
        Best_new_large_pic = new_large_pic;
    end
end
%%  Show small pictures   and   Big picture 
figure
subplot(2,1,1);
imshow(uint8(small_img_matrix)),title('Origin Small Image');
subplot(2,1,2);
imshow(uint8(Best_new_small_pic)),title('Best New Small Image')
imwrite(uint8(Best_new_small_pic),'Best New Small Image.png')

figure
subplot(2,1,1);
imshow(uint8(large_img_matrix)),title('Origin Large Image');
subplot(2,1,2);
imshow(uint8(Best_new_large_pic)),title('Best New Large Image')
imwrite(uint8(Best_new_large_pic),'Best New Large Image.png')
%% K-Means Fuction Return ?
function [centroids,Clusters] = kmeans(K,datas,max_iter_times,convergence_condition)
    %  Return to the center point  K * dimz
    centroids = [];
    [dimy,dimx,dimz] = size(datas);
    %  Initialize the center point randomly  K*3
    centroids = sample_random(K,datas,dimx,dimy);
    for it = 1 : max_iter_times
        %  Calculation its nearest mean
        Clusters = calculate_Ci(centroids,dimx,dimy,datas,K);
        %  Update the center point 
        new_centroids = updatemeans(Clusters,dimx,dimy,dimz,datas,K);
        %  See if it converges 
        convergence = judge(new_centroids,centroids,convergence_condition);
        centroids = new_centroids;
        if convergence
            break
        end
    end
end

%%  Select the initial point randomly  return K * 3 
function sample_Datas = sample_random(num,datas,dimx,dimy)
    % datas  For raw data  num  Is the target number  
    sample_Datas = zeros(num,3);
    a = rand(dimx,1);
    b = rand(dimy,1);
    [vala,posa] = sort(a);
    [valb,posb] = sort(b);
    chose_x = posa(1:num,1);
    chose_y = posb(1:num,1);
    for i=1:num
        sample_Datas(i,:) = datas(chose_x(i),chose_y(i),:);
    end
    
end

%%  Calculate each pixel  Of   Category  return Clusters: N * N val  by   Category  
function Clusters = calculate_Ci(centroids,dimx,dimy,datas,K)
    Clusters = [];
    %  Through each pixel  Calculate a z(i,j)
    for i= 1:dimy
        for j = 1:dimx
            %  obtain xi
            pixel_rgb = [datas(i,j,1),datas(i,j,2),datas(i,j,3)];
            diff = ones(K,1)*pixel_rgb - centroids;
            distance = sum(diff.^2,2);
            [~,class] = min(distance); %  Get the smallest corresponding category index
            Clusters(i,j) = class;
        end
    end
end

%%  Update the center point  return 16 * 3
function newcenters = updatemeans(Clusters,dimx,dimy,dimz,datas,K)
    newcenters = zeros(K,dimz);
    nums =  zeros(K);
    sumred = zeros(K,1);
    sumgreen = zeros(K,1);
    sumblue = zeros(K,1);
    for i=1:dimx
        for j=1:dimy
            class = Clusters(i,j);
            nums(class) = nums(class) + 1;
            sumred(class) = sumred(class) + datas(i,j,1);
            sumgreen(class) = sumgreen(class) + datas(i,j,2);
            sumblue(class) = sumblue(class) + datas(i,j,3);
        end
    end
    for i=1:K
        if nums(i) ~= 0
            sumred(i) = sumred(i) / nums(i);
            sumgreen(i) = sumgreen(i) / nums(i);
            sumblue(i) = sumblue(i) / nums(i);
        end
    end
    newcenters = [sumred,sumgreen,sumblue];
    newcenters = round(newcenters);
end

%%  Judge whether convergence or not 
function convergence = judge(newcenter,oldcenter,condition)
    convergence = 0;
     d = sum(sqrt(sum((newcenter - oldcenter).^2, 2)));
    if d < condition
        convergence = 1;
    end
end

%%  Calculation K-Means  After the picture loss： Calculate all pixel Of rgb sum  then  avg
function Loss1 = calculate_Kmeans_Loss(ori_img,new_img,dimx,dimy)
    Loss = zeros(3,1);
    div_num = dimx*dimy;
    for i=1:dimx
        for j = 1 :dimy
            Loss(1) = Loss(1) + (ori_img(i,j,1)-new_img(i,j,1)).^2 / div_num;
            Loss(2) = Loss(2) + (ori_img(i,j,2)-new_img(i,j,2)).^2 / div_num;
            Loss(3) = Loss(3) + (ori_img(i,j,3)-new_img(i,j,3)).^2 / div_num;
        end
    end
    Loss1 = sum(Loss) / 3;
end

n

    convergence = 1;
end

end

%% Calculation K-Means After the picture loss： Calculate all pixel Of rgb sum then avg
function Loss1 = calculate_Kmeans_Loss(ori_img,new_img,dimx,dimy)
Loss = zeros(3,1);
div_num = dimx*dimy;
for i=1:dimx
for j = 1 :dimy
Loss(1) = Loss(1) + (ori_img(i,j,1)-new_img(i,j,1)).^2 / div_num;
Loss(2) = Loss(2) + (ori_img(i,j,2)-new_img(i,j,2)).^2 / div_num;
Loss(3) = Loss(3) + (ori_img(i,j,3)-new_img(i,j,3)).^2 / div_num;
end
end
Loss1 = sum(Loss) / 3;
end