clear
clc
fg_load=csvread('fg-load.csv');
data1=[];
data=[];
%%对风能进行归一化处理%%
[f_shuju,f_shujus]=mapminmax(fg_load(:,2)',0,1);
fg_load(:,2)=f_shuju';
for i=1:8760
data1=[data1 fg_load(i,:)];
end
clear i;
for i=1:365
data(i,:)=data1(1,i*72-71:i*72);%%data为365*72
end
N=4; %%k数目
[m n]=size(data); %m是数据个数,n是数据维数
ma=zeros(n); %每一维最大的数
mi=zeros(n); %每一维最小的数
u=zeros(N,n); %随机初始化,最终迭代到每一类的中心位置
%%k-means++算法确定初始聚类%%
u(1,:)=data(round(rand(1)*m),:);
dis=inf(m,N);
for i=2:N
for j=1:m
dis(j,i-1)=norm(data(j,:)-u(i-1,:));
diss(j)=min(dis(j,:));
end
rand_diss=rand(1)*sum(diss);
for k=1:m
rand_diss=rand_diss-diss(k);
if rand_diss<=0
u(i,:)=data(k,:);
break;
end
end
end
jj=0;
while 1
jj=jj+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)=norm(tmp{index}(i,:));
end
for i=1:m
for j=1:N
if quan(i,j)>0
quan(i,j)=1;
end
end
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
% if jj>100
% break;
% end
end
re=[];
for i=1:m
tmp=[];
for j=1:N
tmp=[tmp norm(data(i,:)-u(j,:))];
end
[junk index]=min(tmp);
re=[re;data(i,:) index];
end
%%反归一化处理%%
for i=1:24
re(:,i*3-1)=mapminmax('reverse',re(:,i*3-1)',f_shujus)';
u(:,i*3-1)=mapminmax('reverse',u(:,i*3-1)',f_shujus)';%%聚类结果
end
% re1=re(:,1:72)';
% re2=re(:,73)';
%
% [Sw,Sb] = scatter_mat(re1,re2);
%
% PFS=(trace(Sb)/(N-1))/(trace(Sw)/(m-N))
%