基于灰色關(guān)聯(lián)分析法的大氣環(huán)境質(zhì)量評估

1、2011高教社杯全國大學(xué)生數(shù)學(xué)建模競賽承 諾 書我們仔細(xì)閱讀了中國大學(xué)生數(shù)學(xué)建模競賽的競賽規(guī)則.我們完全明白，在競賽開始后參賽隊員不能以任何方式（包括電話、電子郵件、網(wǎng)上咨詢等）與隊外的任何人（包括指導(dǎo)教師）研究、討論與賽題有關(guān)的問題。我們知道，抄襲別人的成果是違反競賽規(guī)則的, 如果引用別人的成果或其他公開的資料（包括網(wǎng)上查到的資料），必須按照規(guī)定的參考文獻(xiàn)的表述方式在正文引用處和參考文獻(xiàn)中明確列出。我們鄭重承諾，嚴(yán)格遵守競賽規(guī)則，以保證競賽的公正、公平性。如有違反競賽規(guī)則的行為，我們將受到嚴(yán)肅處理。我們參賽選擇的題號是（從A/B/C/D中選擇一項填寫）： A 我們的參賽報名號為（如果賽區(qū)設(shè)置

2、報名號的話）： 所屬學(xué)校（請?zhí)顚懲暾娜?參賽隊員 (打印并簽名) ：1. 2. 3. 指導(dǎo)教師或指導(dǎo)教師組負(fù)責(zé)人 (打印并簽名) 日期： 2012年 8 月18日賽區(qū)評閱編號（由賽區(qū)組委會評閱前進(jìn)行編號）：2011高教社杯全國大學(xué)生數(shù)學(xué)建模競賽編 號 專 用 頁賽區(qū)評閱編號（由賽區(qū)組委會評閱前進(jìn)行編號）：賽區(qū)評閱記錄（可供賽區(qū)評閱時使用）：評閱人評分備注全國統(tǒng)一編號（由賽區(qū)組委會送交全國前編號）：全國評閱編號（由全國組委會評閱前進(jìn)行編號）：基于灰色關(guān)聯(lián)分析法的大氣環(huán)境質(zhì)量評估1. 摘要：本文針對全國31個省會的大氣環(huán)境質(zhì)量進(jìn)行建模研究，從分析附錄中數(shù)據(jù)著手，搜集每一年各個城市的污染物濃

3、度。進(jìn)而求得空氣質(zhì)量指標(biāo)。在合理假設(shè)的基礎(chǔ)上，對各個問題進(jìn)行求解。最后根據(jù)預(yù)測和實際情況，提出合理建議。針對對大氣質(zhì)量的綜合評價，根據(jù)所查數(shù)據(jù)，建立空氣質(zhì)量指數(shù)評價模型,實現(xiàn)31個城市空氣質(zhì)量指數(shù)排名；得出?？?，拉薩，南寧，福州空氣質(zhì)量較好，蘭州，烏魯木齊等空氣質(zhì)量較差。針對2012年環(huán)境質(zhì)量的預(yù)測，根據(jù)各個城市地理位置加以分類，進(jìn)行相關(guān)度及殘差檢驗，利用GM（1，1）灰色系統(tǒng)預(yù)測2012年污染物濃度。結(jié)果顯示，顆粒物為主要污染物，華中和華北地區(qū)污染狀況依最為嚴(yán)重，華南地區(qū)較輕。針對不采取節(jié)能減排的預(yù)測，基于GM（1，1）灰色系統(tǒng)預(yù)測分析2006以后十年間的、粉塵、煙塵含量。根據(jù)數(shù)據(jù)顯示：20

4、16年的含量高達(dá)5275萬噸，環(huán)境形勢嚴(yán)峻。在問題二的基礎(chǔ)上加以研究，將2007年以后問題二所得數(shù)據(jù)與附件中的數(shù)據(jù)進(jìn)行差值分析得知：節(jié)能減排之前2007-2010年的總排放量為21838.6萬噸，而節(jié)能減排后的總排放量為9188.8萬噸，在這四年期間的排放總量減少了%，節(jié)能減排效果明顯。對于節(jié)能減排實施前后各省會城市大氣環(huán)境質(zhì)量改善情況的分析，我們?nèi)赃x用問題二中的灰色預(yù)測模型，對20072010年未節(jié)能減排的空氣質(zhì)量指數(shù)進(jìn)行預(yù)測，再與問題一中所求出的節(jié)能減排后的20072010年各省會城市的空氣質(zhì)量指數(shù)做差值分析?？煽闯霰本瑸豸斈君R，呼和浩特，成都改善比高達(dá)79.93%，81.54%，45.

5、50%，38.16%，海口，南寧，福州等環(huán)境較好城市基本無改變。其它城市的改善量都在20%以下。因此節(jié)能減排后各省會城市大氣環(huán)境比減排前環(huán)境狀況有了明顯的改善。關(guān)鍵詞： 指數(shù)評價模型 灰色模型 節(jié)能減排 關(guān)聯(lián)度2. 問題的重述近年來，隨著社會的不斷發(fā)展，人們對周圍環(huán)境重視程度也隨之加深。據(jù)世界衛(wèi)生組織和聯(lián)合國環(huán)境組織發(fā)表的一份報告報道：空氣污染已成為全世界城市居民生活中一個無法逃避的現(xiàn)實。專家們呼吁：“以人為本，建設(shè)文明和諧社會，改善人與自然的環(huán)境，減少污染?！痹诋?dāng)今社會，經(jīng)濟(jì)飛速發(fā)展的同時消耗了能源資源，造成了生態(tài)環(huán)境的日益惡化。我國國情和國力的限制以及生態(tài)環(huán)境十分脆弱的現(xiàn)實情況，我們已經(jīng)沒

6、有“先污染后治理的”資本。因此，我國從2007年8月起，中央財政開始實施節(jié)能減排工作，既是對人類社會發(fā)展規(guī)律認(rèn)識的不斷深化，也是積極應(yīng)對全球氣候變化的迫切需要，取得了顯著的成果。試根據(jù)我國近年污染物總量減排和大氣環(huán)境相關(guān)數(shù)據(jù)，并結(jié)合經(jīng)濟(jì)發(fā)展情況，根據(jù)附錄中的數(shù)據(jù)，結(jié)合你們收集到的相關(guān)資料，建立數(shù)學(xué)模型，完成以下問題：1、建立模型對全國各省會城市的大氣環(huán)境質(zhì)量做出定量的綜合評價，并對2012年各地區(qū)大氣的污染狀況進(jìn)行分析比較。2、假如不采取節(jié)能減排，依照過去幾年的主要統(tǒng)計數(shù)據(jù)，對我國大氣環(huán)境的發(fā)展趨勢做出預(yù)測分析，3、分析討論節(jié)能減排對大氣環(huán)境質(zhì)量改善所起作用。4、對節(jié)能減排實施前后各省會城市大

7、氣環(huán)境質(zhì)量改善情況進(jìn)行科學(xué)分析。5、對下一步實施節(jié)能減排提出建議。3. 問題的分析根據(jù)附表中的含義，分別求出2003-2010年各個城市的空氣質(zhì)量指標(biāo)數(shù)，求其均值。依據(jù)附錄中污染物濃度的不同可對城市的污染狀況做出判斷。用灰色模型預(yù)測發(fā)展趨勢。最后，根據(jù)所得的數(shù)據(jù)，對節(jié)能減排提出合理化建議。要想對各省會的大氣環(huán)境質(zhì)量做出定量的綜合評價，我們需要首先根據(jù)所查得的2003-2013年各個省市的、以及可吸入顆粒物濃度，利用附件中的公式計算，再根據(jù)的計算公式確定的值。對2012年大氣污染，狀況的分析，由于數(shù)據(jù)較多，須首先對其進(jìn)行分類，根據(jù)地理位置的不同，將其分為華北地區(qū)、東北地區(qū)、華中地區(qū)、華東地區(qū)、華

8、南地區(qū)、西南地區(qū)和西北地區(qū)，再采用灰色模型進(jìn)行預(yù)測。我國從2007年8月開始實施節(jié)能減排工作，通過觀察附件中的數(shù)據(jù)，我們利用附件中各污染物排放值，建立灰色理論模型，通過殘差分析對預(yù)測結(jié)果進(jìn)行檢驗，對模型的準(zhǔn)確度加以評價，合格后即可進(jìn)行預(yù)測。對于問題三，需要在問題二的基礎(chǔ)上，將其所得數(shù)據(jù)與附件中已知數(shù)據(jù)進(jìn)行對比分析，使用Excel做出相應(yīng)圖形，形象直觀的反映出節(jié)能減排對大氣環(huán)境質(zhì)量改善所起的作用。對于節(jié)能減排實施前后各省會城市大氣環(huán)境質(zhì)量改善情況的分析。首先，要將31個省會城市在沒有實施節(jié)能減排時的根據(jù)2007年以前的數(shù)據(jù)進(jìn)行預(yù)測，得出20072010的空氣污染指數(shù)，我們利用灰色預(yù)測模型以某一個

9、城市為例，預(yù)測出20072010年節(jié)能減排實施前的城市空氣污染指數(shù)，再和已知的2007年以后數(shù)據(jù)進(jìn)行分析比較。4. 模型的假設(shè)1)假設(shè)從官方獲取的各個省會城市的大氣環(huán)境質(zhì)量的相關(guān)數(shù)據(jù)真實可靠。2)假設(shè)影響空氣質(zhì)量的因素只有，其他影響因素忽略掉。3)假設(shè)相關(guān)數(shù)據(jù)具有獨立性，各個指標(biāo)也不相互影響。 4)假設(shè)城鄉(xiāng)布局不會發(fā)生明顯變化。5)假設(shè)植被、樹木等不會遭到刻意的破壞。6)假設(shè)各個工業(yè)污染源正常排放污染物。5. 符號說明污染物項目的空氣質(zhì)量分指數(shù)的質(zhì)量濃度值與相近的污染物濃度限值的高位值與相近的污染物濃度限值的低位值與對應(yīng)的空氣質(zhì)量指數(shù)與對應(yīng)的空氣質(zhì)量指數(shù)空氣質(zhì)量指數(shù)次累加關(guān)聯(lián)系數(shù)關(guān)聯(lián)度 6.

10、模型的建立與求解6.1. 問題一的建立與求解6.1.1. 空氣質(zhì)量指數(shù)評價模型對空氣質(zhì)量指數(shù)進(jìn)行編程處理大量數(shù)據(jù)，利用國家標(biāo)準(zhǔn)HJ633-2012環(huán)境空氣質(zhì)量指數(shù)(AQI)技術(shù)規(guī)定（試行）中的計算公式算出各個省會城市各個年份的各個指標(biāo)的LAQI值?？諝赓|(zhì)量指數(shù)（AQI）技術(shù)規(guī)定中的公式根據(jù)這些公式建立所需的模型由此，通過查找資料1中華人民共和國國家統(tǒng)計局1，可獲得2003-2010年的各個省會城市的空氣質(zhì)量指標(biāo)數(shù)據(jù)。1.空氣質(zhì)量分指法：空氣質(zhì)量分指法：進(jìn)行單項污染物評價所采用的方法；計算公式：本文在此用空氣質(zhì)量分指數(shù)法來計算各個污染物空氣質(zhì)量指數(shù)。根據(jù)所給附件中表1：和根據(jù)所查得的200320

11、10各省會城市的空氣質(zhì)量指標(biāo)數(shù)據(jù)，用matlab軟件編程（附錄1）進(jìn)行數(shù)據(jù)處理得到結(jié)果見（附錄2）。2.環(huán)境空氣質(zhì)量指數(shù)法3.環(huán)境空氣質(zhì)量指數(shù)：定量描述空氣質(zhì)量狀況的無量綱指數(shù)。它能夠突出高濃度污染物對大氣環(huán)境的影響作用，本文根據(jù)該指數(shù)來作為評判大氣污染的程度。4.環(huán)境空氣質(zhì)量指數(shù)按式計算：由公式得到中國各省會城市年平均的空氣質(zhì)量指數(shù)，進(jìn)行排序表 1結(jié)果如下表所示 表 1 中國各省會城市年平均的空氣質(zhì)量指數(shù)城 市城 市海 口39南 京80拉 薩53合 肥80南 寧60長 沙80福 州61成 都82昆 明67重 慶84廣 州67武 漢85呼和浩特69沈 陽85上 海69西 寧88南 昌69濟(jì) 南

12、89貴 陽72西 安89長 春72石 家 莊94銀 川75北 京100鄭 州78烏魯木齊101哈 爾 濱78太 原103天 津79蘭 州131杭 州79 圖 1 2003-2010各個省會城市的均值表 2 空氣環(huán)境指數(shù)等級標(biāo)準(zhǔn)空氣質(zhì)量指數(shù)空氣質(zhì)量指數(shù)級別空氣質(zhì)量指數(shù)類別050一級優(yōu)51100二級良101150三級輕度污染151200四級中度污染201300五級重度污染300六級嚴(yán)重污染由表 1和表 2得到中國省會城市的污染程度結(jié)果如表 3表 3 中國省會城市的污染程度空氣質(zhì)量指數(shù)級別包括的城市空氣質(zhì)量指數(shù)類別一級?？趦?yōu)二級拉薩 南寧 福州 廣州 昆明 南昌 上海 呼和浩特 長春 貴州 銀川 哈

13、爾濱 鄭州 杭州 天津 長沙 合肥 南京 成都 重慶 沈陽 武漢 西寧 西安 濟(jì)南 石家莊良三級北京 烏魯木齊 太原 蘭州輕度污染表 3可以看出我國大部分城市屬于二級標(biāo)準(zhǔn)?？诘目諝猸h(huán)境質(zhì)量良好，有表 1和表 2可以看出?？凇⒗_、南寧、福州空氣質(zhì)量指數(shù)較低，環(huán)境相對較好。北京、烏魯木齊、太原、蘭州、石家莊空氣環(huán)境質(zhì)量指數(shù)較高，環(huán)境質(zhì)量較差，沒有嚴(yán)重污染的城市。并且根據(jù)附錄1程序結(jié)果數(shù)據(jù)顯示，在所有省會城市中空氣環(huán)境指數(shù)大于50或大于100的環(huán)境污染指標(biāo)均為顆粒物，即顆粒物是這些城市的首要污染物或超標(biāo)污染物。6.1.2. 2012年各地區(qū)大氣的污染狀況進(jìn)行分析比較本題所研究的對象數(shù)據(jù)是2003-

14、2010年空氣質(zhì)量指標(biāo)為小樣本。而神經(jīng)網(wǎng)絡(luò)和多元線性回歸可對大量數(shù)據(jù)進(jìn)行預(yù)測，當(dāng)數(shù)據(jù)不充分時，無法進(jìn)行正常工作，灰色系統(tǒng)的研究對象是“部分信息已知、部分信息未知”的“小樣本”、“貧信息”不確定性系統(tǒng)，它通過對“部分”已知信息的生成、開發(fā)實現(xiàn)對現(xiàn)實世界的確切描述和認(rèn)識。信息不完全是“灰”的基本含義。而本題已知數(shù)據(jù)相對較少，因而采用灰色系統(tǒng)預(yù)測較為合適?；疑到y(tǒng)理論2灰色系統(tǒng)分析方法講解教程2是20世紀(jì)80年代，由華中理工大學(xué)鄧聚龍教授首先提出并創(chuàng)立的一門新興學(xué)科。它是基于數(shù)學(xué)理論的系統(tǒng)工程學(xué)科?；疑到y(tǒng)理論提出了一種新的分析方法關(guān)聯(lián)度3灰色關(guān)聯(lián)度分析 3分析方法，即根據(jù)因素之間發(fā)展態(tài)勢的相似或相

15、異程度來衡量因素間關(guān)聯(lián)的程度，它揭示了事物動態(tài)關(guān)聯(lián)的特征與程度。由于以發(fā)展態(tài)勢為立足點，因此對樣本量的多少沒有過分的要求，也不需要典型的分布規(guī)律，計算量少到甚至可用手算，且不致出現(xiàn)關(guān)聯(lián)度的量化結(jié)果與定性分析不一致的情況。用灰色系統(tǒng)理論建立的微分方程模型稱為灰色模型，即GM模型。用于預(yù)測的模型主要是GM(1，1)模型，它是一階單個變量的預(yù)測模型，其建模過程中僅利用預(yù)測對象本身數(shù)據(jù)的一個時間數(shù)列，而不考慮影響預(yù)測對象的其他各種因素?？陀^世界的很多實際問題，其內(nèi)部的結(jié)構(gòu)、參數(shù)以及特征并未全部被人們了解，人們不可能像研究白箱問題那樣將其內(nèi)部機(jī)理研究清楚，只能依據(jù)某種思維邏輯與推斷來構(gòu)造模型。對這類部分

16、信息已知而部分信息未知的系統(tǒng)，我們稱之為灰色系統(tǒng)。關(guān)聯(lián)度是用來定量描述各變化過程之間的差別，通過對時間序列數(shù)據(jù)進(jìn)行幾何關(guān)系比較，用關(guān)聯(lián)系數(shù)反映兩個被比較序列在某一時刻的緊密程度。關(guān)聯(lián)系數(shù)越大，說明預(yù)測值和實際值越靠近。6.1.3. 模型的建立與求解建立GM模型，實際就是將原始數(shù)列經(jīng)過累加生成后，建立具有微分、差分近似指數(shù)規(guī)律兼容的方程，成為灰色建模，所建模型稱為灰色模型，簡記為GM（Grey Model）。如GM（m,n）稱為m階n個變量的灰色模型，其中GM（1，1）模型是GM（1，n）模型的特例，是灰色系統(tǒng)最基本的模型，也是常用的預(yù)測模型，因此本章重點介紹幾種GM（1，1）模型的建模過程和計

17、算方法，并簡單介紹GM（1，n）建模過程。GM（1，1）模型是GM（1，N）模型的特例，其簡單的微分方程形式（白化形式的微分方程）是利用常數(shù)變易法解得，通解為若初始條件為，則可得到微分方程的特解為或時間響應(yīng)函數(shù)其中白化微分方程中的項中的為的背景值，也稱為初始值; 為常數(shù)（有時也將寫成）。按白化導(dǎo)數(shù)定義有差分形式的微分方程，即顯然，當(dāng)時間密化值定義為1，即當(dāng)時，上式可記為記為離散形式這顯然表明是一次累計生成，因此上述方程可改寫為這實際也表明，模型是以生成數(shù)（是以的一次累加）為基礎(chǔ)的。當(dāng)足夠小時，到不會發(fā)生突變，因此可取與的平均值作為時的背景值，因此，背景值便可記為或于是白化的微分方程可改寫為或K

18、值可從1取到n，構(gòu)成非齊次線性方程組即因此，上述方程可以改寫為矩陣方程形式，即引入下列符號，設(shè) 于是便有令 則解得將求解得到的代入微分方程的解式（也稱時間響應(yīng)函數(shù)），則由于，因此求導(dǎo)還原得上述兩式便為GM（1，1）的時間響應(yīng)式，及灰色系統(tǒng)預(yù)測模型的基本算式，當(dāng)然上述兩式計算結(jié)果只是近似計算值。為簡記，一般可以將GM（1，1）的建模過程記為利用模型進(jìn)行求解的求解：要對2012年各地區(qū)大氣的污染狀況進(jìn)行預(yù)測，需要分別預(yù)測、以及可吸入顆粒物對各地的污染狀況。然而，由于數(shù)據(jù)較多，較為繁瑣，同一地區(qū)的環(huán)境狀況呈現(xiàn)相關(guān)性。因此，根據(jù)各個省會的地理位置不同，可將31個省會分為7個類別，分別為：華北地區(qū)：北京

19、、天津、石家莊、太原、呼和浩特東北地區(qū)：沈陽、長春、哈爾濱華東地區(qū)：上海、南京、杭州、合肥、福州、南昌、濟(jì)南華中地區(qū)：鄭州、武漢、長沙華南地區(qū)：廣州、南寧、海口西南地區(qū)：重慶、成都、貴陽、昆明、拉薩西北地區(qū)：西安、蘭州、西寧、銀川、烏魯木齊根據(jù)中國統(tǒng)計年鑒，查得2003-2010主要城市空氣質(zhì)量指標(biāo)，分別求得各年的污染物濃度。進(jìn)而根據(jù)各地區(qū)的不同，對各污染物的濃度求均值，從而獲得原始數(shù)據(jù)。通過程序?qū)ζ溥M(jìn)行關(guān)聯(lián)度R和殘差分析，從而判斷預(yù)測結(jié)果的好壞。下面以華北地區(qū)為例，對2012年的顆粒物進(jìn)行預(yù)測。設(shè)原始數(shù)據(jù)為通過計算，可以得到1.關(guān)聯(lián)度檢驗2.殘差檢驗：相應(yīng)的預(yù)測模型的模擬序列為殘差序列為相對

20、誤差序列為：平均相對誤差為 檢驗該模型的方法是通過對關(guān)聯(lián)度和殘差的考察來判斷模型的精度，其中平均相對誤差和相對誤差均要求越小越好，關(guān)聯(lián)度要求越大越好，結(jié)合兩者的大小，即可對2012年的大氣環(huán)境污染狀況做出分析（程序見附錄3）。表 4 2012年全國各地區(qū)污染物濃度污染物華北東北華東華中華南西南西北顆粒物92.83079.63675.94484.60955.72672.06292.15368.12642.90147.62567.19648.31463.21761.35342.35242.00941.42442.51736.8826400038.934圖 2 2012年全國各個地區(qū)污染物濃度由上圖

21、可以看出：華北和華中地區(qū)各個污染物濃度均較大，原因可能是該區(qū)經(jīng)濟(jì)比較發(fā)達(dá)，人口比較密集，工業(yè)排氣量大，車輛尾氣排放量大，進(jìn)而導(dǎo)致環(huán)境污染較差；華南地區(qū)污染較輕，原因可能是因為該區(qū)大部分城市沿海，自然條件優(yōu)越，同時輕工業(yè)發(fā)達(dá)，污染較少，工業(yè)廢氣排放量少。6.2. 問題二的建立與求解6.2.1. 對已知條件的分析由于灰色系統(tǒng)的研究對象是“部分信息已知、部分信息未知”的“小樣本”、“貧信息”不確定性系統(tǒng)，可通過對“部分”已知信息的生成、開發(fā)實現(xiàn)對現(xiàn)實世界的確切描述和認(rèn)識。根據(jù)具體灰色系統(tǒng)的行為特征數(shù)據(jù)，充分開發(fā)并利用不多的數(shù)據(jù)中的顯信息和隱信息，尋找因素間或因素本身的數(shù)學(xué)關(guān)系。,建立一個按時間作逐段

22、分析的模型。因此，在已知2000-2006年部分信息已知的前提下，此問題仍然采用了灰色理論系統(tǒng)來對未來十年的排放量進(jìn)行預(yù)測。6.2.2. 模型的建立與求解模型在第一問中已經(jīng)建立，所以此問題中直接進(jìn)行模型的求解。首先，對未來十年的SO2排放總量進(jìn)行預(yù)測：表 5 2000-2006年全國每年的排放總量(萬噸)年份2000200120022003200420052006 1995.11947.81926.62158.72254.92549.32588.8以下即為以排放量來建立GM（1，1）模型的求解過程。第一步，設(shè)原始數(shù)據(jù)為第二步，對原始數(shù)據(jù)進(jìn)行累加生成，即因此累加生成數(shù)據(jù)為第三步，構(gòu)造矩陣 第四步

23、，計算a和u.先求，即根據(jù)逆矩陣的求解方法，得再求的值，即進(jìn)而求得的值為，計算a,u的程序（附錄4）第五步，將的值代入微分方程的時間響應(yīng)函數(shù)，令，得第六步，求導(dǎo)還原得6.2.3. 對上述模型進(jìn)行精度檢驗。常用的方法是回代檢驗，即分別用模型求出各時刻值，然后求相對誤差。先利用時間響應(yīng)函數(shù)求各時刻值（），并計算相對誤差，結(jié)果表 6所示. 表 6 精度檢驗表（程序見附錄5）計算值實測值殘差相對殘差1945.32078.82221.52374.12537.02711.21947.81926.62158.72254.92549.32588.82.5-152.2-62.8-119.212.3-122.40

24、.001280.001290.001150.001110.000980.00097用matlab程序（附錄6）將上表畫成圖形4岡薩雷斯著，阮秋琦譯.數(shù)字圖像處理M，北京，電子工業(yè)出版社，2008。4圖 3 2001-2006年SO2計算值與預(yù)測值再利用時間響應(yīng)函數(shù)模型求各時刻值（），并計算相對誤差，結(jié)果表 7所示.表 7 計算值與實驗原始數(shù)據(jù)值對照表（程序見附表7）計算值實測值殘差相對殘1945.32078.82221.52374.12537.02711.21947.81926.62158.72254.92549.32588.82.5-152.2-62.8-119.212.3-122.40.0

25、01280.001290.001150.001110.000980.00097從殘差檢驗結(jié)果看，數(shù)列曲線擬合較好，相對誤差在0.01即1%左右。從而用（，）模型可以預(yù)測出（2007-2010）全國每年的空氣質(zhì)量指標(biāo)AQI。6.2.4. 進(jìn)一步的預(yù)測利用時間響應(yīng)函數(shù)模型求各時刻值程序（附錄8）表 8 2007-2016年的預(yù)測值年份2007200820092010201120122013201420152016289730963308.3536377840384321461849355275用同樣的方法可以預(yù)測出煙塵排放量、工業(yè)粉末排放量（程序見附錄9）表 9 煙塵排放量、工業(yè)粉末排放量年份20

26、07200820092010201120122013201420152016煙塵11491168.51187.81207.41227.41247.71268.31289.31310.61332.3粉末811784.9759.0734.0709.8686.4663.8641.9620.8600.4用excel對、煙塵、粉末排放量進(jìn)行圖像表示：圖 4 2007-2016年的預(yù)測值由上圖可以看出：在不節(jié)能減排的情況下，含量會急劇上升，到2016年年排放量高達(dá)5275萬噸。煙塵含量也會略有增加，而粉末會減少。會帶來嚴(yán)重的酸雨問題，同時也會消耗大量的煤礦資源和化學(xué)需氧量，空氣質(zhì)量嚴(yán)重下降，一級城市越來越

27、少，也會影響國家年均增長。6.3. 問題三的求解要想討論節(jié)能減排及安排對大氣環(huán)境質(zhì)量改善所起的作用，需要在問題二的基礎(chǔ)上對數(shù)據(jù)加以分析，由于我國在2007年以后開始進(jìn)行節(jié)能減排，因此根據(jù)已知數(shù)據(jù)和預(yù)測數(shù)據(jù)，可以建立下表并用matlab作圖（附錄10），進(jìn)而直觀的反映出節(jié)能減排后使得各個省會城市大氣環(huán)境質(zhì)量的改善情況。表 10 2007-2010實際和預(yù)測值年份2007200820092010預(yù)測2897.43096.33308.93536實際值2468.12321.22214.42185.1煙塵預(yù)測1149.51168.51187.81207.4實際值986.1901.6847.2829.1粉

28、末預(yù)測811.6784.9759.0734.0實際值698.7584.9523.6448.7圖 5 2007-2010年三種污染物實際與預(yù)測值從上圖可以明顯看出，節(jié)能減排的政策能夠有效的改善環(huán)境質(zhì)量。在節(jié)能減排的號召下，煙塵和粉塵的排放量均得到了有效的控制。但得到最有效控制的當(dāng)屬。可以看出是大氣環(huán)境的主要污染物，在節(jié)能減排之前2007-2010年的總排放量為21838.6萬噸，而節(jié)能減排后的總排放量為9188.8萬噸，在這四年期間的排放總量減少了%。滿足了國民經(jīng)濟(jì)和社會發(fā)展第十一個五年規(guī)劃綱要提出了“十一五”期間主要污染物排放總量減少10%的約束性指標(biāo)。從而可以看出“十一五”提出的這項指標(biāo)是行

29、之有效的，酸雨污染減輕，節(jié)約能源，重點城市的環(huán)境質(zhì)量得以保持穩(wěn)定，達(dá)標(biāo)城市的數(shù)量增加，使得“兩控區(qū)”的污染狀況得以減輕，改善了大氣環(huán)境質(zhì)量。6.4. 問題四的求解根據(jù)問題一中所求出的各省會城市空氣質(zhì)量指數(shù)，我們以北京為例，利用北京20032006年空氣質(zhì)量指標(biāo)數(shù)據(jù)來預(yù)測沒有實施節(jié)能減排時北京20072010年的空氣質(zhì)量指數(shù)。再根據(jù)此模型求出每各省會城市空氣質(zhì)量指數(shù)據(jù)，再進(jìn)行分析比較。原始數(shù)據(jù)：累加生成數(shù)據(jù)：模型求解及結(jié)果分析利用GM（1，1）的求解程序得時間響應(yīng)式為先利用時間響應(yīng)函數(shù)，求出k=1，2，3時刻的值，結(jié)果如表 11所示：表 11 精度檢驗實測值、殘差值表 k=1 ,2, 3計算值實

30、測值殘差相對殘差184.1623 296.6607 440.213119629245411.8377 -4.6607 13.78690.0604 -0.0160 0.0304 利用時間響應(yīng)函數(shù)模型 求各時刻值（k=1 2 3），并計算相對誤差，結(jié)果表 12所示.表 12 計算值與實驗原始數(shù)據(jù)值對照表 計算值實測值殘差相對殘差99.3014 126.7048 161.6705100961520.6986 -30.7048 -9.67050.0070 -0.3198 -0.0636從殘差檢驗結(jié)果看，原始累加數(shù)據(jù)的相對殘差均值在3.56%左右；其值小于5%，因此模型通過檢驗。而還原數(shù)列的相對誤差較大

31、，其原因是累加生成數(shù)據(jù)將原始數(shù)據(jù)的隨機(jī)性弱化，正負(fù)誤差可能有相互抵消，當(dāng)數(shù)據(jù)再被還原回來時便表現(xiàn)出來。另一個原因可能是原始數(shù)據(jù)的隨機(jī)性比較大，導(dǎo)致個別數(shù)據(jù)誤差相對較大?？傮w來講該模型還是可以用來預(yù)測未來空氣質(zhì)量指標(biāo)的（模型預(yù)測見附錄11）。對20072010年北京空氣污染指數(shù)預(yù)測所得結(jié)果表 13 20072010年北京空氣污染指數(shù)預(yù)測值年份2007200820092010北京206.3263.212335.8428.5根據(jù)程序計算得到各省會城市空氣質(zhì)量指數(shù)（附件12）如下表根據(jù)所預(yù)測的各省會城市空氣質(zhì)量指標(biāo)數(shù)據(jù)以及問題一中所求的原始數(shù)據(jù)，用excel繪制折線圖6王正東，數(shù)學(xué)軟件與數(shù)學(xué)實驗（第二

32、版），科學(xué)出版社，2010年8月。 10 附錄：附錄1計算各個省會城市每年的空氣質(zhì)量指標(biāo)程序clearclcA=load(A.txt)A=A.*1000;s=A(1:31,1);m=A(1:31,2);n=A(1:31,3);z=length(s);for k=1:zif (s(k)=0)&(s(k)50)&(s(k)150)&(s(k)=0)&(m(k)40)&(m(k)80)&(m(k)180)&(m(k)280)&(m(k)=0)&(n(k)50)&(n(k)150)&(n(k)250)&(n(k)=350)c(k)=1/2*(n(k)-250)+250;endendendendendf

33、or k=1:zAQI(k)=max(max(a(k),b(k),c(k);endV=AQI附錄2各個省會城市每年的空氣質(zhì)量指標(biāo)AQI城市20032004200520062007200820092010平均海 口303340414343384039拉 薩585160565451504853南 寧616663646453506060福 州656261615861576261昆 明687666717559596167廣 州758669636461606067呼和浩特836574767360625969上 海747569686967666569南 昌757570686767656969貴 陽8587

34、69716870636372長 春746875757573687072銀 川918670747167707275鄭 州798180817872758178哈爾濱868277777676767678天 津928185827269767379杭 州858081817980747479南 京858680807974758280合 肥758073758392818380長 沙939586817774716780成 都848388878181817782重 慶999785837978787684武 漢929085868782787985沈 陽939484848584807685西 寧9589829383

35、84968788濟(jì) 南9910089828488878489西 安939690929382828889石家莊1548791968983777494北 京961009615299878686100烏魯木齊89919815093989592101太 原153154959687818475103蘭 州1541531511569091100151131附表3：關(guān)聯(lián)度程序function GM1_1(X0)%format long ;m,n=size(X0); X1=cumsum(X0); %累加 X2=;for i=1:n-1 X2(i,:)=X1(i)+X1(i+1);endB=-0.5.*X2 ;

36、t=ones(n-1,1);B=B,t ; % 求B矩陣YN=X0(2:end) ;P_t=YN./X1(1:(length(X0)-1) %對原始數(shù)據(jù)序列X0進(jìn)行準(zhǔn)光滑性檢驗， %序列x0的光滑比P(t)=X0(t)/X1(t-1)A=inv(B.*B)*B.*YN. ;a=A(1) u=A(2) c=u/a ;b=X0(1)-c ; X=num2str(b),exp,(,num2str(-a),k,),num2str(c); strcat(X(k+1)=,X) %syms k; for t=1:length(X0) k(1,t)=t-1; end kY_k_1=b*exp(-a*k)+c;

37、for j=1:length(k)-1 Y(1,j)=Y_k_1(j+1)-Y_k_1(j);endXY=Y_k_1(1),Y %預(yù)測值CA=abs(XY-X0) ; %殘差數(shù)列Theta=CA %殘差檢驗 絕對誤差序列XD_Theta= CA ./ X0 %殘差檢驗 相對誤差序列AV=mean(CA); % 殘差數(shù)列平均值 R_k=(min(Theta)+0.5*max(Theta)./(Theta+0.5*max(Theta) ;% P=0.5R=sum(R_k)/length(R_k) %關(guān)聯(lián)度Temp0=(CA-AV).2 ;Temp1=sum(Temp0)/length(CA);S2

38、=sqrt(Temp1) ; %絕對誤差序列的標(biāo)準(zhǔn)差%-AV_0=mean(X0); % 原始序列平均值Temp_0=(X0-AV_0).2 ;Temp_1=sum(Temp_0)/length(CA);S1=sqrt(Temp_1) ; %原始序列的標(biāo)準(zhǔn)差TempC=S2/S1*100; %方差比C=strcat(num2str(TempC),%) %后驗差檢驗 %方差比 %-SS=0.675*S1 ; Delta=abs(CA-AV) ;TempN=find(Delta=SS);N1=length(TempN);N2=length(CA);TempP=N1/N2*100;P=strcat(

39、num2str(TempP),%) %后驗差檢驗 %計算小誤差概率附表4：計算au的程序：X0=1995.1 1947.8 1926.6 2158.7 2254.9 2549.3 2588.8;m,n=size(X0);%大小。m、n代表行、列X1=cumsum(X0);%累加生成 X2=;for i=1:n-1 X2(i,:)=X1(i)+X1(i+1); endB=-0.5.*X2;t=ones(n-1,1);B=B,t; YN=X0(2:end);P_t=YN./X1(1:(length(X0)-1)A=inv(B.*B)*B.*YN.;a=A(1) u=A(2) Bb1=B.*Bb2=

40、inv(B.*B)b3=B.*YN.b4=u/ab5=X1(1)-b4b6=-a*b5附錄5：計算值程序for k=1:1:6A=27400*exp(0.066435*k)-25405end附錄6：做圖程序x1=linspace(2001,2006,1);y1=linspace(0,1400,1);y1=3877.1,5888.6,8038.2,10335,12791,15414;x1=2001,2002,2003,2004,2005,2006;plot(x1,y1,r)hold onx2=2001,2002,2003,2004,2005,2006;y2=3942.9,5869.5,8028.

41、2,10283,12832,15421;plot(x2,y2,g)附錄7：計算值程序clcfor k=1:1:6d=1820.3*exp(0.0665*k)end附錄：8clcfor k=7:1:16d=1820.3*exp(0.0665*k)end附錄9：預(yù)測煙塵量：X0=1165.4 1069.8 1012.7 1048.7 1059 1182.5 1088.8;m,n=size(X0);%大小。m、n代表行、列X1=cumsum(X0);%累加生成 X2=;for i=1:n-1 X2(i,:)=X1(i)+X1(i+1); endB=-0.5.*X2;t=ones(n-1,1);B=B

42、,t; YN=X0(2:end);P_t=YN./X1(1:(length(X0)-1)A=inv(B.*B)*B.*YN.;a=A(1) u=A(2) Bb1=B.*Bb2=inv(B.*B)b3=B.*YN.b4=u/ab5=X1(1)-b4b6=-a*b5for k=1:1:16d=1024.8*exp(0.0164*k)end預(yù)測粉末：X0=1092 990.6 941 1021 904.8 911.2 808.4;m,n=size(X0);%大小。m、n代表行、列X1=cumsum(X0);%累加生成 X2=;for i=1:n-1 X2(i,:)=X1(i)+X1(i+1); en

43、dB=-0.5.*X2;t=ones(n-1,1);B=B,t; YN=X0(2:end);P_t=YN./X1(1:(length(X0)-1)A=inv(B.*B)*B.*YN.;a=A(1) u=A(2) Bb1=B.*Bb2=inv(B.*B)b3=B.*YN.b4=u/ab5=X1(1)-b4b6=-a*b5for k=1:1:16 d=1026.1*exp(-0.0335*k)end附錄10：x1=linspace(2007,2010,1);y1=linspace(0,1400,1);x1=2007,2008,2009,2010;y1=2897.4,3096.3,3308.9,35

44、36;plot(x1,y1,r)hold onx2=2007,2008,2009,2010;y2=2468.1,2321.2,2214.4,2185.1;plot(x2,y2,g)x3=2007,2008,2009,2010;y3=1149.5,1168.5,1187.8,1207.4;plot(x3,y3,b)x4=2007,2008,2009,2010;y4=986.1,901.6,847.2,829.1;plot(x4,y4,c)x5=2007,2008,2009,2010;y5=811.6,784.9,759.0,734.0;plot(x5,y5,y)x6=2007,2008,2009,2010;y6=689.7,584.9,523.6,448.7;plot(x6,y6,k)附錄11：利用matlab程序求解得出20072010北京空氣質(zhì)量指數(shù)值：X0=m,n=size(X0);X1=cumsum(X0); X2=;for i=1:n-1 X2(i,:)=X1(i)+X1(i+1);endB=-0.5.*X2;t=ones(n-1,1);B=B,t; YN=X0(2:end);P_t=YN./X1(1:(length(X0)-1)A=inv(B.*B)*B.*YN.;a=A(1) u=A(2) Bb1=B.*Bb2=inv(B.*B)b3=B.*Y