引言

陆地与大气紧密相连,地震前大气的异常情况受到地震学者们的重视。前人进行了大量探索性研究工作,以寻求地震前兆(徐国钧等,1993李贵福等,1996曹新来等,1997张铁宝等,2013)。在这些地震研究中,长期连续完整且均一性较好的气象资料是研究地震前气温异常的数据基础。但是,由于各种原因(如环境干扰、硬件故障等),历史气温资料缺测现象时有发生,导致气温观测资料缺测,造成历史资料的不连续(王建国等,20102013姚会琴等,2012)。

中国许多学者开展了对日、月、年时间尺度的气象资料缺测插补研究,并利用一维车贝雪夫多项式展开、线性回归、标准序列法、基于SVD的迭代等方法对中国部分地区的气象日、月、年值资料进行了恢复性实验(张秀芝等, 1996a, 1996b涂诗玉等,2001张永领等,2006王海军等,2008余予等,2012),但对整点气温值进行缺测插补鲜见文献报道。在国外,Huth等(1995)建立回归模型来插补缺测的日气温数据,Eischeid等(2000)采用空间差值法,插补后建立美国西部40年逐日气温和降水数据集,但这些插补模型只用于1个或数个缺测日数据的插补,不合适用于连续几个月数据缺测的情况。DeGaetano等(1995)引用改进的标准序列法,对美国东北部近400个站的日最高、最低气温缺测值进行了插补。标准序列法和线性回归法解决了插补长期连续缺测数据的问题。整点气温值的缺测插补方法可借鉴日平均、最高、最低气温值的缺测插补方法。由于线性回归法具有更好的统计性能和稳健性,本文对线性回归法进行改进,考虑了距离因素,采用线性回归模型参数求解法,解决了连续数日甚至数月造成的气温缺测问题,为地震前后气温变化特征的研究提供长期连续完整可靠的数据资料,也为今后开展卫星遥感红外亮温与卫星过境时刻气温的对比研究奠定了数据基础。

1 研究区及数据资料
1.1 研究区介绍

本文选取的研究区范围为37°—42°N、113°—119°E,在研究区内收集到15个地震观测站(气温观测站)的气温整点值数据,观测站分布情况如图 1所示。


图 1 气温观测站分布 Fig. 1 Distribution of temperature observatory sites
1.2 气温资料

本文收集了15个地震观测站气温数据,气温指地面以上1.5m处百叶窗测得的空气温度。太阳的热能被地面吸收后,地面再通过辐射、传导和对流把热传给空气,这是空气中热量的主要来源。气温的观测范围-30℃—70℃,精度0.1℃。由于各观测站安装时间、停测时间不同,导致观测时间的长度不等,多数观测数据起止时间为2007年1月1日,截止时间为2014年12月31日。

1.3 气温数据的缺测情况

由于观测环境、仪器设备故障等原因,部分观测站的观测数据缺测,包括仅缺测1个值、缺测1天的值(即24个观测值)、连续缺测几天甚至1个多月导致长达数千个观测值连续缺测等情况。基于震例研究,本文仅讨论唐山站气温数据的缺测情况(表 1),可以看出该站缺测情况较严重,如2008年8月18日—9月15日连续缺测696个整点值,2009年4月2日—5月13日连续缺测1000个整点值,2012年6月2日—24日连续缺测552个整点值,2013年1月21日—2月5日连续缺测384个整点值。唐山站数据连续缺测时间较长,连续缺测几十个值的频率较高,同时该站还存在错误值,如2010年11月22日18时的观测气温整点值是59.701,类似的错误值在每年均有出现。2010、2011年的观测数据较完整。

表 1 唐山站整点气温数据的缺测统计 Table 1 Missing data in integral point temperatures from the Tangshan site
1.4 气温数据的日变情况

图 2为唐山站2012年12月17日—22日气温整点值变化曲线,可以看出气温日变明显,最高温一般出现在14时左右,最低温出现在8时左右,符合日变规律。


图 2 唐山站气温变化曲线 Fig. 2 Daily variation of interpolated integral point temperatures for Tangshan site
1.5 气温数据的月变情况

对唐山站2012年1、4、7、10月的日均值进行分析,绘制了相关曲线,如图 3所示。由图可以看出,气温变化整体表现出明显的夏高冬低的年变规律;1月、7月气温变化平稳,4月气温表现为升温过程,10月气温表现为明显的降温过程。


图 3 唐山站气温日均值变化曲线 Fig. 3 Monthly variation of interpolated daily mean temperatures for Tangshan site
1.6 气温数据的年变情况

唐山站2008—2013年连续的气温整点值数据变化曲线如图 4所示。由图可以看出,气温表现出明显的夏高冬低的年变特征,夏季最高日平均温度接近31℃,冬季最低日平均温度接近-14℃;气温在冬季和夏季处于稳定的状态,气温从3月开始上升,1月—7月处于升温过程,气温上升快,7月、8月气温达到最高,9月后气温开始迅速降低,进入降温阶段,符合季节变化规律;此外,数据缺测明显。


图 4 唐山站气温整点值年变曲线 Fig. 4 Annual variation of integral point temperature of Tangshan site
2 研究方法

应用线性回归法解决连续几日甚至数月的整点气温值缺测问题,并采用交叉验证方法对插补结果进行误差分析。

2.1 线性回归模型的建立

本文对线性回归方法进行改进,考虑了参考站和缺测站之间的距离。改进后该方法更科学,可以更好地去除距离因素的影响。

利用邻近站资料对距离进行加权,建立回归模型,插补缺测站资料的方程式为:

$ {\hat y_i} = \frac{{{a_{1i}}{x_{1i}}}}{{{d_{1i}}}} + \frac{{{a_{2i}}{x_{2i}}}}{{{d_{2i}}}} + \cdots + \frac{{a{}_{mi}{x_{mi}}}}{{{d_{mi}}}} + {a_{m + 1}} $ (1)

其中,${\hat y_i}$为插补值,${x_{mi}}$为临近站数据,${a_{mi}}$为回归模型的参数,${d_{mi}}$为邻近站与缺测站之间的距离,m为临近站站数。利用最小二乘原理求解回归模型的参数,即使观测值和插补值之间差值的平方和(Q)最小:

$ \min Q = {\rm{min}}\sum\limits_{i = 1}^n {{{({y_i} - \hat y{}_i)}^2}} $ (2)

其中${y_i}$为插补站观测值。

因历史同期各要素时空变化规律通常比较相似,选择缺测整点值前后若干整点值的历史同期(不包括缺测值所在的年份)数据,作为拟合回归模型的样本数据,建立线性回归模型,并利用附近站资料,计算缺测记录插补值。

2.2 误差检验方法

本文采用交叉验证的方法对缺测记录的插补结果进行分析,即假设某个站的记录缺测,首先利用插补模型插补整点气温数据,然后对插补值与实际观测资料进行对比和误差分析,并用平均绝对误差(Mean Absolute Error,MAE)代表插补精度(王海军等,2008),其表达式为:

$ MAE = \frac{1}{N}\sum\limits_{i = 1}^N {\left| {{x_{io}} - {x_{ei}}} \right|} $ (3)

其中,${x_{oi}}$为第i个实际整点观测值,${x_{ei}}$为第i个插补整点值,N为插补整点值的个数。

3 插补结果分析
3.1 邻近参考站及时间窗选择

参考站的选择不仅与观测站密度有关,也与插补站及其邻近站所处的地理环境有关(如平原、丘陵、山区等)。同时,时间窗大小也对缺测数据的插补精度有影响。本文采用滑动优选法确定时间窗,时间窗的宽度为气温整点值个数,高度为年数。以选择缺测整点值所在的年份为中心,其前后若干点值历史同期若干年的数据作为样本数据,对于前后无资料的年份,则使用靠近插补年份的资料。

唐山站位于华北平原,周围地势平坦,气温变化相近,故采用最短距离的原则选取临近参考站。根据距离及地形因素,选取了该站周边的北京、昌黎、蓟县、宁河、青光和徐庄子6个观测站,年数为7年。采用15个整点值作为样本资料,建立线性回归模型,插补唐山站的缺测值(包括连续和不连续的单点缺测值)。插补站及其邻近站的基本信息见表 2

表 2 唐山插补站及其邻近台站信息 Table 2 Information of interpolation site and its neighboring sites
3.2 误差检验

采用交叉验证的方法评估上述插补方法,统计了2010年3月的31天实际观测整点值与相应插补值的相关系数,统计结果见表 3。从表中可以看出,3月8日、14日、15日、19日和20日的相关系数相对较低。王海军等(2008)经过对比研究,在平原地区选取了4个参考站,选取年数为8年、天数为15天,插补误差最小。唐山站也位于平原地区,考虑与唐山站的距离及地形因素,选取昌黎、宁河、青光、徐庄子4个参考站,并选取年数为7年、15个整点值的优化模型。通过对比发现,4个参考站的相关系数偏高(表 3)。

表 3 观测值与插补值相关系数 Table 3 Correlation coefficients between interpolated and observed data

利用式(3)对优化后的模型得到的插补结果计算平均绝对误差,并统计其误差的比例分布(表 4)。从表 4可以看出,插补误差在±0.5℃范围内的比例为60.5%,在±0.8℃范围内的比例为80.6%,其误差绝对值大于1.0℃的为9.4%,平均绝对误差为0.82℃。

表 4 唐山站整点气温缺测插补误差比例 Table 4 Proportion of error for interpolated integral point temperature
3.3 插补结果

针对唐山站2008年1月1日—2013年12月18日的气温整点值缺测数据及错误数据,利用唐山邻近站的同期数据和线性回归模型,对缺测数据插补完整,并修正错误数据,绘制气温整点值的年变曲线,如图 5所示。从图中可以看出,在长时间序列缺测的部位,插补值与前后正确的数据衔接吻合,没有出现突升或突降变化。插补后完整连续的数据符合夏高冬低的年变规律,气温6年的变化形态一致。


图 5 插补后的唐山台气温整点值年变曲线 Fig. 5 Annual variation of interpolated integral point temperature for Tangshan site
4 震例应用

据中国地震台网中心测定,2012年5月28日10时22分在河北省唐山市辖区、滦县交界处发生4.8级地震,震源深度8km。

利用插补完整的连续数据,分析2012年3—5月震前气温数据的变化情况。首先,选取无震年份(2008—2011年)同期(3—5月)整点值气温,并计算历年同期气温日均值,以此作为3—5月正常的背景值;其次,将2012年3—5月日均值与历年同期背景值做差值,得到2012年3—5月份数据与历年同期均值的偏移程度;最后,以无震年份同期气温波动范围作为基准,即2008—2011年3—5月所有气温值的标准差作为判断标准。根据以往经验将差值大于2倍标准差视为气温前兆异常。

从2012年3—5月当年日均值与历年同期(2008—2011年)多年日均值的差值及标准差(图 6)中可以看出,从3月27日开始,唐山站数据出现大幅度增温异常现象,4月增温天数也较多,5月1日—11日的气温日均值仍然高于历年同期,特别是5月10日(即震前2天)增幅达到约8℃,且大于2倍标准差;震后差值开始变小,并逐步恢复至平静。


图 6 观测值与同期均值的差值及标准差 Fig. 6 Standard deviation and difference between mean and observed value
5 结论

长期连续完整的历史气象资料是震前气温异常判别研究的重要数据基础,但由于观测环境、仪器故障等原因,造成气温观测数据缺测或错误数据,且部分数据缺测的时间较长。为此,本文利用线性回归模型,插补缺测和错误的气温整点值数据,较好地解决了长期连续缺测的情况。

通过对唐山观测站2008年1月1日—2013年12月18日的气温整点值缺测数据及错误数据进行插补,使得数据完整连续,并应用插补完整的气温整点值数据,分析研究了2012年5月28日唐山4.8级地震的气温前兆异常现象,主要得出以下结论:

(1)唐山观测站的插补值与其前后的观测数据衔接吻合,插补后完整连续的数据符合夏高冬低的年变规律。

(2)插补误差在±0.5℃范围内的比例为60.2%,在±0.8℃范围内的比例为80.3%,其误差绝对值大于1.0℃的为9.6%,平均绝对误差为0.84℃。插补值与观测值的相关系数大部分在0.9以上,可见插补结果真实可靠。

(3)从3月27日起,唐山观测站数据出现增温异常,震前2天增温幅度约8℃。

参考文献
曹新来, 张子广, 李福来, 1997. 1976年唐山7.8级地震前气象要素的异常变化[J]. 华北地震科学, 15(2): 50-59.
李贵福, 解明恩, 1996. 本世纪云南强震(MS ≥ 6.0)的气象特征研究[J]. 地震研究, 19(2): 154-161.
涂诗玉, 陈正洪, 2001. 武汉和宜昌缺测气温资料的插补方法[J]. 湖北气象, (3): 11-13.
王海军, 涂诗玉, 陈正洪, 2008. 日气温数据缺测的插补方法试验与误差分析[J]. 气象, 34(7): 83-91.
王建国, 姚会琴, 高逊, 等, 2010. 天津市地震前兆台网的运行监控与维护管理[J]. 大地测量与地球动力学, 30(S1): 111-115.
王建国, 刘春国, 王伟, 等, 2013. 地震前兆数据库综合管理系统[J]. 大地测量与地球动力学, 33(S1): 114-116.
徐国钧, 谭天明, 兰红林, 等, 1993. 云南破坏性地震与气象要素的关系[J]. 地震研究, 16(2): 148-155.
姚会琴, 李悦, 高逊, 等, 2012. NagVis等开源监控软件在天津地震前兆台网的应用研究[J]. 震灾防御技术, 7(3): 329-333. DOI:10.3969/j.issn.1673-5722.2012.03.012
余予, 李俊, 任芝花, 等, 2012. 标准序列法在日平均气温缺测数据插补中的应用[J]. 气象, 38(9): 1135-1139.
张铁宝, 路茜, 辛华, 等, 2013. 2010年玉树7.1级地震前后气温变化分析[J]. 四川地震, (4): 7-11. DOI:10.3969/j.issn.1001-8115.2013.04.002
张秀芝, 孙安健, 1996a. 利用车贝雪夫多项式进行资料缺测插补的研究[J]. 应用气象学报, 7(3): 344-352.
张秀芝, 孙安健, 1996b. 气候资料缺测插补方法的对比研究[J]. 气象学报, 54(5): 625-632.
张永领, 丁裕国, 高全洲, 等, 2006. 一种基于SVD的迭代方法及其用于气候资料场的插补试验[J]. 大气科学, 30(3): 526-532. DOI:10.3878/j.issn.1006-9895.2006.03.15
DeGaetano A. T., Eggleston K. L., Knapp W. W., 1995. A method to estimate missing daily maximum and minimum temperature observations[J]. Journal of Applied Meteorology, 34(2): 371-380. DOI:10.1175/1520-0450-34.2.371
Eischeid J. K., Pasteris P. A., Diaz H. F., et al, 2000. Creating a serially complete, national daily time series of temperature and precipitation for the western United States[J]. Journal of Applied Meteorology, 39(9): 1580-1591. DOI:10.1175/1520-0450(2000)039<1580:CASCND>2.0.CO;2
Huth R., Nemešová I., 1995. Estimation of missing daily temperatures:can a weather categorization improve its accuracy?[J]. Journal of Climate, 8(7): 1901-1916. DOI:10.1175/1520-0442(1995)008<1901:EOMDTC>2.0.CO;2