0001 function [A,todss]=nt_cluster_jd(x,dsr,flags)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012 if nargin<3; flags=[]; end
0013 if nargin<2; error('!'); end
0014
0015 if ~exist('vl_kmeans');
0016 disp('vl_kmeans() not found, download from http://www.vlfeat.org');
0017 end
0018 if ndims(x)>2;
0019 error('x should be time*channels');
0020 end
0021
0022 disp('entering nt_cluster_jd...');
0023
0024
0025 [C0,C1,A,todss]=nt_bias_cluster(x,dsr,flags);
0026 todss2=todss(:,[1 end]);
0027 z=nt_mmat(x,todss2);
0028
0029 figure(2);clf
0030 subplot 511; plot(x); title('data');
0031 subplot 512; plot(A,'.-'); title('initial cluster map');
0032 subplot 513; plot(z(:,1)); title('DSS1');
0033 subplot 514; plot(z(:,2)); title('DSS2');
0034 drawnow;
0035
0036
0037 old_A=A;
0038 for k=1:10
0039 [~,~,A]=nt_bias_cluster(nt_normcol(z),dsr,flags);
0040 disp(['cluster sizes: ', num2str([sum(A==1), sum(A==2)])]);
0041 c1=nt_cov(x(find(A==1),:));
0042 c0=c1+nt_cov(x(find(A==2),:));
0043 [todss,pwr0,pwr1]=nt_dss0(c0,c1);
0044 todss2=todss(:,[1 end]);
0045 z=nt_mmat(x,todss2);
0046
0047 subplot 515; plot(A,'.-'); title('final cluster map');
0048 if old_A==A; break; end
0049 old_A=A;
0050 end
0051
0052 if nargout==0;
0053
0054 disp(['cluster1: ',num2str(100*numel(find(A==1))/numel(A)), '%']);
0055
0056 z1=nt_mmat(x,todss(:,1));
0057 z2=nt_mmat(x,todss(:,end));
0058 z=nt_mmat(x,todss);
0059 z=nt_normcol(z);
0060 e1=mean(z(find(A==1),:).^2);
0061 e2=mean(z(find(A==2),:).^2);
0062
0063 figure(100); clf
0064 plot(x); hold on
0065 xx=x;
0066 xx(find(A==2),:)=nan;
0067 plot(xx,'k');
0068 axis tight
0069 title('black: cluster2');
0070
0071 figure(101); clf
0072 subplot 121;
0073 plot(pwr1./pwr0,'.-'); xlabel('component'); ylabel('score'); title('DSS cluster A vs B');
0074 subplot 122;
0075 nt_spect_plot(z1,1024,[],[],1);
0076 hold on
0077 nt_spect_plot(z2,1024,[],[],1);
0078 xlim([0 .5])
0079 nt_colorlines([],[1 3]);
0080 legend('first','last'); legend boxoff
0081 hold off
0082
0083
0084 figure(102); clf
0085 subplot 211;
0086 plot(z1); axis tight
0087 title('first DSS component')
0088 subplot 212;
0089 plot(z2); axis tight
0090 title('last DSS component');
0091
0092 figure(103); clf
0093 plot([e1',e2'], '.-'); legend('cluster A', 'cluster B'); title ('power per component');
0094
0095 figure(104);
0096 subplot 121; imagescc(nt_cov(z(find(A==1),:))); title('cluster A');
0097 subplot 122; imagescc(nt_cov(z(find(A==2),:))); title('cluster B');
0098
0099 if 0
0100 figure(105); clf
0101 subplot 211;
0102 nt_sgram(z1,1024,32,[],1);
0103 title('first');
0104 subplot 212;
0105 nt_sgram(z2,1024,32,[],1);
0106 title('last');
0107 end
0108 clear c0 c1 A todss pwr0 pwr1
0109
0110 end
0111
0112 function y=norm2(x,n,ind)
0113 [I,J]=ind2sub([n,n],ind);
0114 for k=1:size(x,1)
0115 a=x(k,1:n);
0116 b=sqrt(a(I).*a(J));
0117 y(k,:)=x(k,:)./b;
0118 end
0119
0120
0121
0122