PAGE10北京市水資源短缺風險綜合評價西南交通大學熊竣熙、周成龍、王義目錄北京市水資源短缺風險綜合評價……………………1摘要……………………21問題提出………………………31.1問題背景………………31.2問題的提出……………32問題分析………………………32.1概念的解釋……………32.2問題的具體分析………33模型假設………………………44符號約定………………………45模型的建立與求解……………55.0數據的預處理…………55.1問題一：確定北京水資源短缺的主要影響因子……55.1.1指標的選取……………………55.1.2使用主成分分析法確定主要因子……………65.1.3主成分分析結果………………75.2問題二：水資源短缺風險等級的劃分………………85.2.1等級的劃分……………………85.2.2構建風險潛在函數……………95.2.3歷年風險等級歸類……………95.2.4結果的合理性檢驗…………115.2.5相關的調控措施及建議……115.3問題三：北京市水資源短缺的預測………………125.3.1利用灰色模型預測北京市水資源短缺風險………………125.3.2預測結果分析………………145.3.3利用BP神經網絡檢驗灰色模型的預測正確性……………145.4問題四：對北京市水行政主管部門的建議………166模型的推廣與評價…………176.1模型的缺陷…………176.2模型的優點…………176.3模型的推廣…………177參考文獻……………………178附錄…………17摘要：本文采用了主成分分析法、模糊數學、灰色模型預測等方法對北京市的水資源短缺風險進行了評價與預測。問題一中，為了找出北京市水資源短缺的主要風險因子，本文首先找出了11個影響水資源短缺的因子，構建了包括來水風險、用水風險、管理制度調節在內的指標體系。利用插值、均值替換法等方法對數據進行修補與校正，接著本文采用主成分分析法，求得累計貢獻率達到86.69%的三個主成分，在每一個主成分中找出影響力較大的因素為：蓄水量、降雨量、地下水埋深和人口總量。問題二中，本文先將風險暫且分為5級，分別為{風險極限，高風險，中度風險，低風險，無風險}，相應的評分為{5,4,3,2,1}。接著選取偏大型柯西分布作為隸屬函數，計算出各級風險的隸屬度，然后本文構建了潛在風險值函數，運用查找到的數據代入風險潛在函數計算出每一年的風險潛在值；將數據極差歸一化后得到每一年的風險隸屬度，從而評判出每一年的風險等級。本文劃分的結果為，前十年為低風險甚至是無風險年，中間十年在向中度風險過度，后十年為中度風險年。為了驗證等級劃分的合理型，本文利用matlab對30年的11項指標進行聚類，30年的數據大致被分為3類，從1979—1994年的數據為一類，1995—1999年的數據為一類，2000—2008年的數據為一類，與等級劃分結果很相似，從而驗證了本文方法的合理型。針對本小問的分析結果，本文提出了一些控制水資源短缺的建議（見正文）。問題三中，本文構建了灰色模型GM（1,1）對11項指標的2009和2010年數據進行預測。由于較早的數據不具有太大的代表性，于是本文取近20年的數據進行預測，將預測的數據標準化后代入風險潛在函數，得到2009和2010年的風險潛在值為，；將其極差歸一化后得到其風險隸屬度，其隸屬風險等級為中度風險。為了檢驗預測的合理性，本文又采用BP神經網絡預測，利用近20年潛在風險值，對2009和2010年度的水資源短缺的潛在風險值進行預測，預測結果與灰色系統很相似，并且風險程度同樣為中度風險，從而驗證了本文預測的合理性。問題四中，本文針對主要風險因子提出一些實用建議，以通過這些方法改善目前北京市的水資源短缺狀況。關鍵詞：主成分分析風險潛在函數模糊數學灰色模型GM（1,1）BP神經網絡一問題提出1.1問題背景水資源，是指可供人類直接利用，能夠不斷更新的天然水體。主要包括陸地上的地表水和地下水。風險，是指某一特定危險情況發生的可能性和后果的組合。水資源短缺風險，泛指在特定的時空環境條件下，由于來水和用水兩方面存在不確定性，使區域水資源系統發生供水短缺的可能性以及由此產生的損失。近年，我國、特別是北方水資源短缺問題日趨嚴重，水資源成為焦點話題。以北京市為例，北京是世界上水資源嚴重缺乏的大都市之一，其人均水資源占有量不足300m31.2問題的提出《北京2009統計年鑒》及市政統計資料提供了北京市水資源的有關信息。利用這些資料和你自己可獲得的其他資料，討論以下問題：評價判定北京市水資源短缺風險的主要風險因子是什么？影響水資源的因素很多,例如：氣候條件、水利工程設施、工業污染、農業用水、管理制度，人口規模等。2建立一個數學模型對北京市水資源短缺風險進行綜合評價，作出風險等級劃分并陳述理由。對主要風險因子,如何進行調控，使得風險降低？3對北京市未來兩年水資源的短缺風險進行預測，并提出應對措施。4以北京市水行政主管部門為報告對象，寫一份建議報告。二問題分析2.1概念解釋2.1.1風險潛在函數由于風險大小是一個模糊的概念，于是本文定義風險潛在函數將風險的大小進行量化。風險潛在函數通過風險的指標綜合求得。2.2問題的具體分析2.2.1問題一的分析問題一要求評價判定北京市水資源短缺風險的主要風險因子。在本小問中無疑會有很多指標和大量的數據，若逐步分析每一個因子，工作量勢必很大。其次，指標之間有一定的信息重疊程度，采用的方法必須避開這一重疊現象。于是本文決定采用主成份分析法進行分析，確定各個因子的貢獻率，以累計貢獻率達到85%-90%來確定主成份。2.2.2問題2的分析問題二要求對北京的水資源短缺做出綜合評價并別劃分等級，由于等級的劃分是一個很模糊的概念，于是本文決定采用模糊數學中的隸屬函數。不妨先將風險等級劃分為5級，即｛風險極限，高風險，中度風險，低風險，無風險｝，相對應的數值為｛5，4，3，2，1｝。然后根據實際情況，本文采用偏大型柯西分布確定每一級風險的隸屬度，再通過計算每年的風險潛在函數，由此將每一年的風險等級歸類。由于等級劃分的合理性未知，本文采用模糊聚類進行驗證。利用matlab編程將11個指標的30年數據進行聚類，一共可以分為4大類。再將風險等級劃分出的4類年份與模糊聚類確定的4類年份進行對比即可。2.2.3問題3的分析問題三要求對北京市未來兩年的水資源短缺進行預測，本文使用灰色系統模型，將各項指標的未來兩年發展情況進行了預測，然后利用預測出的2009和2010的數值帶入風險潛在函數，計算出2009和2010年的水資源短缺風險。為了驗證結果的合理性，本文又采用神經網絡模型進行預測，比較兩種預測的結果從而確定方法的合理性。三模型假設1.所獲得的數據真實可信2.不會出現導致指標產生巨大變動的突發狀況3.水資源短缺風險不受除與之相關統計指標以外其他因素的影響4.30年中大氣中的水循環正常進行四符號約定——一致化以后第個指標第年的數據——標準化以后第個指標第年的數據——第個主成分——第個主成分的貢獻率——第個主成分的累計貢獻率——風險潛在函數——發展灰度——內生控制灰度——灰色模型級比——灰色模型級比界限——隱層神經元數目五模型的建立與求解5.0數據的預處理5.0.1異常數據的處理在不考慮自然災害等突發狀況的影響條件下，查找到的數據應該是平緩變化的，如果發現某個數據異常，則使用均值替換法將其替換：（1）5.0.2缺失數據的處理對于缺失的數據，本文直接調用matlab中的一維插值函數，選用樣條函數進行修補。5.0.3數據的一致化處理在選取的11個指標中，部分指標是極大型指標，部分指標是極小型指標。對于風險函數，某些數據越大對風險函數的促進作用越大，而某些數據則相反。所以將極小型指標極大化處理（2）5.1問題一：確定北京水資源短缺的主要影響因子5.1.1指標的選取對水資源短缺有影響的指標有很多，大致有三類：1.自然狀況；2.產業用水和生活用水；3.管理體制影響。在查閱相關知識后，最終選取了11個指標。潛在風險函數體系如圖5-1所示：圖5-1潛在風險函數體系5.1.2使用主成份分析法確定主要因子Step1：對原始數據的標準化處理問題一涉及到的評價指標一共有8個，評價對象為1979年——2008年。在數據的預處理完成后，將第年的第個指標記為，按照標準化方法進行標準化處理。（3）其中，，得到標準化矩陣，記為Step2：計算相關系數矩陣相關系數矩陣：（4）記第個指標與第個指標的相關系數為，則相關系數矩陣為，其中，，Step3：計算特征值與特征向量計算相關系數矩陣的特征值，及其對應的特征向量，其中，有特征向量組成11個新的指標變量：（5）其中，為第主成分，Step4：選擇個主成分計算特征值的信息貢獻率和累計貢獻率。稱（6）為主成分的信息貢獻率；稱（7）為的累積貢獻率。當時，可以選取前個指標變量作為主成分，代替原來的11個指標進行綜合分析。5.1.3主成分分析結果根據以上步驟，利用matlab編程，首先求出各指標的相關性系數表，從表中可以發現，某些指標具有很強的相關性，如果直接用這些指標對影響力進行綜合評估，必然造成信息的重疊，影響評價的客觀性。主成分分析可以把多個指標轉化成少數幾個不相關的綜合指標。相關系數矩陣的前幾個特征值及其累計貢獻率如表5-1所示5-1主成份分析結果序號特征值分項貢獻率累計貢獻率17.16590.65140.651421.33930.12180.773231.0310.09370.8669可以看出，前三個特征值的累計貢獻率達到86.69%，效果很好，故本文中選取前三個主成份進行分析，其對應的特征向量如表5-2所示表5-2前三個特征值對應的特征向量第一第二第三-0.2520.3490.2920.5220.490.3290.1950.586-0.102-0.570.1810.370.3250.229-0.071-0.1990.377-0.527第一-0.210.236-0.002-0.0220.077第二-0.044-0.207-0.254-0.0520.009第三-0.1330.3650.4720.018-0.04由此，得出第一、二、三主成分分別為：（8）由主成分系數可以看出，第一主成分中，指標4占的比重最大；第二主成份中，指標2、4占的比重最大；第三主成份中，指標6、9占的比重最大。于是本文得出，影響水資源短缺風險的主要因子為：指標2、4、6、9。實際指標為蓄水量、降雨量、地下水埋深、人口總量。蓄水量與降雨量的的一再的減少，增加了水資源短缺的風險；人口總量的一再增加，使得用水量大量增加，在水資源來源減少的情況下，無疑大量增加了水資源供給的負擔；而地下水位的一再下降，使得地下水的使用越來越難。這些都是使得潛在風險函數增加的因素。5.2問題二：水資源短缺風險等級的劃分5.2.1等級的劃分本文暫且將風險等級分為五類，即｛風險極限，高風險，中度風險，低風險，無風險｝，對應的數值為5，4，3，2，1.表5-3風險等級劃分風險度風險極限高風險中度風險低風險無風險風險度函數（）54321注：風險極限是指達缺水程度到了城市的負荷極限根據實際情況，本文選取偏大型柯西分布隸屬函數作為隸屬函數其中、、、為待定系數當風險度為無風險時，取當風險度為中度風險時，隸屬度為0.8，即當風險度達到極限時，隸屬度為1，即計算出系數、、、分別為=1.1086、=0.8942、=0.3915、=0.3699于是得到風險度的隸屬函數為（9）利用風險隸屬函數計算出風險等級的隸屬度如表5-4表5-4風險等級隸屬度風險度風險極限高風險中度風險低風險無風險風險度函數（）10.91260.80.52450.015.2.2在第一問中，本文通過分析得到了主成份。將主成分分別代入下式，得到風險潛在函數（10）其中，為第個主成分的信息貢獻率，風險潛在函數就是對水資源短缺風險的一個綜合評估分別以主成分貢獻量為權重，構建風險潛在函數為（11）將、、代入，得到風險潛在函數的表達式為：（12）5.2.3歷年風險歸類將預處理后的數據標準化處理后，帶入風險潛在函數進行計算，得到歷年的風險潛在值如表5-4所示：表5-4歷年風險潛在值年份潛在風險年份潛在風險年份潛在風險1979-1.227921989-0.3889419991.2813251980-0.731941990-0.7361520000.9569551981-0.742271991-0.7758820011.0018131982-1.0490519920.86814520022.4767661983-0.8808219931.3310120030.6948391984-0.884431994-0.7673920040.5047531985-1.039381995-0.5553820050.8615781986-0.588681996-0.8746620061.3039731987-0.6071119970.03631420071.1701391988-0.926651998-0.5292820080.818082由表5-4所示，本文得出了歷年的風險潛在值，但是僅僅這樣無法將風險歸類，于是本文將風險潛在值進行平移—極差化變換，得到歷年風險隸屬度如表5-5：表5-5歷年風險隸屬度年份風險隸屬度年份風險隸屬度年份風險隸屬度1979019890.22646519990.6773161980013274120000.5897619810.1310919910.12201820010.60186819820.04828119920.5657872002119830.09369119930.69072820030.51900719840.09271819940.1243120040.46769719850.050892199505640151986009535520060.6834319870.16757419970.34125220070.64730419880.08132119980552274由表5-5可以看出，風險最小的一年是1979年，風險度為“無風險”，而風險最大的一年是2002年，該年度的風險隸屬度達到了風險極限，并且概念的風險隸屬度與每一年都相差較大。而本文所研究的水資源短缺風險是在不考慮突發狀況的情況假設下的，數據應該是平緩變換的，于是本文采用均值替換法求該年的風險隸屬度均值替換公式如下：替換后，2002年的風險隸屬度為0.560438。使用matlab編程將歷年劃分等級，得到歷年的潛在風險等級如表5-6：表5-6歷年的風險等級年份風險等級年份風險等級年份風險等級197911989219993198021990220003198121991220013198221992320023198321993320032198421994220042198521995220053198621996220063198721997220073198821998220083由表5-6可以看出風險變化趨勢為，前十年風險度都比較低，而中間十年開始逐漸出現中度風險的現象，到了后面十年風險度幾乎都是中度風險再將風險隸屬度變化做散點圖得到圖5-1，也可以得到相同的結論：圖5-1歷年風險隸屬度變化趨勢由此分析可以看出，如果不立即采取保護水資源的措施，水資源缺乏風險將會愈演愈烈5.2.4結果的合理性檢驗為了檢驗本文對水資源短缺劃分的合理程度，本文采用聚類的方法將11個指標30年數據的變化進行聚類分析。暫不制定聚類數目，使用matlab編程進行聚類分析，得到的結果如圖5-2：圖5-211個指標30年數據的聚類結果從聚類圖可以看出，30年的數據大致被分為3類，從1979—1994年的數據為一類，1995—1999年的數據為一類，2000—2008年的數據為一類。與本文5.2.3的分析結果相差不大，即隨著城市的發展，水資源短缺的風險度由無風險、低風險逐漸在想中度風險變化。如不立即進行調控，后果將難以想象。5.2.5相關的調控措施5.2.5在問題一中，本文通過主成份分析的方法，找出對水資源短缺風險影響較大的因素為：1.蓄水量2.降雨量3.地下水位埋深4.人口總量。通過作圖分析各個指標的變化趨勢如圖5-3：圖5-3主要指標變化趨勢由圖5-3可以看出，北京市的蓄水量在逐年減少，降雨量變動較大，但是總體趨勢也在遞減。地下水位埋深不斷增加為地下水的開采加大了難度。人口總量在逐漸增加，導致了生活用水的增加。有這些不利的趨勢導致了水資源短缺的風險。1.蓄水量：蓄水量的大小由水利工程的數量與規模決定，在短期內無法做出有效的調整，目前北京只有官廳、密云水庫。于是從長期來說，可以增加水庫的容量，甚至是考慮增加水庫的數量，以此增加蓄水量；其次是加強管理制度，科學合理的防水補源。2.地下水位：地下水位主要接受大氣降水補給，其次為側向徑流補給、河渠滲漏補給、地表水灌溉和井灌回滲以及城市工業和生用水的排泄。隨著城市規模擴大和人類活動的加劇，大量的工業廢水、生活污水等通過河渠排放滲入地下，使大部分淺層地下水遭受污染。于是，采取的相應的措施為：適當的采取人工降雨，增加降雨量；提高廢水處理率，減少上游污水的排放；3.降雨量：可以通過人工降雨的方法增加降雨量，增大綠化面積，維持一個正常的水循環。4.人口總量：人口總量對生活用水會產生直接的影響。由圖5-3可以很直觀的看到，總體上人口總量呈增長趨勢，生活用水勢必會大幅度增加，所以限制常住人口的數量也是一個重要的措施。5.3問題三北京市水資源短缺的預測為了對北京市未來兩年水資源短缺的風險進行預測與研究，本文決定采用灰色模型GM(1，1)預測各項指標2009和2010年的數值，在將數據處理之后帶入風險潛在函數計算風險潛在值，得到2009與2010年的水資源短缺情況，具體步驟如下：5.3.1灰色模型的建立分析處理后的數據，共有11個指標，30年的指標數據。由于過早的數據不具有代表性，于是本文取前后二十年的數據進行預測。將數據的數量級化為一致后，將11個指標的平均值記為（13）平均價格的級比規定范圍為（14）經驗證，級比均落在內。對數列進行一次累加，得到數列則記一次累加生成數列（15）取的加權均值，則，，記（16）于是可以得到GM（1,1）的白化微分方程為（17）取為灰導數，為背景值，則方程（12）所對應的灰微分方程為（18）則可得其中，用最小二乘法求得參數的估計值為（19）解得白化微分方程的特解為則可以解得（20）使用matlab編程求解，得到2009和2010的預測值如表5-7：表5-72009與2010各項指標預測值指標農業用水蓄水量入境水量降水量(毫米)污水處理能力200911.98148.52224.5662370.5357440.6512201011.58428.34954.2916358.9641494.4715指標地下水埋深工業用水第三產業人口總量生活水價工業水價200924.65236.11517.564117135.21143.2909201025.95425.826418.2521757.66.16293.73485.3.2預測結果分析將上述灰色模型預測值處理后，代入模型（12），即風險潛在函數表達式，得到2009年和2010年北京市的風險潛在值分別為：，；數據標準化后，其相應的風險等級均為3，風險程度為中度風險。由于預測的未知性，沒有足夠的依據證明本文預測的準確性，于是本文采用第二種預測方法，神經網絡進行對比檢驗。.3.1神經網絡的設計人工神經網絡具有自組織、自適應和自學習能力,以及具有非線性、非局域性、非定性和非凸性等特點，并且在理論上可以逼近任何非線性函數。使用神經網絡進行預測，不需要輸入、輸出之間明確的函數關系，主要通過對數據訓練、學習完成模擬過程，并利用訓練好的網絡對新輸入的數據進行預測。其中BP(BackPropagation)網絡是反向傳播的多層前饋式網絡，是目前使用最為廣泛的一種人工神經網絡。國內學者利用多種預測方法進行預測研究，對比預測模型的預測精度，發現BP神經網絡預測結果與實際值非常接近，效果很好，完全可以滿足對預測的要求。故本文采取BP神經網絡進行進一步的價格預測，然后將其預測結果與灰色模型的預測結果進行比較，從中選出預測精度較高的模型作為分析基礎。對于BP神經網絡,有一個非常重要的定理,即對于任何在閉區間內的一個連續函數以用單隱層的BP神經網絡逼近。因而可以用一個三層BP神經網絡就可以完成任意的m維到n維的映射,即一般的預測問題都可以通過單隱層的BP網絡實現。所以,本文采用三層BP網絡設計預測模型。首先，進行輸入數據的選取。為了和灰色系統預測區分開來，本文決定采用神經網絡直接對風險潛在值進行預測。由于過早的數據沒有代表性，于是選取近20年的風險潛在值進行預測，網絡輸入層的神經元個數為20。輸出層的神經元數為1。數據的預處理在進行網絡預測之前，要對數據進行歸一化處理，將數據處理為區間［0，1］之間的數據。于是本文采用極差化方法對數據進行歸一處理。網絡預測完成后，對預測結果進行分析時，應對預測數據進行反歸一化，然后才能夠與真實值進行比較。依本文的歸一化式子，反歸一化式子為（21）隱層的設計隱層神經元數目的確定是一個很復雜的問題，目前還沒有科學的統一計算式。但有一些經驗公式可用來大致確定隱層神經元的數目范圍。本文采用下式：（22）其中：——輸入層節點數——輸出層節點數——常數。取值為1,2,3…10由上式知：其最佳隱層神經元數目可能在20到30之間。為了得到最佳的隱層神經元數，所以本文增大搜索范圍,以10為初始值逐步增長，即先從一個簡單的神經網絡開始，若不符合要求，逐步增加隱層神經元數目，直到合適時為止。最終確定的最佳神經元數目以實驗調試結果為準，上式確定的范圍僅是一個供參考的初始范圍。BP神經網絡在matlab上的實現現在用matlab中的神經網絡工具箱對2009年和2010年風險潛在值進行預測。先對BP神經網絡結構進行初步設計。由于本文是一個用于預測，較為復雜的神經網絡，輸入量與輸出量之間滿足非線性關系，故輸入層、輸出層傳遞函數分別選取非線性函數S型正切函數tansig、對數函數logsig。再進行對訓練函數的選擇。由于BP神經網絡有眾多的訓練函數，而訓練函數的選取直接關系到預測所需步長和預測精度。為選取最合適的訓練函數，本文取1989年至2007年的數據對2008年的數據進行預測。2008年的潛在風險值數據作為目標輸出，1989年至2007年的數據作為輸入，對神經網絡進行訓練。將各種訓練函數得到的預測結果與2008年的實際風險值進行比較，得出的各組平均誤差如表5-8：表5-8訓練函數精度表訓練函數相對誤差trainrp36.95%trainlm2.22%traingd5.72%trainbfg27.25%traingdx24.61%由上表中結果容易看出，trainlm訓練函數的預測結果最為精確。網絡神經結構初步設計好后，用逐步增長法確定最佳隱層神經數。最后本文發現最佳隱層神經元數22。至此，BP神經網絡的構建完畢。利用神經網絡進行預測分析利用matlab編程實現預測，得到2009年度與2010年度的風險潛在值分別為，；該預測值與灰色系統的預測值相近，且風險等級都為3級，中度風險。由此可見，本文對于北京市2009年度與2010年度的水資源短缺風險預測具有一定的合理性與正確性。5.3.4北京市水資源短缺應對措施在本文的研究中，北京市歷年水資源短缺的潛在風險值是由各個指標所決定的，所以要進行有效的調控，則必須注意11個指標的變化情況。1.北京市的人口總量有穩定增加的趨勢，人口的增加必然會導致生活用水的增加，于是控制北京市常住人口的數量成為了一個主要因素。2.生活用水在逐步增加，而工業用水和農業用水在逐步減少，所以可以降低工業用水的水費，提高生活用水的水費，這樣可以起到較好的限制作用而又可以不影響工業的發展。3.污水處理能力對可使用的地下水量有很大的影響，本文的預測中，污水處理能力在逐步增強，這也許是由于科技的增強所致。如果有關部門能盡力到達本文對污水處理能力的預測值，那么水資源的短缺得到有效的調控。5.4問題四：對北京市水行政主管部門的建議北京作為中國的首都，一個國際化的大都市，經濟、文化等都相當發達。然而今天卻面臨著水資源短缺這一問題，水資源短缺成為了限制北京進一步發展的短板因素，如何調控相關影響因子，是北京的水資源短缺風險減輕至關重要。由題目中給定的數據和2009年北京統計年鑒上的數據可以看到，北京市每年的水資源總量在遞減，且遞減趨勢較大，而用水總量卻沒有顯著減少，如不進行調控，則水資源短缺問題將會越來越嚴重。根據本文對北京市水資源短缺問題的分析，找出了11個與水資源短缺聯系的指標，這11個指標可以歸結為3個方面：用水風險因子，來水風險因子，管理調節因子；11個指標分別為農業用水、工業用水、人口總量、第三產業及生活等用水；地下水位埋深、入境水量、降雨量；蓄水量、污水處理能力、工業水價、生活水價。而根據本文在第三文中做的預測，如果不及時調控水資源短缺的峰下因子，那么水資源短缺風險將進一步加劇。從本文的分析中，11個指標中的主要風險因子為地下水位埋深、人口總量、蓄水量、降雨量。這四個指標在影響水資源短缺風險中占的比重較大，所以應該首先考慮這4個因子的調控。本文經過分析后提出以下調控措施：1.節約用水應當從水費調控。北京市擁有將近1700萬常住人口，生活用水量很大，適當提高水價是控制用水量的一個重要措施。而水價和用水量之間的關系有待進一步研究；2.北京市的蓄水量在總體趨勢上有一定的減小。增加蓄水量是一個較為重要的措施。增加蓄水量可以通過跨流域調水等實現；3.地下水是一個相當重要的因素。隨著每一年地下水位的逐漸下降，北京市水資源短缺的風險在逐步增加。提高污水處理能力是一個保護地下水的實用措施，這樣可以讓北京市擁有更多的可用地下水資源。4.在降雨量整體減小的情況下，人力對其的控制能力十分微弱，可以通過改善環境進行一些微弱的調控。除了以上針對主要風險因子采取的措施外，還可以通過采取其它方法對非主要因子進行調控，本文在此就不再贅述。六模型的推廣與評價6.1模型的缺陷本文運用合理的方法按照題目的要求給出了模型與解答過程，但是也存在一些缺點。1.在2009和2010年的水資源短缺風險預測中，對人口的預測應該使用logstic模型，本文使用灰色系統進行預測，會造成一些誤差，但是考慮到預測的年限較少，所以誤差可以忽略。在其他的指標預測中，用于預測的年份數選取沒有較好的方法。2.本文的研究中忽略了突發事件的影響，是模型過于理想化，這在以后的研究中應該進一步的完善。6.2模型的優點本文存在缺點的同時，也有一些優點。1.本文定義了風險潛在函數，對研究問題有很大的幫助。2.第一問采用主成分分析定義的權重無主觀性，有較強的說服力。3.第二問中將風險等級劃分后，本文并未草草了事，而是運用聚類分析進行驗證；4.第三問中運用灰色模型對水資源缺乏的風險進行預測后，還運用精度較高的神經網絡模型進行檢驗。這樣的檢驗使得本文的分析較為嚴謹。同時，本文也針對分析結果對相關部門提了一些控制水資源短缺的實用建議。6.3模型的推廣本文運用了科學的方法對水資源的短缺進行了評價和預測，有較強的推廣性。本文的模型不僅可以用于水資源短缺的風險研究，還可以用于其它的很多項目，例如工程建設中的風險研究、工程造價中的成本研究、供應鏈中的影響因素的研究等。七參考文獻[1]韓中庚，數學建模方法及其應用，北京：高等教育出版社，2005年6月第1版[2]姜啟源、謝金星，數學建模案例精選，北京：高等教育出版社，2006年7月第1版[3]薛長虹、于凱，大學數學實驗，成都：西南交通大學出版社，2003年10月第1版[4]北京市2009年度統計年鑒，/tjnj/2009-tjnj/，訪問日期（2011-5）八附錄8.1主成分分析源代碼（matlab）clear,clcX=[7.2 0.0518 0.0499 0.0014 0.0431 24.18 14.37 4.37 897.1 8.3333 49.2 0.0518 0.0499 0.0026 0.0431 31.83 13.77 4.94 904.3 8.3333 49.2 0.0518 0.0499 0.0025 0.0397 31.6 12.21 4.3 919.2 8.3333 49.2 0.0518 0.0499 0.0018 0.0397 28.81 13.89 4.52 935 8.3333 49.2 0.0518 0.0499 0.002 0.0397 31.6 11.24 4.72 950 8.3333 49.2 0.0518 0.0499 0.002 0.0397 21.84 14.376 4.017 965 8.3333 49.2 0.0518 0.0499 0.0014 0.0397 10.12 17.2 4.39 981 8.3333 410 0.0518 0.0935 0.0015 0.0379 19.46 9.91 7.18 1028 8.3333 410 0.0518 0.0935 0.0015 0.0379 9.68 14.01 7.26 1047 8.3333 410 0.0518 0.0619 0.0015 0.0379 21.99 14.04 6.4 1061 8.3333 411 0.0518 0.158 0.0023 0.0379 24.42 13.77 6.45 1075 8.3333 410.7 0.0769 0.0853 0.0014 0.0329 21.74 12.34 7.04 1086 8.3333 410.7 0.0719 0.0779 0.0013 0.0329 22.7 11.9 7.43 1094 8.3333 4.347810.5 0.0939 0.0974 0.0018 0.2222 19.94 15.51 10.98 1102 3.5714 2.631612.8 0.1639 0.1522 0.002 0.2222 20.35 15.28 9.59 1112 4.1667 3.030311.7 0.0526 0.0718 0.0012 0.0408 20.93 14.57 10.37 1125 5 3.571411.6 0.0767 0.0579 0.0017 0.0171 19.33 13.78 11.77 1251.1 5.5556 410.2 0.0511 0.0398 0.0014 0.0171 18.95 11.76 9.3 1259.4 3.5714 2.777812.1 0.1153 0.1176 0.0023 0.0171 18.12 11.1 11.1 1240 3.7037 2.564112 0.0684 0.0692 0.0014 0.0171 17.39 10.84 12.2 1245.6 2.6316 1.470614.2 0.2564 0.1727 0.0037 0.0171 18.45 10.56 12.7 1257.2 1.8182 1.08715.2 0.2907 0.1406 0.0027 0.0078 16.49 10.52 13.39 1363.6 1.1111 0.854716.2 0.216 0.189 0.003 0.007 17.4 9.2 12.3 1385.1 0.8929 0.833317.4 0.578 0.3846 0.0027 0.0055 15.5 7.5 11.6 1423.2 0.7299 0.740718.4 0.1071 0.2392 0.0022 0.0047 13.8 8.4 13.6 1456.4 0.6061 0.662319.04 0.093 0.1582 0.0021 0.0039 13.5 7.7 13.4 1492.7 0.5405 0.602420.21 0.0837 0.2179 0.0024 0.0031 13.2 6.8 14.5 1538 0.4785 0.549521.52 0.0814 0.2353 0.0031 0.003 12.8 6.2 15.3 1581 0.4292 0.507622.79 0.0904 0.2899 0.0021 0.0028 12.4 5.8 16.6 1633 0.3906 0.469522.92 0.0773 0.1869 0.0016 0.003 12 5.2 17.9 1695 0.3571 0.4386];Z=zscore(X);[pc,score,variance,t2]=princomp(Z)percent=variance/sum(variance)leiji=cumsum(percent)8.2聚類分析源代碼（matlab）clc,cleara=[7.2 0.0518 0.0499 0.0014 0.0431 24.18 14.37 4.37 897.1 8.3333 49.2 0.0518 0.0499 0.0026 0.0431 31.83 13.77 4.94 904.3 8.3333 49.2 0.0518 0.0499 0.0025 0.0397 31.6 12.21 4.3 919.2 8.3333 49.2 0.0518 0.0499 0.0018 0.0397 28.81 13.89 4.52 935 8.3333 49.2 0.0518 0.0499 0.002 0.0397 31.6 11.24 4.72 950 8.3333 49.2 0.0518 0.0499 0.002 0.0397 21.84 14.376 4.017 965 8.3333 49.2 0.0518 0.0499 0.0014 0.0397 10.12 17.2 4.39 981 8.3333 410 0.0518 0.0935 0.0015 0.0379 19.46 9.91 7.18 1028 8.3333 410 0.0518 0.0935 0.0015 0.0379 9.68 14.01 7.26 1047 8.3333 410 0.0518 0.0619 0.0015 0.0379 21.99 14.04 6.4 1061 8.3333 411 0.0518 0.158 0.0023 0.0379 24.42 13.77 6.45 1075 8.3333 410.7 0.0769 0.0853 0.0014 0.0329 21.74 12.34 7.04 1086 8.3333 410.7 0.0719 0.0779 0.0013 0.0329 22.7 11.9 7.43 1094 8.3333 4.347810.5 0.0939 0.0974 0.0018 0.2222 19.94 15.51 10.98 1102 3.5714 2.631612.8 0.1639 0.1522 0.002 0.2222 20.35 15.28 9.59 1112 4.1667 3.030311.7 0.0526 0.0718 0.0012 0.0408 20.93 14.57 10.37 1125 5 3.571411.6 0.0767 0.0579 0.0017 0.0171 19.33 13.78 11.77 1251.1 5.5556 410.2 0.0511 0.0398 0.0014 0.0171 18.95 11.76 9.3 1259.4 3.5714 2.777812.1 0.1153 0.1176 0.0023 0.0171 18.12 11.1 11.1 1240 3.7037 2.564112 0.0684 0.0692 0.0014 0.0171 17.39 10.84 12.2 1245.6 2.6316 1.470614.2 0.2564 0.1727 0.0037 0.0171 18.45 10.56 12.7 1257.2 1.8182 1.08715.2 0.2907 0.1406 0.0027 0.0078 16.49 10.52 13.39 1363.6 1.1111 0.854716.2 0.216 0.189 0.003 0.007 17.4 9.2 12.3 1385.1 0.8929 0.833317.4 0.578 0.3846 0.0027 0.0055 15.5 7.5 11.6 1423.2 0.7299 0.740718.4 0.1071 0.2392 0.0022 0.0047 13.8 8.4 13.6 1456.4 0.6061 0.662319.04 0.093 0.1582 0.0021 0.0039 13.5 7.7 13.4 1492.7 0.5405 0.602420.21 0.0837 0.2179 0.0024 0.0031 13.2 6.8 14.5 1538 0.4785 0.549521.52 0.0814 0.2353 0.0031 0.003 12.8 6.2 15.3 1581 0.4292 0.507622.79 0.0904 0.2899 0.0021 0.0028 12.4 5.8 16.6 1633 0.3906 0.469522.92 0.0773 0.1869 0.0016 0.003 12 5.2 17.9 1695 0.3571 0.4386];[m,n]=size(a);d=zeros(m,m);fori=1:mforj=i+1:md(i,j)=mandist(a(i,:),a(j,:)');endenddnd=nonzeros(d);nd=union(nd,nd)fori=1:m-1nd_min=min(nd);[row,col]=find(d==nd_min);tm=union(row,col);tm=reshape(tm,1,length(tm));s(i)={char(['第',int2str(i),'次合成，平臺高度為',num2str(nd_min),'時的分類結果為：',int2str(tm)])};%上面大括號{}代表建立數組nd(find(nd==nd_min))=[];iflength(nd)==0breakendends(:)y=pdist(a,'cityblock');yc=squareform(y)z=linkage(y)[h,t]=dendrogram(z)8.3灰色系統預測源代碼（matlab）clc,clear

x0=[1110.710.710.512.811.711.610.212.11214.215.216.217.418.4

19.0420.2121.5222.7922.92];

n=length(x0);xx=1:20;

plot(xx,x0,'o-')

holdon

disp('級比')

lamda=x0(1:n-1)./x0(2:n)

range=minmax(lamda)

x1=cumsum(x0)

fori=2:n

z(i)=0.4*x1(i)+0.6*x1(i-1);

end

B=[-z(2:n)',ones(n-1,1)];

Y=x0(2:n)';

u=B\Y

x=dsolve('Dx+a*x=b','x(0)=x0');

x=subs(x,{'a','b','x0'},{u(1),u(2),x1(1)});

yuce1=subs(x,'t',[0:n-1]);

digits(6),y=vpa(x)%為提高預測精度，先計算預測值，再顯示微分方程的解

yuce=[x0(1),diff(yuce1)]

disp('殘差')

epsilon=x0-yuce%計算殘差

disp('相對誤差')

delta=abs(epsilon./x0)%計算相對誤差

disp('級比偏差值')

rho=1-(1-0.5*u(1))/(1+0.5*u(1))*lamda%計算級比偏差值

%若以上各項檢驗符合標準，則以下可預測將后的所需的數據

yuce2=subs(x,'t',[0:n-1+2]);%這里的n-1為驗證模型時的數據量的上限，可以加上或減去一定的數值進行相應的預測

digits(2),y=vpa(x)%為提高預測精度，先計算預測值，再顯示微分方程的解

yuce3=[x0(1),diff(yuce2)]

yy=1:22;

plot(yy,yuce3,'+-')

holdoff;8.4神經網絡預測源代碼（matlab）%通用BP神經網絡P=[0.1449293520.0413286760.0294740210.5200189410.6581288080.0320072710.095266910.010.2718168330.103054640.6433037760.5465181090.5599028590.5141053390.4683078190.4115898720.5180594770.6500614960.6101280470.505081121]';t=[0.6794];H=[0.0413286760.0294740210.5200189410.6581288080.0320072710.095266910.010.2718168330.103054640.6433037760.5465181090.5599028590.5141053390.4683078190.4115898720.5180594770.6500614960.6101280470.5050811210.6794]';net=newff([01;01;01;01;01;01;01;01;01;01;01;01;01;01;01;01;01;01;01;01],[22,1],{'tansig','logsig'},'trainlm');%輸入參數依次為:'樣本P范圍',[各層神經元數目],{各層傳遞函數},'訓練函數'%訓練函數traingd--梯度下降法，有7個訓練參數.%訓練函數traingdm--有動量的梯度下降法,附加1個訓練參數mc(動量因子，缺省為0.9)%訓練函數traingda--有自適應lr的梯度下降法,附加3個訓練參數:lr_inc(學習率增長比，缺省為1.05;%lr_dec(學習率下降比，缺省為0.7);max_perf_inc(表現函數增加最大比，缺省為1.04)%訓練函數traingdx--有動量的梯度下降法中賦以自適應lr的方法，附加traingdm和traingda的4個附加參數%訓練函數trainrp--彈性梯度下降法,可以消除輸入數值很大或很小時的誤差，附加4個訓練參數:%delt_inc(權值變化增加量，缺省為1.2);delt_dec(權值變化減小量，缺省為0.5);%delta0(初始權值變化，缺省為0.07);deltamax(權值變化最大值,缺省為50.0)%適合大型網絡%訓練函數traincgf--Fletcher-Reeves共軛梯度法;訓練函數traincgp--Polak-Ribiere共軛梯度法;%訓練函數traincgb--Powell-Beale共軛梯度法%共軛梯度法占用存儲空間小,附加1訓練參數searchFcn(一維線性搜索方法,缺省為srchcha);缺少1個訓練參數lr%訓練函數trainscg--量化共軛梯度法,與其他共軛梯度法相比，節約時間.適合大型網絡%附加2個訓練參數