國賽數(shù)學(xué)建模真題城市污染正式版

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

2、報(bào)名號的話）： 所屬學(xué)校（請?zhí)顚懲暾娜?參賽隊(duì)員 (打印并簽名) ：1. 2. 3. 指導(dǎo)教師或指導(dǎo)教師組負(fù)責(zé)人 (打印并簽名)： 日期： 2013 年 8 月26 日賽區(qū)評閱編號（由賽區(qū)組委會評閱前進(jìn)行編號）：2011高教社杯全國大學(xué)生數(shù)學(xué)建模競賽編號專用頁賽區(qū)評閱編號（由賽區(qū)組委會評閱前進(jìn)行編號）：賽區(qū)評閱記錄（可供賽區(qū)評閱時使用）：評閱人評分備注全國統(tǒng)一編號（由賽區(qū)組委會送交全國前編號）：全國評閱編號（由全國組委會評閱前進(jìn)行編號）：城市表層土壤重金屬污染分析 摘 要本文主要是對重金屬污染濃度的分析，畫出濃度的空間分布圖，畫出濃度的空間分布圖，建立相關(guān)的模型，得出該地區(qū)的重金屬污染

3、程度、污染原因以及污染源的所在位置。問題一要求我們對重金屬污染濃度進(jìn)行具體分析。首先采用Cubic插值法，利用Matlab軟件畫出五個功能區(qū)在城區(qū)的位置分布圖及8種重金屬元素空間分布圖，以顏色深淺對重金屬污染程度進(jìn)行直觀表現(xiàn)，顏色淺則表示濃度大并將兩種分布圖進(jìn)行對比，得到8種重金屬的具體分布區(qū)域。其次建立雙權(quán)重因子模糊綜合評判模型，考慮濃度和毒性雙重因子作為權(quán)重判定各區(qū)域內(nèi)不同元素的污染程度，最終得出：工業(yè)區(qū)普遍受到較重的重金屬污染，其次為主干道路區(qū)，再次是生活區(qū)和公園綠地區(qū)，山林區(qū)基本未受到重金屬的污染。問題二要求通過數(shù)據(jù)分析重金屬污染的主要原因，我們發(fā)現(xiàn)各種重金屬不是單一產(chǎn)生的，多種重金屬

4、在產(chǎn)生的過程中往往具有伴隨性。針對這種現(xiàn)象，我們使用因子分析法來確定各種重金屬元素之間的相關(guān)性，然后利用它們的相關(guān)性定量分析各因子對8種重金屬元素產(chǎn)生的貢獻(xiàn)值。最終通過已知數(shù)據(jù)可以得出： Ni、Cr、Cu的污染主要是由主干道路區(qū)汽車尾氣產(chǎn)生， Hg的污染主要是由工業(yè)“三廢”產(chǎn)生， As的污染主要是由生活區(qū)殺蟲劑等藥劑產(chǎn)生， Zn的污染主要是由生活區(qū)垃圾焚燒和主干道路區(qū)汽車輪胎磨損產(chǎn)生。問題三在分析重金屬污染的傳播特征時，可以結(jié)合大氣擴(kuò)散模型和水體對流擴(kuò)散模型的偏微分方程建立重金屬元素的對流擴(kuò)散模型。為求具體污染源的位置，在對擴(kuò)散模型簡化的基礎(chǔ)上，利用反演擴(kuò)散模型方程來求污染源。反演方程可以用遺

5、傳算法來實(shí)現(xiàn)，把問題轉(zhuǎn)化為建立一個非線性最優(yōu)化問題。以為例，其污染源坐標(biāo)為問題四要求我們分析已用模型的優(yōu)缺點(diǎn)，并改進(jìn)模型以研究更好的城市地質(zhì)環(huán)境的演變模式。通過搜集不同歷史時間不同位置的數(shù)據(jù)，利用時間序列分段線性擬合各個指標(biāo)值，較準(zhǔn)確地恢復(fù)任一歷史時間的土壤質(zhì)量指標(biāo)，根據(jù)重金屬元素濃度變化的大致趨勢，采用三次指數(shù)平滑法對未來各指標(biāo)變化進(jìn)行預(yù)測。關(guān)鍵詞：雙權(quán)重因子模糊綜合評判模型 因子分析 遺傳算法 反演 擴(kuò)散方程一、 問題提出隨著城市經(jīng)濟(jì)的快速發(fā)展和城市人口的不斷增加，人類活動對城市環(huán)境質(zhì)量的影響日顯突出。對城市土壤地質(zhì)環(huán)境異常的查證，以及如何應(yīng)用查證獲得的海量數(shù)據(jù)資料開展城市環(huán)境質(zhì)量評價，研

6、究人類活動影響下城市地質(zhì)環(huán)境的演變模式，日益成為人們關(guān)注的焦點(diǎn)。按照功能劃分，城區(qū)一般可分為生活區(qū)、工業(yè)區(qū)、山區(qū)、交通區(qū)及公園綠地區(qū)等，分別記為1類區(qū)、2類區(qū)、5類區(qū)，不同的區(qū)域環(huán)境受人類活動影響的程度不同?，F(xiàn)對某城市城區(qū)土壤地質(zhì)環(huán)境進(jìn)行調(diào)查。為此，將所考察的城區(qū)劃分為間距1公里左右的網(wǎng)格子區(qū)域，按照每平方公里1個采樣點(diǎn)對表層土（0-10 厘米深度）進(jìn)行取樣、編號，并用GPS記錄采樣點(diǎn)的位置。應(yīng)用專門儀器測試分析，獲得了每個樣本所含的多種化學(xué)元素的濃度數(shù)據(jù)。另一方面，按照2公里的間距在那些遠(yuǎn)離人群及工業(yè)活動的自然區(qū)取樣，將其作為該城區(qū)表層土壤中元素的背景值。附件1列出了采樣點(diǎn)的位置、海拔高度及

7、其所屬功能區(qū)等信息，附件2列出了8種主要重金屬元素在采樣點(diǎn)處的濃度，附件3列出了8種主要重金屬元素的背景值。現(xiàn)要求你們通過數(shù)學(xué)建模來完成以下任務(wù)：(1) 給出8種主要重金屬元素在該城區(qū)的空間分布，并分析該城區(qū)內(nèi)不同區(qū)域重金屬的污染程度。(2) 通過數(shù)據(jù)分析，說明重金屬污染的主要原因。(3) 分析重金屬污染物的傳播特征，由此建立模型，確定污染源的位置。(4) 分析你所建立模型的優(yōu)缺點(diǎn)，為更好地研究城市地質(zhì)環(huán)境的演變模式，還應(yīng)收集什么信息？有了這些信息，如何建立模型解決問題？二、 問題分析2.1對問題一的分析 題目要求統(tǒng)計(jì)出八種重金屬元素在城區(qū)的空間分布，需要考慮重金屬元素在X,Y,Z三維空間的濃

8、度分布特征。由附件給出的數(shù)據(jù)分析得，當(dāng)X,Y的坐標(biāo)確定之后，該點(diǎn)的海拔高度Z相應(yīng)確定下來。所以在統(tǒng)計(jì)重金屬的空間分布時，可以先考慮元素濃度和平面位置的關(guān)系，再通過海拔分布圖分析元素分布與海拔之間的關(guān)系。對這八種重金屬元素分布單獨(dú)考慮，利用附件中各元素在城區(qū)的樣本濃度，采用Cubic插值法采用Matlab軟件對各元素含量分布進(jìn)行空間分析，可以分別繪制出8種重金屬元素的空間等高線分布圖。當(dāng)考慮不同區(qū)域內(nèi)各元素的污染程度時，現(xiàn)行重金屬污染評價方法一般采用污染物濃度超標(biāo)賦權(quán)法。其中，地累積指數(shù)評價模型法1 ( Forstner，1989,Forstner et al.1993 ) 是Muller 提出

9、的一種能夠較好地表征土壤中重金屬富集程度的定量指標(biāo)，但其評價法只考慮了濃度而未考慮到重金屬元素本身的生物毒性。對于不同重金屬，因污染物個體的毒性級別不同，污染物濃度超標(biāo)賦權(quán)法有可能掩蓋某些低濃度有機(jī)組分的毒性作用。因此本文采取雙權(quán)重因子模糊綜合評判數(shù)學(xué)模型，將污染物濃度和毒性級別指數(shù)加權(quán)疊加，并作歸一化處理，得到某污染組分的權(quán)重公式考慮引入模糊綜合評價來進(jìn)行基于各元素濃度、毒性響應(yīng)的多因素綜合評價問題。結(jié)合環(huán)境土壤法給出的標(biāo)準(zhǔn)，對各區(qū)的綜合污染進(jìn)行進(jìn)一步精確的分析。2.2對問題二的分析 對于問題（2），分析重金屬污染的主要原因，我們主要要依據(jù)問題（1）中得出的各個功能區(qū)的不同重金屬元素的污染狀

10、況，得出某一重金屬元素主要是由哪一種污染物所帶來的。首先進(jìn)行每個區(qū)域內(nèi)八種重金屬元素的相關(guān)性分析，相關(guān)性比較顯著的重金屬元素可以認(rèn)為來自同一個污染源，由此可以大致確認(rèn)每個區(qū)域污染源的種類。因?yàn)槊總€區(qū)域中重金屬的種類眾多，相互之間又具有比較明顯的相關(guān)性，所以我們可以采用KMO法和Bartlett法對原始數(shù)據(jù)進(jìn)行主成分分析適宜性檢驗(yàn)的基礎(chǔ)上，進(jìn)一步考慮使用主成分分析法，由此分析出每個區(qū)域營銷比較顯著的若干種金屬，從而確認(rèn)出重金屬元素的主要原因。2.3對問題三的分析 基于一個基本假設(shè)，我們發(fā)現(xiàn)重金屬污染物的擴(kuò)散符合對流擴(kuò)散模型，則通過求解擴(kuò)散微分方程，可以得到一個重金屬污染物濃度的函數(shù)表達(dá)式，我們發(fā)

11、現(xiàn)，對于特定的重金屬污染物，如果假設(shè)出幾個污染源的位置以及相應(yīng)的排放強(qiáng)度，即可求出整個區(qū)域的污染物濃度。根據(jù)附件已經(jīng)給出的某一區(qū)域的污染物濃度實(shí)測值，我們可以檢驗(yàn)假設(shè)的污染源位置是否正確，通過不斷驗(yàn)證，一定可以找到幾個污染源的位置，它們可以使得求出的這種重金屬污染物的濃度分布與實(shí)際的分布最為接近，那么這幾個污染源的位置即就是我們需要的污染源的位置。那么原問題便轉(zhuǎn)化為了一個非線性最優(yōu)化的問題，在matlab軟件中，我們借助遺傳算法，可以通過不斷迭代，最終找出這樣幾個位置，滿足要求，進(jìn)而可以得出不同重金屬的污染源位置。 2.4對問題四的分析 前三個問題的求解模型可以比較合理地解決問題，但是也在一定

12、程度上存在局限性。由此我們可以進(jìn)一步收集該城區(qū)突然懲罰的信息，建立更優(yōu)化的模型。三、 模型的假設(shè)（1）假設(shè)所有數(shù)據(jù)都是真實(shí)可靠的；（2）假設(shè)污染源為重金屬含量最高的地方；（3）假設(shè)污染源的位置大概固定，不隨時間變化。（4）假設(shè)污染元素在水平面的擴(kuò)散在各個方向的擴(kuò)散速度相同。（5）重金屬污染物的擴(kuò)散方式為簡單的對流擴(kuò)散，符合對流擴(kuò)散方程；（6）污染物的濃度在海拔高度上大致都分布在0100以內(nèi)，我們只考慮方向上的擴(kuò)散；（7）由于擴(kuò)散到一定程度，擴(kuò)散的影響因素將會增多，我們假設(shè)存在理想最大擴(kuò)散距離；四、 符號說明表1元素As Cd Cr Cu Hg Ni Pb Zn 標(biāo)號/j12345678表2污染

13、因子實(shí)測值污染因子對應(yīng)于各級土壤重金屬環(huán)境質(zhì)量狀況的標(biāo)準(zhǔn)值因子毒性指數(shù)因子權(quán)重值 的相關(guān)系數(shù)重金屬濃度在三個方向的擴(kuò)散系數(shù) 重金屬濃度在三個方向的對流系數(shù) 污染源的排放強(qiáng)度污染源的理想最大擴(kuò)散距離函數(shù)與解原微分方程等效的新函數(shù)五、 模型建立與求解5.1問題一的模型建立與求解5.1.1重金屬空間分布的求解5.1.1.1模型建立與求解 根據(jù)采樣點(diǎn)各重金屬元素的含量分布，采用Cubic插值法，應(yīng)用Matlab軟件對海拔高度和8種重金屬的含量分布進(jìn)行空間分析，可匯出以下空間分布圖。 圖1：海拔分布圖圖2：As空間分布圖圖3：Cd空間分布圖圖4：Cr空間分布圖圖5：Cu空間分布圖圖6：Hg空間分布圖圖7

14、：Ni空間分布圖圖8：Zn空間分布圖 通過比較8種重金屬元素分布圖與功能分布圖，可以得出8種重金屬元素的分布情況，如下表3： 表3: 8種重金屬元素的分布元素主要分布區(qū)域工業(yè)區(qū)、西南方的公園綠地區(qū)、西南方和東南方的交通區(qū)工業(yè)區(qū)、西南方的交通區(qū)、公園綠地區(qū)、西南方和西北方的山區(qū)城區(qū)西南的工業(yè)區(qū)、少量分布在城區(qū)東北的工業(yè)區(qū)城區(qū)南部以及中部的工業(yè)區(qū)城區(qū)西南部的工業(yè)區(qū)和生活區(qū)的交接處城區(qū)西南部的交通區(qū)和工業(yè)區(qū)以及山區(qū)密集的地方城區(qū)南部的生活區(qū)以及工業(yè)區(qū)城區(qū)南部的生活區(qū)、交通區(qū)和工業(yè)區(qū)5.1.1.2模型評價運(yùn)用Cubic插值法，畫圖精準(zhǔn)，有助于空間濃度分析。5.1.2 不同區(qū)域內(nèi)重金屬的污染程度的求解雙

15、權(quán)重因子模糊綜合評判模型是基于各元素濃度、毒性響應(yīng)的雙因素綜合評價問題2。將污染物濃度和毒性級別指數(shù)加權(quán)疊加，并作歸一化處理，得到某污染組分的權(quán)重公式考慮結(jié)合環(huán)境土壤法給出的標(biāo)準(zhǔn)，對各區(qū)的綜合污染進(jìn)行進(jìn)一步精確的分析。5.1.2.1建立模型雙權(quán)重因子模糊綜合評判決策的方法與步驟如下：（1）設(shè)評判因子，代表8種金屬，評價集，代表國家標(biāo)準(zhǔn)中5個土壤等級，分別為清潔、尚清潔、輕污染、中污染和重污染。（2）各因素的隸屬度確定采用呈梯形分布的隸屬函數(shù)。評價單項(xiàng)污染因子 為隸屬度函數(shù)的拐點(diǎn)，用分段直線函數(shù)模擬，。其中為污染因子對應(yīng)于各級土壤重金屬環(huán)境質(zhì)量狀況的標(biāo)準(zhǔn)值， 表示污染因子實(shí)測值。當(dāng)時， 當(dāng)時，

16、公式（1） 當(dāng)時， （3）得到由評價（各金屬）在評價集上的隸屬度所組成的8*5的關(guān)系模糊矩陣。（4）為評價集對因子的權(quán)重分配構(gòu)成的權(quán)重向量， ，為第種污染因子所對應(yīng)的權(quán)重。確定權(quán)重時，以， 表示因子權(quán)重值，為因子毒性指數(shù)，將污染物濃度和毒性指數(shù)加權(quán)疊加，作歸一化處理，得權(quán)重公式： 公式（2） 公式（3）(5)據(jù)矩陣運(yùn)算法則，確定了一個模糊映射，把評價因子集上的權(quán)重向量 映射到評價集，得模糊子集，即評價等級，據(jù)最大隸屬度原則確定評價等級。 5.1.2.2模型求解（1） 先對數(shù)據(jù)進(jìn)行預(yù)處理。將公式（1）中的用每個區(qū)8種重金屬元素含的均值替代，對五個區(qū)8種重金屬元素分別進(jìn)行處理，結(jié)果如下表4 ：表

17、4：各區(qū)各個重金屬含量的均值區(qū)域生活區(qū)工業(yè)區(qū)山區(qū)主干道路區(qū)公園綠地區(qū)As (g/g)6.2704557.2513894.0440915.7080434786.263714286Cd (ng/g)0.2899610.3931110.152320.3600140.280543Cr (g/g)69.0184153.4091738.959758.0539130443.636Cu (g/g)49.40318127.535817.3172762.2149275430.19171429Hg (ng/g)0.0930410.6423550.0409560.4468230.114992Ni (g/g)18.34

18、22719.8116715.4537917.6171014515.28971429Pb (g/g)69.1063693.0408336.5559163.534202960.70857143Zn (g/g)237.0086277.927573.29424242.8548551154.2422857（2）搜集資料，確定等級標(biāo)準(zhǔn)。利用國家土壤環(huán)境質(zhì)量標(biāo)準(zhǔn)(GB 156181995)確定的區(qū)域土壤評價，選取甘肅省土壤重金屬元素背景值和臨界含量作為的標(biāo)準(zhǔn)參照；考慮重金屬生物體的毒性作用時選取Hakanson的重金屬毒性響應(yīng)系數(shù)所確定的毒性指數(shù)作為的標(biāo)準(zhǔn)參照，最終確立土壤重金屬污染程度分級標(biāo)準(zhǔn)及生物毒性指

19、數(shù)，如下表5：表5：土壤重金屬污染程度分級標(biāo)準(zhǔn)及生物毒性指數(shù)(mg·kg ) 3元素/等級1清潔2尚清潔3輕污染4中污染5重污染毒性指數(shù)As (g/g)82.5695.43124.26164.28201.383Cd (ng/g)0.12040. 25230.61.422Cr (g/g)74.8899.541503505005Cu (g/g)28.3740.631202804004Hg (ng/g)0.0920.25920.451.051.51Ni (g/g)77.5182.5694.21102.56154.253Pb (g/g)23.3536.091503505004Zn (g/g)

20、83.68116.752405608006 （3）根據(jù)公式（1）運(yùn)用Matlab軟件計(jì)算各元素隸屬值，得到各區(qū)域的關(guān)系模糊矩陣,以生活區(qū)為例，其他區(qū)域詳細(xì)見附件。=（4）根據(jù)公式（2）、（3）運(yùn)用Matlab軟件，將污染物濃度和毒性級別指數(shù)加權(quán)疊加并作歸一化處理，得到某污染組分的權(quán)重如下表6：表:6：該城市各類區(qū)各參評因子的權(quán)重值區(qū)域生活區(qū)工業(yè)區(qū)山區(qū)主干道路區(qū)公園綠地區(qū)As (g/g)0.01290.01450.01070.01120.0125Cd (ng/g)0.36940.24850.19930.12860.4113Cr (g/g)0.09020.09220.08230.1410.0652

21、Cu (g/g)0.0530.05920.06590.07410.0603Hg (ng/g)0.37170.48010.54380.53710.3701Ni (g/g)0.05430.06840.04350.04590.0436Pb (g/g)0.03780.04150.04520.03920.0362Zn (g/g)0.07770.07850.06350.08010.0568（5）模糊矩陣復(fù)合運(yùn)算，繼續(xù)用Matlab軟件處理得到各評價樣品對評價等級的隸屬度,又根據(jù)最大隸屬度原則，確定各樣品的的污染程度，此即為土壤環(huán)境質(zhì)量分級，最終結(jié)果該城區(qū)內(nèi)不同區(qū)域8種主要重金屬元素的污染程度分析結(jié)果如下表

22、：表 7：模糊綜合評價結(jié)果區(qū)域隸屬度污染等級生活區(qū)0.5267860.4075730.132641001（清潔）工業(yè)區(qū)0.17510.1686140.5731770.16600903（污染）山區(qū)0.960770.0932460.000185001（清潔）交通區(qū)0.19810.1814090.6769760.00071503（污染）公園綠地區(qū)0.4918470.5056420.058511002（輕污染） （5）綜合分析：根據(jù)表7可以看出，生活區(qū)和山區(qū)污染度最低，屬于標(biāo)準(zhǔn)里的清潔范疇；公園綠地區(qū)輕微污染，而工業(yè)區(qū)和交通區(qū)污染最嚴(yán)重，其中又以交通區(qū)的污染程度為最；而縱觀表6各區(qū)的隸屬度值的大小，我

23、們可以得出以下結(jié)論：對生活區(qū)的分析：污染排序Zn>Cu>Hg>Pb>Cd>Cr>As,其中Ni元素尚清潔對工業(yè)區(qū)的分析：污染排序Hg>Cu>Zn>Pb>Cd>As>Cr,其中Ni元素尚清潔對山區(qū)的分析：山區(qū)沒有重金屬污染，尚清潔對交通區(qū)的分析：污染排序Hg>Cu>Zn>Cd>Pb>Cr>As,其中Ni元素尚清潔對公園綠地區(qū)的分析：污染排序Hg>Zn>Cu>Cd>Pb>As,其中Cr、Ni元素尚清潔。5.1.2.3模型評價優(yōu)點(diǎn)：雙權(quán)重因子模糊綜合評判模型考慮濃

24、度和毒性雙重因子作為權(quán)重，查找資料引入毒性這一影響污染程度的指標(biāo)，使不同評價對象間可比性增強(qiáng)，采用評價分值，即利用評價向量的分量形成權(quán)重，通過確定各等級對應(yīng)的分值，對各等級的分值進(jìn)行加權(quán)平均，得到評價分值。從評價分值可以直觀的看出各采樣點(diǎn)土壤重金屬環(huán)境質(zhì)量的優(yōu)劣。缺點(diǎn)：雙權(quán)重因子模糊綜合評判模型計(jì)算冗長，耗時較長，所以對于采樣樣本數(shù)目較多的情況并不可行，并且沒有考慮人為因素的影響。適用：建議對于少量樣本測量時進(jìn)行使用。5.2問題二的模型建立與求解 城市土壤重金屬來源于成土母質(zhì)和人類活動，同一來源的重金屬之間存在著相關(guān)性，根據(jù)相關(guān)性可以判斷土壤重金屬污染來源是否相同。如果重金屬之間存在顯著的正相

25、關(guān)，則其來源可能相同，否則來源可能不止一個。5.2.1建立模型通過查閱資料45,我們采用多元統(tǒng)計(jì)數(shù)學(xué)方法之一的因子分析法，它根據(jù)多個實(shí)測變量（重金屬元素）之間的相關(guān)性，運(yùn)用數(shù)學(xué)變換，將多個變量轉(zhuǎn)換為少數(shù)幾個線性不相關(guān)的綜合指標(biāo)（不相關(guān)的幾個主要原因），從而簡化數(shù)據(jù)處理，其目的在于對大量觀測數(shù)據(jù)，用較少的有代表性的因子（目標(biāo)）來說明眾多變量所提取的主要信息，提示出多個變量的因果關(guān)系。具體方法如下：因子分析從變量的相關(guān)矩陣出發(fā)將一個m維隨機(jī)向量X（數(shù)據(jù)總庫）分解成低于m個且有代表性的公因子和一個特殊的m維向量，是公因子數(shù)取得最佳的個數(shù)，從而使對m維隨機(jī)向量的研究轉(zhuǎn)換成較小個數(shù)的公因子研究。設(shè)一共有

26、n（n=319）個樣本，n個指標(biāo)構(gòu)成樣本空間X： ；。之后我們對數(shù)據(jù)做如下處理步驟：(1)原始數(shù)據(jù)的標(biāo)準(zhǔn)化，標(biāo)準(zhǔn)化公式為：其中為第i個樣本第j個指標(biāo)值，而和分別為j指標(biāo)的均值和標(biāo)準(zhǔn)差。標(biāo)準(zhǔn)化目的在于消除不同變量量綱的影響，而且標(biāo)準(zhǔn)化轉(zhuǎn)化不會轉(zhuǎn)變變量的相關(guān)系數(shù)。（2）KMO檢驗(yàn)，判斷其是否符合相關(guān)性分析條件。KMO值定義為其判斷標(biāo)準(zhǔn)是：表8：KMO檢驗(yàn)法判斷標(biāo)準(zhǔn)0.9<KMO非常適合0.8<KMO<0.9適合0.7<KMO<0.8一般0.6<KMO<0.7不適合KMO<0.6非常不適合若符合上表中的標(biāo)準(zhǔn)，則計(jì)算標(biāo)準(zhǔn)化數(shù)據(jù)的相關(guān)系數(shù)陣，相關(guān)性公式如下

27、：其中X,Y表示不同的重金屬元素。 （3）進(jìn)行正交變換，使用方差最大法，其目的是使因子載荷兩極分化，且旋轉(zhuǎn)后因子仍正交。（4）確定因子個數(shù)，計(jì)算因子得分，進(jìn)行統(tǒng)計(jì)分析。5.2.2模型求解依據(jù)題中給出的數(shù)據(jù)，對此城區(qū)土壤單點(diǎn)樣重金屬元素含量的數(shù)據(jù)標(biāo)準(zhǔn)化處理。然后分步得出結(jié)論：首先，算出KMO=0.778>0.7（符合因子分析法的條件），再給出8種重金屬元素原始含量數(shù)據(jù)的相關(guān)系數(shù)矩陣。各元素相關(guān)系數(shù)矩陣如表5所示，可見As和Ni、Pb的相關(guān)性分別為0.317、0.290（一般），而As和Cd、Zn的相關(guān)性更低，As和Hg幾乎無相關(guān)性，同樣可以看出其他兩兩之間的相關(guān)性，Cd和Pb的相關(guān)性較好達(dá)

28、到0.660，Cr和Ni、Cu和Pb、Cu和Cr相關(guān)性分別達(dá)到0.716、0.520和0.532，而Hg只與Cu的相關(guān)性為0.417（一般），Hg與其它元素的相關(guān)性很弱，剩下的元素之間的相關(guān)性基本上都一般。從成因來分析，相關(guān)性較好的元素可能在成因和來源上有一定的關(guān)聯(lián)。表9：8種重金屬元素相關(guān)性系數(shù)矩陣相關(guān)系數(shù)AsCdCrCuHgNiPbZnAs1.0000.2550.1890.1600.0640.3170.2900.247Cd0.2551.0000.3520.3970.2650.3290.6600.431Cr0.1890.3521.0000.5320.1030.7160.3830.424Cu0

29、.1600.3970.5321.0000.4170.4950.5200.387Hg0.0640.2650.1030.4171.0000.1030.2980.196Ni0.3170.3290.7160.4950.1031.0000.3070.436Pb0.2900.6600.3830.5200.2980.3071.0000.494Zn0.2470.4310.4240.3870.1960.4360.4941.000然后，利用相關(guān)系數(shù)矩陣在SPSS13.0統(tǒng)計(jì)軟件中計(jì)算出相應(yīng)的因子的特征值和累值貢獻(xiàn)率，得到表10，如下：表10：旋轉(zhuǎn)前后因子的特征值和累值貢獻(xiàn)率因子旋轉(zhuǎn)前旋轉(zhuǎn)后總的特征值占總百分量的百

30、分率累計(jì)貢獻(xiàn)率總的特征值占總百分量的百分率累計(jì)貢獻(xiàn)率13.56044.50044.5002.05925.73325.73321.15014.37758.8771.76722.08747.82030.96012.06370.9411.21715.21563.03540.7689.59680.5371.02612.82475.85850.5787.22087.7560.95211.89887.756在積累方差為87.756%（>80%）符合因子分析法的前提下，分析得到5個主因子，這5個因子共占了源資源的87.756%，這5個因子也就是我們要分析得到的重金屬污染的五個主要原因。從表中可以看出主

31、因子1和2占的貢獻(xiàn)率達(dá)到25.733%和22.087%，顯然是重金屬污染的最主要的原因，因子3、4、5（貢獻(xiàn)率在11至16之間）為重要原因。最后，因子分析的主要目的是將具有相近的因子載荷的各個變量置于一個公因子之下，正交方差最多旋轉(zhuǎn)使每一個主因子只與最少個數(shù)的變量有關(guān)系，而使足夠多的因子負(fù)荷均很小，以便對因子的意義作出合理的解釋，輸出見表11。表11：各元素與各因子的關(guān)聯(lián)度元素因子1因子2因子3因子4因子5As0.1310.1540.0230.9700.076Cd0.1700.8770.0820.1020.132Cr0.8820.2090.0020.0140.146Cu0.6140.3620.

32、505-0.0220.029Hg0.0130.1340.9530.0300.083Ni0.8640.0890.0190.2220.196Pb0.1950.8320.1910.1210.208Zn0.2710.2850.1010.0910.904從表7數(shù)據(jù)可得，元素Ni、Cr和Cu的對因子1的關(guān)聯(lián)度分別達(dá)到0.882、0.864和0.614，是此因子產(chǎn)生的主要污染物，從相關(guān)性表里也可以看出Ni和Cr的相關(guān)性好（相關(guān)性系數(shù)達(dá)到0.716），查找資料，可發(fā)現(xiàn)汽車燃燒的石化燃料所產(chǎn)生的尾氣中含有Ni、Cr和Cu，那么因子1極有可能是汽車燃燒的石化燃料所產(chǎn)生的尾氣。元素Cd和Pb對于因子2的關(guān)聯(lián)度分別達(dá)

33、到0.877和0.832，是此因子產(chǎn)生的主要污染物，從相關(guān)性表里也可以看出Cd和Pb的相關(guān)性好（相關(guān)性系數(shù)達(dá)到0.660），查找資料6,可發(fā)現(xiàn)工業(yè)中電鍍、采礦、冶煉、燃料、電池和化學(xué)工業(yè)中排放的廢氣、廢水和廢渣以及汽車排放的尾氣中都含有比較多的Cd和Pb，由于因子1是汽車燃燒的石化燃料所產(chǎn)生的尾氣，所以可以判斷因子2主要是工業(yè)中電鍍、采礦、冶煉、燃料、電池和化學(xué)工業(yè)中排放的廢氣、廢水和廢渣。元素Hg對于因子3關(guān)聯(lián)度達(dá)到0.953，是此因子產(chǎn)生的主要污染物，查找資料6,可發(fā)現(xiàn)工業(yè)中氯堿、塑料、電子等工業(yè)排放的汞廢水中含有大量的Hg，所以可以判斷因子3為工業(yè)中氯堿、塑料、電子等工業(yè)排放的汞廢水。元

34、素As對于因子4的關(guān)聯(lián)度達(dá)到0.970，查找資料6,可發(fā)現(xiàn)生活中使用的殺蟲劑、殺菌劑、滅鼠劑和除草劑等用品產(chǎn)生的廢氣煙塵中含有As，而且其它元素對于此因子的關(guān)聯(lián)度偏低，也正好驗(yàn)證了因子4很有可能是生活中使用的殺蟲劑、殺菌劑、滅鼠劑和除草劑等用品產(chǎn)生的廢氣煙塵。元素Zn的對于因子5的關(guān)聯(lián)度達(dá)到0.904，其它的都比較低，查找資料6,可發(fā)現(xiàn)汽車輪胎磨損產(chǎn)物以及生活中燃燒煤制品和垃圾焚燒產(chǎn)生的廢氣、煙塵，所以因子5為汽車輪胎磨損產(chǎn)物以及生活中燃燒煤制品和垃圾焚燒產(chǎn)生的廢氣、煙塵。綜上所述，通過分析數(shù)據(jù)可以得出此城區(qū)重金屬污染的主要原因有以下五點(diǎn)：a) 汽車燃燒的石化燃料所產(chǎn)生的尾氣，主要引起重金屬N

35、i、Cr、Cu的污染；b) 工業(yè)中電鍍、采礦、冶煉、燃料、電池和化學(xué)工業(yè)中排放的廢氣、廢水和廢渣，主要引起重金屬Cd、Pb的污染；c) 工業(yè)中氯堿、塑料、電子等工業(yè)排放的汞廢水，主要引起重金屬Hg的污染；d) 生活中使用的殺蟲劑、殺菌劑、滅鼠劑和除草劑等用品產(chǎn)生的廢氣煙塵,主要引起重金屬As的污染；e) 汽車輪胎磨損產(chǎn)物以及生活中燃燒煤制品和垃圾焚燒產(chǎn)生的廢氣、煙塵，主要引起重金屬Zn的污染。5.3問題三的模型建立與求解5.3.1.污染物傳播特征的前提假設(shè)重金屬元素傳播途徑主要有水流傳播、土壤滲透、大氣擴(kuò)散等。在城區(qū)中我們主要考慮污染以大氣擴(kuò)散。為了研究方便，作如下假設(shè)（重新申明）：（1）重金

36、屬污染物的擴(kuò)散方式為簡單的對流擴(kuò)散，符合對流擴(kuò)散方程；（2）污染物的濃度在海拔高度上大致都分布在0100以內(nèi)，我們只考慮方向上的擴(kuò)散；（3）由于擴(kuò)散到一定程度，擴(kuò)散的影響因素將會增多，我們假設(shè)存在理想最大擴(kuò)散距離；5.3.2. 模型建立我們根據(jù)對流擴(kuò)散方程建立模型7：其中為重金屬濃度，在三個方向的擴(kuò)散系數(shù),為對流系數(shù)。 模型一：由于附件中給的數(shù)據(jù)都是離散化的點(diǎn)樣本，我們將空間的平面看成有限個點(diǎn)的網(wǎng)格組成，于是我們將原來的方程作離散化處理，用差商代替原方程的微商，通過求解離散點(diǎn)上的解來近似原方程的解，轉(zhuǎn)化如下：根據(jù)對流擴(kuò)散假設(shè)，污染源中心處重金屬濃度會隨著時間慢慢減小，而且其濃度隨時間的變化率與

37、周圍相比比較大，那么原問題就轉(zhuǎn)化為求等步長時間內(nèi)濃度向方向的變化率問題。5.3.3.模型求解上式c為重金屬濃度，Ex，Ey，Ez為三個方向上的擴(kuò)散系數(shù)，ux，uy，uz為對流系數(shù)，Mi為污染源的排放強(qiáng)度。R=f（x，y，z），L為污染源的理想最大擴(kuò)散距離。簡化上式可得求解上式微分方程，對解析解進(jìn)一步變換，結(jié)果如下：則上式求c的問題可以轉(zhuǎn)化為求v函數(shù)的問題對結(jié)果再一步變化可得：由邊界值可解出：其中：最后求得v函數(shù)的解析解為：我們已經(jīng)有了各個樣本點(diǎn)的濃度數(shù)據(jù)，由上式可知，可由， ，我們先研究一種特定的金屬污染物，則其為定值，只要我們假設(shè)出它的污染源的位置，即，那么我們就可以得到所有樣本點(diǎn)的濃度數(shù)據(jù)

38、，如果這些樣本點(diǎn)的濃度數(shù)據(jù)和實(shí)際的符合度很好，那么就說明我們的假設(shè)是正確的，那么原問題就轉(zhuǎn)化成了下面的非線性最優(yōu)化問題：則滿足上式的即就是我們需要的污染源的位置，再轉(zhuǎn)化為原始三維坐標(biāo)即可通過matlab的工具箱，我們采用遺傳算法可以大大提升算法效率和準(zhǔn)確度，得到不同金屬污染的污染源位置坐標(biāo)如下表：表12：各重金屬污染的污染源位置坐標(biāo)表坐標(biāo)區(qū)域編號濃度70568348371427.4112644149434322755.4112734103443232517.3423664979046392114.320261758629320918814405180321523280398.426416650

39、814110736.7318413117218838059.612017717642276328823.15537586431512511.821956973487032109.9522624481827321928.97932843112413617018738109215338258.931426215129663282199212113055223127.9115007553570313831.5259811805117333019.713855334579118354.9817008477582319438.35209831586293328727.541262516259661276

40、149.52169477487412197104.9226244818273219105.315.4問題四的解答5.4.1優(yōu)缺點(diǎn)分析 本文比較有特點(diǎn)的幾個模型中，優(yōu)點(diǎn)明顯的有：（1） 雙權(quán)重因子模糊綜合評判模型考慮了環(huán)境質(zhì)量界限模糊性的特性，打破用一個確定性指標(biāo)來評價環(huán)境質(zhì)量的缺陷，具有一定的合理性和優(yōu)越性，并且結(jié)合了不同元素的毒性響應(yīng)系數(shù)，對污染程度的判斷更合理。（2） 反演擴(kuò)散模型比較合理地模擬重金屬在該城區(qū)內(nèi)的擴(kuò)散。通過對模型的簡化，比較簡單的將模型轉(zhuǎn)換為二維偏微分方程。由于遺傳算法有很好的的全局優(yōu)化性能，減少了陷入局部極值的風(fēng)險(xiǎn)，特別是對非線性優(yōu)化問題有很強(qiáng)的適應(yīng)性。傳統(tǒng)的反演方法在源

41、強(qiáng)和濃度的關(guān)系上，只能給出一種虛擬的聯(lián)系，而遺傳算法可以給出較為直接的關(guān)系，從而為污染預(yù)報(bào)和污染的控制提供了基礎(chǔ)。用遺傳算法來反演，方法簡單且通用性強(qiáng)。而缺點(diǎn)則有： 由于不考慮模型里的對流項(xiàng)以及海拔上面的影響，導(dǎo)致其自身的確定誤差較大，其次遺產(chǎn)算法本身就是一個枚舉算法，算法在運(yùn)行中具有隨機(jī)性，導(dǎo)致其求出的污染點(diǎn)的分布不是完全符合，但具有穩(wěn)定的收斂性，是可以接受的。5.4.2 城市地質(zhì)環(huán)境的演變模式模型5.4.2.1 所需要搜集的數(shù)據(jù)查閱相關(guān)的數(shù)據(jù)8，得到了有關(guān)地層的幾個重要性指標(biāo)：土壤重金屬元素含量，pH值，有機(jī)質(zhì)及養(yǎng)分含量，不同種類無機(jī)礦物的含量，N，P等營養(yǎng)元素的含量，農(nóng)藥、有機(jī)污染物殘留

42、含量等多方面影響土壤地質(zhì)環(huán)境的因素。要能研究地質(zhì)演化的信息，必須搜集不同歷史時間的數(shù)據(jù)。對于數(shù)據(jù)的要求為，將所要研究的區(qū)域數(shù)據(jù)的取樣點(diǎn)越密集，研究結(jié)果越準(zhǔn)確。5.2.2.2 時間序列分析根據(jù)已有研究8，在城市發(fā)展對土壤質(zhì)量影響的分析、規(guī)律統(tǒng)計(jì)和建模方面，主要有三種模型：基于基本統(tǒng)計(jì)理論的回歸分析和自回歸趨勢分析，灰色關(guān)聯(lián)分析和灰色建模、基于GIS的圖形空間分析模型。我們可以利用時間序列分段線性擬合各個指標(biāo)值，較準(zhǔn)確地恢復(fù)任一歷史時間的土壤質(zhì)量指標(biāo)，根據(jù)重金屬元素濃度變化的大致趨勢，采用三次指數(shù)平滑法對未來各指標(biāo)變化進(jìn)行預(yù)測。六、 參考文獻(xiàn)1.張江華等，內(nèi)梅羅指數(shù)和地質(zhì)累積指數(shù)在土壤重金屬評價中

43、的差異探討，環(huán)保與分析，第31 卷第8期：P43-44，2010 2. 張志斌，王軻，蒲瑞豐，基于雙權(quán)重因子模糊綜合模型的金昌市土壤重金屬污染評價，西北師范大學(xué)學(xué)報(bào)，2009 3. 國家土壤環(huán)境質(zhì)量標(biāo)準(zhǔn)(GB 156181995)4.王雄軍、賴健清、魯艷紅、李德勝、周繼華、王建武，基于因子分析法研究太原市土壤重金屬污染的主要來源，生態(tài)環(huán)境，第17卷第2期：671-676,2008年；5.史貴濤、陳振樓、許世遠(yuǎn)、王利、張菊、李海雯、李麗娜，上海城市公園土壤及灰塵中重金屬污染特征，環(huán)境科學(xué)，第28卷第2期：238-242，2007年；6.百度百科“砷” 百度百科“鎘污染” 百度百科“鉻污染” 百度

44、百科“銅污染” 百度百科“汞污染” 百度百科“鎳污染” 百度百科“鉛污染” 百度百科“鋅污染” 2011/9/10；7.閔濤等，二維穩(wěn)態(tài)對流-擴(kuò)散方程反演的迭代算法，水動力學(xué)研究與進(jìn)展 2007.118.吳新民等，長江三角洲土壤質(zhì)量演變趨勢分析，地理與地理信息科學(xué)，第51 卷第2 期：P88-89，2008七、 附件7.1作圖程序clcclearclose allx=xlsread('2011a.xls','B4:B322');y=xlsread('2011a.xls','C4:C322');z=xlsread('2011

45、a.xls','D4:D322');c=xlsread('2011a.xls','E4:E322');M=xlsread('2011a.xls',2,'B4:I322');xi=linspace(min(x),max(x),150);yi=linspace(min(y),max(y),150);xi,yi=meshgrid(xi,yi);for i=1:8Mi(:,:,i)=griddata(x,y,M(:,i),xi,yi);endci=griddata(x,y,c,xi,yi);zi=griddata(

46、x,y,z,xi,yi,'cubic');marker='*','o','s','','p'color='k','r','y','c','b'mat='As','Cd','Cr','Cu','Hg','Ni','Pb','Zn'str='等高線','生活區(qū)',&#

47、39;工業(yè)區(qū)','山林區(qū)','交通區(qū)','綠地區(qū)'for j=1:8%等高線圖形 a=M(:,j); figure contourf(xi,yi,Mi(:,:,j); % set(h,'ShowText','on','TextStep',get(h,'LevelStep')*2) % clabel(C,h,0:10:50,50:50:300,300:100:500); title('金屬 ',matj, ' 二維等高線分布圖') xlabel(

48、'X') ylabel('Y') colormap summer colorbar grid on hold on for i=1:5 loc=c=i; plot(x(loc),y(loc),markeri,'markerfacecolor',colori,'MarkerEdgeColor',colori); end legend(str,'location','best') end%三維體圖figuresurf(xi,yi,zi);colormap hsvtitle('三維圖立體（顏色條表

49、示分類）')xlabel('X')ylabel('Y')colorbarhidden offhold onfor i=1:5 loc=c=i; plot3(x(loc),y(loc),z(loc),markeri,'markerfacecolor',colori);endstr1='三維圖'legend(str,'location','best')7.2模糊關(guān)系矩陣程序clcclearx=xlsread('x.xlsx');s=xlsread('s.xlsx');for dtt=1:5 for i=1:8 for stn=1:5 if stn=1 if x(i,dtt)<=s(i,stn) A(i,stn,dtt)=1; elseif x(i,dtt)>=s(i,stn)&&x(i,dtt)<=s(i,stn+1) A(i,stn,dtt)=(s(i,2)-x(i,dtt)/(s(i,2)-s(i,1); elseif x