(19)国家知识产权局 (12)发明 专利申请 (10)申请公布号 (43)申请公布日 (21)申请 号 202210650692.9 (22)申请日 2022.06.09 (71)申请人 中国科学院地理科 学与资源研究所 地址 100101 北京市朝阳区大屯路甲1 1号 (72)发明人 桑燕芳　李鑫鑫　任智慧　 (74)专利代理 机构 江苏圣典律师事务所 32 237 专利代理师 贺翔 (51)Int.Cl. G06F 30/20(2020.01) G06F 111/08(2020.01) (54)发明名称 一种实现水文气候时间序列趋势无偏估计 的方法 (57)摘要 本发明公开了一种实现水文气候时间序列 趋势无偏估计的方法， 包括： 利用蒙特卡罗试验 获取五种自然演变类型的时间序列差分后各自 对应的一阶和二阶自相关系数的95%置信区间； 剔除待分析时间序列的突变成分和季节成分， 对 剩余成分进行差分后求解一阶和二阶自相关系 数， 与上述各种类型对应的95%置信区间对比， 确 定其自然演变类型： 若为白噪声， 利用广义最小 二乘法直接评估趋势显著性； 若为AR过程， 利用 Newey‑West方法消除异方差和自相 关性后评估 趋势显著性； 若为长持续过程， 拟合并去除长持 续成分后再评估趋势显著性。 本发 明针对不同自 然演变类型采取针对性的方法检测趋势类型并 评估显著性， 消除自然演变特征的影响， 实现对 水文气候时间序列趋势的无偏估计 。 权利要求书2页 说明书7页 附图5页 CN 114970187 A 2022.08.30 CN 114970187 A 1.一种实现水文气候时间序列趋势无偏估计的方法， 其特 征在于， 步骤如下： 1)分别生成与待分析 时间序列S(t)长度相同的白噪声、 AR(1)过程、 AR(2)过程、 单位根 过程、 长持续过程五种类型 的时间序列， 对生成的各时间序列进行差分处理后求解各自对 应的一阶自相关系数和二阶自相关系数； 2)重复上述步骤1)， 直至各类型的时间序列差分处理后的一阶自相关系数和二阶自相 关系数的统计特征趋于稳定， 进而获取各类型的时间序列差分处理后的一阶自相关系数和 二阶自相关系数对应的95％置信区间； 3)识别出时间序列S(t)中的突变成分B0， 求解多年平均的季节成分S0， 剔除时间序列S (t)的突变成分B0和季节成分S0， 将剩余成分作为 新时间序列S ’(t)＝S(t) ‑B0‑S0； 4)对新时间序列S ’(t)做差分处理后， 求解其一阶自相关系数AC_diff(1)和二阶自相 关系数AC_dif f(2)； 5)将一阶自相关系数AC_diff(1)和二阶自相关系数AC_diff(2)与步骤2)中得到的各 类型的时间序列差分处理后的一阶自相关系 数和二阶自相关系 数对应的95％置信区间进 行对比， 来判定时间序列S(t)的具体自然演变 类型； 6)当时间序列S(t)的自然演变类型为白噪声过程时， 利用广义最小二乘方法估计时间 序列S(t)的线性趋势斜率β1， 并利用t检验方法评估其显著性： 若通过显著性检验， 则判定 时间序列S(t)存在确定性趋势， 趋势斜 率的无偏估计结果β ＝β1； 7)当时间序列S(t)的自然演变类型为AR(1)或AR(2)过程时， 利用Newey ‑West方法处理 该时间序列的异方差和自相关性后求得趋势斜率为β2； 若通过显著性检验， 则判定时间序 列S(t)存在确定性趋势， 趋势斜 率的无偏估计结果β ＝β2； 8)当时间序列S(t)的自然演变类型为单位根过程时， 则时间序列S(t)呈现明显的随机 性趋势， 不存在确定性线性趋势， β ＝0； 9)当时间序列 S(t)的自然演变类型为长持续过程时， 依次取不同的长持续特性di值， 再 利用长持续特性模型拟合新时间序列S ’(t)的长持续成分Mi(t)， 得到剩余成分Si(t)＝S’ (t)‑Mi(t)； 分别估计剩余成分Si(t)的趋势斜率βi与对应的残差Ri(t)； 当残差最小时得到 最终的趋势斜 率的无偏估计结果β ＝β3； 10)得到时间序列S(t)的趋势斜率无偏估计结果β， 利用广义最小二乘方法估计时间序 列S(t)的线性趋势斜率β0； 对比β0与β 的差别， 得到时间序列S(t)的自然演变特征对其趋势 无偏估计结果的影响程度 λ＝|β0‑β |。 2.根据权利要求1所述的实现水文气候 时间序列趋势无偏估计的方法， 其特征在于， 所 述的步骤1)具体包括： 11)利用蒙特卡罗方法生成白噪声的时间序列y1(t)； 12)利用一阶自回归 模型生成AR(1)过程的时间序列y2(t)如下： y2(t)＝ρ×y2(t‑1)+u(t) 式中， t表示时序； ρ 为一阶自相关系数， 且|ρ |<1， u(t)是均值为0 的符合独立同分布的 白噪声序列； 13)利用二阶自回归 模型生成AR(2)过程的时间序列y3(t)如下： y3(t)＝ρ1×y3(t‑1)+ρ2×y3(t‑2)+u(t) 式中， ρ1和ρ2分别为一阶和二阶自相关系数， ρ1+ρ2<1， ρ2‑ρ1<1，‑1<ρ2<1；权　利　要　求　书 1/2 页 2 CN 114970187 A 214)生成单位 根过程的时间序列y4(t)如下： y4(t)＝y4(t‑1)+u(t)； 15)利用ARFIMA模型生成长持续过程的时间序列y5(t)。 3.根据权利要求1所述的实现水文气候 时间序列趋势无偏估计的方法， 其特征在于， 所 述步骤5)具体包括： 51)当AC_diff(1)和AC_diff(2)属于白噪声的95％置信区间内， 则时间序列S(t)的自 然演变类型判定为白噪声过程； 52)当AC_diff(1)和AC_diff(2)属于单位根过程的95％置信区间内， 则时间序列S(t) 的自然演变 类型判定为单位 根过程； 53)当AC_diff(1)和AC_diff(2)属于AR(2)过程的95％置信区间内， 则时间序列S(t)的 自然演变 类型判定为AR(2)过程； 54)当AC_diff(1)和AC_diff(2)属于长持续过程或AR(1)过程的95％置信区间内， 利用 DFA方法求 解新时间序列S ’(t)的标度指数α； 进一 步判别S(t)的自然演变 类型。 4.根据权利要求3所述的实现水文气候 时间序列趋势无偏估计的方法， 其特征在于， 所 述步骤54)具体包括： 541)利用DFA方法获取新时间序列S ’(t)的波动函数F(s)和时间尺度s的双对数散点图 (ln(F(s))、 ln(s))； 542)识别双对数散点图的结构突变点B1； 543)利用最小二乘法对区间B1<s<L/4的ln(F(s))和ln(s)进行线性拟合， 线性趋势为标 度指数α， L 为时间序列S(t)的序列长度； 544)若α ＝0.5， 则时间序列S(t)的自然演变 类型判定为AR(1)过程； 545)若α >0.5， 则时间序列S(t)的自然演变 类型判定为长持续过程。 5.根据权利要求1所述的实现水文气候 时间序列趋势无偏估计的方法， 其特征在于， 所 述步骤9)具体包括： 91)设置不同的长持续特性di值， 即di＝[‑0.5:0.01:0.5]； 92)依次取不同的di值， 利用长持续特性模型拟合新时间序列S ’(t)的长持续成分Mi (t)： Mi(t)＝u(t)/(1 ‑N)di 式中， N为自回归 模型AR(n)的滞后算子， n 为模型阶数； 93)得到剩余成分Si(t)＝S’(t)‑Mi(t)； 94)利用广义最小二乘方法估计剩余成分Si(t)的线性趋势斜率βi； 若通过显著性检验， 则βi＝βi； 若未通过显著性检验， 则βi＝0； 95)得到趋势成分Ti(t)＝ci+βi×t， ci为趋势成分的截距项； 96)新时间序列S ’(t)减去长持续成分Mi(t)和趋势成分Ti(t)后得到残差Ri(t)： Ri(t)＝S’(t)‑Mi(t)‑Ti(t)； 97)当残差Ri(t)取最小 值时， 对应的斜率值为时间序列S(t)趋势斜率的无偏 估计结果β ＝β3。权　利　要　求　书 2/2 页 3 CN 114970187 A 3