clc%清除命令窗口变量
clear%清除工作空间变量
close all%关闭FIGURE图像
num = xlsread('片片追溯形变量数据（研磨前后）.xls');%读入数据
num(:,10:12) = num(:,10:12)*100;%预处理
n = randperm(85);
input_train=num(n(1:80),1:4)';%训练数据的输入数据
output_train=num(n(1:80),9:12)';%训练数据的输出数据
input_test=num(n(81:85),1:4)';%测试数据的输入数据
output_test=num(n(81:85),9:12)'; %测试数据的输出数据
%选连样本输入输出数据归一化
[inputn,inputps]=mapminmax(input_train);%训练数据的输入数据的归一化
[outputn,outputps]=mapminmax(output_train);%训练数据的输出数据的归一化de
%% BP网络训练
% %初始化网络结构
net=newff(minmax(inputn),[15,4],{'logsig','tansig'},'traincgp');%大括号里面修改传递函数 ，可以选择（purelin,logsig,tansig）,最后一个是训练函数，可以选
%择（traingd traingdx,trainlm,trainc,traincgb,traincgp,trainda,trainoss,trainlm,trainr,trainp,trains,trainscg）
net.trainParam.epochs=5000;%最大迭代次数
net.trainParam.lr=0.01;%学习率
net.trainParam.goal=0.001;%学习目标
%网络训练
net=train(net,inputn,outputn);
%% BP网络预测
%预测数据归一化
inputn_test=mapminmax('apply',input_test,inputps);
output_test1=mapminmax('apply',output_test,outputps);
%网络预测输出
an=sim(net,inputn_test);
%网络输出反归一化
BPoutput=(mapminmax('reverse',an,outputps));
BPoutput(2:4,:) = BPoutput(2:4,:)/100 ;
output_test(2:4,:) = output_test(2:4,:)/100 ;
% [nummax,mx1] = max(BPoutput,[],2);
% [nummin,mx2] = min(BPoutput,[],2);
%% 结果分析
figure(1)
subplot(4,1,1)
plot(BPoutput(1,:),'r-o')%预测的结果数据画图-代表虚线，O代表圆圈标识，r代表红色
hold on
plot(output_test(1,:),'k-*');%期望数据，即真实的数据画图，-代表实现，*就是代表*的标识,k代表黑色
hold off
legend('预测输出','期望输出')%标签
subplot(4,1,2)
plot(BPoutput(2,:),'r-o')%预测的结果数据画图-代表虚线，O代表圆圈标识，r代表红色
hold on
plot(output_test(2,:),'k-*');%期望数据，即真实的数据画图，-代表实现，*就是代表*的标识,k代表黑色
hold off
legend('预测输出','期望输出')%标签
subplot(4,1,3)
plot(BPoutput(3,:),'r-o')%预测的结果数据画图-代表虚线，O代表圆圈标识，r代表红色
hold on
plot(output_test(3,:),'k-*');%期望数据，即真实的数据画图，-代表实现，*就是代表*的标识,k代表黑色
hold off
legend('预测输出','期望输出')%标签
subplot(4,1,4)
plot(BPoutput(4,:),'r-o')%预测的结果数据画图-代表虚线，O代表圆圈标识，r代表红色
hold on
plot(output_test(4,:),'k-*');%期望数据，即真实的数据画图，-代表实现，*就是代表*的标识,k代表黑色
hold off
legend('预测输出','期望输出')%标签
error=BPoutput-output_test;%误差
mse1 = sum(sum(error.^2))/20;%均方根差
map = sum(sum(abs(error./output_test)))/20;%平均误差百分比
[nummax,mx1] = max(abs(error),[],2);%最大误差
[nummin,mx2] = min(abs(error),[],2);%最小误差
figure(2)
plot(error,'-*')
title('BP网络预测误差','fontsize',12)
ylabel('误差','fontsize',12)
figure(3)
plot((output_test'-BPoutput')./BPoutput','-*');
title('神经网络预测误差百分比')
legend('单位面积制冷能耗','单位面积制热能耗')
