2011年高教社杯全国大学生数学建模竞赛优秀论文.doc

2011高教社杯全国大学生数学建模竞赛承 诺 书我们仔细阅读了中国大学生数学建模竞赛的竞赛规则.我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。我们知道，抄袭别人的成果是违反竞赛规则的, 如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。我们参赛选择的题号是（从A/B/C/D中选择一项填写）： 我们的参赛报名号为（如果赛区设置报名号的话）： 所属学校（请填写完整的全名）： 参赛队员 (打印并签名) ：1. 2. 3. 指导教师或指导教师组负责人 (打印并签名)： 日期： 年月 日赛区评阅编号（由赛区组委会评阅前进行编号）：2011高教社杯全国大学生数学建模竞赛编 号 专 用 页赛区评阅编号（由赛区组委会评阅前进行编号）：赛区评阅记录（可供赛区评阅时使用）：评阅人评分备注全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：城市表层土壤重金属污染分析摘要本文利用克里格插值法，多因子综合指数评价法，相关分析，因子分析，聚类分析，Grid data插值模型分别进行了空间分析，污染程度分析，污染原因的分析以及确定污染源的位置。针对问题一，本文利用Sufer软件，通过克里格插值法画出8种重金属元素的等值图，并利用单污染指数评价法和多因子综合指数评价法对分别对同功能区不同采集点和不同功能区重金属综合污染程度进行分析，得出1. 不同的人类活动造成城市土壤中不同类型的重金属积累2. 工业区，交通区，生活区，公园绿地区属于重污染区，山区属于轻度污染区3.不同功能区污染指数的排序为：工业区交通区生活区公园绿地区山区。针对问题二，本文首先对8种重金属元素进行相关分析，得出Cd与Pb，Cr与Cu，Cr与Ni，Cu与Pb之间相关性显著，可能具有同源关系或复合关系。然后对8种重金属元素进行因子分析，得出Cr Ni Cu在因子1有较高的载荷，As Cd Pb Zn在因子2上有较高的载荷，可能主要来自类似的一个源，Hg在因子3上远远高于在其他元素在该因子的载荷，可能由单一污染源产生。 最后通过对两个模型综合分析，得出Cr Cu Ni在两种分析方法中都具有很强的关联性，分析污染产生的原因可能与电镀，合金制造等重工业有关，As Cd Pb Zn Cu在两种分析方法中呈现出一些差异，因此其污染来源情况比较复杂，可能与尾气排放，人类活动丢弃物以及农用废物有关，Hg在两种分析方法中均表现出独立性，因此伴随Hg产生的其他重金属元素则相对较少，主要来源可能是科学仪器制造，医学仪器制造工厂等对Hg有特定需求的工厂。针对问题三，本文通过对污染传播特征的研究，发现污染源的重金属浓度最高，并向周围呈递减趋势扩散。本文利用聚类分析找出8种重金属元素浓度的所有极大异常点，通过其传播特性可以预见这些极大异常点距离污染源是非常接近，从而对污染源的确定提供依据。然后我们又采用Grid data插值模型，运用局部搜索方法，成功地搜索出了每种重金属元素的污染源位置并与极大异常点坐标接近，证明了该方法的可行性。最后统一各元素污染源坐标得出综合污染源有8处，单一污染源有9处，并用逐步多元回归模型对部分点进行了检验，验证了解的准确性。针对问题四，本文采用的Grid data插值法优点在于能够快速定位污染源，但是由于在搜索过程中考虑到计算量巨大，所以搜索步长为100，由此产生的误差较大。为更好地研究城市地质环境的演变模式，还应考虑土壤类型、土地利用、道路等因素对重金属的影响。关键词：克里格插值；综合指数评价；相关分析；因子分析；GridData插值目录一、问题的重述与分析2二、模型的假设3三、符号说明4四、模型的建立与求解44.1重金属元素的空间分布及污染程度44.1.1 重金属元素的空间分布44.1.2 不同区域重金属的污染程度的分析54.2 重金属污染的原因分析74.2.1 Pearson相关系数74.2.2 相关系数的求解与分析84.2.3 因子分析84.2.4 因子分析的求解与结论94.3 污染源位置的确定104.3.1 聚类分析104.3.2 Grid data插值模型114.3.3 插值模型的建立与求解114.4 模型的检验134.5 模型优缺点分析及讨论14五、结论14六、参考文献15一、问题的重述与分析1.问题的重述随着城市经济的快速发展和城市人口的不断增加，人类活动对城市环境质量的影响日显突出。对城市土壤地质环境异常的查证，以及如何应用查证获得的海量数据资料开展城市环境质量评价，研究人类活动影响下城市地质环境的演变模式，日益成为人们关注的焦点。按照功能划分，城区一般可分为生活区、工业区、山区、主干道路区及公园绿地区等，分别记为1类区、2类区、5类区，不同的区域环境受人类活动影响的程度不同。现对某城市城区土壤地质环境进行调查。为此，将所考察的城区划分为间距1公里左右的网格子区域，按照每平方公里1个采样点对表层土（010 厘米深度）进行取样、编号，并用GPS记录采样点的位置。应用专门仪器测试分析，获得了每个样本所含的多种化学元素的浓度数据。另一方面，按照2公里的间距在那些远离人群及工业活动的自然区取样，将其作为该城区表层土壤中元素的背景值。附件1列出了采样点的位置、海拔高度及其所属功能区等信息，附件2列出了8种主要重金属元素在采样点处的浓度，附件3列出了8种主要重金属元素的背景值。现要求你们通过数学建模来完成以下任务：(1) 给出8种主要重金属元素在该城区的空间分布，并分析该城区内不同区域重金属的污染程度。(2) 通过数据分析，说明重金属污染的主要原因。(3) 分析重金属污染物的传播特征，由此建立模型，确定污染源的位置。(4) 分析你所建立模型的优缺点，为更好地研究城市地质环境的演变模式，还应收集什么信息？有了这些信息，如何建立模型解决问题？2.问题的分析本文将所有问题细化为一下几个步骤：(1)利用软件画出8种重金属污染元素的空间分布图；(2)建立综合评价指标对污染程度进行评价；(3)对数据进行统计分析，研究各元素之间的相关性，从而分析污染原因；(4)运用插值法解出污染源坐标，并用逐步回归进行检验；二、模型的假设1.提供的数据真实有效；2.8种重金属元素能够全面反映某一区域污染程度，不受其他元素影响；3.插值方法所产生的模型误差在可承受范围之内。三、符号说明符号意义单位第种污染物的污染分指数-实测浓度-评价标准-土壤污染综合指数-土壤污染物中最大的污染分指数-各污染分指数的算术平均数-其余未说明的符号会在文章中进行说明四、模型的建立与求解4.1重金属元素的空间分布及污染程度4.1.1 重金属元素的空间分布本文利用Surfer软件，通过克里格插值法画出8种重金属元素的标准化数据的等值图。 （a）As浓度等值图 （b） Cd 浓度等值图 (c)Cr浓度等值图 (d)Cu浓度等值图 (e) Hg浓度等值图 (f) Ni浓度等值图 (g) Pb浓度等值图 (h) Zn浓度等值图图1 8种重金属元素浓度等值图4.1.2不同区域重金属的污染程度的分析4.1.2.1评价方法与分级标准1.单污染指数评价法 以土壤单项污染物的实测值与评价标准相比，用以表示土壤中该污染物的污染程度。(1)式中，为第种污染物的污染分指数；为其实测浓度；为其评价标准。2.多因子综合指数评价法各类土壤一般为多种重金属所污染，因而土壤污染评价多应用综合指数法进行污染综合评价。(2)式中，为土壤污染综合指数；为土壤污染物中最大的污染分指数；为各污染分指数的算术平均数。3.分级标准下表是国内现行的指数分级标准表1 土壤综合污染指数分级标准污染指数分级污染等级污染水平PN0.71安全清洁0.7PN12警戒线尚清洁1PN3轻污染土壤作物已受污染235重污染土壤作物已受严重污染4.1.2.2结论与分析（1） 同一功能区不同元素污染分析将数据代入公式(1)得到各个功能区的8种污染物的污染分指数和综合指数。表2 8种元素污染分指数和综合指数功能区11.7417922.2304722.2264073.7426772.658291.4912412.2292293.434922.0142723.0239371.7228839.66177818.352911.6107023.0013194.02796931.123361.171691.2567631.3119181.1701721.2564061.1792251.06223641.5855662.7693461.8727094.71327712.766351.4322852.0494923.51964351.7399192.1580341.4076152.2872613.2854671.2430751.9583372.235389结论：1.从污染程度上说，山区的污染指数明显低于其他4个区，说明人类活动对城市土壤重金属分布有重要影响2.从污染类型上说，工业区以Hg Cu的积累为特征, 而生活区和公园绿地区则以Cu Zn的积累为特征, 交通区以Hg Cu Zn 的积累为特征，即不同的人类活动造成城市土壤中不同类型的重金属积累（2） 不同功能区的综合污染分析将数据代入公式(2)算得5种功能区的综合污染指数得到表3表3 5种功能区的综合污染指数功能区生活区工业区山区交通区公园绿地区PN4.55862114.681591.65285811.874113.744166将计算得到的5种功能区的综合污染指数用SIGMAPLOT软件绘出，以形象直观的反应各污染区的综合污染指数大小。图2 5种功能区的综合污染指数通过对比表1，不同功能区的污染程度具有一下特点：1工业区，交通区，生活区，公园绿地区属于重污染区，其中，工业区和交通区的污染非常严重，PN值分别高达14.68159和11.87411，而山区的人类活动相对较弱，PN值为1.652858，属于轻度污染。2对不同区域PN值对比发现，其不同区域污染程度呈现如下特点：工业区交通区生活区公园绿地区山区 4.2 重金属污染的原因分析为了研究重金属污染的原因，首先要分析8种元素之间的相关性，对其进行分类，找到污染产生的共同点，进而确定污染源的类型，不同污染物的产生就是因为不同类型的污染源的存在。4.2.1 Pearson相关系数当两个变量都是正态连续变量，而且两者之间呈线性关系时，表现这两个变量之间相关程度用积差相关系数，主要有Pearson简单相关系数。计算公式为：(3) 或 (4)其中，表示相关系数，表示样本数，和分别表示两变量的变量值，和分别为表示变量的平均值。 还可用下式表示： (5) 其中，和分别表示两变量的样本标准差4.2.2 相关系数的求解与分析为了研究8种元素可能的同源关系或是复合关系，本文通过SPSS对各元素浓度进行相关分析。表4 8种重金属元素的相关性元素AsCdCrCuHgNiPbZnAs1Cd0.225*1Cr0.189*0.352*1Cu0.160*0.397*0.532*1Hg0.0640.265*0.1030.417*1Ni0.317*0.329*0.716*0.495*0.1031Pb0.290*0.660*0.383*0.520*0.298*0.307*1Zn0.247*0.431*0.424*0.387*0.196*0.436*0.436*1由上表可知，8种元素中，Cd与Pb，Cr与Cu，Cr与Ni，Cu与Pb之间相关性显著，说明这些元素可能具有同源关系或复合关系。4.2.3因子分析因子分析是指研究从变量群中提取共性因子的统计技术，本文为了进一步研究8种重金属元素之间的同源关系，利用因子分析对公因子的各个元素载荷进行分析，找出同类污染源产生的重金属。其基本步骤可分为：1.建立原始变量间的相关系数矩阵。2.求该系数矩阵的特征值、累计贡献率，并确定公因子个数。3.建立因子载荷阵。4.将因子载荷阵实行方差最大正交旋转，使得能对公因子给出合理的解释。4.2.4 因子分析的求解与结论利用SPSS软件对8种重金属元素浓度进行因子分析，获得表5和表6结果。表5 方差贡献率成 分初始特征值提取平方和载入旋转平方和载入合计方差的%累积%合计方差的%累积%合计方差的%累积%13.56044.50044.5003.56044.50044.5002.22827.84927.84921.15014.37758.8771.15014.37758.8771.77622.20450.05330.96512.06370.9410.96512.06370.9411.67120.88870.9440.7689.59680.53750.5787.22087.75660.4325.39993.15670.3013.76996.92480.2463.076100.000表6 因子载荷元素旋转前旋转后因子1因子2因子3因子1因子2因子3As0.426-0.2000.6810.1150.794-0.202Cd0.7110.2810.2820.1940.6420.463Cr0.735-0.444-0.3030.8940.1420.102Cu0.7560.125-0.3650.5970.1250.591Hg0.4080.673-0.2970.0050.0080.842Ni0.723-0.515-0.1900.8790.227-0.002Pb0.7640.3140.2370.2300.6360.531Zn0.699-0.0370.1230.4440.4910.258因子分析方法可以用来判别土壤中重金属的来源，在本文的因子分析中Bartlett球度检验相伴概率为0，小于显著性水平0.05，因此本研究中的数据适合于作因子分析。本研究中提取前3个主成分可以解释总变量的70.941%，即对前3个主成分进行分析已经能够反映全部数据的大部分信息。Cr Ni Cu在因子1有较高的载荷，As Cd Pb Zn在因子2上有较高的载荷，这说明Cr Ni Cu和As Cd Pb可能主要分别来自类似的一个源。Hg在因子3上远远高于在其他元素在该因子的载荷，这说Hg可能由单一污染源产生。 综合以上两种方法分析污染原因：Cr Cu Ni在两种分析方法中都具有很强的关联性，根据其共同特性，分析产生的原因可能与电镀，合金制造等重工业有关，As Cd Pb Zn Cu在两种分析方法中呈现出一些差异，因此其污染来源情况比较复杂，可能与尾气排放，人类活动丢弃物以及农用废物有关，Hg在两种分析方法中均表现出独立性，因此伴随Hg产生的其他重金属元素则相对较少，可能主要来源可能是科学仪器制造，医学仪器制造工厂等对Hg有特定需求的工厂。4.3 污染源位置的确定通过对污染传播特征的研究，发现污染源的重金属浓度最高，并向周围呈递减趋势扩散，在计算污染源时，只要找到重金属浓度的极大值点，便能确定污染源的位置。4.3.1 聚类分析聚类分析是研究分类问题的一种多元分析方法，根据一组样品的多个观测值定出能度量样品相似程度的统计量，把相似程度大的归为一类。由于现实的分类往往伴随着模糊性，聚类问题采用模糊数学语言描述有其方便之处，这便是模糊聚类方法。根据实际问题的需要模糊聚类过程可以常采用两种方式：一种是通过建立模糊等价关系进行聚类分析的方法，另一种是基于模糊划分的模糊聚类方法，即模糊C一划分。本文采用前一种方式。本文利用聚类分析找出8种重金属元素浓度的所有极大异常点，通过其传播特性可以预见这些极大异常点距离污染源是比较近的，这样便能为之后污染源的确定提供依据。以8种重金属元素浓度为基础，进行聚类，找出每种重金属元素分类较少点，即为极大异常点。得到各种元素异常点如下：表7 各元素异常点元素取样点编号As178,84,41,29,30,6Cd95,22Cr20,22Cu8,22Hg257,182,8,9Ni135,22Pb16,8,6Zn61,36,224.3.2 Grid data插值模型分析图1，我们可以看出每种重金属的污染源分布较散，在此我们采用Grid data插值模型，运用局部搜索方法，搜索出每种重金属元素的污染源位置。4.3.3 插值模型的建立与求解将所有取样点位置的横坐标（x）、纵坐标（y）、以及8重元素的浓度作为变量，对横坐标（x）、纵坐标（y）进行插值搜索，搜索出每种元素的污染源的横坐标（x）、纵坐标（y），和污染源浓度。运用Grid data插值法，选择步长100进行搜索，求出每种元素的污染源位置如下：表8 As元素污染源位置横坐标x(m)纵坐标y(m)污染程度3300750062.25275100990091.45905600990024.605812600300023.924779001510050.9202184001010032.2364表9 Cu元素污染源位置横坐标x(m)纵坐标y(m)污染程度1700470013268.843300750018603.92 表10 Hg元素污染源位置横坐标x(m)纵坐标y(m)污染程度2300230021328.7713800230015142.4315100930018473.18表11 Ni元素污染源位置横坐标x(m)纵坐标y(m)污染程度7002900268.795330075001769.5213900960040.155015200590043.8626表12 Pb元素污染源位置横坐标x(m)纵坐标y(m)污染程度19003900866.766330075002024.4951005100425.021表13 Zn元素污染源位置横坐标x(m)纵坐标y(m)污染程度950045004131.83119002900976.1871390096003767.79表14 Cr元素污染源位置横坐标x(m)纵坐标y(m)污染程度2900680011295.063300750011279.07表15 Cd元素污染源位置横坐标x(m)纵坐标y(m)污染程度1700470010200.43360056003616.7953005100794.38160001160013560.0860007900779.26381003200884.649910015900863.521125001900900.4071770027001006.161790010000728.601结果分析：1.通过聚类分析与Grid data插值模型搜索出的结果可以看出，异常点的取样位置与各污染源的污染程度最大值的位置相接近，且都是从污染程度都是从最大值的位置向四周逐步衰减，从而可以看出Grid data插值模型搜索出的各种元素污染源的位置是可靠的。2.根据插值模拟出的污染源坐标，可以找到坐标非常相近的不同污染源，考虑到模型误差，本文将相距2Km之内的污染源看作是同一污染源，坐标取其平均值。表16 综合污染源含污染源元素种类平均值坐标As Cu Ni Pb Cr(3300,7500)As Cd(5800,10750)As Hg Zn Cd(12700,2525)As Cd(18150,10050)Cu Pb Cd(1767,4433)Hg Zn Ni(14300,9500)Pb Cd(5200,5100)Cr Cd(3250,6200)经过上述方法，仍旧有一部分污染源呈现独立性。表17 单一污染源污染源元素坐标As(7900,15100)Hg(2300,2300)Ni(700,2900)Ni(15200,5900)Zn(9500,4500)Cd(6000,7900)Cd(8100,3200)Cd(9100,15900)Cd(17700,2700)4.4 模型的检验为了进一步检验插值法的有效性，本文通过逐步多元回归对其中一些点进行检验。1、功能区属性变量的量化处理如下表：表18 五种功能区量化处理表功能区生活区工业区山区交通区公园绿化区量化处理1 0 0 0 00 1 0 0 00 0 1 0 00 0 0 1 00 0 0 0 12、以8种重金属浓度分别为因变量，每个取样点的坐标，以及功能区为自变量，建立多元回归模型，我们以As元素为例，用SAS编程，求的结果如下表：表19 回归方程检验VariableParamete