数学建模数据分析——趋势性检验和平稳性检验

在数学建模比赛中，经常需要对数据进行分析和预处理，常见的比如趋势分析（上升/下降/无明显趋势）和突变分析，很多时候靠人的经验观察得出结论，但这是不够严谨的。于是我们通常会采用一些更科学的方法，下面我们就来详细的捋一遍数据分析检验方法：

文章目录

数学建模数据分析——趋势性检验和平稳性检验
- 时间序列趋势性检验方法
- - 斜率法
  - Cox-Stuart检验法
  - Mann-Kendall检验法
- 时间序列平稳性检验方法
- - 平稳时间序列定义与性质
  - 检验方法
  - - 时序图检验法
    - 自相关图检验法
    - ADF检验法
    - KPSS检验
- 非常有用的相关链接

时间序列趋势性检验方法

斜率法

原理

斜率法就是使用最小二乘法对时序数据进行拟合，根据拟合的直线的斜率K来判断序列的数据走势，当K>0时，则代表上升趋势；当k<0时，则代表下降趋势；

优缺点

优点为使用简单；缺点是要求数据是线性的，当数据波动较大时无法准确拟合，导致丧失了评价意义。

Cox-Stuart检验法

原理

直接考虑数据的变化趋势，若数据呈现上升趋势，则排在后面的数据的值要比排在前面的数据的值显著的大；反之，若数据呈现下降趋势，则排在后面的数据的值要比排在前面的数据的值显著的小。该方法利用前后两个时期不同数据的差值的正负来判断数据的总的变化趋势。

算法步骤

取xix_ixi和xi+cx_i+cxi+c组成一对(xi,xi+c)(xi,xi+c)(xi,xi+c)。这里如果nnn为偶数，则c=n/2c=n/2c=n/2，如果nnn是奇数，则c=(n+1)/2c=(n+1)/2c=(n+1)/2。当n为偶数时，共有n’=cn’=cn’=c对，而nnn是奇数时，共有n’=c−1n’=c-1n’=c−1对。
用每一对的两元素差Di=xi−xi+cD_i=x_i−x_{i+c}Di=xi−xi+c的符号来衡量增减。令S+S_+S+为正的DiD_iDi的数目，S−S_−S−为负的DiD_iDi的数目。显然当正号太对时有下降趋势，反之有增长趋势。在没有趋势的零假设下他们因服从二项分布B(n’,0.5)B(n’,0.5)B(n’,0.5)。
用p(+)p(+)p(+)表示取到正数的概率，用p(−)p(-)p(−)表示取到负数的概率，这样就得到符号检验方法来检验序列是否存在趋势性。

优缺点

优点是不依赖趋势结构，可以快速判断趋势是否存在；缺点是未考虑数据的时序性，仅从符号检验来判断

实现代码

假设有下面一个时序数据，波动性较大，无法根据经验判断数据呈上升还是下降趋势。

import matplotlib.pyplot as plta_array = [206,223,235,264,229,217,188,204,182,230,223,227,242,238,207,208,216,233,233,274,234,227,221,214,226,228,235,237,243,240,231,210]plt.plot(range(1, len(a_array) + 1), a_array)

统计正的DiD_iDi出现次数为14次，负的DiD_iDi出现次数为2次，接下来判断是否存在显著的上升趋势

做出原假设H0H_0H0：该数据不存在上升趋势；

做出替代假设H1H_1H1：该数据存在上升趋势。

运行CS趋势检验方法，可以发现，使用二项式离散分布检验得到p-value值为0.004181，小于0.05，因此有理由认为原假设H0H_0H0不成立，从而支持假设H1H_1H1：该数据存在上升趋势

import numpy as np
import pandas as pd
import scipy.stats as statsdef cox_staut(list_c):lst=list_c.copy() ## 副本n0 = len(lst) ##列表长度## 判断奇偶性if n0%2 == 1:del lst[int((n0-1)/2)] #奇数删除中间数值## 配对前后数据c=int(len(lst)/2)## 统计正Di和负Di出现次数pos=neg=0for i in range(c):diff=lst[i+c]-lst[i]if diff>0:pos+=1elif diff<0:neg+=1else:continuen1=pos+negk=min(pos,neg)## 使用二项式离散分布检验趋势是否显著p=2*stats.binom.cdf(k,n1,0.5)print('fall:%i  rise:%i  p-value:%f'%(neg, pos, p))cox_staut(a_array)

fall:2  rise:14  p-value:0.004181

Mann-Kendall检验法

原理
Mann-Kendall检验不需要样本遵循一定的分布，也不受少数异常值的干扰。在Mann-Kendall检验中，原假设H0H_0H0为时间序列数据(X1，…,Xn)(X1，…,Xn)(X1，…,Xn),是nnn个独立的、随机变量同分布的样本；备择假设H1H_1H1 是双边检验，对于所有的k，j≤nk，j≤nk，j≤n，且k≠jk≠jk=j，XkX_kXk和XjX_jXj的分布是不相同的。若原假设是不可接受的，即在α置信水平上，时间序列数据存在明显的上升或下降趋势。对于统计量ZZZ，大于0时是上升趋势；小于0时是下降趋势。

算法步骤

将数据按采集时间列出：x1,x2,…,xnx_1,x_2,…,x_nx1,x2,…,xn，即分别在时间1,2,…,n1,2,…,n1,2,…,n得到的数据。
确定所有n(n−1)/2n(n-1)/2n(n−1)/2个xj−xkxj−xkxj−xk差值的符号，其中j>kj > kj>k，，这些差值是：x2−x1，x3−x1,…,xn−x1，x3−x2，x4−x2，…，xn−xn−2，xn−xn−1x_2 - x_1，x_3 - x_1, … , x_n - x_1，x_3 - x_2，x_4 - x_2，…，x_n - x_{n-2}，x_n - x_{n-1}x2−x1，x3−x1,…,xn−x1，x3−x2，x4−x2，…，xn−xn−2，xn−xn−1
令sgn(xj−xk)sgn(x_j−x_k)sgn(xj−xk)作为指示函数，依据xj−xkx_j−x_kxj−xk的正负号取值为1,0或-1，即
sgn(xj−xk)={1,xj−xk>0−1,xj−xk<00.xj−xk=0sgn(x_j-x_k)= \begin{cases} 1, x_j − x_k > 0\\ -1, x_j − x_k < 0\\ 0. x_j − x_k = 0 \end{cases} sgn(xj−xk)=⎩⎪⎨⎪⎧1,xj−xk>0−1,xj−xk<00.xj−xk=0
计算S=∑k−1n−1∑j−k+1nsgn(xj−xk)S=∑^{n−1}_{k−1}∑^n_{j−k+1}sgn(x_j−x_k)S=∑k−1n−1∑j−k+1nsgn(xj−xk)。即差值为正的数量减去差值为负的数量。如果SSS是一个正数，那么后一部分的观测值相比之前的观测值会趋向于变大；如果S是一个负数，那么后一部分的观测值相比之前的观测值会趋向于变小。
如果n≤10n \leq 10n≤10，依据Gilbert (1987, page 209, Section 16.4.1)中所描述的程序，接下来要在概率表 (Gilbert 1987, Table A18, page 272) 中查找SSS。如果此概率小于α\alphaα(认为没有趋势时的截止概率)，那就拒绝零假设，认为趋势存在。如果在概率表中找不到n(存在结数据——tied data values——会发生此情况)，就用表中远离0的下一个值。比如S=12S=12S=12，如果概率表中没有S=12S=12S=12，那么就用S=13S=13S=13来处理也是一样的。
如果n>10n > 10n>10，则依以下步骤6-10来判断有无趋势。这里遵循的是Gilbert (1987, page 211, Section 16.4.2)中的程序。
计算S的方差如下：VAR(S)=118[n(n−1)(2n+5)−∑p−1gtp(tp−1)(2tp+5)]VAR(S)=\frac{1}{18}[n(n−1)(2n+5)−∑^g_{p−1}t_p(t_p−1)(2t_p+5)]VAR(S)=181[n(n−1)(2n+5)−∑p−1gtp(tp−1)(2tp+5)]。其中ggg是结组（tied groups）的数量，tpt_ptp是第p组的观测值的数量。例如：在观测值的时间序列23,24,29,6,29,24,24,29,23{23, 24, 29, 6, 29, 24, 24, 29, 23}23,24,29,6,29,24,24,29,23中有g=3g = 3g=3个结组，相应地，对于结值(tiied value)23有t1=2t_1=2t1=2、结值24有t2=3t_2=3t2=3、结值29有t3=3t_3=3t3=3。当因为有相等值或未检测到而出现结时，VAR(S)VAR(S)VAR(S)可以通过Helsel (2005, p. 191)中的结修正方法来调整。
计算MK检验统计量Z_{MK}:
ZMK={S−1VAR(S),S>00,S=0S+1VAR(S),S>0Z_{MK}=\begin{cases} \frac{S-1}{\sqrt{VAR(S)}},& S>0\\ \qquad0\quad,& S=0\\ \frac{S+1}{\sqrt{VAR(S)}},& S>0 \end{cases} ZMK=⎩⎪⎪⎨⎪⎪⎧VAR(S)S−1,0,VAR(S)S+1,S>0S=0S>0
设想要测试零假设。H0H_0H0(没有单调趋势)对比替代假设Ha(有单调增趋势)，其1型错误率为α，0<α<0.50（注意α是MK检验错误地拒绝了零假设时可容忍的概率——即MK检验拒绝了零假设是错误地，但这个事情发生概率是α，我们可以容忍这个错误）。如果ZMK≥Z1−αZ_{MK}≥Z_1−αZMK≥Z1−α，就拒绝零假设H0，接受替代假设HaH_aHa，其中Z1−αZ_1−αZ1−α是标准正态分布的100(1−α)th百分位。
测试上面的H0与Ha（有单调递减趋势），其1型错误率为alpha，0<α<0.5，如果ZZMK≤–Z1−α，就拒绝零假设H0，接受替代假设Ha
测试上面的H0与Ha（有单调递增或递减趋势），其1型错误率为alpha，0<α<0.5，如果|ZMK|≥Z1−α2，就拒绝零假设H0，接受替代假设Ha，其中竖线代表绝对值。

实现代码

笔者自己找了一个时序数据，长度为2000，由于数据较长，就不再展示了。由时序图可以看出该数据波动性较大，无法根据经验判断数据呈上升还是下降趋势。

plt.plot(range(1, len(b_array) + 1), b_array)

根据上述提供的思路，做出假设，并判断数据的变化趋势。

运行MK趋势检验方法，可以发现，P-value值为4.10e-12，小于给定的显著性水平，因此有理由认为该时序数据呈现上升趋势。

import math
from scipy.stats import norm, mstats
import numpy as npdef mk_test(x, alpha=0.05):n = len(x)# calculate Ss = 0for k in range(n-1):for j in range(k+1, n):s += np.sign(x[j] - x[k])# calculate the unique dataunique_x, tp = np.unique(x, return_counts=True)g = len(unique_x)# calculate the var(s)if n == g:  # there is no tievar_s = (n*(n-1)*(2*n+5))/18else:  # there are some ties in datavar_s = (n*(n-1)*(2*n+5) - np.sum(tp*(tp-1)*(2*tp+5)))/18if s > 0:z = (s - 1)/np.sqrt(var_s)elif s < 0:z = (s + 1)/np.sqrt(var_s)else: # s == 0:z = 0# calculate the p_valuep = 2*(1-norm.cdf(abs(z)))  # two tail testh = abs(z) > norm.ppf(1-alpha/2)if (z < 0) and h:trend = 'decreasing'elif (z > 0) and h:trend = 'increasing'else:trend = 'no trend'return trend, p, hmk_test(b_array)

('increasing', 4.100275674545628e-12, True)

时间序列平稳性检验方法

平稳时间序列定义与性质

平稳时间序列按照限定条件的严格程度可以分为以下两种类型：

严平稳时间序列：严平稳顾名思义，是一种条件非常苛刻的平稳性，它要求序列随着时间的推移，其统计性质保持不变。

设{Xt}\{X_t\}{Xt}为一时间序列，对任意的正整数mmm，任取t1,t2,⋯,tm∈Tt_1, t_2, \cdots, t_m ∈ Tt1,t2,⋯,tm∈T，对任意整数τ\tauτ，其联合概率密度函数满足：

Ft1,t2,⋯,tm(x1,x2,⋯,xm)=Ft1+τ,t2+τ,⋯,tm+τ(x1,x2,⋯,xm)F_{t_1, t_2,\cdots,t_m} (x_1, x_2, \cdots, x_m)= F_{t_{1+\tau}, t_{2+\tau}, \cdots, t_{m+\tau}}(x_1, x_2, \cdots, x_m) Ft1,t2,⋯,tm(x1,x2,⋯,xm)=Ft1+τ,t2+τ,⋯,tm+τ(x1,x2,⋯,xm)
则称时间序列{Xt}\{X_t\}{Xt}为严平稳时间序列。严平稳的条件只是理论上的存在，现实中用得比较多的是宽平稳的条件。

宽平稳时间序列：宽平稳也叫弱平稳或者二阶平稳（均值和方差平稳），满足：

均值为常数，即
E(Xt)=u,∀t∈TE(X_t)=u,\quad \forall t \in T E(Xt)=u,∀t∈T
方差为常数
D(Xt)=γ(t,t)=γ(0),∀t∈TD(X_t)=\gamma(t,t)=\gamma(0), \quad \forall t \in T D(Xt)=γ(t,t)=γ(0),∀t∈T
自协方差为常数，即自协方差函数和自相关系数只依赖于时间的平移长度，而与时间的起点无关。
γ(t,s)=γ(k,k+s−t),∀t,s,k∈T\gamma(t,s)=\gamma(k, k+s-t),\quad \forall t,s,k \in T γ(t,s)=γ(k,k+s−t),∀t,s,k∈T
因此，可以记γ(k)\gamma(k)γ(k)为时间序列{Xt}\{X_t\}{Xt}的延迟kkk自协方差函数。

在现实生活中，时间序列是很难满足严平稳时间序列的要求的，因此，一般所讲的平稳时间序列在默认情况下都是指宽平稳时间序列。

由于平稳时间序列具有这些优良性质，因此，对于一个平稳时间序列来说，其待估计的参数量就变得少了很多，因为他们的均值、方差都是一样的，因此，可以利用全部的样本来估计总体的均值和方差，即：
μ^=xˉ=∑i=1nxnγ^(0)=∑t−1n(xt−xˉ)2n−1\hat \mu = \bar x = \frac{\sum_{i=1}^nx}{n}\\ \hat \gamma(0) = \frac{\sum_{t-1}^{n}(x_t - \bar x)^2}{n-1} μ^=xˉ=n∑i=1nxγ^(0)=n−1∑t−1n(xt−xˉ)2
所以当拿到一个时间序列后，需要对其进行平稳性检验。

检验方法

平稳性检验的主要方法是看时序图、ACF图和单位根检验，其中单位根检验方法有ADF、PP、KPSS等。其中ADF和PP检验都拒绝原假设，即认为序列平稳；KPSS拒绝原假设，即认为序列非平稳（KPSS零假设和ADF/PP恰好相反）。当出现不一致时，需要根据属于下面哪种平稳过程来判断，有一种观点认为KPSS检验结果更强更鲁棒，因为ADF和PP检验将差分平稳模型作为零假设，它检验随机游走或带漂移的随机游走效果奇好，但它们都需要假设是否包含常数均值和时间趋势，因此拒绝零假设的功效（low power）较低，而KPSS则完全不需要选择假设类型。

时序图检验法

那么，当拿到一个时间序列后，应该如何对其进行平稳性的检验呢？目前，对时间序列的平稳性检验主要有两种方法，一种是图检法，即根据时序图和自相关图进行直观判断，另一种是构造检验统计量的方法，目前主要有单位根检验法。
对于图检法，我们一般绘制时间序列的时序图，如下图所示，如果时间序列是平稳的，那么序列应该是围绕某一个均值上下随机波动，而下图中的序列明显具有一定的增长趋势，因此，可以断定该序列肯定不是平稳时间序列。

自相关图检验法

对于平稳时间序列，其自相关图一般随着阶数的递增，自相关系统会迅速衰减至0附近，而非平稳时间序列则可能存在先减后增或者周期性波动等变动。如下图所示，该时间序列随着阶数的递增，自相关系数始终在[0.8,1][0.8,1][0.8,1]之间，衰减微弱，因此，可以判断该时间序列不是平稳时间序列。

from statsmodels.graphics.tsaplots import plot_acf, plot_pacfplot_acf(b_array)
plot_pacf(b_array)

除了上述两种方法，便是单位根检验，下面这张图对单位根检验做了很好地说明。

ADF检验法

在使用很多时间序列模型的时候，如 ARMA、ARIMA，都会要求时间序列是平稳的，所以一般在研究一段时间序列的时候，第一步都需要进行平稳性检验，除了用肉眼检测的方法，另外比较常用的严格的统计检验方法就是ADF检验，也叫做单位根检验。

ADF检验全称是 Augmented Dickey-Fuller test，顾名思义，ADF是 Dickey-Fuller检验的增广形式。DF检验只能应用于一阶情况，当序列存在高阶的滞后相关时，可以使用ADF检验，所以说ADF是对DF检验的扩展。
原理

ADF检验就是判断序列是否存在单位根：如果序列平稳，就不存在单位根；否则，就会存在单位根。

所以，ADF检验做出的假设为：

**原假设H0H_0H0

数学建模时序数据分析——趋势性检验和平稳性检验相关推荐

数学建模与数据分析 || 1. 数学建模简介
数学建模简介文章目录数学建模简介 1. 数学建模比赛的理解 2. 一般数据分析的流程 3. 机器学习与统计数据分析 4. 各种编程软件仅仅是工具,对问题的观察视角和解决问题的策略才是关键 4.1 ...
2020暑期数学建模（数据分析）学习笔记
总算忙完所有课程论文,购买了视频课程. 第1讲层次分析法综合评价课已学,B站视频也看了,略过. 第2讲 TOPSIS法(优劣解距离法) 综合评价课已学,B站视频也看了,略过. 第3讲插值算法针 ...
数学建模之Shapiro-wilk夏皮洛-威尔克检验
Shapiro-wilk夏皮洛-威尔克检验作用:检验数据是否满足正态分布一.说明直接利用SPSS即可,得到的结果的最后一项为P值,一般样本数在 3 到 50 之间宜使用二.操作步骤就到这里啦 ...
数学建模与数据分析中的时间序列分析
时间序列分析将某种现象的指标数值按照时间顺序排列而成的数值序列. 文章目录 (1) 时间序列的基本概念 1. 组成要素 2. 时间序列的分类 3. 时间序列的分解 ① 长期变动趋势 T ② 季节趋势 ...
数学建模与数据分析 || 3. 面向数据的特征提取方法: 探索性数据分析
面向数据的特征提取方法: 探索性数据分析文章目录面向数据的特征提取方法: 探索性数据分析 1. 原始数据的准备 1.1 导入 python 模块 1.2 导入数据集并进行宏观认识 1.3 数据集描 ...
数学建模（数据分析C题）-建模思路
前言: 参考E038 的2019薄利多销优秀论文的模型建立与求解 http://dxs.moe.gov.cn/zx/a/hd_sxjm_sxjmlw_2019qgdxssxjmjslwzs/19102 ...
数学建模及数据分析上的插值处理——第三部分实践插值实战
2.二维网格节点插值例已知平面区域0<=x<=1400,0<=y<=1200的高程数据见下表(单位:m).求该区域地表面积的近似值,并用插值数据画出该区域的等高线图和三 ...
数学建模与数据分析中的灰色关联分析
灰色关联分析可以在数据量比较少的情况下,分析出主要因素.次要因素等文章目录 (1) 数理统计传统方法的问题 (2) 灰色关联分析 1. 基本思想 2. 进行系统分析 3. 用于综合评价模型 (1) ...
R语言——单位根检验/平稳性检验
导入数据: 得到数据集x: 这里RM1是指货币供应量的同比增长率,想要检验RM1的平稳性. 先安装程序包tseries: 调用tseries: 进行单位根检验: P值为0.2843,大于0.05,不能 ...
python平稳性检验_Python数据分析0.3 用statsmodels进行ADF平稳性检验
#statsmodels用于数据的统计建模分析 #此例为ADF平稳性检验的例子 from statsmodels.tsa.stattools import adfuller as ADF import ...

数学建模时序数据分析——趋势性检验和平稳性检验