灰色预测模型GM

灰色预测模型GM(1,1)§1预备知识灰色预测是就灰色系统所做的预测。所谓灰色系统是介于白色系统和黑箱系统之间的过渡系统，其具体的含义是：如果某一系统的全部信息已知为白色系统，全部信息未知为黑箱系统，部分信息已知，部分信息未知，那么这一系统就是灰箱系统。一般地说，社会系统、经济系统、生态系统都是灰色系统。例如物价系统，导致物价上涨的因素很多，但已知的却不多，因此对物价这一灰色系统的预测可以用灰色预测方法。灰色系统理论认为对既含有已知信息又含有未知或非确定信息的系统进行预测，就是对在一定方位内变化的、与时间有关的灰色过程的预测。尽管过程中所显示的现象是随机的、杂乱无章的，但毕竟是有序的、有界的，因此这一数据集合具备潜在的规律，灰色预测就是利用这种规律建立灰色模型对灰色系统进行预测。平面上有数据序列nnyxyxyx,,,,,,2211，大致分布在一条直线上。设回归直线为：baxy，要使所有点到直线的距离之和最小（最小二乘），即使误差平方和niiibaxyJ12最小。J是关于a,b的二元函数。由0120211niiiiniiiiibxaybJxbxayaJ00112niiiniiiiibaybxaxyx则得使J取极小的必要条件为：iiiiniiiynbxayxxbxa12（*）22222iiiiiiiiiiiiixxnyxxxybxxnyxyxna(1)以上是我们熟悉的最小二乘计算过程。下面提一种观点，上述算法，本质上是用实际观测数据ix、iy去表示a与b，使得误差平方和J取最小值，即从近似方程bbbxxxayyynn2121中形式上解出a与b。把上式写成矩阵方程。yixxiiyx,jjyx,令nyyyY21，baxxxYn11121令11121nxxxB，则baBY左乘TB得baBBYBTT注意到BTB是二阶方阵，且其行列式不为零，故其逆阵(BTB)-1存在，所以上式左乘1BBT得YBBBbaTT1（2）可以具体验算按最小二乘法求得的结果（1）与（2）式完全相同，下面把两种算法统一一下：由最小二乘得结果：方程（*）iiiiniiiynbxayxxbxa12方程组改写为：nniiiyyyxxxbanxxx21212111令：11121nxxxB，nyyyY21，baaˆ（*）化为YBaBBTTˆ所以YBBBaTT1ˆ以后，只要数据列njyxjj,,2,1,大致成直线，既有近似表达式nibaxyii,,2,1当令：nyyyY21，11121nxxxB，baaˆ则有aBYˆyBBBaTT1ˆ（2）（2）式就是最小二乘结果，即按最小二乘法求出的回归直线baxy的回归系数a与b。推广：多元线性回归设有m个变量mxxx,,,21，每个自变量有n个值，因变量y有n个值mnmnnnmmmmxbxbxbayxbxbxbayxbxbxbayn2211222212121212111121（1）如n个人，每人有m个指标。女生：人：1x（体重）公斤2x（胸围）厘米3x（呼吸差）厘米ky（肺活量）毫升111x=3521x=6931x=0.71600212x=4022x=7432x=2.52600313x=4023x=6433x=2.02100414x=4224x=7434x=32650515x=3725x=7235x=1012400616x=4526x=6836x=1052200717x=4327x=7837x=4032750818x=3728x=6638x=21600919x=4429x=7039x=30227501010x=4220x=6530x=32500方程组（1）是n个方程m个数据mmnnnmmbbbaxxxxxxxxxY21212221211211111用X表示增广矩阵：n行，m+1列baXYˆ，mbbbb21ˆ，baXXYXTTˆYXXXbaTT1ˆ其中XXT为11mm阶矩阵。由此可解出：mbbba,,,,21注意：方程组中mbbba,,,,21不知，意思是：如果线性关系成立mmxbxbxbay2211当mbbba,,,,21为多少时，iy到mmxbxbxba2211的距离之和为最小。或说，当所有iy到(mmxbxbxba2211)距离之和为最小时的mbbba,,,,21就是我们要求的最佳系数。§2前言为什么要讲GM(1,1)模型？80年代初，华中理工大学邓聚龙教授提出了灰色系统理论，先后发表过灰色控制、灰色预测、灰色决策、灰色系统理论等多部专著，较详细在阐述了灰色系统理论的产生、理论、方法与应用。在80年代中后期到90年代初，举行了十数次国际、国内有关灰色系统理论的研讨会，在全国形成一股灰色系统理论研究与应用热潮。邓聚龙先生因灰色系统理论方面的供献，获得国家科技进步一等奖。什么叫灰？用邓先生自己的话来讲：“完全已知的系统称作白系统；完全未知的系统称作黑系统或黑箱；部分已知、部分未知的系统称作灰色系统。”在此，已知或未知到什么程度没有具体说明。所以，“灰”的内涵不是很清楚。举个例子讲，已知某量的真值x在闭区间[a,b]上，不可能落在[a,b]之外，但具体落到区间[a,b]的什么位置则是完全不知道的。那么，这个量称作灰量，可具体表示为[a,b]，称其为区间灰数。显然，区间灰数是客观实际中存在的，除了知道真值x在[a,b]上，而不在[a,b]之外，不再有任何已知信息，这就是灰量的最基本原型。由于灰色系统理论从一开始就没有建立在严格的集合论基础之上，使之缺乏必要的数学支撑，这大大限制了灰色系统理论和应用的发展。虽然灰色系统理论在控制、预测、决策等领域有着广泛的应用；但就其精华而言，还在于GM(1,1)模型。即便是现在，在特定情况下，GM(1,1)还有用，还在被应用，并且预测效果很好。其使用限制条件是：原始数据单调，预测背景呈现稳定发展趋势；其优势是：适用于原始观测数据较少的预测问题，由于数据量很小，无法应用概率统计方法寻找统计规律，而GM(1,1)模型恰恰弥补了这个空白，由于GM(1,1)算法简单易行，预测精度相对较高，所以在一些特定问题中，GM(1,1)仍然是决策者乐于选择的预测模型。上面讲到的背景稳定的发展趋势是指下述情况：如化工设备的腐蚀量，随着使用时间的推移腐蚀不断增加，呈现出稳定的发展趋势，并且腐蚀量的测量通常比较困难(如停产才能测量)，所以实际观测数据较少。这类问题很适合GM(1,1)模型预测。§3GM(1,1)预备知识(二)3.1回忆一阶线性常系数微分方程uaxdtdx(1)其解为：~aueauxxat)0((2)其中a，u为给定的常数。一阶线性常系数微分方程(1)的解(2)是指数型曲线，如下图所示atex图象ateauxx0图象aueauxxat0图象3.2在预备知识中，讲述了最小二乘法：若数据点)(iiyx，ni,,2,1近似落在一条直线上，设这条直线为y＝ax+b，a,b为参数。理想的直线要求：每个数据点)(iiyx，ni,,2,1，到该直线的距离平方和最小――即最小二乘。用最小二乘法求出参数a与b，这相当于形式上的解线性方程组：baxyiini,,2,1(3)当令nyyyy21，11121nxxxB，baaˆ则（3）化为aBYˆ，aBBYBTTˆYBBBaTT1ˆ（4）由此求出abaˆ，可得回归直线baxy（5）上述形式上的求解结果，本质上是用最小二乘法求解回归参数的过程，故有下面结论。结论：一组数据点（n个），且近似线性关系10txa0a0x(0)–u/a0txa0a0x(0)0txa0a0baxyii则下述表达式可求出回归系数a与b。nnTTyyyYxxxBYBBBba21211111，：上述形式上的计算，本质是使点),(iiyx到直线y＝ax＋b的距离平方和最小，即是最小二乘法得来的结果。§4GM(1，1)模型G表示Grey(灰)，M表示Model(模型)，前一个“1”表示一阶，后一个“1”表示一个变量，GM(1,1)则是一阶，一个变量的微分方程模型。给定等时间间隔的数据列，且设数据列单调：),n(),2(),1()(,n21xxxkxk，k表示时刻，kxkx＝)(表示t＝k时刻某量的观测值，不妨设1kkxx，1,,2,1nk，将数据列记成：0n030201)0(,,xxxxx)0(x表示原始数据序列。比如：697.3,390.3,337.3,278.3,874.2)0(x。)1(kx={2.874,6.152,9.489,12.879,16.558}对原始数据作一次累加生成：即令(1)(0)1()()1,2,,kixkxkkn得一次累加生成数序列为：)1()1(21(1)1(,,,nxxxx）2.874,3.278,3.337,3.390,3.697在此，)1(kx={2.874,6.152,9.489,12.879,16.558}给定的原始数据序列)0(kx已经是单增序列，经一次累加后生成的累加数序列具有更强烈的单调性。我们知道指数序列是单调的，但是，单调序列却不一定是指数型的，不过强烈的单调序列可近似看做是指数的，即可用指数型曲线进行弥合。如果用指数曲线来弥合一次累加生成序列，那么，这条指数曲线一定是某个一阶线性常系数微分方程uaxdtdx)1()1((6)的满足某个初始条件的一条积分曲线：aueauxxat)1()1()1(即aueauxxat)1()0()1(（7）其中a，u是待确定的未知参数，该微分方程中的导数dtdx)1(可用差商近似表示。tkxtkxdtdxt)1()1(0)1(limt为时间间隔，将时间间隔t看做是单位时间间隔，并且认为时间被充分细化(秒，毫秒。微秒……事实上只要单位时间内函数的增量相对很小，这个单位时间间隔也可以是日，月，年等。)此时有kxkxdtdx)1()1()1(1注意到一次累加生成数)()1(tx在时刻t＝k＋1与t＝k时的差为：11)0()1()1(kxkxkx而)1(dtdx是在[k,k＋1]上某一点取值，既然是近似，索性将dtdx)1(的值取在点k+1，即11)0()1()1(1)1(