摘要：正确和及时的水文预报，有利于防汛、抗旱及发挥水利设施的作用。本 文通过MATLAB及EXCEL对给定数据进行处理，按照假设的抽样方式采样本点。 针对问题一：对数据进行处理，抽取其中具有代表性的数据，按时段利用样本点的预测雨量与相应站点的观测雨量的“误差的期望”、“误差的
标准
excel标准偏差excel标准偏差函数exl标准差函数国标检验抽样标准表免费下载红头文件格式标准下载
差”来评价方法一、方法二的准确性。本文论的两种方法总体上都较好。通过“误差的3阶中心矩”来比较方法一和方法二的优良性。 针对问题二：根据气象部门提供的雨量范围，将雨量划分为等级，通过预测的雨量的等级与观测的雨量的等级差的概率来考虑公众的感受。在评价方法中适当扩大等级之间定义预报的准确性。若结果是等级差为“0”及“1”时认为公众是满意的。此时方法一：公众满意 95%的预报，方法二：公众满意94.17%的预报。 【关键字】：误差 标准差 3阶中心矩 等级 1. 问题重述 C题 雨量预报方法的评价 雨量预报对农业生产和城市工作和生活有重要作用，但准确、及时地对雨量作出预报是一个十分困难的问题，广受世界各国关注。我国某地气象台和气象研究所正在研究6小时雨量预报方法，即每天晚上20点预报从21点开始的4个时段（21点至次日3点，次日3点至9点，9点至15点，15点至21点）在某些位置的雨量，这些位置位于东经120度、北纬32度附近的53×47的等距网格点上。同时设立91个观测站点实测这些时段的实际雨量，由于各种条件的限制，站点的设置是不均匀的。 气象部门希望建立一种科学评价预报方法好坏的数学模型与方法。气象部门提供了41天的用两种不同方法的预报数据和相应的实测数据。预报数据在文件夹FORECAST中，实测数据在文件夹MEASURING中，其中的文件都可以用Windows系统的“写字板”程序打开阅读。 FORECAST中的文件lon.dat和lat.dat分别包含网格点的经纬度，其余文件名为
_dis1和_dis2，例如f6181_dis1中包含2002年6月18日晚上20点采用第一种方法预报的第一时段数据（其2491个数据为该时段各网格点的雨量），而f6183_dis2中包含2002年6月18日晚上20点采用第二种方法预报的第三时段数据。 MEASURING中包含了41个名为<日期>.SIX的文件，如020618.SIX表示2002年6月18日晚上21点开始的连续4个时段各站点的实测数据（雨量），这些文件的数据格式是： 站号 纬度 经度 第1段 第2段 第3段 第4段 58138 32.9833 118.5167 0.0000 0.2000 10.1000 3.1000 58139 33.3000 118.8500 0.0000 0.0000 4.6000 7.4000 58141 33.6667 119.2667 0.0000 0.0000 1.1000 1.4000 58143 33.8000 119.8000 0.0000 0.0000 0.0000 1.8000 58146 33.4833 119.8167 0.0000 0.0000 1.5000 1.9000 …… 雨量用毫米做单位，小于0.1毫米视为无雨。 (1) 请建立数学模型来评价两种6小时雨量预报方法的准确性； (2) 气象部门将6小时降雨量分为6等：0.1—2.5毫米为小雨，2.6—6毫米为中雨，6.1—12毫米为大雨，12.1—25毫米为暴雨，25.1—60毫米为大暴雨，大于60.1毫米为特大暴雨。若按此分级向公众预报，如何在评价方法中考虑公众的感受？ 2. 问题的背景与分析 虽然我国幅员辽阔，地形复杂，但降水量的空间分布仍有一定规律。由于受季风影响，我国降水量的地域分布总的趋势是:从东南沿海向西北内陆逐渐减少。台湾、福建、海南、广东、浙江等省年降水量一般在1500～2000毫米(台湾省部分地方可超过4000毫米)，长江中下游地区在1200～1400毫米，淮河流域和秦岭山地在750～1000毫米，黄河下游、渭河流域在500～750毫米，其他一些地方大部分在200～500毫米，局部地方达750～1000毫米。在塔里木盆地等沙漠干旱地区，年降水量仅有几十毫米，其中新疆的托克逊年最少，降水量只有几毫米。而台湾的最多降水量达8409毫米。 我国降水量季节分配的特点是夏季最多，冬季最少，春、秋季介于两者之间。大的雨带一般5月中旬出现在华南地区，6月中旬北推到长江流域，7月中旬北推到淮河流域以北，至8月下旬到9月上旬开始雨带又逐步往南推移。 随着我国经济的快速发展，气象信息越来越受到广泛的关注，无论是人民的生活，还是经济建设，从城市保障到商业活动都与气象信息息息相关。向社会提供准确及时的气象信息是我们的宗旨；满足人民对气象信息的多种需求是我们的目标。因此，准确的对雨量预报，有着十分重要的意义。 我们建立的这个数学模型就是本着对雨量预报的准确度出发，从误差，预报的准确率等因素来全面评价雨量预报的准确性，而比较两种预报方法哪一种更准确。 3. 模型的假设与符号说明 模型的假设： 1．各观测点的降雨量可用以该观测点为中心，r为半径的圆内的各网格点的降雨量的平均值做为该站的预报降雨量。 2．各观测站在整个观测系统中所占的权重相同，不受地理位置分布的影响。 3．假定公众对降雨量的感受仅取决于雨量的大小与个人因素及时间无关。 4. 观测点的设置是不均匀的。 5. 利用坐标轴，x轴，y轴分别表示纬度和经度的坐标，通过把点的纬度和经度分别看作横坐标和纵坐标，用欧式距离计算公式r= 来作为两点之间的距离。 6. 雨量用毫米来做单位，小于0.1毫米视为无雨。 7．假设所选站点的分布在区域内具有随机性和普遍性 符号说明： r：实测站点与预测站点之间的距离。 X0ij(Pi,tjk) 为观测点的降雨量 i站j日k时段 X1ij(Pi,tjk) 为对应观测点第一种方法经平均后的降雨量。 X2ij(Pi,tjk) 为对应观测点第二种方法经平均后的降雨量。 其中：i:1—u;u=15 j:1—v;v=8 k:1—q;q=2 p为观测点的位置； t为时段 ：表示i站观测点j日k时段用第一种方法的测量的误差。 ：表示i站观测点j日k时段用第二种方法的测量的误差。 和 ：第一种方法测量误差的期望和标准差。 和 ：第二种方法测量误差的期望和标准差。 h：降雨等级空报或漏报是，公众满意度的折扣因子。范围在[0,1]内。 4．问题分析 1、本题中两种预报方法所包含的数据都比较繁多，虽然可以编程序实现所有数据的录入和调用，但是如果只是评价这两种方法的优劣，我们可以使用采集样本点的方法，去评价这两种方法的准确性。采集样本的方法考虑了： （1）. 天数的选择每隔五天取一天，均匀采样。 （2） 站点选取较密度区域的15个站点。 2、同时我们还要在评价的方法中去考虑公众的感受。也就是说，无论取哪种方法，公众对降雨预报评价应该是与预报准确度成正比的。预报结果和实际降雨越接近，公众对预报水平的评价就会越高。表一将降雨的等级分为，小雨，中雨，大雨，暴雨，大暴雨，特大暴雨7个等级。预报的降雨等级和实际的降雨等级相差越大，公众对预报的评价就会越低，如果预测的数值与实测的数值都在一个等级上时，公众的信服力就会增加。所以可根据预报的降雨等级和实际结果所在的等级差距不同给出不同的评分。 5．模型的建立 由题中所给数据可以看出各个站点的位置分布是不规则的，如果应用常规的方法去做，工作量会很大，而且应用MATLAB软件画雨量图来后可以发现有很多的站点的值上一样的（或者是非常相近的），所以我们不妨在91个实测站点抽样出15个出来（见表1），还有，并非每个时段的数据都是公众所关心的，比如，第一时段（21点至次日3点）时人们都处于睡梦中，对天气的要求不高；第三时段（9点至15点）时人们都处于工作时间，也对天气要求不高。而在第二时段（次日3点至9点）和第四时段（15点至21点）处于人流密集且白天范围较大的时段。在一个月中的，并不是每天的天气都是变化无常的，应此不用考虑全部41天的情况，我们不妨抽随机样取8天作为这41天的样本进行分析。 分析所给定的观测点和预测网格节点的坐标图（见图）。可知观测点是分布是不均匀的。因此为了以同样的权重对数据进行分析。按图分成等分区域进行抽 样。每一抽样在中间取一个观测点做为采样点数据，一共15个观测点。 用电脑分析查找15个观测点中，分别以该观测点为中心，r为半径的圆内的各网格点的降雨量的平均值做为该站的预报降雨量，取平均后即为该站点相对应的该种预报方法的降雨量。 （表一：选取的15个观测站点） 表2：（蓝色，红色分别表示预测和实测站点的经纬度分布） 为了对这两种预报方法进行客观、准确的评估，我们将建立以下两种模型。 首先，找到用计算机编程算得的对应观测站点相邻节点的各预报降雨量值(2种方法)再求出其平均值。 另：我们只取八个预报日。 设：15个站点和相对应的预报均值的数学表达如下： X0ij(Pi,tjk) 为观测点的降雨量 i站j日k时段 X1ij(Pi,tjk) 为对应观测点第一种方法经平均后的降雨量。 X2ij(Pi,tjk) 为对应观测点第二种方法经平均后的降雨量。 其中：i:1—u;u=15 j:1—v;v=8 k:1—q;q=2 p为观测点的位置； t为时段 所得结果如表一所示：（如6月22日） 误差分析表 第二时段 第四时段 法一 法二 法一 法二 F622 -3.1587 -3.1638 -0.5020 -0.4879 F627 -1.2594 -1.2573 0.0009 -0.0064 F704 -1.6767 -1.6402 -0.6030 -0.5999 F709 -3.3648 -3.3656 -0.6978 -0.6981 F714 -3.2918 -3.2747 -0.6801 -0.6794 F719 -3.3451 -3.3454 -0.0636 -0.0765 F724 -3.3660 -3.3652 -0.2636 -0.2607 F729 -3.3679 -3.3698 -0.6876 -0.6873 一阶比较 0 1.0e-015*0.3331 1.0e-015*0.0069 1.0e-015*-0.0208 二阶比较 0.6552 0.6657 0.0729 0.0706 三阶比较 0.6311 0.6415 0.0111 0.0100 6月22日实测值图像如下图所示： 表4：6月22日降雨量的实测值 （方法1第2时段） 表5：6月22日降雨量的实测值 （方法1第4时段） 表6：6月24日降雨量的实测值 （方法2第2时段） 表7：6月24日降雨量的实测值 （方法2第4时段） 先求出预报降雨量相对于观测站实际降雨量的误差（2种方法） ： 现认为其降雨量的误差构成一组数。现我们评价二种方法就是要评价二种方法所构成误差数列的统计量。用于比较二种预测方法的优劣。 分别对两种方程进行数据统计运算，对第一种预报方法：先求 的期望和标准差。 期望 : 标准差 : 再对日期统计求得期望和标准差 同理，对第二种方法求得期望和标准差 和 所求得的标准差如下表所示：（6月22日） 标准差分析表 第二时段 第四时段 法一 法二 法一 法二 F622 10.5853 10.5902 1.4656 1.4870 F627 10.6794 10.6380 1.6612 1.6194 F704 10.6555 10.7291 1.4160 1.4158 F709 10.6321 10.6319 1.3973 1.3953 F714 10.6351 10.6438 1.4020 1.4022 F719 10.6344 10.6352 0.9097 0.9779 F724 10.6295 10.6293 1.4560 1.4551 F729 10.6305 10.6297 1.3983 1.3984 一阶比较 0 1.0e-014 * -0.1110 1.0e-014 * 0.0180 1.0e-014 * 0.0208 二阶比较 0.0006 0.0013 0.0394 0.0296 三阶比较 -0.0000 0.0001 -0.0111 -0.0074 由题中第二问可知根据降雨量将降雨等级分为： 降雨量（毫米） 0.1—2.5 2.6—6 6.1—12 12.1—25 25.1—60 大于60.1 等级 小雨 中雨 大雨 暴雨 大暴雨 特大暴雨 由于对公众的降雨预报是分等级来预报的，而公众对各等级的感知度不同，例如，公众能十分清楚的判断有雨和无雨的差别，但是在有雨的情况下，对相差只有一个等级差距的雨量的差别却不一定能够清楚的判断。另外，公众对空报和漏报的感受又是不同的，对漏报的厌恶度明显会高于空报的厌恶度。所以我们在评分的时候也要考虑这个问题。则我们在原正确率的评分的的基础上，制定如下的评分细则：（1）无雨报有雨，或者有雨报成无雨时，在所得分数上所乘的折扣因子h（h在[0,1]内）比其他情况下的折扣因子要小一点,而且有雨报成无雨比无雨报成有雨的折扣因子更小.(2)有雨的情况下,空报的折扣因子大于漏报的折扣因子。 表C:：方法一的预测值与实际值等级差的误差的绝对值的总个数分别为 等级 0 1 2 3 4 个数 163 65 6 4 2 表D:：方法二的预测值与实际值等级差的误差的绝对值的总个数分别为 等级 0 1 2 3 4 个数 150 76 7 4 3 6．模型的结果分析与评价 用matlab做简单的计算，可以算出8天内两种方法分别在第二和第四时段的期望(误差)，然后算出期望的标准差，就可以初步比较两种预测方法的优劣。从计算的数据可以看出是第一种的方法比较好。如果结合第二问的把降雨量划分的等级，求出预测值与实际值等级差，然后划分级别（0，1为满意；2为较满意；3为不满意；4为极不满意；），就可以求出第一种方法的满意度为S1=234/240；第二种方法的满意度为S2=233/240；这样一来，就可以比较出第一种预测方法较好。 第二问：从等级的划分去考虑公众的感受。由计算得，用第一种方法预测公众的满意程度：为满意的概率P11=228/240；较满意的为概率P12=6/240；不满意的概率P13=4/240；极不满意的概率P14=2/240；用第二种方法预测公众的满意程度：为满意的概率P21=226/240；较满意的为概率P22=7/240；不满意的概率P23=4/240；极不满意的概率P24=3/240。 结合所有的数据来看，两种天气预报的准确率都比较高，并且公众的满意度也是比较高的，从公众平时的出行，日常生活以及工作，对这两种天气预报的预报还是信赖的过的。 模型的评价及改进方向： 1． 通过误差的期望、误差的标准差进行评价，验证了两种方法的准确性是合理的。 2． 通过误差的3阶中心矩，预报等级与观测等级差的概率可以认为方法一略优。当数字特征期望、方差无法对所描述的问题进行评价时，3阶中心矩是一个重要的评价指标。 3． 文中采集样本点的方式是科学的、合理的、值得推广。 4． 预报正确，公众满意，两种方法预报的正确率为何95%、94.17%.. 5． 由于时间的限制，本文没有考虑从个观测站的角度讨论文中所涉及的数字特征。但是文中所采用的从各时段的讨论的数字特征结果显示是合理的。 6． 文中忽略了地域的特殊气候的情况。 7． 系统抽样的数据不多，对这两种方法的评价只是在某方面的数值上的。 7．参考文献 范正森 谢兆鸿等，数学建模技术，北京，中国水利水电出版社，2003 王敏生 王庚, 现代数学建模方法, 北京,科学出版社 2006 冯杰 黄力伟等，数学建模原理与案例，北京，科学出版社，2007 李志林 欧宜贵 数学建模及典型案例分析 北京，化学工业出版社 2006 附录：Matlab程序 %15个站点误差的q期望和标准差及一，二，三阶比较 a=xlsread('Book1.xls','Sheet1','A1:G91');%导入数据 plot(a(:,2),a(:,3),'r-*')%画出91个观测点的分布图 grid on; hold on b=xlsread('Book2.xls','Sheet1','A1:AU53');%导入预测点的经度 c=xlsread('Book3.xls','Sheet1','A1:AU53');%导入预测点的纬度 plot(b,c,'bo')%画出预测点的分布图 d=xlsread('Book4.xls','Sheet1','A1:G15'); %选取8天为研究对象，分别计算第一个观测点（k=1）每天2，4时段两种预测方法的预测值 k=1; %第二时段误差 e=xlsread('f6222-1.xls','Sheet1','A1:AU53');%方法一的预测数据 f=xlsread('f6222-2.xls','Sheet1','A1:AU53');%方法二的预测数据 [E1,E2]=forecast(e,f,k);%E1,E2分别为第k个点方法一，二的预测值 E1-d(k,5); %方法一 E2-d(k,5); %方法二 %第四时段误差 e=xlsread('f6224-1.xls','Sheet1','A1:AU53');%方法一的预测数据 f=xlsread('f6224-2.xls','Sheet1','A1:AU53');%方法二的预测数据 [E1,E2]=forecast(e,f,k);%E1,E2分别为第k个点方法一，二的预测值 E1-d(k,7); %方法一 E2-d(k,7); %方法二 W1=xlsread('1.xls','Sheet1','A1:E8'); W2=xlsread('2.xls','Sheet1','A1:E8'); W3=xlsread('3.xls','Sheet1','A1:E8'); W4=xlsread('4.xls','Sheet1','A1:E8'); W5=xlsread('5.xls','Sheet1','A1:E8'); W6=xlsread('6.xls','Sheet1','A1:E8'); W7=xlsread('7.xls','Sheet1','A1:E8'); W8=xlsread('8.xls','Sheet1','A1:E8'); W9=xlsread('9.xls','Sheet1','A1:E8'); W10=xlsread('10.xls','Sheet1','A1:E8'); W11=xlsread('11.xls','Sheet1','A1:E8'); W12=xlsread('12.xls','Sheet1','A1:E8'); W13=xlsread('13.xls','Sheet1','A1:E8'); W14=xlsread('14.xls','Sheet1','A1:E8'); W15=xlsread('15.xls','Sheet1','A1:E8'); W(:,:,1)=W1; W(:,:,2)=W2; W(:,:,3)=W3; W(:,:,4)=W4; W(:,:,5)=W5; W(:,:,6)=W6; W(:,:,7)=W7; W(:,:,8)=W8; W(:,:,9)=W9; W(:,:,10)=W10; W(:,:,11)=W11; W(:,:,12)=W12; W(:,:,13)=W13; W(:,:,14)=W14; W(:,:,15)=W15; Q=zeros(8,4); for i=1:8 for j=1:4 w=0; for k=1:15 w=W(i,j,k)+w; Q(i,j)=w/15; end end end Q %对15个站点误差的期望 for i=1:8 for j=1:4 w=0; for k=1:15 w=(W(i,j,k)-Q(i,j))^2+w; R(i,j)=sqrt(w/15); end end end R %某日某时段15个站点误差的标准差 L=sum(R)/8 %对日期统计求得标准差的期望 T=sum(Q)/8 %对日期统计求得误差的期望 for j=1:4 w=0; for i=1:8 w=Q(i,j)-T(j)+w; S(j)=sqrt(w/8); end end S %对日期统计求得标准差 for j=1:4 w=0; w1=0; w2=0; for i=1:8 w=Q(i,j)-T(j)+w; w1=(Q(i,j)-T(j))^2+w1; w2=(Q(i,j)-T(j))^3+w2; end Z1(j)=w/8; %一介中心距 Z2(j)=w1/8; %二介中心距 Z3(j)=w2/8; %三介中心距 end Z1 Z2 Z3 for j=1:4 w=0; w1=0; w2=0; for i=1:8 w=R(i,j)-L(j)+w; w1=(R(i,j)-L(j))^2+w1; w2=(R(i,j)-L(j))^3+w2; end Z5(j)=w/8; %一介中心距 Z6(j)=w1/8; %二介中心距 Z7(j)=w2/8; %三介中心距 end Z5 Z6 Z7 %选取8天为研究对象，分别计算第一个观测点（k=1）每天2，4时段两种预测方法的预测值 %15个观测点实测数据分等级 a=xlsread('Book1.xls','Sheet1','A1:G91');%导入数据 plot(a(:,2),a(:,3),'r-*')%画出91个观测点的分布图 grid on; hold on b=xlsread('Book2.xls','Sheet1','A1:AU53');%导入预测点的经度 c=xlsread('Book3.xls','Sheet1','A1:AU53');%导入预测点的纬度 plot(b,c,'bo')%画出预测点的分布图 d=xlsread('f622.xls','Sheet1','A1:G15'); for i=1:15 for j=5:2:7 if d(i,j)<0.1 p(i,j)=0; elseif 0.1<=d(i,j)<2.5 p(i,j)=1; elseif 2.5<=d(i,j)<6 p(i,j)=2; elseif 6<=d(i,j)<12 p(i,j)=3; elseif 12<=d(i,j)<25 p(i,j)=4; elseif 25<=d(i,j)<60 p(i,j)=5; elseif 60<=d(i,j) p(i,j)=6; end end end p(:,5:2:7); %实测数据雨量分级 %q求预测值 function [E1,E2]=forecast(e,f,k) d=xlsread('f622.xls','Sheet1','A1:G15'); x=d(:,2); y=d(:,3);%导入15个观测点位置坐标 b=xlsread('Book2.xls','Sheet1','A1:AU53');%导入预测点的经度 c=xlsread('Book3.xls','Sheet1','A1:AU53');%导入预测点的纬度 m=1; e=e'; f=f'; %计算第k个点两种方法的预测值 for i=1:53 for j=1:47 D(m)=sqrt((b(i,j)-x(k))^2+(c(i,j)-y(k))^2); m=m+1; end end [sortE,turnV]=sort(D); D=0; m=1; for i=1:4 E=0; F=0; jj=mod(turnV(i),47); ii=((turnV(i)-jj)/47)+1; E=e(ii,jj)+E; F=f(ii,jj)+F; end E1=E/4; %方法一所得观测点的预测值 E2=F/4; %方法二所得观测点的预测值 %画1图 a=xlsread('B1.xls','Sheet1','A1:G91'); plot(a(:,2),a(:,3),'ro'); grid on gtext('1');gtext('2');gtext('3');gtext('4');gtext('5');gtext('6');gtext('7');gtext('8'); gtext('9');gtext('10');gtext('11');gtext('12');gtext('13');gtext('14');gtext('15'); figure(1); %6月22日方法一第二时段的降雨量 x=xlsread('B11.xls','Sheet1','A1:G91'); a=x(:,2);%经度 b=x(:,3);%纬度 x1=a([1:47],:); y1=b([1:53],:); z=xlsread('6222-1.xls','Sheet1','A1:AU53'); stem3(x1,y1,z,'r.'); figure(2); %6月22日方法二第一时段的降雨量 x=xlsread('B11.xls','Sheet1','A1:G91'); a=x(:,2);%经度 b=x(:,3);%纬度 x1=a([1:47],:); y1=b([1:53],:); z=xlsread('6222-2.xls','Sheet2','A1:AU53') stem3(x1,y1,z,'g.'); figure(3); %6月22日方法一第四时段的降雨量 x=xlsread('B11.xls','Sheet1','A1:G91'); a=x(:,2);%经度 b=x(:,3);%纬度 x1=a([1:47],:); y1=b([1:53],:); z=xlsread('6224-1.xls','Sheet3','A1:AU53'); stem3(x1,y1,z,'m.'); figure(4); %6月22日方法二第四时段的降雨量 x=xlsread('B11.xls','Sheet1','A1:G91'); a=x(:,2);%经度 b=x(:,3);%纬度 x1=a([1:47],:); y1=b([1:53],:); z=xlsread('6224-2.xls','Sheet1','A1:AU53'); stem3(x1,y1,z,'b.');
浙公网安备 33021202002483