一种判别水文气候过程自然演变类型的方法

技术领域

本发明属于水文气候科学技术领域，具体指代一种判别水文气候过程自然演变类型的方法。

背景技术

准确揭示水文气候过程的演变特征、掌握其未来演变情势，是科学评估气候变化以及合理应对气候变化影响的基本依据和必要前提。受多种复杂因素(包括随机因素)的共同作用和影响，实际水文气候过程十分复杂，特别是近几十年来受全球变化影响的日益加剧，水文气候过程更为复杂多变。

当前水文气候过程的演变特征研究主要关注趋势和自然演变特征，其中自然演变特征主要关注短持续特性和长持续特性。当自相关函数C(s)随着时间间隔的增大快速衰减为0或以指数函数衰减时，该序列呈现短持续特性；反之，当自相关函数C(s)衰减的很慢或以幂函数形式衰减时，该序列呈现长持续特性。针对短持续特性的研究较早，研究者相继提出了AR模型、MA模型和ARMA模型等用以描述时间序列的短持续特性。Hurst对尼罗河的长期观测数据进行分析，发现序列的自相关函数呈现缓慢衰减状态，这种现象被后人称为“Hurst现象”。由于具有“Hurst现象”的序列取值在较长时间时仍存在联系，故又被称为长持续特性、长相关特性或长距离依赖特性，并以参数d表示序列的长持续特性大小。当d>0时，序列呈现长持续特性；当d<0时，序列呈现反持续特性；当d＝0时，序列呈现短持续特性或不存在相关性。特别地，当d＝1时，时间序列转化为单位根过程，其不仅产生明显的随机性趋势，时间序列之间还会出现“伪回归”现象。因此，判别和区分水文气候过程的不同自然演变类型，是准确揭示水文气候过程演变特征、科学评估气候变化和合理应对气候变化影响的重要前提。

目前，水文气候过程自然演变特征的研究主要基于长持续特性评估方法。长持续特性评估方法主要分为三类：非参数估计法、参数估计法和半参数估计法；非参数估计方法主要包括自相关系数分析法、R/S方法、重标方差法、去趋势波动分析(DetrendedFluctuation Analysis，DFA)等；参数估计方法主要有均值模型、波动率模型、周期图法、功率谱密度法和小波分析法；半参数方法主要包括GPH方法、局部Whittle(Local Whittle，LW)方法、以及基于GPH估计量和LW估计量的修正方法。然而，上述长持续特性评估方法均基于水文气候过程存在长持续特性的基本前提和假设，而忽略了该假设在实际中是否真正成立。已有研究表明，长持续特性评估方法易将短持续特性误检测为长持续特性，且随着数据样本量的增加，AR(1)过程的d值逐渐趋于0，此时可准确检测AR(1)过程并不属于长持续特性。但受限于观测数据长度，无法对长持续特性和短持续特性进行准确评估和区分，为此树轮等代用数据被用于长持续特性研究。研究表明，代用数据会低估长持续特性，而低分辨率代用数据则会高估长持续特性，即代用数据也无法用于区分水文气候过程的短持续特性和长持续特性。实际时间序列中往往存在噪声、趋势等成分，导致多数长持续特性估计方法的结果出现偏差或错误。相比上述长持续特性评估方法，DFA方法可以有效地滤去时间序列中的各阶趋势成分，适合于非平稳时间序列的长持续特性分析，其原理是消除时间序列中外部因素产生的趋势成分，然后对剩余成分进行长持续特性估计，因此被广泛地应用于水文气候过程长持续特性研究。

近年来，有学者提出利用DFA方法可对AR(1)过程和长持续过程进行准确区分。AR(1)过程在整个时间尺度上的d>0，在大时间尺度上d＝0，而长持续过程在整个时间尺度和大时间尺度上的长持续特性均为d>0。基于此，DFA方法可避免AR(1)过程被误判为长持续过程的结果。然而，该方法只适用于AR(1)过程，实际观测时间序列很难利用AR(1)过程进行准确描述，需要考虑更高阶的自相关特性。但对于更高阶的AR(p)(p>1)过程，DFA方法仍无法对短持续过程和长持续过程进行准确区分。通常自相关函数用于确定AR过程的阶数，但由于受确定成分的影响，自相关系数存在较大偏差。

发明内容

针对上述现有技术的不足，本发明的目的在于提供一种判别水文气候过程自然演变类型的方法，以解决现有技术中忽略长持续特性评估方法的适用性，导致无法准确判别和区分水文气候过程不同自然演变类型的缺陷。

为达到上述目的，本发明采用的技术方案如下：

本发明的一种判别水文气候过程自然演变类型的方法，步骤如下：

1)分别生成与待分析时间序列TS(t)长度相同的白噪声、AR(1)过程、AR(2)过程、单位根过程、长持续过程五种类型的时间序列，对生成的各时间序列进行差分处理后求解各自对应的一阶自相关系数和二阶自相关系数；

2)重复上述步骤1)，直至各类型的时间序列差分处理后的一阶自相关系数和二阶自相关系数的统计特征趋于稳定，进而获取各类型时间序列差分处理后的一阶自相关系数和二阶自相关系数对应的95％置信区间；

3)识别出时间序列TS(t)中的突变成分B

4)对新时间序列TS’(t)做差分处理后，求解其一阶自相关系数AC_diff(1)和二阶自相关系数AC_diff(2)；

5)将一阶自相关系数AC_diff(1)和二阶自相关系数AC_diff(2)与步骤2)中得到的白噪声、AR(1)过程、AR(2)过程、单位根过程、长持续过程的时间序列差分处理后的一阶自相关系数和二阶自相关系数对应的95％置信区间进行对比，来确定时间序列TS(t)的具体自然演变类型。

进一步地，所述的步骤1)具体包括：

11)利用蒙特卡罗方法生成白噪声的时间序列y

12)利用一阶自回归模型生成AR(1)过程的时间序列y

y

式中，t表示时序；ρ为一阶自相关系数，且|ρ|<1，u(t)是均值为0的符合独立同分布的白噪声序列；

13)利用二阶自回归模型生成AR(2)过程的时间序列y

y

式中，ρ

14)生成单位根过程的时间序列y

y

15)利用ARFIMA模型生成长持续过程的时间序列y

进一步地，所述步骤5)具体包括：

51)当AC_diff(1)和AC_diff(2)属于白噪声的95％置信区间内，则时间序列TS(t)的自然演变类型判定为白噪声过程；

52)当AC_diff(1)和AC_diff(2)属于单位根过程的95％置信区间内，则时间序列TS(t)的自然演变类型判定为单位根过程；

53)当AC_diff(1)和AC_diff(2)属于AR(2)过程的95％置信区间内，则时间序列TS(t)的自然演变类型判定为AR(2)过程；

54)若AC_diff(1)和AC_diff(2)属于长持续过程或AR(1)过程的95％置信区间内，则利用DFA方法求解新时间序列TS’(t)的标度指数α，进一步判别TS(t)的自然演变类型。

进一步地，所述步骤54)具体包括：

541)利用DFA方法获取新时间序列TS’(t)的波动函数F(s)和时间尺度s的双对数散点图(ln(F(s))、ln(s))；

542)识别双对数散点图的结构突变点B

543)利用最小二乘法对区间B

544)若α＝0.5，则时间序列TS(t)的自然演变类型判定为AR(1)过程；

545)若α>0.5，则时间序列TS(t)的自然演变类型判定为长持续过程。

本发明的有益效果：

本发明利用蒙特卡罗试验确定各种自然演变类型统计特征的置信区间，以此为依据准确区分白噪声、单位根过程、AR(2)过程，可避免AR(2)过程被误判为长持续过程的错误结果；

本发明利用DFA方法进一步区分AR(1)过程和长持续过程，避免直接采用LocalWhittle方法而将AR(1)过程误判为长持续过程的错误操作，从而准确判别和区分水文气候过程自然演变的五种主要类型，较常规方法判别结果更为合理可靠，可为科学评估气候变化以及合理应对气候变化影响提供科学依据。

附图说明

图1为本发明方法的流程图；

图2a为AR(2)过程的人工生成序列S11的示意图；

图2b为AR(2)过程的人工生成序列S12的示意图；

图2c为AR(2)过程的人工生成序列S13的示意图；

图2d为AR(2)过程的人工生成序列S14的示意图；

图3a为AR(1)过程的人工生成序列S21的示意图；

图3b为AR(1)过程的人工生成序列S22的示意图；

图3c为长持续过程的人工生成序列S23的示意图；

图3d为长持续过程的人工生成序列S24的示意图。

具体实施方式

为了便于本领域技术人员的理解，下面结合实施例与附图对本发明作进一步的说明，实施方式提及的内容并非对本发明的限定。

参照图1所示，本发明的一种判别水文气候过程自然演变类型的方法，步骤如下：

1)随机设置各参数的取值，分别生成与待分析时间序列TS(t)长度相同的白噪声、AR(1)过程、AR(2)过程、单位根过程、长持续过程五种类型的时间序列，对生成的各时间序列进行差分处理来消除随机性趋势的影响，然后求解差分处理后的各时间序列的一阶自相关系数和二阶自相关系数；

其中，所述的步骤1)具体包括：

11)利用蒙特卡罗方法生成白噪声的时间序列y

12)利用一阶自回归模型生成AR(1)过程的时间序列y

y

式中，t表示时序；ρ为一阶自相关系数，且|ρ|<1，u(t)是均值为0的符合独立同分布的白噪声序列；

13)利用二阶自回归模型生成AR(2)过程的时间序列y

y

式中，ρ

14)生成单位根过程的时间序列y

y

15)利用ARFIMA模型生成长持续过程的时间序列y

2)设置各参数的不同取值，重复上述步骤1)生成各类型的不同时间序列，进行蒙特卡罗试验(Monte-Carlo试验)，直至各类型的时间序列差分处理后的一阶自相关系数和二阶自相关系数的统计特征(均值和标准差)趋于稳定；

3)获取各类型时间序列差分处理后的一阶自相关系数和二阶自相关系数的均值±2倍标准差的区间，作为其对应的95％置信区间，实现对五种自然演变类型(白噪声、AR(1)过程、AR(2)过程、单位根过程、长持续过程)之间差异的准确量化和区分；

4)采用Mann-Kendall检验方法识别出时间序列TS(t)中的突变成分B

5)对新时间序列TS’(t)做差分处理后，求解其一阶自相关系数AC_diff(1)和二阶自相关系数AC_diff(2)；

6)将AC_diff(1)和AC_diff(2)与步骤3)中得到的白噪声、AR(1)过程、AR(2)过程、单位根过程、长持续过程的时间序列差分处理后的一阶自相关系数和二阶自相关系数对应的95％置信区间进行对比，来确定时间序列TS(t)的具体自然演变类型；具体包括：

61)当AC_diff(1)和AC_diff(2)属于白噪声的95％置信区间内，则时间序列TS(t)的自然演变类型判定为白噪声过程；

62)当AC_diff(1)和AC_diff(2)属于单位根过程的95％置信区间内，则时间序列TS(t)的自然演变类型判定为单位根过程；

63)当AC_diff(1)和AC_diff(2)属于AR(2)过程的95％置信区间内，则时间序列TS(t)的自然演变类型判定为AR(2)过程；

64)当AC_diff(1)和AC_diff(2)属于长持续过程或AR(1)过程的95％置信区间内，则利用DFA方法求解新时间序列TS’(t)的标度指数α，进一步判别TS(t)的自然演变类型。

优选地，所述步骤64)具体包括：

641)利用DFA方法获取新时间序列TS’(t)的波动函数F(s)和时间尺度s的双对数散点图(ln(F(s))、ln(s))；

642)识别双对数散点图的结构突变点B

643)利用最小二乘法对区间B

644)若α＝0.5，则时间序列TS(t)的自然演变类型判定为AR(1)过程；

645)若α>0.5，则时间序列TS(t)的自然演变类型判定为长持续过程。

示例：

由于人工生成序列的自然演变类型等情况已知，利用人工生成序列，有利于检验本发明方法的有效性，而实测水文时间序列的自然演变类型往往未知，无法准确判断本发明方法判别水文气候过程自然演变类型结果的准确性。为证明该发明方法判别时间序列自然演变类型结果的准确性，设计方案时生成两类人工序列，分别用于验证本发明方法对区分AR(2)过程和长持续过程的有效性，以及区分AR(1)过程和长持续过程的有效性。第一类时间序列的自然演变特征为AR(2)过程，序列长度相同，但时间序列的自回归系数不同，分别记为S11、S12、S13和S14(图2a，图2b，图2c，图2d)。第二类序列具有相同的序列长度，但自然演变类型不同，其中S21、S22为AR(1)过程(图3a，图3b)，S23和S24为长持续过程(图3c，图3d)。选用不同方法对上述人工生成的时间序列的自然演变类型分别进行判定，两类时间序列的判别结果分别见表1和表2(不同方法对人工生成序列的自然演变类型判别结果)：

表1

表2

自然演变类型判别结果显示：DFA方法无法准确识别AR(2)过程，导致对时间序列自然演变类型的误判，而本发明可以准确区分AR(2)过程和长持续过程。对于AR(1)过程和长持续过程，Local Whittle方法无法对其进行准确区分。相比于直接评估时间序列长持续特性的Local Whittle方法，本发明方法首先对时间序列进行差分处理，利用差分时间序列一阶和二阶自相关系数识别AR(2)过程，可以克服DFA方法的局限性，避免了AR(2)过程被误判为长持续过程；其次，利用DFA方法区分AR(1)过程和长持续过程，进一步克服由于忽略Local Whittle方法的局限性导致对时间序列自然演变类型结果的误判，最终得到准确的时间序列自然演变类型结果。

对比上述时间序列自然演变类型判别结果，可以得到以下几点重要结论：

(1)DFA方法可以识别AR(1)过程和长持续过程，但无法准确区分AR(2)过程和长持续过程；

(2)相比DFA方法，Local Whittle方法是基于时间序列存在长持续特性的基本假设，会将AR(1)过程误判为长持续过程；

(3)本发明方法首先识别AR(2)过程，避免了AR(2)过程被误判为长持续过程的错误，然后利用DFA方法进一步区分AR(1)过程和长持续过程，克服了其他常规方法的局限性，因此对时间序列自然演变类型的判定结果更加合理可靠，可为揭示水文气候过程演变特征和科学评估气候变化等提供科学依据。

本发明具体应用途径很多，以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以作出若干改进，这些改进也应视为本发明的保护范围。