附录
(1)BP神经网络模型程序:
p = [91 65 106.6 4838.9 48.8;
92 70 99.6 5160.3 46.6; 110 212 98.8 5425.1 44.7; 124 263 96.8 5854.02 42.1; 136 414 98.5 6280 39.4; 142 421 98.4 6859.6 38.2; 140 666 100.2 7702.8 37.7; 140 1040 100.1 8472.2 37.1; 145 1425 103 9421.6 37.7; 147 1559 100.9 10493 36.7; 186 1761 101.2 11759.5 35.8; 190 2206 104 13785.8 36.3; 200 2316 106.1 15780.76 37.9; 215 2433 99.1 17174.65 36.5; 243 3088 102.1 19109 35.7]'; t = [1309 1614 1620 2094 2537 2900 3270 3391 4089 5058 6000 6944 8100 9399.7 10815]';
[pn,minp,maxp,tn,mint,maxt]=premnmx(p,t);
net = newff(minmax(pn),[5,14,1],{'tansig' 'tansig' 'purelin'},'trainlm'); net.trainParam.epochs=20000;%训练次数设置 net.trainParam.goal=1e-6;%训练目标设置
net.trainParam.lr=0.001;%学习率设置,应设置为较少值,太大虽然会在开始加快收敛速度,但临近最佳点时,会产生动荡,而致使无法收敛 net=train(net,pn,tn); an=sim(net,pn);
y=postmnmx(an,mint,maxt)
20
[m,b,r]=postreg(y,t); %计算误差 All_error=[]; for i=1:15
m=(t(i)-y(i))/t(i);
All_error=[All_error,m];
disp(['百分相对误差为:',num2str(m)]); end figure
xx=1:length(All_error); %plot(xx,All_error); %title('误差变化图'); %计算仿真误差 E = t- y;
MSE=mse(E)
%对BP网络进行仿真 echo off figure
plot((1996:2010),t,'-*',(1996:2010),y,'-o') xlabel('年份')
ylabel('旅游流量(万)') title('仿真图')
plot(p,t,'*r',p,y,':b')
title('*为真实值,:为预测值');
(2)灰色理论GM(1,1)模型程序:
y=input('请输入数据 ');%输入数据请用如例所示形式:[1309 1614 1620 2094 2537 2900 3270 3391 4089 5058 6000 6944 8100 9399.7 10815]or [50.15 79.35 81.64 111.29 134.6 161.39 191.1 197.47 240.81 320.02 390.89 463.67 559.38 675.61 818] n=length(y); yy=ones(n,1); yy(1)=y(1); for i=2:n
yy(i)=yy(i-1)+y(i);%对原始灰色数据序列作一次累加 end
B=ones(n-1,2); for i=1:(n-1)
B(i,1)=-(yy(i)+yy(i+1))/2;%B矩阵 B(i,2)=1; end BT=B'; for j=1:n-1
YN(j)=y(j+1);
21
end
YN=YN';
A=inv(BT*B)*BT*YN; a=A(1);%求解a u=A(2);%求解U t=u/a;
t_test=input('请输入需要预测个数:'); i=1:t_test+n;
yys(i+1)=(y(1)-t).*exp(-a.*i)+t;%预测方程 时间响应式 yys(1)=y(1);
for j=n+t_test:-1:2
ys(j)=yys(j)-yys(j-1);%还原后的预测结果及后几年的预测值 end x=1:n;
xs=2:n+t_test; yn=ys(2:n+t_test);
plot(x,y,'^r',xs,yn,'*-b');%原数据与预测数据图 disp(['预测值为: ',num2str(ys(2:n+t_test))]); e0=[];
%计算百分相对误差 for i=2:n
e(i)=y(i)-ys(i); m=e(i)/y(i); e0=[e0,e(i)];
disp(['百分相对误差为:',num2str(m)]); end
%计算关联度
max1=max(abs(e0)); r=1;
for k=1:n-1
r=r+0.5*max1/(abs(e0(k))+0.5*max1); end
r=r/n; % r 表示关联度
disp(['关联度为:',num2str(r)]); %计算百分相对误差 for i=2:n
det=abs(ys(i)-y(i));
disp(['百分绝对误差为:',num2str(det)]); end
(3)多元线性回归模型程序: ①开始的程序: y=[1309 1614 1620
2094 2537 2900
3270 3391 4089
22
5058
6000 6944 8100 9399.7 10815]';
x1=[50.15 79.35 81.64 111.29 134.6 161.39 191.1 197.47 240.81 320.02 390.89 463.67 559.38 675.61 818]'; x2=[91 92 110 124 136 142 140 140 145 147 186 190 200 215 243]';
x3=[65 70 212 263 414 421 666 1040 1425 1559 1761 2206 2316 2433 3088]'; x4=[106.6 99.6 98.8 96.8 98.5 98.4 100.2 100.1 103 100.9 101.2 104 106.1 99.1 102.1]';
x5=[4838.9 5160.3 5425.1 5854.02 6280 6859.6 7702.8 8472.2 9421.6 10493 11759.5 13785.8 15780.76 17174.65 19109]';
x6=[48.8 46.6 44.7 42.1 39.4 38.2 37.7 37.1 37.7 36.7 35.8 36.3 37.9 36.5 35.7 ]'; x=[ones(15,1) x1 x2 x3 x4 x5 x6]; [b,bint,r,rint,stats]=regress(y,x) rcoplot(r,rint) ②修正程序:
y=[1309 1614 1620 2094 2537 2900 3270 3391 4089 5058
6000 6944 8100 9399.7 10815]';
x1=[50.15 79.35 81.64 111.29 134.6 161.39 191.1 197.47 240.81 320.02 390.89 463.67 559.38 675.61 818]'; x2=[91 92 110 124 136 142 140 140 145 147 186 190 200 215 243]';
x3=[65 70 212 263 414 421 666 1040 1425 1559 1761 2206 2316 2433 3088]'; x4=[106.6 99.6 98.8 96.8 98.5 98.4 100.2 100.1 103 100.9 101.2 104 106.1 99.1 102.1]';
x5=[4838.9 5160.3 5425.1 5854.02 6280 6859.6 7702.8 8472.2 9421.6 10493 11759.5 13785.8 15780.76 17174.65 19109]';
x6=[48.8 46.6 44.7 42.1 39.4 38.2 37.7 37.1 37.7 36.7 35.8 36.3 37.9 36.5 35.7 ]'; x=[ones(15,1) x1 x2 x3 x4 x5 x6]; [b,bint,r,rint,stats]=regress(y,x) rcoplot(r,rint)
y1=-206.055+14.1299*x1+0.6997*x2+0.2806*x3+27.6024*x4-0.1361*x5-31.5573*x6 %计算重新预测值
s=abs(y-y1)./y %计算相对误差 ③拟合程序及结果:
y=[1309 1614 1620 2094 2537 2900 3270 3391 4089 5058 6000 6944 8100 9399.7 10815]';
x1=[50.15 79.35 81.64 111.29 134.6 161.39 191.1 197.47 240.81 320.02 390.89 463.67 559.38 675.61 818]';
23
x2=[91 92 110 124 136 142 140 140 145 147 186 190 200 215 243]';
x3=[65 70 212 263 414 421 666 1040 1425 1559 1761 2206 2316 2433 3088]'; x4=[106.6 99.6 98.8 96.8 98.5 98.4 100.2 100.1 103 100.9 101.2 104 106.1 99.1 102.1]';
x5=[4838.9 5160.3 5425.1 5854.02 6280 6859.6 7702.8 8472.2 9421.6 10493 11759.5 13785.8 15780.76 17174.65 19109]';
x6=[48.8 46.6 44.7 42.1 39.4 38.2 37.7 37.1 37.7 36.7 35.8 36.3 37.9 36.5 35.7 ]'; plot(x1,y,'o') plot(x2,y,'o') plot(x3,y,'o') plot(x4,y,'o') plot(x5,y,'o') plot(x6,y,'o')
④预测2011-2015的旅游人数: x=16:20 for i=1:5
x1=12.7.*x(i)+806;
x2=-0.0027.*x(i)^3+1.5.*x(i)^2-200.*x(i)+9200;
x3=-5.6.*10^(-8).*x(i)^3+0.00076.*x(i)^2+1.3.*x(i)+1600; x4=[101.2 103.4 98.6 100.5 106.7]; x5=0.64.*x(i)-1800;
x6=[35.8 36.3 37.9 36.5 35.7 ]; end
y=-206.055+14.1299.*x1+0.6997.*x2+0.2806.*x3+27.6024.*x4-0.1361.*x5-31.5573.*x6
(4)时间序列程序:
clc,clear
y=[1309 1614 1620 2094 2537 2900 6000 6944 8100 9399.7 10815]; m1=length(y);
n=2; %n 为移动平均的项数 for i=1:m1-n+1
yhat1(i)=sum(y(i:i+n-1))/n; end yhat1
m2=length(yhat1); for i=1:m2-n+1
yhat2(i)=sum(yhat1(i:i+n-1))/n;
24
3270 3391 4089 5058