南方医科大学 2012级生物医学工程学院
摘要
要评价两种预测方法的准确性,就要在相同地点、相同时段对预报值和实测值进行比较,本文首先针对这一问题,我们通过对监测站位于网格中的位置分类处理的方法对数据进行压缩,将2491个网格点的预报值变为91个观测站点的预报值;在问题一中,我们引入切比雪夫距离,构造对比矩阵,判断哪种方法误差小的方法,来评价两种预报方法准确性的好坏;在问题二中我们通过公众对预报结果满意度好坏给予相应的权重再求和(每个时段的满意度之和)并且应用了模糊数学法评价两种预报方法的准确性。由问题一和问题二所得出结果可知,两种预测方法各有各的优点和缺点,我们可以根据实际情况,选择较佳的预测方法。
我们应用MATLAB软件编程,处理数据并求解结果。检验结果表明第二种预测方法所得的结果相对准确。
一.问题的重述
雨量预报对农业生产和城市工作和生活有重要作用,但准确、及时地对雨量作出预报是一个十分困难的问题,广受世界各国关注。现有一个地方的气象台和气象研究所正在研究6小时雨量预报方法,每天晚上20点预报,预报从21点开始的 4个时段(21点至次日3点,次日3点至9点,9点至15点,15点至21点)的降雨量,预报位置位于东经120度、北纬32度附近的53×47的等距网格点上。同时,设立了91个观测站点,这些站点实测了41天各时段的实际雨量(站点的设置是不均匀的),气象部门还提供了两种不同预测方法41天的预报雨量。我们要利用提供的数据,解决下面的问题:
(1) 请建立数学模型来评价两种6小时雨量预报方法的准确性;
(2) 气象部门将6小时降雨量分为6等:0.1—2.5毫米为小雨,2.6—6毫米为中雨,6.1—12毫米为大雨,12.1—25毫米为暴雨,25.1—60毫米为大暴雨,大于60.1毫米为特大暴雨。若按此分级向公众预报,如何在评价方法中考虑公众的感受?
二.符号和名词的说明
----------给天数的编号,从预测的第1天到第41天;
----------给观察站的编号,从1到91号;
----------给时段的编号,从1到4 时段;
----------观察站点;
---------网格上的点;
---------观察站点与网格点的距离;
-------第i天第j站第k时段观测值与实际值偏差;
-------第i天第j站第k时段方法一得到的预测值;
-------第i天第j站第k时段方法二得到的预测值;
----------第i天第j站第k时段实际测量值;
--------比较矩阵;
-------对比函数;
三.问题分析:
题目中第一问是建立数学模型来评价两种6小时雨量预报方法的准确性。由于给的数据较多不能全部考虑,而又不能主观的随机抽取若干点处理,故对监测点分布的规律作分类考虑,用最能代替监测点的一个或几个网格点来代替监测点。尽可能的让监测点检测到的实际值利用率更高,大大的减少了数据量,而又不失准确性。而后在评价两种方法的问题中,我们通过绘制雨量折线图简单的表现并评判两种预测的贴切程度,引用了切比雪夫距离用最大误差代替整体误差的方法,建立误差矩阵,并通过方法一方法二的误差矩阵,得到一个对比矩阵。我们建立一个评价函数,通过评价函数求和的方式直观的解决了两种方法谁更优越的问题。
在第二问中,由于前一问已得出了方法一较为优越,故我们选用方法一。
第二问是:气象部门将6小时降雨量分为6等:0.1—2.5毫米为小雨,2.6—6毫米为中雨,6.1—12毫米为大雨,12.1—25毫米为暴雨,25.1—60毫米为大暴雨,大于60.1毫米为特大暴雨。若按此分级向公众预报,如何在评价方法中考虑公众的感受?
由于公众的感受是主观心理感受,不好分析。我们将其量化,并分等级,对于每次的预测误差评估归类以后,每个等级都有一个满意度权值,最后求和便得结果。
五.模型解决:
问题一:
雨量预报的方法选取十分重要,对于不同的时间段,不同的区域,使用的预测方法也因而相异。在不同的时段,选用更好的预报方法可以减小预报误差,为生活和生产提供更好的指导。
根据问题提供数据,我们得到了47*53个网格点的不同经度和纬度。并且在不同的网格点上,数据给出了预测1和预测2的一天四个时段(以六小时为间隔)的雨量预报值。同时,还建立了91个实际观察站,对一天内相同的四个时段雨量采样、处理并统计。最终给出了6月18日至28日、7月1日至30日在内41天的预测值和实测值。
首先,因为预测值远远多于实际测量值,并且由于地理原因和人力限制我们不能将观察站精确的建立在我们的预测点(47*53的网格格点)上。这样我们的预测点的值并不能直接的代表我们在这91个观察点上进行预测的值。在与实测值对比之前,我们需要对这47 x 53个点上的预测值进行一个处理。
由于我们的预测点相对密集,并且经纬度的变化情况不会太大,对于预测点我们采用最邻近插值的方法。对于某个观察站S (x,y),它只可能落在网格点上、落在网格线上、落入网格中三种情况:
1)落在网格点上的观察站正好可以直接将该站的预测值与实际值比较。
2)落在网格线上的观察站,在横向或者纵向上一定有与之相邻近的两个点 d1 (x1,x2)、d2 (x1,y2)。S点与两点间的某个更邻近时,我们可以取这一点与S相合,将这一点的预测数据作为S的预测数据。若与两者的距离相等,那么取两点上预测值的均值作为S点的预测值。
3)落在网格中的观测站,图示如下,那么周围临近点有四个d1(x1,y1)、d2(x2,y2)、d3(x3,y3)、d4(x4,y4)。观测站与网格点的距离(红线)为:。取最近的点,即,将上的预测值作为该观测站预测值。
图1 图2
使用最邻近插值,用原本的47 x 53个网格点去贴和这91个观察站点,将原本的比较复杂且庞大的数据2491(网格点数) x 4 x 41个简化为91 x 4 x 41个。
接着,我们总体的评价一下方法一和方法二。
通过绘制雨量折线图简单的表现并评判两种预测的贴切程度。下图是在第1站和第21站的第一个时段三种不同预测对于这41天预测的变化值。
对于某天我们一共进行了四次预测,再次我们引入切比雪夫距离作为判断一次降雨量的,每次预测与实际值比较获得误差值,并从这四个误差中提取出最大的误差作为第i天j站点预测的整体误差:。
这样我们就可以得到误差矩阵——方法一: 方法二
通过上述的方法得到的是一个41 x 91矩阵,矩阵中第i行第j列的元素代表第i天第j个站点在四个时段中预测与实际偏差的最大值即与实际的切比雪夫距离。
通过方法一方法二的误差矩阵,我们可以得到一个对比矩阵:
其中元素为。
我们建立一个评价函数
也就是说每次评价中,我们将第一种对于第二种方案的优于评价为1,即误差第一种相对于第二种更小,误差相等则平等均不具有优势,评价为0,若劣势于则评价为-1。
即通过判断矩阵中元素的值是否大于来评价这次的预测哪个方案更具有优势。
我们可以得到站点优势值.
下表为方案一对于方案二的评价矩阵的列向量累加,即方案一中每个站点预测的优势值。
方案一 | 站点1 | 站点2 | 站点3 | 站点4 | 站点5 | 站点6 | 站点7 | 站点8 | ```` |
优势值 | 9 | 9 | 12 | -15 | -9 | 0 | 5 | -3 | ```` |
通过计算累计的优势值R累加,即.若优势值大于零,则表明整体上方案一的误差较方案二小,方案一更加优秀。
通过计算得到R=194.
说明在对比每个站点四个时段预测的最大偏差中,方案一的最大偏差相对小一些,整体预测的误差应该比方案二要小。
同时,我们比较每个时段的误差值优势值。方法与上述类似不再赘述给出计算公式。
得出第一时段。通过观察每个时段的R的值发现,每个时段预测1的精确值的实际误差程度都比预测二要好,直观的体现了在插值过后比较点预测的结果方案一预测的相合度比方案二更加优秀。
问题二:
由于问题一得到方案一的预测更为精确,那么我们选择方案一作为研究目标。
等级 | 雨级 | 雨量(毫米) |
0 | 无雨 | —0.1 |
1 | 小雨 | 0.1—2.5 |
2 | 中雨 | 2.6—6 |
3 | 大雨 | 6.1—12 |
4 | 暴雨 | 12.1—25 |
5 | 大暴雨 | 25.1—60 |
6 | 特大暴雨 | 60.1— |
降雨量等级分类表
通过上述雨量分类,我们可以将方案一中的雨量数据表转化为雨量等级预报表。将矩阵转化为预报矩阵,P中元素最大为6,最小为0,i代表天数,j代表观测站点,k代表时段。
同样的根据上表将实际的天气预报矩阵给出,即通过上述分段将转换为。
由于群众对于天气预报的误差,心态各不相同,较为主观,我们简单的将预报错误带来的不满足感分为6级,也就是报错等级的最大误差到报错等级的最小误差。并得到满意度矩阵。且将群众的心情划分为下表的等级。
误差等级 | 心情 | 满意度积累权值 |
0 | 满意 | 1 |
1 | 勉强接受 | 0.9 |
2 | 有些不满 | 0.8 |
3 | 不满 | 0.4 |
4 | 十分不满 | 0.2 |
5 | 愤怒 | 0 |
6 | 愤怒 | 0 |
预报在较小误差内可以让人接受,也不至于影响人的情绪,一旦预报偏差过大则会影响生活生产、交通运输、商铺运营、城市排水等.所以在预报误差只有2以内时满意度积累值变化较小,当误差等级变大满意度变化递减变快。
通过统计计算出,发现,也就是说最大误差等级不会到达六级。我们采用模糊评估法,给出一组评估向量。
对于每次的预测误差评估归类以后,每个等级都有一个满意度权值。
误差等级 | 次数 | 占全部预测的百分比 |
0 | 12603 | |
1 | 2239 | |
2 | 58 | |
3 | 18 | |
4 | 3 | |
5 | 3 | |
MATLAB程序:
%导入测量的数据存入A中
for i=18:28
k=i-17;
A(:,:,k)=importdata(strcat('C:\Users\computer\Desktop\插值拟合作业\MEASURING','\0206',num2str(i),'.SIX'));
end;
for i=1:9
k=k+1;
A(:,:,k)=importdata(strcat('C:\Users\computer\Desktop\插值拟合作业\MEASURING','\02070',num2str(i),'.SIX'));
end;
for i=10:30
k=k+1;
A(:,:,k)=importdata(strcat('C:\Users\computer\Desktop\插值拟合作业\MEASURING','\0207',num2str(i),'.SIX'));
end;
%导入网格点的经度,纬度存入Lon,Lat中并转化为行向量
Lon=importdata('C:\Users\computer\Desktop\插值拟合作业\FORECAST\lon.dat');
Lat=importdata('C:\Users\computer\Desktop\插值拟合作业\FORECAST\lat.dat');
lon=reshape(Lon,2491,1)';
lat=reshape(Lat,2491,1)';
mlat(1:91)=A(:,2,1);%实测站点的经度
mlon(1:91)=A(:,3,1);%实测站点的纬度
%导入第一种预测的各个时间段的所有数据存入dis1中
K=0,I=1;
for I=1:4
K=0;
for j=618:628
K=K+1;
dis1(K,:,I)=reshape(importdata(strcat('C:\Users\computer\Desktop\插值拟合作业\FORECAST','\f',num2str(j),num2str(I),'_dis1')),2491,1)';
end;
for j=701:730
K=K+1;
dis1(K,:,I)=reshape(importdata(strcat('C:\Users\computer\Desktop\插值拟合作业\FORECAST','\f',num2str(j),num2str(I),'_dis1')),2491,1)';
end;
end;
%[LAT LON]=meshgrid(mlat,mlon);
%mrain=A(:,i,j)';
%对网格点上的预测降雨量进行插值使其与实测点的91个站的降雨量逼近,并存入rain中
for j=4:7
for i=1:41
Rain(1,:)=dis1(i,:,j-3);
r0=griddata(lat,lon,Rain,mlat,mlon,'nearest');
rain(i,:,j-3)=r0;
end;
end;
B=A(:,4:7,:);
%将实际测量值以41x91x4的格式存入me中
for i=1:4
for j=1:41
tem=B(:,i,j);
me(j,:,i)=tem';
end;
end;
%第二种预测得到的降雨量误差K=0,I=1;
for I=1:4
K=0;
for j=618:628
K=K+1;
dis2(K,:,I)=reshape(importdata(strcat('C:\Users\computer\Desktop\插值拟合作业\FORECAST','\f',num2str(j),num2str(I),'_dis2')),2491,1)';
end;
for j=701:730
K=K+1;
dis2(K,:,I)=reshape(importdata(strcat('C:\Users\computer\Desktop\插值拟合作业\FORECAST','\f',num2str(j),num2str(I),'_dis2')),2491,1)';
end;
end;
%对网格点上的预测降雨量进行插值使其与实测点的91个站的降雨量逼近,并存入rain中
for j=4:7
for i=1:41
Rain2(1,:)=dis2(i,:,j-3);
r0=griddata(lat,lon,Rain2,mlat,mlon,'nearest');
rain2(i,:,j-3)=r0;
end;
end;
%降雨量与实测值的误差
delta1=abs(rain-me);
delta2=abs(rain2-me);
%降雨量误差值的切比雪夫距离
for i=1:41
for j=1:91;
del1(i,j)=max(delta1(i,j,:));
del2(i,j)=max(delta2(i,j,:));
end;
end;
%对比矩阵cpa1,如果预测1误差比2小赋值1,比2大赋值-1
cpa1=del1-del2;
cpa1(cpa1<0)=-1;
cpa1(cpa1>0)=1;
rain(rain<0.1)=0;
rain(rain>=0.1&rain<2.6)=1;
rain(rain>=2.6&rain<6.1)=2;
rain(rain>=6.1&rain<12.1)=3;
rain(rain>=12.1&rain<25.1)=4;
rain(rain>=25.1&rain<60.1)=5;
rain(rain>=60.1)=6;
本文来源:https://www.2haoxitong.net/k/doc/f020e733dcccda38376baf1ffc4ffe473368fd06.html
文档为doc格式