热带气旋的运动规律研究——数学建模

1、热带气旋的运动规律研究摘要近年来，随着热带气旋对我国东南沿海地区的影响日益增大，由热带气旋造成的自然灾害对我国沿海地区造成的巨大的经济损失甚至人员伤亡也逐年增加，对热带气旋运动规律的了解以及对其发展方向的正确预测越来越受到大家的关注以及重视。本文围绕西太平洋海域的热带气旋进行了讨论，对热带气旋各指标之间的关系、热带气旋发源地的分布规律及其运动特征、季节变化对气旋强度的影响、热带气旋最佳路径的预测等问题分别建立了模型，并对求解结果进行了详细的分析。针对问题一，由附录给出的数据说明可知热带气旋的强度是根据风速进行划分的，所以包含的隐性条件即是强度与风速的关系可以不用求解。于是我们采用线性回归的方式

2、分别对强度与气压以及风速与气压进行建模计算，得出两组数值确实存在线性回归的关系，并得出了两两之间的函数关系式。针对问题二，为了求得热带气旋的发源地分布情况，先通过c+程序筛选出所有发源地的热带季风的数据。首先利用spass软件绘制出发源地的经纬度散点图，接着利用地统计学中kriging插值法，对海域上未采样点的热带气旋的出现概率、气压以及风速进行插值，绘制出热带气旋发源地的空间分布图，结果显示发源地的分布规律与气象局发布的热带气旋主要发源地分布相符。针对问题三，为了研究季节变化对热带气旋强度的影响，首先根据国家气象局的规定按月份划分出春夏秋冬并筛选出不同季节的热带气旋的强度数据，然后采用单因素

3、方差分析的方法，比较春夏秋冬四组的组内差异和组间差异。根据结果显示，季节变化对热带季风强度变化的影响是显著的。其次考虑到1949年到2011年时间跨度较长，为了检验时间变化是否会对测试结果产生影响，我们又分别对春夏秋冬四组数据分别进行单样本t检验，从而优化季节变化对热带气旋强度影响的讨论。根据结果显示，春季，夏季，秋季，冬季各个组间差异不显著，所以，单因素方差分析方法找出的季节间影响关系成立。关键词：线性回归分析 krging插值法 单样本t检验 单因素方差分析一、问题重述热带气旋(Tropical Cyclone)是发生在热带或副热带洋面上的低压涡旋，是一种强大而深厚的热带天气系统。热带气旋

4、通常在热带地区离赤道平均3-5个纬度外的海面(如西北太平洋，北大西洋，印度洋)上形成，其移动主要受到科氏力及其它大尺度天气系统所影响，最终在海上消散、或者变性为温带气旋，或在登陆陆地后消散。登陆陆地的热带气旋会带来严重的财产和人员伤亡，是自然灾害的一种。我国的东南沿海地区也常年遭受热带气旋的影响，当地居民的生产生活也遭受了严重的影响，因此掌握热带气旋的运动规律对我们进一步监控热带气旋，做好相关防范措施等都有积极的作用。附件提供了自1949年至2011年西北太平洋(含南海, 赤道以北, 东经180°以西) 海域热带气旋每6小时的位置和强度, 按年份分别放在单独的文本文件中。请利用这些数

5、据建立数学模型，回答以下问题：（1）探讨热带气旋强度、气压与风速之间的关系；（2）给出热带气旋发源地的空间分布，并评价不同发源地产生的热带气旋的运动特征；（3）研究季节变化对热带气旋强度的影响；（4）提出一种能预测热带气旋最佳路径的方法，并验证其有效性。二、问题分析随着极端天气对人类日常生活的影响越来越大，对各种极端天气的科学研究正越来越得到人们的重视。西北太平洋上的热带气旋对我国一直有着较大的影响，对该片海域上的气旋的研究也成为我国气象科学中不可缺少的一部分。热带气旋是发生在热带或副热带洋面上的低压涡旋，是一种强大而深厚的热带天气系统，登陆陆地的热带气旋会带来严重的财产和人员伤亡，是自然灾害

6、的一种。本文旨在对建国以来出现在我国东南海域的热带气旋的指数及运动规律进行研究，主要包含了四个相关性问题。这四个问题分别为热带气旋指数之间的关系问题、发源地分布规律及运动特征问题、季节变化对热带气旋强度影响问题和对热带气旋最佳路径的预测问题。2.1 问题一：探讨热带气旋强度、气压与风速之间的关系问题一要求探讨出热带气旋多个指标之间的内在关系。能够猜测，热带气旋的指标之间应该存在一定的可以量化的函数关系，为了达到这样一个目的，本文考虑采取线性回归的方式进行研究。但应当注意到，在附录所给的数据解释中，热带气旋的强度是根据热带气旋的风速进行分级的，即本题的隐藏条件便是热带气旋的强度与风速之间的关系是

7、不需要进行讨论的，所以本文着重讨论的是热带气旋强度与气压以及气压与风速之间的关系。对于数据的选取问题，附件所给的数据表中给出了两种风速，即近中心最大风速和平均风速，通过对附件所给的数据的初步浏览可知，只有近中心最大风速测量得比较全面，平均风速缺测比较严重，所以采用的是近中心最大风速；对于强度数据，由于我们要预测以及预防的是热带气旋是否有可能演变为台风进而对沿海地区造成影响，而对一个气旋危害程度研究最需要的即是其最大强度，一旦最大强度超过了5，即达到了台风的强度，我们就有必要对其进行预防以及检测，这也是我们气象工作的核心目的所在，所以本文采用的是每一个热带气旋的最大强度进行研究。要想探讨两个因素

8、直接的关系，本文采用了线性回归的研究方式，通过spss中的线性回归方程工具对附件所给的数据进行分析，从而得到本问题所需要的结果。2.2 问题二：热带气旋的发源地分布问题二要求确定热带气旋发源地的分布，传统统计学方法不考虑测定参数和测定位置之间的关系，只根据空间不同位置的测定数据，无法估计未测定区域所考察的热带气旋的相关特性的空间分布。所以本文考虑采用地统计学方法。Kringing插值法是地统计学中一种重要的数据处理方法。本文通过Kriging空间插值法，对未采样的地区进行模拟，利用MATLAB工具绘制出土壤重金属的空间分布规律图，得出不同区域出现热带季风概率的分布图，并分析了其空间分布特征及其

9、变异规律。2.3 问题三：研究季节变化对热带气旋的影响首先，我们按照国家气象局的规定，按月份划分成春夏秋冬。春季：3月5月，夏季：6月8月，秋季：9月11月，冬季：12月2月。对春夏秋冬四组数据进行单因素方差分析，比较春夏秋冬四组的组内差异和组间差异。从而比较出季节变化对热带气旋强度的影响。其次，我们考虑到附件数据是从1949年2011年的，时间跨度较长，气候变化对热带气旋也会产生影响，所以，我们采用单样本t检验的方法，检验春夏秋冬组内是否有显著性差异，从而优化季节变化对热带气旋强度影响的讨论。三、条件假设与符号说明3.1 条件假设1、 不考虑外部其他因素对热带气旋各指标的相互关系的影响；2、

10、 假设只有季节变化会对热带气旋的强度产生影响；3、 假设附件所给的统计数据能够反映西北太平洋海域热带气旋的产生及运动状况。3.2 符号说明符号 解释与说明r(h) 变异函数 C0 块金效应 C1 基台值减去C0的值 Ck 第个k个点的扩散函数四、模型准备4.1 数据的筛选首先对数据进行预处理，由附件可知，当WND=0时，表示数据缺测。用spss删除WND=0的相关数据。模型假设中，风速取两分钟平均近中心最大风速而不取两分钟的平均风速，是因为对建筑物或者是船只产生影响时，要考虑在最大风速下产生的最坏的结果，而不能用平均风速来衡量。筛选条件见附件对筛选完成的数据进行线性回归分析4.2 热带季风发源

11、地分布概率及各项指标检验4.2.1起源地空间图利用C+C程序提取1949到2011年每一种热带气旋的发源地的数据再对所得数据进行坐标投影变换，用MATLAB软件分别绘制出散点分布图,如图所示： 空间分布散点图结合PS把散点图和地图进行处理得到热带气旋发源地的经纬度在地图上空间分布的。 有地图的空间分布散点图 台风源地分布在西北太平洋广阔的洋低纬洋面上。西北太平洋热带扰动加强发展为台风的初始位置，在经度和纬度方面都存在着相对集中的地带。在东西方向上，热带扰动发展成台风相对集中在4个海区：南海中北部的海面、菲律宾群岛以东和琉球群岛附近海面、 马里亚纳群岛附近海面、马绍尔群岛附近海面。所以我们按照经

12、度的分布给热带气旋的发源地进行划分。所以100°E119°E、126°E140°E、144°E146°E、150°E180°E。 4.2.2热带气旋各元素概率分布检验 Kriging空间插值法是以变异函数为基础的，其参数设置和变异函数模型的选择对内插值效果影响很大3。变异函数的计算一般要求数据符合正态分布，否则可能存在比例效应，比例效应的存在会使实验变异函数产生畸变，抬高基台值和块金值，增大估计误差，变异函数点的波动大，甚至会掩盖其固有的结构，消除比例效应主要是通过原始数据取对数来消除，所以本文先对5个影响因素进行

13、分布检验。 检验数据的正态分布性有多种方法：频率分布直方图法、2法检验法、Q-Q图、P-P图、Kolmogorov-Smirnov(K-S)检验法等。本文利用SPSSl9.0统计软件中的P-P图对数据集5个影响因素进行了正态分布检验，若样点在正态分布P-P图上呈直线散布，则被检验数据基本上成一条直线。由上图可见，对5中因素进行正态分布检验后，5中因素的正态分布P-P图大致是一条直线，所以其均服从正态分布五、模型的建立与求解5.1 问题一：探讨热带气旋强度、气压与风速之间的关系回归分析方法，是研究要素之间具体的数量关系的一种强有力的工具，能够建立反映地理要素之间具体的数量关系的数学模型，即回归模

14、型。本小结利用线性回归的分析方法，对热带气旋的强度与气压以及风速与气压之间的关系进行探究，从而得出待求的两组数据之间的相互关系。1、模型的建立回归分析方法，是研究要素之间具体的数量关系的一种强有力的工具，能够建立反映地理要素之间具体的数量关系的数学模型，即回归模型。一元线性回归模型的基本结构形式：假设有两个地理要素（变量）x和y，x为自变量，y为因变量。则一元线性回归模型的基本结构形式：a和b为待定参数；=1,2,n为各组观测数据的下标； a为随机变量。如果记a和b 分别为参数a与b的拟合值，则得到一元线性回归模型。2、模型的求解 本小节利用spass软件对上述两组数据进行线性回归分析，结果如

15、下：系数a模型非标准化系数标准系数tSig.B标准 误差试用版1(常量)63.993.315203.399.000中心最低气压-.062.000-.629-194.466.000a. 因变量: 强度强度：y，中心气压：xY=-0.062x + 63.993由强度的标准化残差直方图可以看出，残差大致符合正态分布。系数a模型非标准化系数标准系数tSig.B标准 误差试用版1(常量)632.86853.93411.734.000中心最低气压-.615.055-.047-11.233.000a. 因变量: 风速风速：y，中心气压：x Y=-0.615x + 632.868 由上述结果可以看出，所选取的

16、两组数据确实存在相应的函数关系。5.2 问题二：热带气旋的发源地分布热带季风的发源地具有空间分布的特点，有效地对热带季风发源地出现概率进行统计能帮助我们更科学地预测以及预防热带季风所造成的灾害。热带季风发源地的空间变异性研究是对西北大西洋海域热带季风的变化特征、趋势走向及不同方向差异性等空间结构特征的探讨。空间插值方法是分析热带季风发源地的空间分布特征及其变异规律的有效工具。因此，本文先对热带季风发源地的空间变异性进行讨论。1 模型的建立 利用Kriging插值法进行空间插值 时，均要求变量具 有符合本征假设的 规定。所谓本征假设课表示为以下公式：（1） 假设 x 表示采样点的位置坐标。同一重

17、金属元素浓度均值存在并且不取决于x，任意距离为 h 的两点间的差值的数学期望为0： EZ(x)= m，对于任意x；EZ(x)- Z(x+h)=0 (1)(2 ) 式中，E表示数学期望；Z(x)和Z(x+ h)表示空间区域内距离为h的两个位置的观测值。（2） 对于任何距离 h ，变量 Z ( x) - Z ( x + h) 具有有限方差，此方差不取决于 x ，即对于任何 x 和 h ，下式成立： VarZ ( x + h) - Z ( x)=EZ ( x + h ) - Z ( x) - EZ ( x + h) - Z ( x) =E Z x+h-Z ( x)2=2g (h)式中，g (h) 称

18、为变异函数或半方差函数，变异函数可用于表征随机变量的空间变异结构，或空间连续性。 Kriging 插值法的一般公式为： Z（x0）*=i=1niZ(xi) (3)式中，Z(xi)为位于区域内xi位置的观测值（i=1,2,，n；n为样本总数）；x0为某个未采样点；Zx0* 表示采样点的观测值，即插值；i 为权，且满足下式： i=1ni=1 (4)选取 li ，使Z*(x0)的估计无偏差，并使其方差 s 小于任意观测值其他任意线性组合的方差。由于 li 值的在普通Kriging插值中是通过样本半方差求得的故选择一个较好的变异函数模型来拟合样本的变异性对于Kriging插值的效果影响较大。常用的变异

19、函数模型有球形模型、指数模型，表达式如下：球状模型 h=c0+c11.5ha-0.5ha3,0<h<ah=c0+c1, h>ah=0,h>0 (5)指数模型 h=0,h=0h=c0+c11-e-ha,h>0 (6) 其中h为分离距离，C0为块金值，C1为部分基台值，C0+C1为基台值，a为变差距离，即曲线到达基台值时所对应的分离距离。对于指数模型和高斯模型而言。变差距离，即曲线到达基台值时所对应的分离距离。对于指数模型和高斯模型而言变差距离即为变异函数曲线达到95时的分离距离。一般认为块金值代表随机变异的量基台值代表变量空间变异的结构性方差而块金系数是块金值与基台

20、值的比值，用于反映变量的空间自相关程度。 2 模型求解 由前面分析可知，热带气旋的5个元素均成对数正态分布状态，适合于Kriging 法插值。变异函数模型的选取要比较两个变异函数模型的优劣，采用如下的一些统计指标。 1） 计算残差 残差值可用于定量估计源数据的 z 值与 Kriging 插值后的值之间的一致性。根据残差结果可以计算平均残差值和标准残差等。 Zns=Zdat-Zgrd (7)式中，Zns为残差值，Zdat为源数据的Z值，Zgrd为Kriging插值后的Z值。2） 均方根预测误差 RMSPE 21nik=1nZi*xk-Zixk2=RMSPE （8）3） 相对均方差RMSE 1s2

21、nik=1niZi*xk-Zi(xk)2=RMSE (9)在计算RMSPE 和 RMSE 时，需要从采样点中随机预留出一定数量的点作为检测点，而采用其他的采样点作为插值的源数据，通过 Kriging 插值计算各监测点的估计值，式（8）（9）中，Zi ( xk ) 和 Zi (xk ) 分别表示检测点的Z观测值和Z估计值。 对球状模型和指数模型用MATLAB软件的Kriging插值工具箱进行Kriging 插值，并对结果进行统计 ,这里使用 RMSPE 作为变异函数模型效果判别空间分布与气压的三维曲面图 空间分布与风速的三维曲图热带气旋发源地空间频率分布图通过以上3图可以看出西北太平洋热带气旋发

22、源地分布规律，中国台风网公布的主要发源地区域相符合，可以作为分析热带气旋发源地分布以及移动规律的数据基础。5.3 问题三：研究季节变化对热带气旋的影响要求研究季节变化对热带气旋强度的影响。首先，我们按照国家气象局的规定，按月份划分成春夏秋冬。春季：3月5月，夏季：6月8月，秋季：9月11月，冬季：12月2月。我们采用单因素方差分析的方法，比较春夏秋冬四组的组间差异。从而比较出季节变化对热带气旋强度的影响。5.3.1 模型一：单因素方差分析单因素方差分析,用于完全随机设计的多个样本均数间的比较，其统计推断是推断各样本所代表的各总体均数是否相等。1. 模型的建立完全随机设计的单因素方差分析是把总变

23、异的离均差平方和SS及自由度分别分解为组间和组内两部分，其计算公式如下。SS组间=离均平方和/组间自由度SS组内=离均平方和/组内自由度SS总=SS组间+SS组内单因素方差分析：核心就是计算组间和组内离均差平方和。两组或两组以上数据，大组全部在一组就是组内，以每一组计算一均数，再进行离均平方和的计算：SS组间=组间离均平方和，MS组间=SS组间/组数-1（注：离均就有差的意思了！）SS组内=组内离均平方和，MS组内=SS组内/全部数据-组数F值=MS组间/MS组内查F值，判断见下面的分析步骤部份2. 模型的求解从表格中，我们看出，季节的变化对气旋强度的影响是显著的。5.3.2 模型二：单样本t

24、检验t检验是用t分布理论来推论差异发生的概率，从而比较两个平均数的差异是否显著。它与Z检验、卡方检验并列。t检验是戈斯特为了观测酿酒质量而发明的。戈斯特在位于都柏林的健力士酿酒厂担任统计学家，基于Claude Guinness聘用从牛津大学和剑桥大学出来的最好的毕业生以将生物化学及统计学应用到健力士工业程序的创新政策。戈斯特于1908年在Biometrika上公布t检验，但因其老板认为其为商业机密而被迫使用笔名（学生）。实际上，跟他合作过的统计学家是知道“学生”的真实身份是戈斯特的。1、模型的建立单总体t检验是检验一个样本平均数与一个已知的总体平均数的差异是否显著。当总体分布是正态分布，如总体

25、标准差未知且样本容量小于30，那么样本平均数与总体平均数的离差统计量呈t分布。单总体t检验统计量为：    独立样本t检验统计量为：  S1 和 S2 为两样本方差；n1 和n2 为两样本容量。（上面的公式是1/n1 + 1/n2 不是减！）配对样本t检验统计量为：2、模型的求解本文采用spss软件含有的t检验函数进行模型的求解，结果如下：单个样本统计量N均值标准差均值的标准误春季49582.722.084.030夏季238362.652.077.013秋季247953.05

26、2.107.013冬季42282.581.773.027单个样本检验检验值 = 0 tdfSig.(双侧)均值差值差分的 95% 置信区间下限上限春季91.7774957.0002.7162.662.77夏季197.13223835.0002.6522.632.68秋季228.27024794.0003.0543.033.08冬季94.4824227.0002.5772.522.63由单样本t检验的结果可知，春季，夏季，秋季，冬季各个组内差异不显著。所以，单因素方差分析方法找出的季节间影响关系成立。六、模型评价本文所建立模型的优缺点主要体现在一下几个方面：对问题二：优点：通过地统计学而非传统统

27、计学进行分析，增加了科学性； 缺点：kriging插值依赖于数据的分布，难以普遍适用。对问题三：优点：采用单样本t检验的方式验证了不同年份的不同春季对热带季风强度的影响几乎不大，进一步证明了是季节的变化很大程度地影响了热带季风的强度； 缺点：单样本t 检验要求数据成正太分布，普遍使用的可能性不大。附录附录1 matlab画二维图clc; clear all; close all;filename = 'C:UsersAdministratorDesktop13Ä꽨ģ³õʼֵ.xls'num1, txt1, raw1 = xlsread(filename, '¸½¼þ1');A = num1(:, 2:3);x = A(:, 2);y = A(:, 1);figure; hold on; box on;plot(x, y, 'r.');view(2)colorbarcolormap(0 0 0; 0.1 0.5 0.4; 0.2 0.