K均值据类C语言和matlab实现

一、均值聚类算法的流程如图1所示：

图1：K-means流程图
算法描述：
Step1:读入数据记为x,输入聚类的点的起始行cow1、结束行位置cow2、分类个数K.
Step2:初始化聚类中心(取数据集合前K个向量)，记为K_position，K_position1.
Step3:遍历数据集合x,判断每一个向量属于哪一类.并将标志矩阵存储在数据集合x的下一行，用数字1:K来表示该向量属于哪一类.
Step4:更新聚类中心，将K_position1赋值给K_position，对K_position1更新为当前聚类中心值。
Step5:判断K_position1==K_position？，满足条件结束while循环,不满足条件，继续执行step3-step4.

二、仿真结果
本代码采用matlab实现，运行结果如图2所示。图中显示了对150个数据分为2类的运行结果，迭代了5次。第一类的聚类中心为[6.30 2.89 4.96 1.70]，第二类的聚类中心为[5.01 3.36 1.56 0.29]。数据集中归为第1类的有97组数值，归为第2类的有53组数值。

图2：运行结果

当取集合中不同的数据进行分类时，结果如表1所示。
序号 起始位置 结束位置 聚类中心1 聚类中心2

表1：聚类结果
使用系统自带的库函数计算结果如表2所示。

表2：库函数聚类结果
对比表1和表2结果可知，二者聚类中心基本相同。但是观察每一个向量的归属点发现，自己编写算法和库函数计算结果会出现细微的差距。
例如取51-100行元素进行测试，每一个元素的归属结果如表3所示，第二行表示自己实现的算法结果，第三行表示库函数实现的结果。

表3：归类结果
由表3可以知道，第6、12、17点出现误判，比较表1和表2第二行结果发现两者聚类中心基本一致。出现误判的原因可能是，本代码，计算出某一个向量距离2个聚类中心距离相等，将这个向量归属于第一类的缘故。index=find(s==min(s)); x(loop,lp1)=index(1)。index的值可能不唯一，所以取第一个索引，归类
matlab代码

clear,clc,close all;
x=textread('IrisData.txt')';
cow1=input('请输入行起始位置:');
cow2=input('请输入列结束位置');
len=size(x);len=len(2);
if(cow1>=1&&cow2<=len&&cow1<cow2)
x=x(1:4,cow1:cow2);
[l1,l2]=size(x);
K=input('请输入将您的数据分为几类');
K_position=x(:,[1:K]);
K_position1=x(:,[1:K]);
flag=1;
loop=l1+1;
while(flag||~isequal(K_position,K_position1))
    %~isequal(roundn(K_position,-3),roundn(K_position1,-3))
    flag=0;
    for lp1=1:l2
        %判断这个点距离哪个聚类中心近
        s=zeros(1,K);
        for lp2=1:K
            for lp3=1:l1
                s(lp2)=s(lp2)+(x(lp3,lp1)-K_position1(lp3,lp2))^2;
            end
        end
        index=find(s==min(s));
        x(loop,lp1)=index(1);
    end
    %然后更新聚类中心
    temp=zeros(l1+1,K);
    for lp4=1:l2
        temp(1:l1,x(loop,lp4))=temp(1:l1,x(loop,lp4))+x(1:l1,lp4);
        temp(l1+1,x(loop,lp4))=temp(l1+1,x(loop,lp4))+1;
    end
    %得到新聚类中心。
    for lp5=1:K
        temp(1:l1,lp5)=temp(1:l1,lp5)/temp(l1+1,lp5);
    end
    K_position=K_position1;
    K_position1=temp(1:l1,:);
    loop=loop+1;
    disp(strcat(num2str(loop-1-l1),'轮计算'));
end
x_index=zeros(K,cow2-cow1+1);
lp11=ones(1,K);
for lp10=1:cow2-cow1+1
    x_index(x(end,lp10),lp11(1,x(end,lp10)))=lp10;
    lp11(1,x(end,lp10))=lp11(1,x(end,lp10))+1;
end
len_K=zeros(1,K);
disp('-----------------');
for lp12=1:K
    disp(['第',num2str(lp12),'类聚类中心为',num2str(K_position1(:,lp12)')]);
    cnt=0;
    for lp13=1:cow2-cow1+1
        if(x_index(lp12,lp13))
            cnt=cnt+1;
        end
    end
    len_K(1,lp12)=cnt;
    disp(x_index(lp12,1:cnt));
    disp('---------------------------');
end
else
    disp('输入有误');
end
disp('自己算的');
K=K_position1'
disp('库函数');
x1=textread('IrisData.txt')';
x1=x1(1:4,cow1:cow2)';
[idx,center ]= kmeans(x1,2);
center=center'
idx=idx';

C语言代码

#include<stdio.h>
#include<stdlib.h>
#include<math.h>
/*
为了以后特定行数的读取，採用宏定义的方法来确定行数
程序中read（）读取该文件的行数，当未知文件行数时可使用返回值做全局变量
*/
#define N 150   //行
#define L 4       //列
#define cow1 0 //行起始位置
#define cow2  49 //行结束位置
#define cow3 50 //总行数cow2-cow1+1;
#define K 2   //2类
//const char file_name[50] = "d:\\dat.txt";
const char file_name[50] = "IrisData.txt";
void read(FILE* fp)
{
    int row = 0;
    char mid;
    while (!feof(fp))
    {
        mid = fgetc(fp); //从txt文本中读取一个字符赋值给mid
        if (mid == '\n') //假设这个字符为换行符
            row++; //记录txt数据行数
    }
    row++; //最后一行没有换行符
   // printf("行数为%d\n", row);
    rewind(fp); //回文件起始位置 
}
//定义函数来比较两个数组是否相等,只能是两个数组
int Comp_Array(double K_position[][L], double K_position1[][L], int len1, int len2)
{
    int flag2;//初始化返回值
    flag2= 1;
    for (int i = 0; i < len1; i++) 
    {
        for (int j = 0; j < len2; j++) 
        {
            if (K_position[i][j] != K_position1[i][j]) 
            {
                flag2 = 0;
                break;
            }
        }
    }
    return flag2;
}
//寻找数组中最小元素的索引
int Find_min_index(double Dis[], int L_end)
{
    int Dis_min_index = 0;
    for (int i = 1; i < L_end; i++) 
    {
        if (Dis[i] < Dis[Dis_min_index]) 
        {
            Dis_min_index = i;
        }
    }
    return Dis_min_index;
}

//替换聚类中心
void K_p_replace(double K_position[][L], double K_position1[][L],int l1, int l2)
{
    for (int i = 0; i < l1; i++) 
    {
        for (int j = 0; j < l2; j++) 
        {
            K_position[i][j] = K_position1[i][j];

        }
    }
}

//计算新的聚类中心
void K_P1_process(double K_position1[][L], int l1, int l2,int K_index[],int l3,double data2[][L])
{
    //初始化聚类中心K_position1的值，将其中的元素全部赋值为0
    for (int i = 0; i < l1; i++)
    {
        for (int j = 0; j < l2; j++)
        {
            //printf("数字%f", K_position1[i][j]);//传递过来的数值没有错误
            K_position1[i][j] = 0;
        }
    }
    //遍历数据集合data2,将其叠加到对应的位置
    int cnt[K] = { 0 };
    for (int i = 0; i < l3; i++) 
    {
        //通过数据集合K_index判断该点属于哪个类
        for (int j = 0; j < l2; j++) 
        {
            K_position1[K_index[i] - 1][j] = K_position1[K_index[i] - 1][j] + data2[i][j] ;//这块直接除以l3导致错误
            //printf("数字%f,%d\t", K_index[i] - 1);//传递过来的数值没有错误
        }
        //printf("\n");
        cnt[K_index[i] - 1] = cnt[K_index[i] - 1] + 1;
    }
    //计算平均值
    for (int i = 0; i < l1; i++) 
    {
        for (int j = 0; j < l2; j++) 
        {
            //printf("数字%f,%d\n", K_position1[i][j],cnt[i]);//传递过来的数值没有错误
            K_position1[i][j] = K_position1[i][j] / cnt[i];
            //printf("数字%f\t", K_position1[i][j]);            
        }
        //printf("\n");
    }
    //printf("\n\n\n");
}
int main(int argc, char* argv[])
{
    //读取文本文件的值存储在data中
    FILE* fp;
    double data[N][L] = { 0.0 };   //二维数组
    int index[N] = { 0 }; //二维数组行下标
    int K_index[cow3] = { 0 };//定义每一个数据对应哪一个类
    double temp;
    int i, j;
    double data2[cow3][L] = { 0.0 };//定义选取的文本数据
    double K_position[K][L] = { 0.0 };//定义聚类中心
    double K_position1[K][L] = { 0.0 };//定义聚类中心
    int count = 0;  //计数器，记录已读出的浮点数
    if ((fp = fopen(file_name, "rb")) == NULL) {
        printf("请确认文件(%s)是否存在!\n", file_name);
        exit(1);
    }
    read(fp);     //读取行数
    while (1 == fscanf(fp, "%le", &temp)) //lf,le都能够。但别的都不能够，%e也不行
    {
        data[(index[count % L])++][count % L] = temp;
        count++;
    }
    fclose(fp);   //关闭句柄
    /******处理数据****************/
    for (i = 0; i < N; i++)
    {
        //printf("第%d行数据：\n", i + 1);
        for (j = 0; j < L; j++)
        {
            //printf("%5.3le ", data[i][j]);//.16f能够,le时以科学计数法显示
        }
        //printf("\n");
    }
    //对数据进行聚类分析
    //对数据进行初始化
    //printf("选择的数据-----------------\n");
    for (int i = 0; i < cow3; i++)
    {
        //printf("第%d行数据：\t", i + 1);
        for (int j = 0; j < L; j++)
        {
            data2[i][j] = data[cow1 + i][j];
            //printf(" % f", data2[i][j]);
        }
        //printf("\n");
    }
    //printf("初始的聚类中心+++++++++++++++++++++\n");
    for (int i = 0; i < K; i++)
    {
        for (int j = 0; j < L; j++)
        {
            K_position[i][j] = data[cow1 + i][j];
            K_position1[i][j] = K_position[i][j];
        }
        //printf("\n");
    }
    int flag = 1;
    while (flag || !Comp_Array(K_position, K_position1, K, L))
    {
        flag = 0;
        //判断当前点data2离哪个聚类中心K_position近
        for (int i = 0; i < cow3; i++)
        {
            double Dis[K] = { 0.0 };///定义点到聚类中心的距离
            for (int j = 0; j < L; j++)
            {
                for (int k = 0; k < K; k++)
                {
                    Dis[k] = Dis[k] + pow((data2[i][j] - K_position1[k][j]), 2);
                }
            }
            //寻找K个值中的最小值，并且返回索引
            int min_index = Find_min_index(Dis, K);
            K_index[i] = min_index + 1;
        }
        //更新聚类中心
        K_p_replace(K_position, K_position1, K, L);
        K_P1_process(K_position1, K, L, K_index, cow3, data2);//计算新的聚类中心
    }
    for (int i = 0; i < K; i++) 
    {
        printf("第%d类聚类中心为:", i + 1);
        for (int j = 0; j < L; j++)
        {
            printf("%f\t", K_position1[i][j]);
        }
        printf("\n有数据:\t");
        for (int j = 0; j < cow3; j++)
        {
            if (K_index[j] == i + 1) 
            {
                printf("%d  ",j+1);
            }
        }
        printf("\n");
    }
}