数学建模社区-数学中国
标题:
matlab 灰色系统预测 GM(1,1) 数学建模
[打印本页]
作者:
佛自业障
时间:
2018-10-31 09:23
标题:
matlab 灰色系统预测 GM(1,1) 数学建模
本文代码主要是基于邓聚龙教授在20实际80年代提出的灰色系统理论。
: h4 }4 E- \) m1 F
GM0.m
0 X4 E' y& c* A4 O' O( l
%该函数为GM（1，1）模型返回还原值
9 V v; a7 h3 |1 V/ M1 e
function f=GM0(x0,t) %数据数列
7 M# a3 W, Z1 H' }
[M,N]=size(x0); %算出数据数列的大小
0 ]( ?. y% G4 E. i d
x1(1)=x0(1); %累加生成数列
% M: a) Y; q* y! E5 V9 V9 r& Z: t5 W
for i=2:N;
8 x! \+ K- C/ O& x
x1(i)=x1(i-1)+x0(i);
7 }9 K2 f8 h( [# p1 o( m8 E& p
end
4 N9 v- B" d! S
x2=[]; %累加生成数列均值生成数列
' m! S# U. U/ |
for j=1
N-1);
: I% D/ ^9 v( L K6 s* t% G
x2(j)=(x1(j)+x1(j+1))/2;
6 Q0 _$ |" d# ?6 d8 Q/ R4 h
end
* e' k4 [3 j/ m# }( U/ O: J3 I
x=x0; %数据数列镜像
3 V" ^. w8 |% s" l8 n- t! u
x(1)=[]; %删除第一个数据
3 R3 C, ^, T; v* I
Y=x'; %数据列向量
# O0 d/ [1 U, S
global a;
4 S; J; A- D* [0 ?
global b;
7 ~* G# {0 _' @9 y$ x* x
B(:,1)=-x2';
( v! i* s5 z, {* c3 ~
B(:,2)=1;
g) Z7 C: {! w. q) L- P* w
A=inv(B'*B)*B'*Y; %求参量a，b组成的参数向量
$ i* Q% N+ L8 [* ~
a=A(1,1); %求参数a
) z" ]3 e) D' k/ v- D6 _
b=A(2,1); %求参数b
5 [3 @" ~$ g, g% Q: R$ j$ H6 J9 r
f=(1-exp(a))*(x0(1)-b/a)*exp(-a*(t-1));
+ X/ }+ [/ y8 C" @3 }) V6 ^, J/ t1 g
f
' n( w+ r8 @; h# M. j1 S
( J1 q' |, ]" g% L' k+ @
GM1.m:
. T7 D; b9 V$ K" G5 i
%该函数为GM（1，1）模型中数据数列进行光滑比检验
7 c# ~1 e( W0 `& s; V z4 X. t* s$ {
function f=GM1(x0) %数据数列
8 W, P5 W8 K2 G& E3 t* B
N=max(size(x0)); %算出数据数列的大小
& n' `# ?1 u F! u+ }
x1=cumsum(x0); %累加生成数列
+ V2 f" C' c/ j' @- z8 |9 @
global J;
! O# E5 _/ x9 r9 f
global J1;
4 Z4 k* ~/ [3 K4 G: f/ T
global J2;
! D/ ]8 S- g: j
x0(1)=[];
) Z" T, l! [2 @ C1 i
x1(N)=[];
: t& t# N* x0 g! v
global r;
. ~6 V/ _, D7 x5 \! d& @* u. C& q9 `2 i9 W
r=x0./x1;
/ J v8 W" @: g) L8 o
for j=2
N-1); %判断数据数列是否满足准光滑条件1
* i V# f9 _ @' V: _2 `+ L( Q
if(r(j)>=0.5||r(j)<0)
4 {1 [4 Z+ E( v. [) N: D( c0 m1 {8 d
J1=0;
8 F- g3 i6 F" E$ O5 Q: Z
break;
{) B/ M* Q, ?# J/ n: P
else
# u4 I" \: \9 a! I
J1=1;
1 E( ] L0 ]+ u5 }4 }
end
( g! }- Q+ H& ^& k/ p
end
1 O3 O" Q$ h: P# T
for l=1
N-2); %判断数据数列是否满足准光滑条件2
9 |0 S. c& e# X4 v
if((r(l+1)/r(l))>=1)
) J9 u5 P: i3 Q
J2=0;
4 k# d% a' t: J; Z H( z8 z
break;
, ]( ~+ d% m& H. W# t% \$ U+ B
else
7 B$ {4 F: Z; N/ R6 u3 d7 F% Y
J2=1;
& s# R" L$ y% E8 Y# c" @
end
. `% a) m9 \2 ^3 z+ r
end
# ^+ h6 L; T3 G0 ]9 D2 x
J=J1+J2;
1 f G* z: M y$ h4 U+ L% N% x
if(J==2) %判断数据数列是否为准光滑数列
- u m! ~% U% K. a6 c
disp('数据为准光滑数列')
. A! {* A- ^3 |) a! Y# u
else
( M3 ?: @& G7 t6 \0 w6 M- j5 O
disp('数据不是准光滑数列')
N1 q2 n: R1 V/ }. M
end
, c4 Z1 B8 p) M7 J- c+ U
) P. z0 k! Y. a8 K
GM2.m
" ]) r. k+ t* w1 C# v
%该函数为GM（1，1）模型还原值参数计算
* B8 f2 @4 m1 G2 v
function f=GM2(x0) %数据数列
; Y/ {& J) Q, j! e- L& _
[M,N]=size(x0); %算出数据数列的大小
8 l" [7 G8 ^; n
x1(1)=x0(1); %累加生成数列
8 u0 h. {; {& i6 A4 u- ~8 h5 C
for i=2:N;
9 W1 k( G/ n5 }: ^
x1(i)=x1(i-1)+x0(i);
3 c. e. `6 A. S, A
end
9 c5 J L& |7 ~. G, ^7 c
x2=[]; %累加生成数列均值生成数列
$ ]4 b8 E5 X& z F
for j=1
N-1);
3 y0 O, s$ K" ^4 j
x2(j)=(x1(j)+x1(j+1))/2;
' V+ _. Y# }1 L1 T0 v
end
- B3 U4 o# A. p0 n- e1 g' P& q
x=x0; %数据数列镜像
1 j6 ]9 n2 V# \4 V% Y
x(1)=[]; %删除第一个数据
W7 ^' R4 g% ]$ ?% i
Y=x'; %数据列向量
8 x; v/ J. t. W( Q; J. m
global a;
. P: l7 e; t' R- \
global b;
: O! {0 @- H% p" t" k# R
B(:,1)=-x2';
5 A7 J( Z {3 b8 V( E5 g; C
B(:,2)=1;
4 O2 _6 ?( A, }8 |3 u5 K2 Z- ]
A=inv(B'*B)*B'*Y; %求参量a，b组成的参数向量
1 `, |2 B) @5 O) C4 a# W1 B+ ~ C
a=A(1,1); %求参数a
) A1 ^: m; n. C0 }" q
disp('参数a为：')
3 R9 `' F1 I$ S4 I+ M) B# J
a
' H2 t4 W" N) n. j
b=A(2,1); %求参数b
9 @- b$ T- ~. I0 Q. D0 ?% [& m
disp('参数b为：')
& @8 M2 D2 ?1 a* \7 w$ e
b
4 i5 d" F: }% { x
- h5 Z7 ?) X1 a% M. L
GM3.m
$ t3 ?+ `* U' x; w
%该程序实现G(1,1)模型的精度检验
( q% |1 X1 r; D- z- ~6 X: P
%包括平均相对误差，绝对关联度，均方差比值，小误差概率检验
+ |7 [4 l9 S3 \
function f=GM3(x0)
; w c1 u3 l e: B1 J i
N=max(size(x0));
+ m: X; U# Q# I7 \, @
x=GM0(x0,1:N); %利用已有程序GM2得出数据列模型估计值
6 T+ ?& C8 V# X9 T2 l; \7 z" w) N* P
x(1)=x0(1); %更正第一个估计值
) P" Q/ V# y! M7 V) W5 F
disp('模型模拟估计值为')
8 R6 c c: t5 s
x
+ A3 _; {: k* {$ O- d* O1 f5 S
A=x-x0; %计算绝对残差序列
" }, C% G' P; }* Z; F
disp('模型估计值绝对残差序列为：')
' ^! \2 W: f, M" Y' J/ |
A
: l7 H9 _, q8 p% {$ l; J2 ]
G=abs(A);
- q5 m D- B( h- f$ ?
Amin=min(G); %计算最小绝对值绝对残差
5 n2 t( K9 C O8 v* u5 C }
Amax=max(G); %计算最大绝对值绝对残差
7 g7 X1 ?/ { F0 f1 w0 d4 b" S
B=A./x0; %计算相对误差序列
& ~, c* J0 z( l9 D
disp('模型估计值相对误差序列为：')
+ e0 M# k# ~) R% @
B
( t& u- r/ }3 U
P=sum(abs(B))/N; %计算平均相对误差
7 x4 E. z- O. k) H
disp('模型估计值平均相对误差为：')
5 ]3 n6 u* ?3 n0 R- V" H$ r' M
P
! B' }& s+ u/ F+ J& D* g
for i=1:1:N %通过循环计算关联系数序列
7 F% ~6 v4 y2 t v+ ^
D(i)=(Amin+0.5*Amax)/(G(i)+0.5*Amax);
: G6 j) V1 c1 s- d" V
end
8 M. K7 R5 n# X& x) m: [; E/ a
R=sum(D)/N;
' t. O7 k7 A* R6 v' ? H) c
disp('关联度为：')
# M/ Z# S, ^. o r8 a; O
R
6 `9 @7 n/ w4 s- t o+ G, R
x_=sum(x0)/N; %计算数据的均值
6 }/ p. B) t+ ~( y
S1=(sum((x0-x_).^2)/(N-1))^0.5; %计算数据序列方均差
" d) d o% Y" I
A_=sum(A)/N; %计算残差平均值
7 z' N1 d. y# u% V7 J% F$ T/ r& h
S2=(sum((A-A_).^2)/(N-1))^0.5; %计算残差序列方均差
! _, v: u4 B l6 E' Z. |4 x
C=S2/S1; %计算方均差比值
0 i2 x! J% }! W0 n1 \' y2 ?1 X
disp('均方差比值为：')
* \0 x m% c8 p* g4 Z
C
( O# |. p( _0 L
S0=0.6745*S1;
! i( D% i/ T& Z* V; a# a
E=A-A_;
2 P. J Y2 d) ~+ \9 v/ h0 ~ [7 I
F=find(E<S0);
( v9 k- ~0 |* U8 ]3 I; D( E% O* `8 b2 o
M=max(size(F)); %计算小残差个数
- ^! Z H" L( K" P% c5 }
p=M/N; %计算小误差概率
: J$ w9 b% j% l9 w: H
disp('小误差概率为：')
* O$ C% W4 r5 n% V% G+ [7 e4 T, i3 u
p
P3 u* w( h: m4 l+ D& V Q# V
- \6 M3 W6 C5 }" O8 h
) B( o! d3 b" M& X* s) k+ M
8 t' G# ^$ B4 l1 z( u
欢迎光临 数学建模社区-数学中国 (http://www.madio.net/)
Powered by Discuz! X2.5
N-1); N-1); %判断数据数列是否满足准光滑条件1 N-2); %判断数据数列是否满足准光滑条件2 9 |0 S. c& e# X4 vN-1);