clear,close all
%无随机取,取1到9个连续整周波扰动,在十个整周波上,(含噪),200个测试集,全部模型
w=100*pi;
%t=1:1/1040:1.2;
t=0:1/6400:0.2;%取1028 2^10=1024
t1=0.071;
t2=0.142;
%pure sin
%figure(1);
% subplot(5,1,1);
u01=sin(w*t);
u0=awgn(u01,25);
% plot(t,u0);
%randperm(10)%随机取10个数,从1-10中间的整数
%swell
for c1=1:1:200;
%subplot(5,1,2);
u1(c1,:)=sin(w*t).*(1+0.2*(heaviside(t-t1)-heaviside(t-t2)));
u(c1,:)=awgn(u1(c1,:),25);
%plot(t,u);
end
%sag
for c1=1:1:200;
v1(c1,:)=sin(w*t).*(1-0.8*(heaviside(t-t1)-heaviside(t-t2)));
v(c1,:)=awgn(v1(c1,:),25);
%subplot(5,1,3);
%plot(t,v);
end
%interupt
for c1=1:1:200;
vi1(c1,:)=sin(w*t).*(1-1*(heaviside(t-t1)-heaviside(t-t2)));
vi=awgn(vi1,25);
%subplot(5,1,4);
%plot(t,vi);
end
%homonics
for c1=1:1:200;
vv1(c1,:)=sin(w*t)+0.1*sin(3*w*t)+0.1*sin(5*w*t)+0.1*sin(7*w*t);
vv(c1,:)=awgn(vv1(c1,:),25);
%subplot(5,1,5);
%plot(t,vv);
end
%pluse
for c1=1:1:200;
[z11,t3] = gensig('pulse',t1,0.2,0.00015625);%1.2 0.0011673152
z1(c1,:)=sin(w*t)+[0.4*z11]';%幅度0.1-0.8
z(c1,:)=awgn(z1(c1,:),25);
%subplot(7,1,6);plot(t,z);
end
%osc
for c1=1:1:200;
zz1(c1,:)=sin(w*t)+0.4*exp(10*(t1-t2))*sin(20*w*t).*(heaviside(t-t1)-heaviside(t-t2));%取20比较象 osc.m
zz(c1,:)=awgn(zz1(c1,:),25);
end
%figure (2);
for c1=1:1:200;
[c0,l0]=wavedec(u0,10,'db4');
[Ea0,Ed0] = wenergy(c0,l0);
[cu,lu]=wavedec(u(c1,:),10,'db4');
[Eau,Edu(c1,:)] = wenergy(cu,lu);
[cv,lv]=wavedec(v(c1,:),10,'db4');
[Eav,Edv(c1,:)] = wenergy(cv,lv);
[cvi,lvi]=wavedec(vi(c1,:),10,'db4');
[Eavi,Edvi(c1,:)] = wenergy(cvi,lvi);
[cvv,lvv]=wavedec(vv(c1,:),10,'db4');
[Eavv,Edvv(c1,:)] = wenergy(cvv,lvv);
[cz,lz]=wavedec(z(c1,:),10,'db4');
[Eaz,Edz(c1,:)] = wenergy(cz,lz);
[czz,lzz]=wavedec(zz(c1,:),10,'db4');
[Eazz,Edzz(c1,:)] = wenergy(czz,lzz);
end
t3=1:10;
%plot(t3,Ed0,'c',t3,Edu,'g',t3,Edv,'k',t3,Edvi,'m',t3,Edvv,'y');
plot(t3,Ed0,'g',t3,Edu,'r',t3,Edv,'k',t3,Edvi,'m',t3,Edvv,'y',t3,Edz,'b',t3,Edzz,'g*');
legend('Ed0','Edu','Edv','Edvi','Edvv','Edz','Edzz');
X1=[Edu(1:100,:)',Edv(1:100,:)',Edvv(1:100,:)'];
X22=ones(200,1);
X2=Edvi;%outage
X44=ones(200,1);
X4=X44*Ed0 ;%pure
%X=[X1,X2(1:90,:)',X4(1:90,:)']
%X=[X4(1:200,:),X1,X2(1:200,:),Edz(1:200,:),Edzz(1:200,:)];
%-----------------四份类例子
%算法一 多类支持向量机(一对多算法)
%训练数据产生(四类)
xapp=[X4(1:100,:);Edu(1:100,:);Edv(1:100,:);Edvv(1:100,:);X2(1:100,:);Edz(1:100,:);Edzz(1:100,:)];
yapp=[1*ones(1,100) 2*ones(1,100) 3*ones(1,100) 4*ones(1,100) 5*ones(1,100) 6*ones(1,100) 7*ones(1,100)]';
%------------------%用svm分类,只要特征与目标的列(类似clouds数据的列都是50000)一致即可,要测试的模式行与X一致,列可以是一个也可以是多个,C4_
%5返回的结果直接是分类后的结果-目的是对列向量进行分类,行代表属性
%试验数据产生
xtest=[X4(101:200,:);Edu(101:200,:);Edv(101:200,:);Edvv(101:200,:);X2(101:200,:);Edz(101:200,:);Edzz(101:200,:)];
ytest=[1*ones(1,100) 2*ones(1,100) 3*ones(1,100) 4*ones(1,100) 5*ones(1,100) 6*ones(1,100) 7*ones(1,100)]';
%一对多算法
% 参数设置
c = 1000;
lambda = 1e-7;
kerneloption= 2;
kernel='gaussian';
verbose = 1;
%---------------------One Against All algorithms----------------
nbclass=7;
%训练
[xsup,w,b,nbsv]=svmmulticlassoneagainstall(xapp,yapp,nbclass,c,lambda,kernel,kerneloption,verbose);
%计算训练样本自身分类误差
[ypred] = svmmultival(xapp,xsup,w,b,nbsv,kernel,kerneloption);
fprintf( '\nRate of correct class in training data : %2.2f \n',100*sum(ypred'==yapp)/length(yapp));
%测试
[ypred,maxi] = svmmultival(xtest,xsup,w,b,nbsv,kernel,kerneloption);
%识别率
fprintf( '\nRate of correct class in test data : %2.2f \n',100*sum(ypred'==ytest)/length(ytest));
%-----------------------------------------------------------------
%XX=[X4(201:400,:)',Edu(201:400,:)',Edv(201:400,:)',Edvv(201:400,:)',X2(201:400,:)',Edz(201:400,:)',Edzz(201:400,:)'];%测试集
%t=treefit(X',Y'); % Create decision tree
%sfit = treeval(t,XX'); % Find assigned class numbers