nntwarn off
clear
%数据输入
p=xlsread('输入因子.xls');
t=xlsread('输出因子.xls');
pn=premnmx(p);
tn=premnmx(t);
[R,Q]=size(pn);
iitst=2:4:Q;
iival=4:4:Q;
iitr=[1:4:Q 3:4:Q];
vval.P=pn(:,iival);vval.T=tn(:,iival);
test.P=pn(:,iitst);test.T=tn(:,iitst);
ptr=pn(:,iitr);ttr=tn(:,iitr);
P=pn;
T=tn;
R=size(P,1);
S2=size(T,1);
S1=12;%隐含层节点数
S=R*S1+S1*S2+S1+S2;%遗传算法编码长度
aa=ones(S,1)*[-1,1 ];
popu=100;%种群规模
initPpp=initializega(popu,aa,'gabpEval');%初始化种群
gen=300;%遗传代数
%下面调用gaot工具箱，其中目标函数定义为gabpEval
[x,endPop,bPop,trace]=ga(aa,'gabpEval',[],initPpp,[1e-6 1 1],'maxGenTerm',gen,...
'normGeomSelect',[0.08],[tarithXover'],[2],'nonUnifMutation',[2 gen 3]);
%绘收敛曲线图
figure(1)
plot(trace(:,1),1./trace(:,3),'r-');
hold on
plot(trace(:,1),1./trace(:,2),'b-');
xlabel('Generation');
ylabel('Sum—Squared Error');
figure(2)
plot(trace(:,1),trace(:,3),'r-');
hold on
plot(trace(:,1),trace(:,2),'b-');
xlabel('Generation');
ylabel('Fittness');
%下面将初步得到的权值矩阵赋给尚未开始训练的BP网络
net=newff(minmax(pn),[12,1],{'tansig','purelin'},'trainlm');
[W1,B1,W2,B2,P,T,A1,A2,SE,val]=gadecod(x);
net.iw{1,1}=W1;
net.lw{2,1}=W2;
net.b{1}=B1;
net.b{2}=B2;
%数据输入2：网络有关参数
net.trainParam.epochs=2000;
net.trainParam.goal=0.001;
tml=cputime;
[net,tr]=train(net,ptr,ttr,[],[],vval,test);
tm2=cputime;
y=sim(net,vval.P);
e=vval.T-y;
error=mse(e,net)%error为网络的误差向量
r=norm(error);%r为网络的整体误差
tm=tm2-tml;%测试时间
%输出训练后的权值和阈值
iw=net.IW{1}
bl=net.b{1}
lw2=net.LW{2}
b2=net.b{2}
save net51 net； %保存最好的网络
figure(3)
plot(tr.epoch,tr.perf,'-',tr.epoch,tr.vperf,':',tr.epoch,tr.tperf,'-.')
legend('Training','Validation','Test',-1);
yhbel('平方差');xlabel('时间');
load net51 net；
yuce=sim(net,test.P);
%加载训练后的bp网络
%结果反归一化
yl=y;
yucel=yuce;
t1=((y1+1)*(max(t)-min(t)))/2+min(t);
yuce2=((yucel+1)*(max(t)-min(t)))/2+min(t);
%计算误差和标准差
wucha=yuce2-t(iitst);
s=std(wucha);
%作图
figure(4)
plot(iitst,t(iitst),'*-',iitst,yuce2,'o：');
figure(5)
plot(iitst,wucha,'r','Linewidth',3);
