谱聚类:Ng算法

来源:互联网 发布:手机版开淘宝店铺 编辑:程序博客网 时间:2024/05/21 16:59

文章: On Spectral Clustering: Analysis and an algorithm

作者: Andrew Y. Ng  




 在算法第二步并没有讲清楚,矩阵A是表示的是各点之间的亲和值?我按照这个算法写了后对于半月形的数据集无法得出好的结果,在参考其他人的资料后,才意识到中途有一个图构造的过程,即取与点xi最邻近的k个点,只保留这一部分信息,来计算图的权重矩阵W。代码如下:


function re = spectral(X,k,KnearNum,sigma)%UNTITLED2 Summary of this function goes here%   Detailed explanation goes here    %Matrix A    [m,n] =size(X);    A = zeros(m,m);    D = zeros(m,m);   %度    %KnearNum = 10;    %构图采用10邻近    d = sigma;    %1-step: A    for i = 1:m        for j =1:m            if i~=j                A(i,j) = exp(-(norm(X(i,:)-X(j,:)).^2)/(2*d.^2));                %D(i,i) = D(i,i)+ A(i,j);            end        end    end        A_s= zeros(m,m);    in_A=zeros(m,m);    for i=1:m        [A_s(i,:) in_A(i,:)] = sort(A(i,:),'descend');    end        %W 权重向量    w= zeros(m,m);    for i=1:m        for j = 1:KnearNum            w(i,in_A(i,j))=1;  %将亲和点之间的权重赋值为1        end    end        for i=1:m        D(i,i)=sum(w(i,:));    end    %lapl      L = D^(-0.5)*w*D^(-0.5);        %3-step    [a b] = eig(L);    bx = eig(L);    [b_s in] =sort(bx,'descend');    a_s = a(:,in);%     x=[]%     for i=1:k%         x=[x;a(i,:)]%     end    x= a(:,1:k);        %4-step    y=zeros(m,k);    for i=1:m        s = (sum(x(i,:).^2)).^(1/2);        for j=1:k            y(i,j) = x(i,j)/s;        end    end        %5-step :kmeans    [u re hi]=kmean(y,k);        %6-step :result                  end


所调用kmean函数,代码:

%N是数据一共分多少类%data是输入的不带分类标号的数据%u是每一类的中心%re是返回的带分类标号的数据%hi  每一类的统计数据function [u re hi]=kmean(data,N)       [m n]=size(data);   %m是数据个数,n是数据维数    ma=zeros(n);        %每一维最大的数    mi=zeros(n);        %每一维最小的数    u=zeros(N,n);       %随机初始化,最终迭代到每一类的中心位置    for i=1:n       ma(i)=max(data(:,i));    %每一维最大的数       mi(i)=min(data(:,i));    %每一维最小的数       for j=1:N            u(j,i)=ma(i)+(mi(i)-ma(i))*rand();  %随机初始化,不过还是在每一维[min max]中初始化好些       end          end    u    while 1        pre_u=u;            %上一次求得的中心位置        for i=1:N            tmp{i}=[];      % 公式一中的x(i)-uj,为公式一实现做准备            for j=1:m                tmp{i}=[tmp{i};data(j,:)-u(i,:)];            end        end                quan=zeros(m,N);        for i=1:m        %公式一的实现            c=[];            for j=1:N                c=[c norm(tmp{j}(i,:))];            end            [junk index]=min(c);            quan(i,index)=1;                   end                for i=1:N            %公式二的实现           for j=1:n                u(i,j)=sum(quan(:,i).*data(:,j))/sum(quan(:,i));           end                   end                if norm(pre_u-u)<0.1  %不断迭代直到位置不再变化            break;        end    end        re=[];    figure(2);hold on;    for i=1:m        tmp=[];        for j=1:N            tmp=[tmp norm(data(i,:)-u(j,:))];        end        [junk index]=min(tmp);        re=[re; index];%         if index==1%             plot(data(i,1), data(i,2), 'r.'); %         elseif index==2%             plot(data(i,1), data(i,2), 'b.'); %          elseif index==3%             plot(data(i,1), data(i,2), 'k.'); %          elseif index==4%             plot(data(i,1), data(i,2), 'g.'); %          elseif index==5%             plot(data(i,1), data(i,2), 'm.'); %          elseif index==6%             plot(data(i,1), data(i,2), 'm.');%         end    hi=[];    for i=1:N         hi=[hi length(find(re(:)==i))];    endend
运行结果:

判断knearNum对精确度的影响:





0 0