数学建模-城市表层土壤重金属污染分析

2011高教社杯全国大学生数学建模竞赛承诺书我们仔细阅读了中国大学生数学建模竞赛的竞赛规则.我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。我们知道，抄袭别人的成果是违反竞赛规则的,如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。我们参赛选择的题号是（从A/B/C/D中选择一项填写）：我们的参赛报名号为（如果赛区设置报名号的话）：所属学校（请填写完整的全名）：参赛队员(打印并签名)：1.2.3.指导教师或指导教师组负责人(打印并签名)：日期：年月日赛区评阅编号（由赛区组委会评阅前进行编号）：2011高教社杯全国大学生数学建模竞赛编号专用页赛区评阅编号（由赛区组委会评阅前进行编号）：赛区评阅记录（可供赛区评阅时使用）：评阅人评分备注全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：城市表层土壤重金属污染分析摘要本文考虑城市表层土壤重金属污染问题。随着人口增长、人口密度增加、经济发展，环境污染问题日益成为人们关注的焦点。问题一，为考察重金属元素的空间分布，本文对已有采样点进行griddata插值，运用matlab软件，给出各重金属元素的空间分布。进一步建立综合评价模型，对8种主要重金属元素进行赋权，通过影响率和加权求出每个功能区污染程度得分。评判标准以背景值1为分界线，分值高者越高污染程度越深。模型求解结果以工业区污染程度最高，主干道区次之。问题二，为分析污染原因，本文建立多维评价模型，先以数据间差异性分析，判断影响因素是主要是人为还是自然形成，再以潜在生态危害系数法对重金属污染、各功能区主要金属污染进行分析，通过对结果数据的分析再结合网络与书籍，基于重金属元素的传播特征，得到产生重金属污染的主要原因。主要污染元素为Hg和Cd，工业区和主干道区的污染最为严重。问题三，为确定污染源，本文建立聚类分析模型和重心定位模型，通过分析传播特征拟定出金属污染物是延污染源点向四周辐射的途径，对组内污染程度高的采样点赋权求重心，所得重心即为污染源位置，用spss进行模型求解得到9个污染源位置。问题四，为进一步研究地质演变模式，本文首先评价了前面所建模型的优缺点，根据前面建模中所欠考虑的几个因素,如时间,考虑在信息可收集的情况下,建立重金属污染的扩散模型。关键词：重金属污染聚类分析熵值法重心定位1问题重述随着城市经济的快速发展和城市人口的不断增加，人类活动对城市环境质量的影响日显突出。对城市土壤地质环境异常的查证，以及如何应用查证获得的海量数据资料开展城市环境质量评价，研究人类活动影响下城市地质环境的演变模式，日益成为人们关注的焦点。按照功能划分，城区一般可分为生活区、工业区、山区、主干道路区及公园绿地区等，分别记为1类区、2类区、……、5类区，不同的区域环境受人类活动影响的程度不同。现对某城市城区土壤地质环境进行调查。为此，将所考察的城区划分为间距1公里左右的网格子区域，按照每平方公里1个采样点对表层土（0~10厘米深度）进行取样、编号，并用GPS记录采样点的位置。应用专门仪器测试分析，获得了每个样本所含的多种化学元素的浓度数据。另一方面，按照2公里的间距在那些远离人群及工业活动的自然区取样，将其作为该城区表层土壤中元素的背景值。附件1列出了采样点的位置、海拔高度及其所属功能区等信息，附件2列出了8种主要重金属元素在采样点处的浓度，附件3列出了8种主要重金属元素的背景值。要求通过数学建模来完成以下任务：(1)给出8种主要重金属元素在该城区的空间分布，并分析该城区内不同区域重金属的污染程度。(2)通过数据分析，说明重金属污染的主要原因。(3)分析重金属污染物的传播特征，由此建立模型，确定污染源的位置。(4)分析所建立模型的优缺点，为更好地研究城市地质环境的演变模式，在可以收集到信息的情况下，建立模型解决问题。2问题分析2.1问题一分析运用matlab的griddata插值计算后得到8个重金属含量的插值图，从8钟重金属元素的空间分布图可以直观的看出8种主要重金属元素在该城区的空间分布情况。并用excel对重金属数据的分析，更加具体的得到重金属在不同区域的比重。对污染程度，建立综合评价模型，把8种主要重金属元素作为污染程度的8个指标，从As,Cd…依次编号为1,2…运用熵值法求出各个指标对污染程度的权重。在计算出不同地区金属元素的影响率,通过加权求和得出各个地区重金属污染程度的得分。以背景值的平均值为1，得分越高说明污染程度越严重。2.2问题二分析对重金属进行差异性分析，以此来判断重金属污染主要归因于自然因素还是人为因素，分清主要来源对象。再运用潜在生态危害指数的计算方法，得出土壤重金属潜在的生态风险，以不同功能区的横向分析和以不同重金属的纵向分析，再加上重金属含量与海拔高度关系的趋势分析，再以这些数据分析为基础，通过网络与书籍，查找资料，详细的得出重金属污染的主要原因。2.3问题三分析基于文献资料和前两问对数据的分析，本文得出重金属污染物在相对理想状态下的传播特征。根据传播特征，运用聚类分析法对采样点进行分类，找到用于定位的采样点，给每个样点加权求质量，建立重心定位模型定位得出的污染源位置。2.4问题四分析主要根据前面所建模型忽略的重金属传播的因素，如土壤PH值。以及考虑不周的地方进行再考虑。基于重金属在表层土壤中的传播特征，根据扩散定律结合物质守恒定律，可以给出扩散方程模型来研究地质环境的演变模式。模型假设1、假设不受风向的影响。2、假设不受河流、高山阻隔的影响。3、假设忽略重金属传播随时间变化的浓度扩散速度不同。4、假设重金属污染物传播不考虑在不同载体中的速度不同。4符号约定符号含义V第i个编号下底j个指标的值P第i编号下第j个指标的比重e第i项编号下第j个指标的熵值g第j项指标的差异系数w第j项指标的权重v重金属背景值得平均值R第j个项指标的影响率Grade各区域污染程度的得分比例C表层重金属浓度E单个重金属的潜在生态危害系数RI多种重金属的潜在生态危害系数C单个重金属的污染系数CV变异系数σ标准差V归一化处理后的数据C污染浓度5模型的建立与求解与检验5.1问题一5.1.1模型建立建立各地区重金属污染程度综合评价模型：Grade=过程如下：熵值法求权重对二级指标进行归一化处理，消除量纲的影响。公式如下：计算第i编号下第j个指标的比重计算第i项编号下第j个指标的熵值，其中计算第j项指标的差异系数。对于第j项指标，指标值的差异越大，对方案评价的作用就越大，墒值就越小。定义差异系数：计算第j项因素的权重加权求各地区重金属污染程度得分以重金属元素的背景值为标准值，影响率为各地区超标重金属与对应标准平均值的比例。其表达式为：R以重金属元素的背景值为1分，通过加权计算各区域污染程度的得分比例，得分越高污染程度越严重。其表达式为：Grade=5.1.2模型求解运用matlab作图可得（编程见附录A）：图1城市表土壤重金属采样区域模拟图图2功能区划分与采样点图2As、Cd、Cr、Cu空间分布图3Hg、Ni、Pb、Zn空间分布用excel作8种重金属元素在不同区域的分布柱形图：图4As、Cd、Cr、Cu区域分布图4Hg、Ni、Pb、Zn区域分布从8种重金属的空间分布图和柱形图结合，可以直观的看到，Cr、Cu、Ni、Zn的分布情况比较相近，在y值较小的范围，含量都比较小，大部分地区含量分布均匀，其中Zn的含量在x为1000~2000范围内有两个高浓度点，多处在工业区的位置。另外，As、Cd、Pb浓度延x轴方向逐渐减低，主要在y值高的地方，浓度较高。Hg的分布差异比较大，大部分地区表层土壤含量都在200以下，但在主干道区和工业区浓度明显增高。三类区土壤重金属含量在五个区域中最小，说明人类活动严重影响土壤重金属含量。运用matlab求解(见附录B)，得出各主要重金属元素的的权重为：w=（0.0233,0.0475,0.0651,0.1831,0.4903,0.0222,0.0446,0.1238）图5主要重金属的污染权重在运用excel的功能计算出以下数据：表1区域重金属污染程度得分表元素1类区2类区3类区4类区5类区平均值影响率平均值影响率平均值影响率平均值影响率平均值影响率As(μg/g)3.60.747.251.014.040.125.710.596.260.74Cd(ng/g)1301.23393.112.02152.320.17360.011.77280.541.16Cr(μg/g)311.2353.410.7238.960.2658.050.8743.640.41Cu(μg/g)13.22.74127.548.6617.320.3162.213.7130.191.29Hg(μg/g)351.66642.3617.3540.960.17446.8211.77114.992.29Ni(μg/g)12.30.4919.810.6115.450.2617.620.4315.290.24Pb(μg/g)311.2393.042.0036.560.1863.531.0560.710.96Zn(μg/g)692.43277.933.0373.290.06242.852.52154.241.24Grade2.8411.741.197.972.66有不同区域重金属污染程度得分情况可以看出，在3类区得分最低，也就是重金属污染程度最低。另外2类区和4类区重金属污染严重，大比例超过背景值。5.2问题二5.2.1模型建立建立潜在生态危害指数模型RI=Er过程如下：1、对重金属含量进行差异性分析。用变异系数消除不同水平变量单位和水平高低的影响，更好的反映单位均值上的离散程度。变异系数计算公式：CV=2、潜在生态危害指数的计算方法。单个重金属的污染系数：C单个重金属的潜在生态危害系数：E多种重金属的潜在生态危害系数：RI=以Hakanson的标准化毒性系数为评价指标。通过网上查资料可知，As,Cd,Cr,Cu,Hg,Ni,Pb,Zn的毒性系数分别为：10,30,2,5,40,5,5,1。其中，Er表2Er指数类型所处范围污染程度指数类型所处范围污染程度潜在生态危害系数(ErE轻微生态危害潜在生态危害指数(RI)RI轻微生态危害40≤E中等生态危害150≤RI中等生态危害80≤E强生态危害300≤RI强生态危害160很强生态危害RI≥600很强生态危害E极强生态危害3、重金属元素含量与海拔关系对附件2的数据进行量化整理，以重金属元素含量与对应背景值的比值减1，消除量纲的影响，再用量化后的值作出随海拔高度的趋势图。不同地区同一海拔的重金属元素含量也可能不同，为了消除区域间的差异，对每个区域分别进行曲线作图。5.2.2模型求解 对附件2的数据进行处理，求出总体重金属元素的平均值、标准差、变异系数。表3总体重金属元素数据分析金属元素平均值标准差变异系数CV/%背景值As(μg/g)5.683.0253.173.6Cd(ng/g)302.40224.9974.40130Cr(μg/g)53.5170.00130.8231Cu(μg/g)55.02162.92296.1113.2Hg(ng/g)299.711629.50543.6935Ni(μg/g)17.269.9457.5912.3Pb(μg/g)61.7450.0681.0831Zn(μg/g)201.20339.23168.6069 从表中，可以得到主要重金属元素的变异系数都超过了50%，说明重金属污染主要是人为因素。在其中，Hg的变异系数高达543.69%，Cr、Cu、Hg、Zn的变异系数也都超过了100%，人为活动对造成重金属污染是十分严重的。 再根据多种重金属的潜在生态危害系数公式计算得城市表层土壤重金属潜在生态风险程度分级值。表4重金属潜在潜在生态危险系数项目ErRI(分级)AsCdCrCuHgNiPbZn最大值83.69373.8059.41957.7618285.7157.9376.2154.5019765.76最小值4.479.230.990.879.791.743.170.4830.74平均值15.7769.783.4520.84342.537.029.962.92472.27由表可得主要重金属在区域上的潜在生态系数大小关系为Hg>Cd>Cu>As>Pb>表5各功能区潜在生态危害系数项目ErRI(分级)AsCdCrCuHgNiPbZn1类区17.4266.914.4518.71106.337.4611.153.43235.862类区20.1490.723.4548.31734.128.0515.014.03923.833类区11.2335.152.516.5646.816.285.901.06115.504类区15.8683.083.7523.57510.657.1610.253.52657.845类区17.4064.742.8211.44131.426.229.792.24246.07从功能区的生态危险系数可以看出，在横向上，2类区和4类区生态危害系数最大，其中2类区RI=923.83，4类区RI=657.84，已经属于很强的生态危害，这说明重金属污染存在地区差异，与2、4类区的功能有关。在8种元素纵向上，Hg的危害系数最大，在2和4类区的危害系数Eri分别达到了734.12和510.65，均属极强生态危害。其次，Cd的 重金属元素量化值随海拔高度的趋势，以1类区为例，作出如下图。图6类区重金属元素量化值随海拔高度的趋势由图6，可以看出重金属含量在不同高度的量化值并没有一定的相关性。所以海拔高度并不是重金属污染的主要原因。通过查阅资料对数据处理结果分析： 在各种人为因素中，主要来自工业区和主干道区人类活动的影响，其中这种污染原因主要是产生了大量的Hg和Cd。 工业是产生重金属污染的主要原因之一。工业过程中广泛使用重金属元素，工矿企业将未经严格处理的废水直接排放，使得它周围土壤容易富集高含量的有毒金属[1]。城市表层土壤Hg的危害系数最大。其中由图5可以看出，Hg含量均值在工业区和主干道区最多。这是因为企业排污、废水处理等都含有Hg，并通过雨淋或自然进入土壤，甚至通过风的传播或废弃物的移动将污染传播到更远的地方，这些都是产生Hg的原因。 城市表层土壤Cd污染危害系数居第二，在图4中，我们可以看出Cd含量在除山区外的所有地区含量都很高。在生活区，由于一些人会从事农业活动，在污水灌溉、垃圾堆肥以及大量使用农药化肥等都会产生Cd。而在主干道区，由于含铅汽油的燃烧，气体轮胎、车体老化磨损等会产生Cd、Zn等重金属污染。在公园和绿地区Cd含量高的主要原因是污泥堆肥的广泛使用[2]。5.3问题三5.3.1模型建立土壤重金属污染物是由于人类活动，土壤重的微量有害元在土壤重的含量超过背景值，过量沉积而引起的含量过高，统称为土壤重金属污染[3]。基于这点，需要根据各采样点主要重金属在土壤中的含量超过背景值的比率有关的Grade来选取一组过高值，以平均值为界。 土壤重金属污染可以通过多种载体进行传播，如随大气、污水、固体废弃物、农用物资等为载体，忽略不同方式扩散的速度和随时间变化扩散速度都会有所不同，总体是以污染源为中心向四周扩散的。污染源可以看成是球形扩散的，建立重心定位模型，以周围扩散点的污染程度作为权值，用聚类分析找到同一污染源的空间污染散点坐标，最后求出重心，即为污染源位置。另外，由于不同金属有不同的传播特征可以推断不同金属的污染源可能来自不同的地方，因此需要根据相关性用聚类的方法，将样本点关于各重金属元素进行分类，确定不同的污染源。建立重心定位模型G(x,y.z)过程如下：1、各采样点污染浓度计算。重金属污染是随着污染源向四周扩散的，其浓度是随距离污染源的位置逐渐降低的。 广义浓度是指单位体积里所含溶质的物质的量。是与体积有关的。但由于所计算的区域是点，故可忽略体积因素。由于重金属物浓度越大的地方，其污染程度得分Grade也越高，其变化趋势是一样的，所以用问题一中所得重金属污染程度作为污染物浓度的替代值。 2、选取数据 用问题一中的方法求各个采样点的污染程度Grade，综合考虑各重金属污染的得分，选取Grade平均值以上的点做聚类分析。最后选取出40个采样点数据，结果见附表A。3、对各指标间进行相关性分析 8种物质不一定来自同一个污染源，对8个重金属指标进行相关性分析，把相关性的高的指标作为同时作为聚类分析的变量，划分出不同的区域。4、K-means聚类分析 聚类分析是直接比较各事物之间的性质，将性质相近的归为一类，将性质相近的归为一类，将性质差别较大的归入不同的类。所以我们根据采样点的指标为分类标准，对采样点进行分类。 首先对原始数据进行标准化处理，消除量纲的影响。ZScores：标准化变换公式为：v 数据选择，根据污染浓度，先对数据进行第一次选择。污染浓度在其平均值以上的采样点说明离污染源比较接近。 K-means算法： 将已经测定n个采样点的m种指标数据列成一个二维矩阵，亦称作观测矩阵：v=把vij其中n-1迭代过程：给出初始聚类中心，s用下列公式更新uu其中l为迭代次数。更新ss如果maxi=sil+1重心定位把分类点进行整理，用Grade的作为散点的质量，求分类采样点的重心位置即为污染源。每个散点对应坐标为xi,yi,zG(x,y.z)则x=y=(z=(5.3.2模型求解对8个指标两两间做相关性分析，分析结果见附表6。表6种重金属两两间的相关系数AsCdCrCuHgNiPbZnAs1.00Cd0.251.00Cr0.190.351.00Cu0.160.400.531.00Hg0.060.260.100.421.00Ni0.320.330.720.490.101.00Pb0.290.660.380.520.290.311.00Zn0.250.430.420.380.190.440.491.00根据相关性，我们把8个指标分成三组（分别为：Cd、Pb、Zn为一组，Cr、Cu、Ni为一组，As、Hg为一组）分别做聚类分析的变量，运用spss做聚类分析，以Cd、Pb、Zn为例的聚类分析结果见附录表A。根据聚类分析，我们共得出9组数据，用excel求重心，最终得到9个污染源的位置：表7污染源坐标指标组别X/mY/mZ/mGradeCd、Pb、Zn1166070822115437112026091222793531253626201287As、Hg3299601846609235518622366532383369277345Cd、Pb、Zn82614229202962181816520211157054365897123468用matlab做图结果如下：图7污染源位置图5.4问题四建立扩散模型：∂C问题一的griddata插值用matlab求解模型直观但不能界定好采样点内的区域。另外用综合评价模型可以客观全面的分析各指标对污染程度的影响，比较具有可行性。问题二的潜在生态危害系数模型，分析了各因素的危害系数和不同功能区的危害性，便于通过Er对于问题三建立的模型，能够考虑到不同重金属污染源可能来自不同地方，分差异性大的采样点进行分类，但计算过程复杂，可操作性比较弱。而且模型忽略了时间、山川阻隔等影响，对污染传播方式有比较大的影响。城市地质环境的演变，还涉及到重金属污染物排放后经过的时间；经过一定的衰竭期后土壤中金属的含量变化；重金属在扩散过程中的沉淀；土壤的PH值等信息。通过这些信息，我们可以建立重金属在表层土壤中的演变模式。根据扩散定律结合物质守恒定律，得出扩散模型[4]。由于重金属的扩散是在一定范围内沿平面分布污染，可假设模型扩散分布沿圆形区域，根据扩散定律，扩散物质在单位时间内沿法线方向n流过单位面积的曲面的质量与物质浓度Cx,y,z,t沿法线方向的方向导数∂C扩散物质在时段dt内沿法线方向流过面积为dS的曲面的质量dm为：dm=-D其中Dx,y,z任取一封闭曲面∝，所围区域为ω，则从时刻t1到时刻tm=t由高斯公式αDm=经过时间t的累积，金属含量产生变化，变化的含量为ω由质量守恒有t综上，得到扩散模型∂C在分析城市地质演变模式中，可以通过此扩散模型，来检查土壤中的重金属含量的变化，以此来探测城市地质的变化。6参考文献[1]郑喜坤，鲁安淮，高翔，赵谨，郑德圣，土壤中重金属元素污染现状与防治方法/view/b41f0a7d5acfa1c7aa00cc51.html，2011年9月11日[2]王应刚，辛晓云，郭翠花，太原市土壤重汞污染及成因研究，2000。[3]土壤重金属/view/1735172.htm，2011年9月11日[4]胡素华，张彤，张世英，正态逆高斯扩散模型的MCMC估计，系统工程理论与方法应用：第2期，133—138,2006。[5]刘卫国，MATLAB程序设计与应用，北京：高等教育出版社，2006。[6]邬学军周凯宋军全，数学建模竞赛辅导教程，杭州：浙江大学出版社，2009。[7]黄润龙，数据统计与分析技术:SPSS软件实用教程，北京：高等教育出版社，2004。[8]肖华勇，实用数学建模与软件应用，西安：西北工业大学出版社，2008。[9]姜启源谢金星叶俊，数学模型，北京：高等教育出版社，2011。[10]徐争启倪师军庹先国张成江，潜在生态危害指数法评价重金属毒性系数计算，环境科学与技术：第31卷，112—115，2008。7附录附表A根据重金属污染程度选取的40个数据编号AsCdCrCuHgNiPbZn97.411123.988.17151.641600025.8172.36926.841824.17357.538.7843.891490117.798.22349.889.621066.2285.582528.481350041.7381.641417.862573.57514.142.3447.671380017.369.96269.894116.5845994.7947.17190019.971.32215.12323.77515.628.7477.47180115.170.18178.57934.29232.929.1740.021714.299.7938.6595.412921.87424.573.459.72152027.883.7175.71614.081092.967.96308.61104028.2434.8966.731456.56202.333.2837.671339.2911.36140.62110.92229.841213.5920.841364.85115142.5181.481818.473018.3863096.68114.8164534.8130.361626.0256.35525.259.35117.5380020.2169.96726.02614.9303.856.0265.8663.2140.0590.693760.82494.69568.6306.0270.4190016.879.67196.7346.56223.940.0825.1795015.432.28117.35363.5396.3138.3758.9717024.291.762893.472345.2417.835.7147.5880116.0394.93204.082258.23420.635.8173.8656122.4960.91241.6317823.72452.137.2243.6116016.3573.271635.92208.5614744.46130.5515632.8228.641013.471438.721024.977.6171.24128.5719.35227.41389.391577.78315.549.7628.0355018.9545.73109.291635.6288078.29121.1229325.61171.14540406.56852.772.59158.6731121.2124.24377.141611.451044.594.78136.9720222.3472.48602.04345.51778.774.6692.4833019.7110.2351.6331011.68870.570.8435.1730229.178.15435.442755.41302.534.2227.640814.868.7218.241616.91568.554.59113.4626423.2282.4399.9454.9635.568.42227.7617619.596.33538.98394.29449.167.2251.6431515.7106.97294.693110.53635.364.03101.3519028.3162.64615.11548.9307.354.3957.2132625.72131.93256.94354.49754.899.8897.9224324.9100.79323.372445.2284.340.4744.7334916119.46184.611747.56580.139.7885.4627017.4966.69315.92155.51516.491.9789.0418919.8121.72494.8375.51687.885.5272.8520119103.2403.27144.49359.5258.15123.277712.9106.47853.98附表B以Cd、Pb、Zn为变量的spss聚类分析结果聚类成员案例号编号聚类距离19.0031.2592182.002.52038.0031.0124257.002.366541.002.4226232.002.443793.0021.220829.002.45296.0031.47910145.0021.2231122.0031.5621230.0021.649135.002.8271461.001.5851549.002.440164.0021.2701736.001.58518234.002.41519225.002.53920178.0021.5902120.0021.58322143.003.77323157.002.95224163.0021.5692540.0021.3062616.0031.9592734.0021.02328310.0021.35929275.002.82030161.002.3073145.002.5153239.002.2633331.002.85234154.002.8003535.002.93536244.002.86237174.002.4663815.002.2563937.002.6854014.002.747ANOVA聚类误差均方df均方dfFSig.Cd12.0092.4053729.657.000Pb11.9152.4103729.060.000Zn14.6252.2633755.504.000F检验应仅用于描述性目的，因为选中的聚类将被用来最大化不同聚类中的案例间的差别。观测到的显著性水平并未据此进行更正，因此无法将其解释为是对聚类均值相等这一假设的检验。每个聚类中的案例数聚类12.000232.00036.000有效40.000缺失.000附录AMATLAB插值法绘图x=data(:,1)';y=data(:,2)';z=data(:,3)';v1=data(:,5);v2=data(:,6);v3=data(:,7);v4=data(:,8);v5=data(:,9);v6=data(:,10);v7=data(:,11);v8=data(:,12);set(gca,'XTick',0:1000:30000);set(gca,'YTick',0:1000:20000);subplot(2,2,1);[X,Y,v1]=griddata(x,y,v1,linspace(0,28654)',linspace(0,18449),'v4');title('As');subplot(2,2,2)