GM模型及其Matlab實(shí)現(xiàn)

1、灰色系統(tǒng)預(yù)測(cè) GM(1,1) 模型及其 Matlab 實(shí)現(xiàn)三天三夜 72 小時(shí)：讀懂題目 -查找文獻(xiàn)資料 - 選擇題目 - 重查找文獻(xiàn)資料 - 精讀其中幾篇 - 查找資料的資料。 。在數(shù)學(xué)建模中常常會(huì)遇到數(shù)據(jù)的預(yù)測(cè)問(wèn)題，有些賽題中，預(yù)測(cè) 占主導(dǎo)地位，例如：2003 年A題SARS的傳播問(wèn)題；2005 年A題 長(zhǎng)江水質(zhì)的評(píng)價(jià)和預(yù)測(cè)問(wèn)題；2006 年B題 艾滋病療法的評(píng)價(jià)及療效的預(yù)測(cè)問(wèn)題；2007 年A題 中國(guó)人口增長(zhǎng)預(yù)測(cè)問(wèn)題。有些問(wèn)題則是需要在求解的過(guò)程中進(jìn)行預(yù)測(cè)，如2009年D題“會(huì)議籌備”對(duì)與會(huì)人數(shù)的確定等。參考資料：灰色系統(tǒng)理論及其應(yīng)用 第五版 作者：劉思峰，黨耀國(guó)等著 出 版時(shí)間： 20

2、10.05 校超星數(shù)字圖書(shū)館可閱讀?；疑Ｐ?Gray Model)有嚴(yán)格的理論基礎(chǔ)，最大優(yōu)點(diǎn)是實(shí)用。 用灰色模型預(yù)測(cè)的結(jié)果比較穩(wěn)定， 不僅適用于大數(shù)據(jù)量的預(yù)測(cè)， 在數(shù) 據(jù)量較少時(shí) (3) 預(yù)測(cè)結(jié)果依然較準(zhǔn)確。預(yù)備知識(shí)(1) 灰色系統(tǒng)白色系統(tǒng)是指系統(tǒng)內(nèi)部特征是完全已知的 , 即人們不僅知道該系 統(tǒng)的輸入輸出關(guān)系， 而且知道實(shí)現(xiàn)輸入輸出關(guān)系的結(jié)構(gòu)與過(guò) 程；黑色系統(tǒng)是指系統(tǒng)內(nèi)部信息完全未知的 , 即人們只知道該系統(tǒng)輸 入輸出關(guān)系， 但不知道實(shí)現(xiàn)輸入輸出關(guān)系的結(jié)構(gòu)與過(guò)程； 而 灰色系統(tǒng)是介于白色系統(tǒng)和黑色系統(tǒng)之間的一種系統(tǒng)， 灰色系統(tǒng)其內(nèi) 部一部分信息已知，另一部分信息未知或不確定。例如，一個(gè)加有電

3、壓的電阻，也是一個(gè)系統(tǒng)，根據(jù)歐姆定律，匸U/R, 當(dāng)電阻的大小知道后， 便可由多大電壓算出能得到多大電流。 電壓與 電流之間有明確的關(guān)系或函數(shù)，這便是白色系統(tǒng)。因此，這樣的系統(tǒng) 要求有明確的作用原理， 一個(gè)有明確作用原理的系統(tǒng)必定是具有確定 結(jié)構(gòu)的，必定是有物理原型的。 然而許多社會(huì)經(jīng)濟(jì)系統(tǒng)都沒(méi)有物理原 型，雖然知道影響系統(tǒng)的某些因素，但很難明確全部因素，更不可能 確定因素之間的映射關(guān)系。 這種沒(méi)有確定的映射關(guān)系 (函數(shù)關(guān)系 ) 的系 統(tǒng)是灰色系統(tǒng)。(2) 灰色預(yù)測(cè)灰色預(yù)測(cè)， 是指對(duì)系統(tǒng)行為特征值的發(fā)展變化進(jìn)行的預(yù)測(cè)， 對(duì)既 含有已知信息又含有不確定信息的系統(tǒng)進(jìn)行的預(yù)測(cè)， 也就是對(duì)在一定 范圍

4、內(nèi)變化的、 與時(shí)間序列有關(guān)的灰過(guò)程進(jìn)行預(yù)測(cè)。 盡管灰過(guò)程中所 顯示的現(xiàn)象是隨機(jī)的、雜亂無(wú)章的，但畢竟是有序的、有界的，因此 得到的數(shù)據(jù)集合具備潛在的規(guī)律。 灰色預(yù)測(cè)是利用這種規(guī)律建立灰色 模型對(duì)灰色系統(tǒng)進(jìn)行預(yù)測(cè)。目前使用最廣泛的灰色預(yù)測(cè)模型就是關(guān)于數(shù)列預(yù)測(cè)的一個(gè)變量、一階微分的GM(1,1)模型。它是基于隨機(jī)的原始時(shí)間序列，經(jīng)按 時(shí)間累加后所形成的新的時(shí)間序列呈現(xiàn)的規(guī)律可用一階線(xiàn)性微分方 程的解來(lái)逼近。 經(jīng)證明， 經(jīng)一階線(xiàn)性微分方程的解逼近所揭示的原始 時(shí)間序列呈 指數(shù)變化規(guī)律。因此，當(dāng)原始時(shí)間序列隱含著指數(shù)變化規(guī) 律時(shí)，灰色模型GM(1,1)的預(yù)測(cè)是非常成功的。目前，灰色模型GM(1,1)己

5、廣泛應(yīng)用于工程技術(shù)、社會(huì)、經(jīng)濟(jì)、 農(nóng)業(yè)、生態(tài)、環(huán)境等各種系統(tǒng)的預(yù)測(cè)中。1 灰色預(yù)測(cè)基礎(chǔ)知識(shí)灰色系統(tǒng)理論認(rèn)為 :系統(tǒng)的行為現(xiàn)象盡管是朦朧的 ,數(shù)據(jù)是復(fù)雜 的,但它畢竟是有序的 ,是有整體功能的。在建立灰色預(yù)測(cè)模型之前 , 需先對(duì)原始時(shí)間序列進(jìn)行數(shù)據(jù)處理 , 經(jīng)過(guò)數(shù)據(jù)預(yù)處理后的數(shù)據(jù)序列稱(chēng) 為生成列。對(duì)原始數(shù)據(jù)進(jìn)行預(yù)處理 , 不是尋求它的統(tǒng)計(jì)規(guī)律和概率分 布, 而是將雜亂無(wú)章的原始數(shù)據(jù)列通過(guò)一定的方法處理, 變成有規(guī)律的時(shí)間序列數(shù)據(jù) , 即以數(shù)找數(shù)的規(guī)律 , 再建立動(dòng)態(tài)模型。 灰色系統(tǒng)常用 的數(shù)據(jù)處理方式有累加和累減兩種 , 通常用累加方法。灰色預(yù)測(cè)通過(guò)鑒別系統(tǒng)因素之間發(fā)展趨勢(shì)的相異程度 ,并對(duì)原始

6、 數(shù)據(jù)進(jìn)行生成處理來(lái)尋找系統(tǒng)變動(dòng)的規(guī)律 , 生成有較強(qiáng)規(guī)律性的數(shù)據(jù) 序列, 然后建立相應(yīng)的微分方程模型 , 從而預(yù)測(cè)事物的未來(lái)發(fā)展趨勢(shì)。 灰色預(yù)測(cè)的數(shù)據(jù)是通過(guò)生成數(shù)據(jù)的模型所得到的預(yù)測(cè)值的逆處理結(jié) 果?；疑A(yù)測(cè)是以灰色模型為基礎(chǔ)的 , 在諸多的灰色模型中 , 以灰色系 統(tǒng)中單序列一階線(xiàn)性微分方程模型 GM(1,1 )模型最為常用。下面簡(jiǎn)要 地介紹GM(1,1)模型。設(shè)有原始數(shù)據(jù)列x(0)=(x(0)(1),x (0)(2),x (0) (n),n 為數(shù)據(jù)個(gè)數(shù)。如果根據(jù)x(0)數(shù)據(jù)列建立GM(1,1)來(lái)實(shí)現(xiàn)預(yù)測(cè)功能，則基本步驟 如下:(1) 原始數(shù)據(jù)累加以便弱化隨機(jī)序列的波動(dòng)性和隨機(jī)性 ,得到新

7、數(shù)據(jù)序列：x(1)= X (1)(1),x (1)(2),x(1)(n)其中,x (1) (t)中各數(shù)據(jù)表示對(duì)應(yīng)前幾項(xiàng)數(shù)據(jù)的累加。tx二-x(0)(k),t=1,2,nk 二 對(duì)x 建立下述一階線(xiàn)性微分方程：即GM(1,1)模型。dx(1)+ ax(1) = u dt16 / 15其中,a,u為待定系數(shù)，分別稱(chēng)為發(fā)展系數(shù)和灰色作用量,a的有效區(qū)間是(-2,2),并記a,u構(gòu)成的矩陣為灰參數(shù)a二勺1只要求出參數(shù) 2丿a,u,就能求出x，進(jìn)而求出x(0)的未來(lái)預(yù)測(cè)值。(3) 對(duì)累加生成數(shù)據(jù)做均值生成 B與常數(shù)項(xiàng)向量Yn,即1(X (1) + X (2)-2(X(1)(2) X(1)(3)匸 2(X

8、 (n -1) +X (n)1x(0)(2)Yn二x(0)(1)(t)(7) 對(duì)建立的灰色模型進(jìn)行檢驗(yàn),步驟如下：計(jì)算x(0) (t)與?(0) (t)之間的殘差e(0) (t)和相對(duì)誤差q(0) (t):e(0) (t)=x (0) (t)- ?(0) (t)(0) (0) (0)q (t)=e(t)/x(t)等等。(8) 利用模型進(jìn)行預(yù)測(cè)：?(0) = x(0)(1), ?(0)(2),?(0)(n ),。0)( n+1),。0)( n+m)原數(shù)列的模擬未來(lái)數(shù)列的預(yù)測(cè)應(yīng)用舉例取某高校2000年2005年的某專(zhuān)業(yè)招生數(shù)據(jù)建模，見(jiàn)表1。表1某高校專(zhuān)業(yè)招生數(shù)據(jù)表年招生人數(shù)200013220019

9、22002118200313020041872005207以表1中的數(shù)據(jù)構(gòu)造原始數(shù)據(jù)列 於0),即乂。)=X(0)(1),X (0) (2),X (0) (3),X (0)(4),X (5),X (0) (6)=132,92,118,130,187,207對(duì)X(0)進(jìn)行一次累加(1 AGO)生成數(shù)列：kZX(1) (k) = v X0) (i)即疋=X(1),X (2),X (1)(3),X (4),X (5),X (6)=132,224,342,472,659,866-178-28311-4071-565 .51B=】一 762 .51經(jīng)計(jì)算可得和數(shù)據(jù)陣B、數(shù)據(jù)列YnTAa =a,u,Yn=

10、 (92,118,130,187,207)1-0.205=|56 .7878184進(jìn)一步得到灰色預(yù)測(cè)模型GM(1,1)為uuX (1) (k) = (X(0)(1) - a) e) + a=(132 + 277.0137483) e0.205(kJ) 277.01374830.205( kJ)=409.0137483 e 277.0137483預(yù)測(cè)值及預(yù)測(cè)精度見(jiàn)表2。表2某高校專(zhuān)業(yè)招生預(yù)測(cè)值及預(yù)測(cè)精度表年GM(1,1)模型計(jì)算值1 AGO還 原 值實(shí) 際 值誤 差擬合相對(duì) 誤差(%)2000132132132132002001225.060879622493921-12002339.29544

11、1834211411843.382003479.52123472140130-10-7.692004651.65196591721871582005862.9466129866211207-4-1.9320061122.316167259252-7-2.78由表2知預(yù)測(cè)精度較高。2006年某專(zhuān)業(yè)招生人數(shù)預(yù)測(cè)值為259人。由于人數(shù)為整數(shù)，所以結(jié)果取整數(shù)部分。GM(1,1)也是一種長(zhǎng)期預(yù)測(cè)模型， 在沒(méi)有大的市場(chǎng)波動(dòng)及政策 性變化的前提下，該預(yù)測(cè)值應(yīng)是可信的。眾所周知，影響招生人 數(shù)的因素很多且難以預(yù)測(cè)。因此，在采用灰色系統(tǒng)理論進(jìn)行定量 預(yù)測(cè)時(shí)，如果存在對(duì)預(yù)測(cè)對(duì)象影響較大的因素，就要在定性分析 的基

12、礎(chǔ)上，尋找原始數(shù)據(jù)信息的突變點(diǎn)的量化值，然后再對(duì)預(yù)測(cè) 值進(jìn)行必要的修正，使預(yù)測(cè)值更接近實(shí)際情況，提高預(yù)測(cè)值的可 信度，為科學(xué)決策提供可靠的數(shù)據(jù)。另外，若作長(zhǎng)期預(yù)測(cè)，要考 慮對(duì)上限值的約束條件。2灰色預(yù)測(cè)的MATLABS序2.1典型程序結(jié)構(gòu)灰色預(yù)測(cè)中有很多關(guān)于矩陣的運(yùn)算，這可是MATLAB勺特長(zhǎng)，所以 用MATLAB是實(shí)現(xiàn)灰色預(yù)測(cè)過(guò)程的首選。用 MATLAB編寫(xiě)灰色預(yù)測(cè)程 序時(shí),可以完全按照預(yù)測(cè)模型的求解步驟，即(1) 對(duì)原始數(shù)據(jù)進(jìn)行累加。(2) 構(gòu)造累加矩陣B與常數(shù)向量Yn。(3) 求解灰參數(shù)。(4) 將參數(shù)帶入預(yù)測(cè)模型進(jìn)行數(shù)據(jù)預(yù)測(cè)。下面以某公司收入預(yù)測(cè)問(wèn)題為例介紹灰色預(yù)測(cè)的MATLAB現(xiàn)過(guò)程

13、。已知某公司 1999 2008年 的 利 潤(rùn) 為 ( 單 位:元/ 年): 89677,99215,109655,120333,135823,159878,182321,209407,246619.0 , 現(xiàn)在要預(yù)測(cè)該公司未來(lái)幾年的利潤(rùn)情況。具體的MATLABS序如下：clearsyms a u;c=a u; % 灰參數(shù) cA=89677,99215,109655,120333,135823,159878,182321,209407,246619.0 ; % 原始數(shù)據(jù)Ago=cumsum(A); %原始數(shù)據(jù)一次累加n=length(A); % 原始數(shù)據(jù)個(gè)數(shù)for i=1:(n-1)C(i)=

14、(Ago(i)+Ago(i+1)/2; % 生成累加矩陣end% 計(jì)算待定參數(shù)的值Yn=A;Yn(1)=;Yn=Yn;E=-C;ones(1,n-1);c=inv(E*E)*E*Yn;c=c; a=c(1);u=c(2);% 預(yù)測(cè)后續(xù)數(shù)據(jù)F=;F(1)=A(1);for i=2:(n+10)F(i)=(A(1)-u/a)/exp(a*(i-1)+u/a ;endG=;G(1)=A(1);for i=2:(n+10)G(i)=F(i)-F(i-1); %得到預(yù)測(cè)出來(lái)的數(shù)據(jù)endt1=1999:2008;t2=1999:2018;Gplot(t1,A,o,t2,G) % 原始數(shù)據(jù)與預(yù)測(cè)數(shù)據(jù)的比較

15、運(yùn)行該程序 , 得到的預(yù)測(cè)數(shù)據(jù)如下 :G =1.0e+006 * Columns1through140.0897 0.0893 0.1034 0.1196 0.1385 0.1602 0.18540.2146 0.2483 0.2873 0.3325 0.3847 0.4452 0.5152 Columns15through200.5962 0.6899 0.7984 0.9239 1.0691 1.2371該程序還顯示了預(yù)測(cè)數(shù)據(jù)與原始數(shù)據(jù)的比較圖514121086420x 10199820002002200420062008201020122014201620182.2灰色預(yù)測(cè)程序說(shuō)明(1)

16、 先熟悉程序中各條命令的功能，以加深對(duì)灰色預(yù)測(cè)理論的理解。(2) 在實(shí)際使用時(shí)，可以直接套用該段程序，把原數(shù)據(jù)和時(shí)間序列數(shù)據(jù) 替換就可以了。(3) 模型的誤差檢驗(yàn)可以靈活處理，圖中給出的是預(yù)測(cè)數(shù)據(jù)與原始數(shù) 據(jù)的比較圖，同樣也可以對(duì)預(yù)測(cè)數(shù)據(jù)進(jìn)行其他方式的精度檢驗(yàn)。3灰色預(yù)測(cè)應(yīng)用實(shí)例3.1實(shí)例一：長(zhǎng)江水質(zhì)的預(yù)測(cè)(CUMCM 2005A)長(zhǎng)江的水質(zhì)問(wèn)題是一個(gè)復(fù)雜的非線(xiàn)性系統(tǒng)，但是由于數(shù)據(jù)樣本少 需要預(yù)測(cè)的時(shí)間長(zhǎng)，直接應(yīng)用神經(jīng)網(wǎng)絡(luò)很難取得理想的效果?？紤]到 污水排放量的變化規(guī)律是一個(gè)不確定的系統(tǒng)，且題中給出污水排放量 數(shù)據(jù)樣本比較少，還要求做出長(zhǎng)達(dá)10年的預(yù)測(cè)，因此采用灰色預(yù)測(cè)方 法來(lái)預(yù)測(cè)未來(lái)的污水排

17、放量。對(duì)原題附件4中的數(shù)據(jù)進(jìn)行整理可以得到10年的長(zhǎng)江污水量排 放數(shù)據(jù),如下表所列。1995 2004年長(zhǎng)江污水排放量年份1995199619971998199920002001200220032004污水 量/億 噸174179183189207234220.5256270285以前述的程序段為基礎(chǔ)，將上表的數(shù)據(jù)代入，并更新時(shí)間軸數(shù)據(jù)，即得 到新程序。clear syms a u;c=a u; % 灰參數(shù) cA二174 179 183 189 207 234 220.5 256 270 285; %原始數(shù)據(jù)Ago二cumsum(A); %原始數(shù)據(jù)一次累加n=length(A); %原始數(shù)據(jù)個(gè)

18、數(shù)for i=1:( n-1)C(i)=(Ago(i)+Ago(i+1)/2; %生成累加矩陣end%計(jì)算待定參數(shù)的值Yn 二A; Yn (1)=;Yn二Yn：E=-C;ones(1,n-1);c=inv(E*E)*E*Yn;c=c;a=c(1);u=c(2);% 預(yù)測(cè)后續(xù)數(shù)據(jù)F=;F(1)=A(1);for i=2:(n+10)F(i)=(A(1)-u/a)/exp(a*(i-1)+u/a ;endG=;G(1)=A(1);for i=2:(n+10)G(i)=F(i)-F(i-1); % 得到預(yù)測(cè)出來(lái)的數(shù)據(jù)endt1=1995:2004;t2=1995:2014;G,a,u %輸出預(yù)測(cè)值、

19、發(fā)展系數(shù)和灰色作用量plot(t1,A,o,t2,G) %原始數(shù)據(jù)與預(yù)測(cè)數(shù)據(jù)的比較550500450400350300250200150 19941996199820002002Ifr2004200620082010201220141995 2014預(yù)計(jì)長(zhǎng)江排放的污水量3.2實(shí)例二：預(yù)測(cè)與會(huì)代表人數(shù)(CUMCM 2009D)1問(wèn)題描述該題要求為會(huì)議籌備組制訂一個(gè)預(yù)訂賓館客房、租借會(huì)議室、租 用客車(chē)的合理方案。為了解決這個(gè)問(wèn)題，需要先預(yù)測(cè)與會(huì)代表的人數(shù)。 預(yù)測(cè)的依據(jù)是代表回執(zhí)數(shù)量及往屆的與會(huì)人員數(shù)據(jù)。已知本屆會(huì)議的回執(zhí)情況表及以往幾屆會(huì)議代表回執(zhí)和與會(huì)情況表。要解決的問(wèn)題是：根據(jù)這些數(shù)據(jù)預(yù)測(cè)本屆與會(huì)代表人數(shù)。本屆會(huì)議代表的回執(zhí)中有關(guān)住房要求的信息單位:人要求合住1合住2合住3獨(dú)住1獨(dú)住2獨(dú)住3男154104321076841女784817592819說(shuō)明：表頭第一行中的數(shù)字 1、2、3分別指每天每間120160元、161200元、201300 元三種不同價(jià)格的房間。合住是指要求兩人合住一間。獨(dú)住是指可安排單人間或一人單獨(dú)住 一個(gè)雙人間。以往幾屆會(huì)議代表回執(zhí)和與會(huì)情況單位：人屆次第一屆第二屆第三屆第四屆發(fā)來(lái)回執(zhí)的代表 數(shù)量315356408711發(fā)來(lái)回執(zhí)但未與 會(huì)的代表數(shù)量891151