基于ASAMC算法的气象数据多变点估计

许 欢, 谭常春

(合肥工业大学 数学学院,安徽 合肥 230009)

基于ASAMC算法的气象数据多变点估计

许 欢, 谭常春

(合肥工业大学 数学学院,安徽 合肥 230009)

文章主要研究了安徽省某气象台站年平均气温、月平均气温的结构性变化情况。该文在对气温数据进行正态性检验的基础上,运用ASAMC(annealing stochastic approximation Monte Carlo)算法对气温数据进行统计分析,估计出平均气温结构性变化的位置,并探索发生结构性变化的气象因素和非气象因素。

ASAMC算法;Bayes多变点模型;M-H算法;平均气温

如今大数据时代,变点问题在各个领域引起越来越多的关注,尤其是在气象和金融方面。有关内容可参考文献[1-2]。文献[3]应用Bayes方法检验正态分布的均值是否改变,判断正态分布中是否产生变点;文献[4]运用似然比方法讨论了多元正态分布的均值变点的检测和估计;文献[5]运用信息论准则讨论了在变点个数已知的条件下,正态分布均值多变点的检测和估计;文献[6]给出了独立正态分布序列中变点个数的估计;文献[7]提出了检验正态分布方差变点的3个统计量;文献[8]讨论了在均值已知的条件下,正态分布方差变点的假设检验问题,运用变点的参数检验、非参数检验和Bayes检验3种方法,并对它们进行比较分析;文献[9]研究了在方差未知的条件下,多元正态分布均值变点的个数和位置的检测及估计。

由于样本空间的维度、参数都不是固定的,这给寻找样本空间的变点数量和变点位置带来了一定的困难。随着Bayes方法在变点检测领域的广泛应用,越来越多的人使用Bayes方法检测和估计序列中的多变点位置,很好地解决了之前的问题。Bayes方法的优点之一是允许在不同的子空间中有不同的参数,并且计算相对简便。大多数Bayes变点问题使用MCMC(Markov chain Monte Carlo)算法或者Gibbs采样器算法去解决已知或未知数量的变点问题。文献[10]提出SAMC(stochastic approximation Mmonte Carlo)算法和ASAMC(annealing stochastic approximation Monte Carlo)算法,2种算法在不同维度的多参数子空间中使用Bayes后验概率函数估计出变点。文献[11]比较了ASAMC算法与SAMC算法,发现ASAMC算法不仅计算量小，而且更容易找到多变点的个数及位置。

在气象数据分析中,一系列无差异气象数据能显示出该地区的气象特征,一般可用变点检测方法判断降水、气压、湿度、气温数据序列是否发生了结构性变化。文献[1]运用统计方法提出标准正态检验(standard normal homogeneity test,SNHT);文献[9]改进了该方法,并应用于气象数据序列的均一性检验,提出参考台站相关系数;文献[12]研究了气温数据序列非均一性检验的3种方法;文献[13]研究了定点观测均一性检验方法;文献[14]通过研究、分析中国气象台站历史资料,认为台站迁移对观测数据序列影响最大;文献[15-16]用Γ分布检测气温序列中的变点;文献[17]通过研究原始气象数据序列得出台站迁移对年平均温度影响较大。ASAMC算法不参考周围站台气候变化趋势的影响,直接研究单站台气象数据序列的均一性,而且一次能检测出多个变点的数量及位置。

本文主要研究1955-2010年间安徽省年度平均气温、月度平均气温变结构问题。在平均气温为正态分布的条件下,运用Bayes方法建立平均气温序列多变点模型,并根据对数后验概率值估计出变点的位置。

1 正态分布多变点模型描述

设z=(z1,z2,…,zn)为一列来自于正态分布的观测数据,满足

其中,γ、δ、λ为超参数;ak=(k+1)(γlnδ-ln Γ(γ))+ln(n-1-k)!-ln(n-1)!+klnλ。

2 ASAMC方法介绍

ASAMC算法的具体步骤如下：

(1) 产生一组长度为n的样本数据z=(z1,z2,…,zn)，然后检验样本z是否服从正态分布。本文利用样本z=(z1,z2,…,zn)画出概率图,观察这些点是否在一条直线附近。若是,则认为z服从正态分布;否则,认为z不服从正态分布,要进行修正。一般作对数变换后,再用画正态概率图的方法检验。

其中,q(X1,x)为proposal分布。设接受概率为:

存在概率μ,若μ≤α,则X1=X1;否则X1=x。

(5) 多次重复步骤(4),找出U2,…,UN及相应的X2,…,XN。

(6) 画出频数直方图。观察频数直方图,若各组频数分布有显著差别,则认为有变点;反之,则认为没有变点。

(7) 若认为有变点,则找出最终最大的Ui(x)及Xi,其中Xi中所有1的位置就是变点的位置。

3 安徽省年平均气温实例分析

本文所使用的数据来源于安徽省气象局,为安徽省58XXX台站从1955—2010年的年平均气温。由于数据年限跨度较长,期间发生了很大的变化,例如设备的改进、监测站地址和监测方法的改变等,都有可能引起监测数据的突变,使用年平均气温避免气温周期变化的自然因素。安徽省1955—2010年的年平均气温如图1所示。

图1 安徽省1955—2010年的年平均气温

对年平均气温数据进行正态性检验,其概率如图2所示。由图2可知数据点近似在一条直线附近,可以认为这组数据近似服从正态分布。

假设在这一组观察数据中变点的个数不超过4个,且变点的位置不出现在尾处(ci∈{1955,1956,…,2009}),受文献[10]启发,选取均匀分布为Proposal分布,超参数分别选定为λ=1,γ=0.05,δ=0.05,概率μ=0.131 4,迭代次数为600。根据样本的对数后验概率函数U(x),为了使相对样本频率大于80%，把空间划分为2个子空间。子空间为(-∞,-1)时，相对样本频率为118.667%；子空间为(-1,+∞)时，相对样本频率为81.333%。

图2 安徽省1955—2010年的年平均气温的正态概率

运用ASAMC算法计算出100组对数后验概率函数值,选10组最大的对数后验概率函数值，见表1所列。将所有的最大对数后验概率函数值按照从大到小的顺序排列,找到排序第1的对数后验概率函数值,即100组中最大的对数后验概率函数值。

表1 10组最大的对数后验概率值

注：有下划线的表示对数后验概率值对应的变点，下同。

结合表1可以看出,大多数最大的后验概率函数值都是带有2个变点的。根据变点可能出现的位置画出频数直方图,如图3所示。

图3 年度变点可能出现位置的频数直方图

由图3可知,各组频数分布有显著差别,因此认为存在2个变点。

根据最终最大的对数后验概率函数值找到变点的具体位置为(1993,2001),如图4所示。

变点位置

导致长期气候数据序列发生变化的因素主要有气候因素和非气候因素，非气候因素主要是观测台站迁移及观测方法的改变。通过查阅该台站的历史资料,发现该台站从1955—2010年期间,先后于1981年、1991年、2002年和2009年发生站点迁移,因此可以认为安徽省58XXX台站的年平均气温发生的结构性变化是由于台站迁移因素造成的,而不是由气候变化本身的因素导致的。

4 安徽省月度平均温度分析

考虑到在计算年平均温度时,使用加权平均方法可能对检测的结果产生影响，因此进一步考虑月度平均气温的变结构问题。为了避免气温周期变化的因素,将相同月份的月平均温度数据组成1组。对这12组数据应用ASAMC算法，分别检测是否存在变点并估计变点的位置。通过计算对数后验概率,并画出频数直方图,发现12组数据中只有2月、3月、4月、9月存在变点,频数直方图如图5～图8所示，对数后验概率值见表2所列。

从图5～图8和表2中可以看出，2月份、4月份和9月份均存在2个变点,变点分别发生在1972年和1990年、1980年和2008年、1997年和2001年;3月份有1个变点,其位置在1996年。7个变点中只有1个变点在1991年附近,所占比例为14.3%,1个变点在2002年附近，所占比例为14.3%。

通过以上分析发现，月度平均气温序列的变结构点位置与年度平均气温的变结构点的位置并不完全一致，其中月度平均气温结构变化发生在1980、1990、2001、2008年,应该是由该台站发生迁移的原因造成的;而导致月度平均气温在1972、1996、1997年发生变化的原因有待进一步分析。7个变点中有4个是由于台站迁移造成的,所占比例为57.2%,表明台站迁移对气温序列的均一性影响很大。

图5 2月份变点可能出现位置的频数直方图

图6 3月份变点可能出现位置的频数直方图

图7 4月份变点可能出现位置的频数直方图

图8 9月份变点可能出现位置的频数直方图

2月份变点对数后验概率值3月份变点对数后验概率值4月份变点对数后验概率值9月份变点对数后验概率值(1972,1990)-63.425(1996)-45.566(1980,2008)-37.493(1997,2002)-24.476(1973,1990)-64.056(1957,1996)-46.877(1991,2008)-37.535(1998,2000)-25.890(1974,1990)-64.948(1957,1993)-47.603(1980)-37.810(1997,2000)-26.086(1972,1990,1997)-66.137(1957,1963,1996)-48.229(1991,2007)-37.961(1995,2001)-26.641(1973,1991,1993)-66.936(1957,1965,1992)-49.313(1980,1991,2008)-40.705(1974)-27.738(1973,1991,1994)-68.037(1957,1965,1994)-49.676(1981,1991,2008)-41.332(1973,1997,2002)-28.104(1972,1976,1991)-68.113(1957,1964,1989)-50.332(1991,2006,2008)-41.399(1971,1997,2001)-28.579(1972)-68.545(1957,1992,1996)-50.646(1980,1991,2006)-41.677(1970,1997,2001)-28.869(1991,1992,1998)-53.712(1965,1975,1991,2005)-46.365(1989,1991,2006,2008)-46.816(1980,1997,2001)-29.106(1985,1995,2001)-29.134

5 结 论

本文运用ASAMC算法研究安徽省58XXX台站在1955-2010年的年度平均气温和月度平均气温的结构性变化问题，得到了年度平均气温序列在1993年和2002年发生了结构性变化。由12组月度平均气温数据序列发现，部分月度平均气温序列存在1个或2个结构变点,部分序列不存在结构性变化。这表明ASAMC算法在多变点的检测和估计方面是有效的,与以往气温序列1次只能估计出1个变点的算法相比,ASAMC算法可以同时检测和估计出温度序列中的多变点的个数和位置。另外,在分析气象因素数据时(比如年平均温度),期望气象数据序列是均一的(不存在变点)。若气象数据不满足均一性,则需要对数据进行适当修正。在检测出气象数据序列存在结构变化后,对其进行修正使之满足均一性条件需要进行更深入的研究。

[1] HAWKINS D M.Testing a sequence of observtions for a shift in location[J].J Amer Statist Assoc,1977,72(357):180-186.

[2] SCHIPPER S,SCHMID W.Sequential methods for detecting changes in the variance of economic time series[J].Sequential Analysis,2001,20(4):235-262.

[3] CHERNOFF H,ZACKS S.Estimating the current mean of a normal distribution which is subjected to changes in time[J].The Annals of Mathematical Statistics,1964,35(3):999-1018.

[4] JAMES B,JAMES K L,SIEGMUND D.Asymptotic approximations for likelihood ratio tests and confidence regions for a change-point in the mean of a multivariate normal distribution[J]. Statistica Sinica,1992,2:69-90.

[5] MIAO B Q,ZHAO L C,KRISHNAIAH P R.On detection of change points usingmean vectors[J].Acta Mathematicae Applicate Sinica (English Series),1993,9(3):193-203.

[6] LEE C B.Estimating the number of change points in a sequence of independent normalrandom variables[J].Statistics and Probability Letters,1995,25(3):241-248.

[7] WANG J L,BHATTI M I.Three tests for a change-point in variance of normal distribution[J].Chinese Journal of Applied Probability and Statistics,1998,14(2):113-121.

[8] 王黎明.均值已知时,正态分布方差变点的假设检验[J].曲阜师范大学学报,1995,21(3):45-49.

[9] 缪柏其,赵林城,谭智平.关于变点个数及位置的检测和估计[J].应用数学学报,2003,26(1):26-39.

[10] LIANG F.Annealing stochastic approximation Monte Carlo for neural network training[J].Machine Learning,2007,68(3):201-233.

[11] KIM J,CHEON S.Bayesian multiple change-point estima tion with annealing stochastic approximation Monte Carlo[J].Computational Statistics,2010,25(2):215-239.

[12] 宋超辉,孙安建.非均一性气温序列订正方法的研究[J].高原气象,1995,14(2):215-220.

[13] 李庆祥,刘小宁,张洪政,等.定点观测气候序列的均一性研究[J].气象科技,2003,31(1):3-10.

[14] 吴增祥.气象台站历史沿革信息及其对观测资料序列均一性影响的初步分析[J].应用气象学报,2005,16(4):461-467.

[15] 吴必文,温华阳,惠军.基于Γ分布的气压序列非均一性检验方法初探[J].应用气象学报,2008,19(4):497-501.

[16] 吴必文,温华阳,曹军.一种新的基于Γ分布气温序列非均一性检验方法[J].合肥工业大学学报(自然科学版),2008,31(9):1521-1524.

[17] 周建平,孙照渤,倪冬鸿,等.中国气象台站迁移对年平均均一性的影响[J].大气科学学报,2013,36(2):139-146.

(责任编辑 朱晓临)

Multiple change-point estimation for meteorological data based on ASAMC algorithm

XU Huan, TAN Changchun

(School of Mathematics, Hefei University of Technology, Hefei 230009, China)

In this paper, the investigation on structural changes of annual and monthly mean temperature based on the statistics of a meteorological station in Anhui Province is conducted. Firstly, the normality test of temperature statistics is conducted. Then the temperature data is analyzed by using the annealing stochastic approximation Monte Carlo(ASAMC) algorithm, and the structural changes of mean temperature are estimated. The meteorological and non-meteorological factors leading to the structural changes are also explored by this method.

annealing stochastic approximation Monte Carlo(ASAMC) algorithm; Bayesian multiple change-point model; M-H algorithm; mean temperature

2015-06-12

国家自然科学基金资助项目(11201108);全国统计科研计划重点资助项目(2012LZ009)和中央高校基本科研业务费专项资金资助项目(JZ2015HGXJ0177)

许 欢(1990-),女,安徽安庆人,合肥工业大学硕士生; 谭常春(1977-),男,安徽合肥人,博士,合肥工业大学副教授,硕士生导师.

10.3969/j.issn.1003-5060.2016.12.024

O212.7

A

1003-5060(2016)12-1719-05