聚类之isodata算法
近期上課學習了isodata的迭代自組織數據分析算法,該算法是一個動態的聚類算法,與K-均值算法有相似之處,即聚類中心同樣是通過樣本均值的迭代運算來決定的,但不同k-均值算法的是,isodata算法是可以通過計算的中間過程,對聚類中心的數量進行調整。對一個算法,要想真正的了解,是能夠在學習完該算法的基礎上,能夠自主地編寫,在編寫過程中你會發現在學習過程中的很多不足。因此,盡管網上已經有了成熟的isodata算法程序,樓主還是不才,自己編寫程序進行該算法的學習。該算法主要參照了蔡元龍老師的《模式識別》一書,程序是用matlab編寫。該算法寫的比較急,必然會出現一些不足,請大家能夠多多指出,共同進步。
首先,是對一系列參數的初始化,能力有限,所以在代碼中pars.L中的參數必須是1,即一次迭代運算中可以合并的聚類中心的最多對數最多是1,這一點會及時更新。
%programed by Lu Qi,my email:qqlu1992@gmail.com clear all clc % train_x=[0 0;0 1;1 0;1 1;7 7;8 8;9 9;8 7;8 9];load data_train %參數的初始化 %step one [pars.num_train,pars.length]=size(train_x);%num_train代表樣本的數量,length代表樣本的維數 pars.num_initial_cluster=5; pars.num_cluster=pars.num_initial_cluster; initial=randperm(pars.num_train); pars.iter_max=20; initial_point=initial(1:pars.num_initial_cluster); pars.z{1}=train_x(initial_point,:); %第一次初始聚類中心 clear initial initial_point pars.K=5; pars.k=1; %k代表分裂處理時的加減 pars.sita_N=3; %每一聚類域中最少的樣本數目,如少于此數就不作為一個獨立的聚類 pars.sita_S=1; %聚類域中樣本距離分布的標準差 pars.sita_C=1; %兩聚類中心之間的最小距離,如小于此數,兩個聚類進行合并 pars.L=1; %一次迭代運算中可以合并的聚類中心的最多對數 pars.iter=1;接著,是通過上一次的聚類中心進行歸類并更新各聚類中心值及類內和類間的平均距離。 for iter=1:pars.iter_maxfor i=1:pars.num_train %計算到聚類中心的分類,樣本的數量for j=1:pars.num_clusterd(i,j)=norm(train_x(i,:)-pars.z{iter}(j,:));end end[min_d,index]=min(d,[],2);clear min_dclear pars.g % pars.g(1:pars.num_cluster,1)=0;for i=1:pars.num_trainfor j=1:pars.num_clusterif (index(i))==jpars.g(j,i)=[;i]; %離每一個聚類最近的樣本點elsepars.g(j,i)=0;endendend % pars.g(1:pars.num_cluster,1)=[]; % for j=1:pars.num_cluster % if ~exsit('pars.g(j,:),) % pars.g(j,1:num_train)=0; % end % end%step threefor j=1:pars.num_clusterN(j)=sum(pars.g(j,:)~=0,2); %求每個聚類中的樣本數目endfor j=1:pars.num_clusterif N(j)<pars.sita_Ndelete_clusters=[;j]; %刪除該行endendif exist('delete_clusters')if delete_clusters~=0pars.g(delete_clusters,:)=[];endclear delete_clustersendpars.num_cluster=size(pars.g,1); %刪除后的聚類數量%step four%x修正各聚類中心值temp_sum_d=0;for j=1:pars.num_clustertemp=find(pars.g(j,:)); %找每一類中不等于0的數pars.z{iter+1}(j,:)=sum(train_x(temp,:))/(length(temp));q_max=length(temp);temp_d=0;for q=1:q_maxtemp_d=temp_d+norm(train_x(temp(q),:)-pars.z{iter+1}(j,:));end%step five,計算各聚類域中諸聚類中心間的平均距離pars.aver_d{iter+1}(j)=temp_d/(q_max); temp_sum_d=temp_sum_d+temp_d;end%step sixpars.aver_sum_d{iter+1}=temp_sum_d/pars.num_train;clear temp_d temp_sum_d q_max qpars.iter=iter;%進行分裂合并的處理工作[pars]=settle_split_merge(train_x,pars);clear pars.gfprintf('iter=%d\n',iter); end在該程序中,有很多編程的小細節要注意,寫完程序之后調了好長時間的BUG。在程序中,有一個分程序是,分裂合并處理 [pars]=settle_split_merge(train_x,pars),該程序如下:function [pars]=settle_split_merge(train_x,pars) % clear all % clc % load temporary % pars.k=0.5; % 進行合并或者分裂的處理 while 1if ~((pars.iter==pars.iter_max)||(rem(pars.iter,2)==0)||(size(pars.g,1)>=2*pars.K)) %跳到第8步%step 8,計算每聚類中樣本距離的標準差向量for j=1:pars.num_clusterfor i=1:pars.length %維數temp=find(pars.g(j,:)); %找每一類中不等于0的數temp_length=length(temp);for q=1:temp_lengthd(q,i)=(train_x(temp(q),i)-pars.z{end}(j,i))^2; %某個樣本的某一維到聚類中心的某一維的距離ends(i,j)=sqrt(sum(d(q,:))/temp_length); %橫坐標代表維數,縱坐標代表聚類數endend%step 9[s_jmax,index]=max(s,[],1);%step 10temp=length(s_jmax); %聚類數tt=0; %用于跳轉的標識符for j=1:tempif s_jmax(j)>pars.sita_Sif ((pars.aver_d{pars.iter+1}(j)>pars.aver_sum_d{pars.iter+1})&&(length(find(pars.g(j,:)))>(2*(pars.sita_N+1))))||(temp<=(pars.K/2))a=zeros(pars.length);dimension=index(j);a(dimension,dimension)=1;pars.z{pars.iter+1}(size(pars.z{pars.iter+1},1)+1,:)=pars.z{pars.iter+1}(j,:)+(pars.k*a*s(:,j))';pars.z{pars.iter+1}(j,:)=pars.z{pars.iter+1}(j,:)-(pars.k*a*s(:,j))';tt=1;pars.num_cluster=pars.num_cluster+1;clear a dimensionbreak; %跳轉到第二步!!!要進行修改endendendif tt==1break;endend% 第11步,合并處理if pars.iter==pars.iter_max %最后一次迭代判斷par.sita_C=0;endtemp=length(pars.z{pars.iter+1});for i=1:tempfor j=1:id_merge(i,j)=norm(pars.z{pars.iter+1}(i,:)-pars.z{pars.iter+1}(j,:));endendpp=max(d_merge(:));for i=1:size(d_merge,1)for j=i:size(d_merge,1)d_merge(i,j)=pp+1;endendclear pp[xsort,index]=sort(d_merge(:));if xsort(1)<pars.sita_Cfor i=1:pars.Lif xsort(i)<pars.sita_C[pp_raw,pp_col]=ind2sub(size(d_merge),index(i));pars.z{pars.iter+1}(min(pp_raw,pp_col),:)=(length(find(pars.g(min(pp_raw,pp_col)))~=0)*pars.z{pars.iter+1}(min(pp_raw,pp_col),:) ...+length(find(pars.g(max(pp_raw,pp_col)))~=0)*pars.z{pars.iter+1}(max(pp_raw,pp_col),:)) .../(length(find(pars.g(min(pp_raw,pp_col)))~=0)+length(find(pars.g(max(pp_raw,pp_col))))~=0);pars.z{pars.iter+1}(max(pp_raw,pp_col),:)=[];pars.num_cluster=pars.num_cluster-1; % temp=find(pars.g(j,:)); %找每一類中不等于0的數 % temp_length=length(temp);endendendbreak; end通過這兩個程序就可以運行出來結果了,例如我一共做了5類的樣本,最終結果如下:
?the 1th cluster center is 5.5000 ?5.5000?
?the 2th cluster center is 5.5000 -5.5000?
?the 3th cluster center is 0.5000 ?0.5000?
?the 4th cluster center is-5.5000 -5.5000?
?the 5th cluster center is-5.5000 ?5.5000
由于在參數初始化中,初始的聚類中心是隨機選的,所以最終的聚類結果和初始聚類中心影響比較大,有時我的會出現四類,而忽略了最后一類,甚至有時會出現三類,如圖:
the 1th cluster center is 5.5000 -5.5000?
the 2th cluster center is 5.5000 ?5.5000?
the 3th cluster center is-5.5000 -5.5000?
the 4th cluster center is-5.5000 ?5.5000
個人感覺isodata算法比較適合高維的樣本,如果是比較簡單的樣本聚類,效果甚至不如K-means。
如果有什么不足,請大家聯系我:qqlu1992@gmail.com
總結
以上是生活随笔為你收集整理的聚类之isodata算法的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: Dogleg“狗腿”最优化算法
- 下一篇: 矩阵之LU分解