时间序列分析基础-描述性统计分析（建模前准备工作，代码+数据）

专题一 时间序列分析基础

本文从股票时间序列数据入手，全面且详细的介绍了时间序列数据分析的基础知识，以帮助读者更好地理解和分析这一关键领域。在文章中提供了清晰的定义和相关示例，以便更好地理解时间序列分析的概念和方法。
本文内容包含了简单趋势分析、白噪声分析、自相关性分析等时间序列的描述性统计分析，是时间序列建模前的准备工作。这些内容可以帮助新手更好地理解和应用时间序列的描述性统计分析，为时间序列建模做好准备。

在接下来的时间序列专题文章中，会继续扩展并深入探讨与时间序列数据分析相关的更多主题，包括：

1. 特殊日期的趋势分析：介绍如何分析和建模节假日、每周、每日和每年等特殊日期对时间序列数据的影响。这包括实时调整和节假日效应的处理。
2. Python时间序列建模：提供更具体的示例和案例，以帮助读者快速上手Python中常用的时间序列建模库，如Prophet、ARIMA、Exponential Smoothing等。为读者提供代码示例以及如何应用这些模型进行分析。
3. 金融大数据应用案例：提供实际的金融时间序列数据应用案例，以展示如何应用时间序列分析在金融领域中解决问题，如股票价格预测、风险管理、投资策略等。
4. 学习资料和资源推荐：提供有关金融大数据和时间序列分析的学习资料、教程、在线课程和相关书籍的推荐，进一步学习和深化知识。

关注gzh“finance 褪黑素”获取案例代码等更多金融大数据等相关学习资料。

1.1 时间序列定义

时间序列是指按照时间顺序排列的一系列数据点的集合。这些数据点通常是按照固定时间间隔或时间点采集的，例如每天、每小时、每分钟或每秒等。

对某一个或者一组变量$x(t)$进行观察测量，将在一系列时刻$t_1,t_2,\cdots ,t_n$所得到的离散数字组成的序列集合，称之为时间序列。

例如: 某股票A从2020年6月1日到2022年6月1日之间各个交易日的收盘价，可以构成一个时间序列；某地每天的最高气温可以构成一个时间序列。

1.2 时间序列的特征

• 趋势（Trend）： 趋势是时间序列数据中长期上升或下降的总体趋势。趋势反映了数据随时间推移呈现出的总体增长或减少的变化模式。

• 季节性（Seasonality）： 季节性是时间序列数据中周期性重复出现的模式。这些重复出现的模式通常是由于某种周期性因素，例如每年的季节性变化或每周的周期性变化。

• 周期性（Cyclic Patterns）： 周期性是时间序列数据中非固定周期重复出现的模式。与季节性不同，周期性的周期长度可能不固定，而且不一定有明显的规律。

• 噪声（Noise）： 噪声是时间序列数据中随机波动或不规律性的部分。噪声通常由随机因素引起，使得时间序列数据不完全符合趋势、季节性和周期性模式。

• 自相关性（Autocorrelation）： 自相关性指时间序列数据中一个数据点与前一时刻（或多个时刻）的数据点之间的相关性。自相关性分析可以帮助发现时间序列中的滞后模式。

1.3 时间序列数据收集和处理（以股票市场为例）

1.3.1 获取股票代码/指数代码和名称列表

• 从东方财富网、同花顺、新浪财经等网站找到相应股票的股票代码。

• 使用Python爬取对应股票的交易信息

pip install tushare

!pip install tushare

# 导入需要的包
import pandas as pd
import tushare as ts #获取股票数据需要安装的库
import numpy as np
import matplotlib.pyplot as plt #绘图
import seaborn as sns     #seaborn画出的图更好看，且代码更简单
sns.set(color_codes=True) #seaborn设置背景

# 根据股票代码和时间范围获取股票数据
df = ts.get_k_data('600519', start='2020-01-01', end='2022-12-31')
# 简单查看一下股票数据：
df.head()

	date	open	close	high	low	volume	code
0	2020-01-02	1022.186	1024.186	1039.246	1010.186	148099.0	600519
1	2020-01-03	1011.186	972.746	1011.186	971.086	130318.0	600519
2	2020-01-06	965.046	972.176	987.086	961.486	63414.0	600519
3	2020-01-07	971.686	988.716	993.186	970.586	47853.0	600519
4	2020-01-08	979.236	982.326	989.686	976.766	25008.0	600519

# 此时如果想要将股票数据获取到Excel文件中,代码如下：
df.to_excel('股价数据.xlsx', index=False)

1.3.2 绘制股价走势图

已经有了股价数据后，我们可以通过可视化的方式将其展示出来，这里我们首先利用set_index()函数将日期设置为行索引，这样方便等会直接用pandas库进行绘图，代码如下：

df.set_index('date', inplace=True)#columns
df.head()

	open	close	high	low	volume	code
date
2020-01-02	1022.186	1024.186	1039.246	1010.186	148099.0	600519
2020-01-03	1011.186	972.746	1011.186	971.086	130318.0	600519
2020-01-06	965.046	972.176	987.086	961.486	63414.0	600519
2020-01-07	971.686	988.716	993.186	970.586	47853.0	600519
2020-01-08	979.236	982.326	989.686	976.766	25008.0	600519

df['close'].plot(figsize=(18,4)) #绘制股价图
# plt.xticks(rotation=45) # 设置X轴标签旋转45度
plt.legend() #设置图例
plt.show()

1.3.3 收益率的计算

# 计算对数收益率
df_lgreturn=np.log(df['close'])-np.log(df['close'].shift(1)) #差分
# df['close'].pct_change().dropna()
df_lgreturn.head()

df_lgreturn=df_lgreturn.dropna()#删除缺失值
df_lgreturn.head()

df_lgreturn.plot(figsize=(18,4))
# plt.xticks(rotation=45) # 设置X轴标签旋转45度
plt.show()

1.4 时间序列图形分析

1.4.1 相关系数

对于两个向量，我们希望定义它们是不是相关。一个很自然的想法，用向量与向量的夹角来作为距离的定义，夹角小，就距离小，夹角大，就距离大。

早在中学数学中，我们就经常使用余弦公式来计算角度：$$cos<\vec{a},\vec{b}>=\frac{\vec{a}\cdot \vec{b}}{|\vec{a}||\vec{b}|}$$

而对于$\vec{a}\cdot \vec{b}$我们叫做内积，例如$$(x_1,y_1)\cdot (x_2,y_2)=x_1{x}_2+y_1{y}_2$$

我们再来看相关系数的定义公式,$X$和$Y$的相关系数为：$${\rho }_{xy}=\frac{Cov(X,Y)}{\sqrt{Var(X)Var(Y)}}$$

而我们的根据样本的估计计算公式为：$${\rho }_{xy}=\frac{{\sum }_{t=1}^T(x_t-\bar{x})(y_t-\bar{y})}{\sqrt{{\sum }_{t=1}^T(x_t-\bar{x}{)}^2{\sum }_{t=1}^T(y_t-\bar{y}{)}^2}}=\frac{\overrightarrow{(X-\bar{x})}\cdot \overrightarrow{(Y-\bar{y})}}{|\overrightarrow{(X-\bar{x})}||\overrightarrow{(Y-\bar{y})}|}$$

我们发现，相关系数实际上就是计算了向量空间中两个向量的夹角！ 协方差是去均值后两个向量的内积！

如果两个向量平行，相关系数等于1或者-1，同向的时候是1，反向的时候就是-1。如果两个向量垂直，则夹角的余弦就等于0，说明二者不相关。两个向量夹角越小，相关系数绝对值越接近1，相关性越高。 只不过这里计算的时候对向量做了去均值处理，即中心化操作。而不是直接用向量$X$,$Y$计算。

对于减去均值操作，并不影响角度计算，是一种“平移”效果，如下图所示：

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
a = pd.Series([9,8,7,5,4,2])
b = a - a.mean() # 去均值
plt.figure(figsize=(10,4))
a.plot(label='a')
b.plot(label='mean removed a')
plt.legend()

​

1.4.2 自相关函数 （Autocorrelation Function, ACF)

相关系数度量了两个向量的线性相关性，而在平稳时间序列{$r_t$}中，我们有时候很想知道，$r_t$与它的过去值$r_{t-i}$的线性相关性。 这时候我们把相关系数的概念推广到自相关系数。

$r_t$与$r_{t-l}$的相关系数称为$r_t$的间隔为$l$的自相关系数，通常记为${\rho }_l$。具体的：
$${\rho }_l=\frac{Cov(r_t,r_{t-l})}{\sqrt{Var(r_t)Var(r_{t-l})}}=\frac{Cov(r_t,r_{t-l})}{Var(r_t)}$$
这里用到了弱平稳序列的性质: $$Var(r_t)=Var(r_{t-l})$$

对一个平稳时间序列的样本{$r_t$}，$1\le t\le T$,则间隔为$l$的样本自相关系数的估计为:
$${\hat{\rho }}_l=\frac{{\sum }_{t=l+1}^T(r_t-\bar{r})(r_{t-l}-\bar{r})}{{\sum }_{t=1}^T(r_t-\bar{r}{)}^2},0⩽l⩽T-1$$

则函数$ \large \hat \rho_1,\hat \rho_2 ,\hat \rho_3…$称为$r_t$的样本自相关函数(ACF)

当自相关函数中所有的值都为0时，我们认为该序列是完全不相关的；因此，我们经常需要检验多个自相关系数是否为0。

1.检验是否存在自相关性

混成检验

原假设 H 0 : ∃ i ∈ { 1 , . . . , m } , ρ i ≠ 0 H0:\exists i \in \{1,...,m\}, \rho_i \ne 0 H0:∃i∈{1,...,m},ρi​=0
备择假设 $$H1:{\rho }_1=...={\rho }_m=0$$

混成检验统计量：$$Q(m)=T(T+2)\sum _{l=1}^m\frac{{\widehat{{\rho }_l}}^2}{T-l}$$
$Q(m)$渐进服从自由度为$m$的$ \chi^2 $分布

决策规则：
$$Q(m)>{\chi }_{\alpha }^2,拒绝H_0$$
即,$Q(m)$的值大于自由度为$m$的卡方分布$100(1-\alpha )$分位点时，我们拒绝$H_0$。

大部分软件会给出$Q(m)$的p-value，则当p-value小于等于显著性水平 $\alpha$时拒绝H0。

下面给出示例：

from scipy import stats
import statsmodels.api as sm  # 统计相关的库

data = df_lgreturn # 使用收益率序列
m = 10 # 我们检验10个自相关系数

acf,q,p = sm.tsa.acf(data,nlags=m,qstat=True)  ## 计算自相关系数 及p-value
out = np.c_[range(1,11), acf[1:], q, p]
output=pd.DataFrame(out, columns=['lag', "AC", "Q", "P-value"])
output = output.set_index('lag')
output

	AC	Q	P-value
lag
1.0	-0.005176	0.019559	0.888775
2.0	0.011402	0.114591	0.944315
3.0	-0.020999	0.437372	0.932419
4.0	0.038581	1.528482	0.821585
5.0	-0.014422	1.681154	0.891265
6.0	-0.051001	3.593148	0.731538
7.0	0.028434	4.188252	0.757857
8.0	0.036009	5.144016	0.742078
9.0	-0.034035	5.999056	0.740013
10.0	0.017096	6.215102	0.796879

我们取显著性水平为0.05，可以看出，所有的p-value都大于0.05；则我们接受原假设$H_0$。

因此，我们认为该序列，即贵州茅台收益率序列，是没有显著的序列相关性的

2.绘制ACF图形

计算原始时间序列 $y(t)$ 的自相关系数 ACF(k)，表示滞后阶数为 $k$ 的自相关系数，即：

$$ACF(k)=\text{Corr}(y(t),y(t-k))$$

其中，$\text{Corr}(\cdot )$ 表示相关系数。

定阶：

ACF图形是描述时间序列数据的自相关结构的图形。在ACF图中，横轴表示滞后阶数（lag），纵轴表示相关系数。通常，ACF图在滞后阶数为0时取值为1（因为每个数据点与自身的相关系数为1），之后随着滞后阶数的增加逐渐减小。截尾点之后的值接近于零，截尾点前的值可能出现振荡。截尾点之后第一个超过显著性水平的相关系数通常用来确定MA模型的阶数。

from statsmodels.graphics.tsaplots import plot_acf as ACF #自相关图

fig = ACF(data,lags = 20)
plt.show()

3.绘制PACF图形

计算 PACF(k)，表示在去除其他滞后阶数的影响后，当前时点与滞后阶数 $k$ 之间的偏相关系数。PACF(k)的计算方法如下：

$$PACF(1)=ACF(1)$$

$$PACF(k)=\left[\text{Corr}(y(t),y(t-k))-\sum _{i=1}^{k-1}PACF(i)\times \text{Corr}(y(t-i),y(t-k))\right],\text{ for }k>1$$

其中，${\sum }_{i=1}^{k-1}PACF(i)\times \text{Corr}(y(t-i),y(t-k))$ 表示对所有滞后阶数 $i<k$ 的 $PACF(i)$ 与 $ACF(y(t-i)$, $y(t-k))$ 的乘积求和。

通过以上递归计算，我们可以得到时间序列数据的 PACF 序列，用来帮助我们确定AR模型的阶数。PACF 在时间序列分析中是一个重要的工具，可以帮助我们理解数据的自相关结构，并选择适当的模型进行预测和分析。

定阶：

PACF图形是描述时间序列数据的偏自相关结构的图形。在PACF图中，横轴也表示滞后阶数（lag），纵轴表示偏相关系数。PACF图的截尾点之后的值也接近于零，截尾点前的值通常显示出截尾点之后的指数衰减。截尾点之后第一个超过显著性水平的偏相关系数通常用来确定AR模型的阶数。

from statsmodels.graphics.tsaplots import plot_pacf as PACF   #偏自相关图

fig = PACF(data,lags = 20)
plt.show()

1.5 时间序列的平稳性

平稳时间序列粗略地讲，一个时间序列，如果均值没有系统的变化（无趋势）、方差没有系统变化，且严格消除了周期性变化，就称之是平稳的。

# AR MA 
dataset=pd.DataFrame()
dataset['close'] = df['close']
dataset['closeDiff_1'] = df['close'].diff(1)  # 1阶差分处理
dataset['closeDiff_2'] = dataset['closeDiff_1'].diff(1)  # 2阶差分处理
dataset.plot(subplots=True,figsize=(18,12))

1.5.1 定义

1. 严平稳：

如果对所有的时刻$t$，任意正整数$k$和任意$k$个正整数$ \large (t_1,t_2…t_k)$, $$(r_{t_1},r_{t_2},......r_{t_k})$$
的联合分布与$$(r_{t_1+t},r_{t_2+t},......r_{t_k+t})$$的联合分布相同，我们称时间序列{$r_t$}是严平稳的。

也就是，$$(r_{t_1},r_{t_2},......r_{t_k})$$的联合分布在时间的平移变换下保持不变，这是个很强的条件。而我们经常假定的是平稳性的一个较弱的方式

2. 弱平稳：

若时间序列{$r_t$}满足下面两个条件:$$E(r_t)=\mu ,\mu 是常数$$ $$Cov(r_t,r_{t-l})={\gamma }_l,{\gamma }_l只依赖于l$$
则时间序列{$r_t$}为弱平稳的。即该序列的均值，$r_t$与$r_{t-l}$的协方差不随时间而改变，$l$为任意整数。

在金融数据中，通常我们所说的平稳序列，是弱平稳的。

3. 差分

差分（这里为前向），就是求时间序列{$r_t$}在$t$时刻的值$r_t$与$t-1$时刻的值$r_{t-1}$的差不妨记做$d_t$,则我们得到了一个新序列{$d_t$}，为一阶差分，对新序列{$d_t$}再做同样的操作，则为二阶差分。

通常非平稳序列可以经过$d$次差分，处理成弱平稳或者近似弱平稳时间序列。回头看上图，我们发现二阶差分得到的序列比一阶差分效果更好。

1.5.2 平稳性检验

上述的操作中对于平稳性我们主要是根据肉眼分辨，有一定误差，在金融时间序列中，我们经常使用ADF检验来检验时间序列的平稳性。

ADF是一种常用的单位根检验方法，它的原假设为序列具有单位根，即非平稳，对于一个平稳的时序数据，就需要在给定的置信水平上显著，拒绝原假设。

data2 = dataset['close'] # 贵州茅台收盘价
data2.plot(figsize=(15,4))

看图形，这里显然是非平稳的。接着我们进行ADF单位根检验。

#statsmodels.org
from statsmodels.tsa.stattools import adfuller as ADF  #平稳性检测

temp = np.array(data2)
t = ADF(temp)  # ADF检验

output=pd.DataFrame(index=['Test Statistic Value', "p-value", "Lags Used", "Number of Observations Used","Critical Value(1%)","Critical Value(5%)","Critical Value(10%)"],columns=['value'])
output['value']['Test Statistic Value'] = t[0]
output['value']['p-value'] = t[1]
output['value']['Lags Used'] = t[2]
output['value']['Number of Observations Used'] = t[3]
output['value']['Critical Value(1%)'] = t[4]['1%']
output['value']['Critical Value(5%)'] = t[4]['5%']
output['value']['Critical Value(10%)'] = t[4]['10%']
output

	value
Test Statistic Value	-2.30762
p-value	0.169522
Lags Used	0
Number of Observations Used	727
Critical Value(1%)	-3.439377
Critical Value(5%)	-2.865524
Critical Value(10%)	-2.568891

t

(-2.3076196321727642,
 0.16952226399336345,
 0,
 727,
 {'1%': -3.439376877165393,
  '5%': -2.865523768488869,
  '10%': -2.5688914082860164},
 7104.773703721578)

可以看出，p-value为0.169522，大于显著性水平0.05。原假设：序列具有单位根即非平稳。不能被拒绝。因此收盘价序列为非平稳的。

我们将序列进行1次差分后再次检验！

data2Diff = data2.diff()  # 差分
data2Diff.plot(figsize=(15,4))

temp = np.array(data2Diff)[1:] # 差分后第一个值为NaN,舍去
t = ADF(temp)  # ADF检验
print("p-value:   ",t[1])

p-value:    6.887197234767805e-24

可以看出，p-value非常接近于0，拒绝原假设，因此，该序列为平稳的。

可见，经过1次差分后的序列是平稳的，对于原序列，d的取值为1即可。

1.6 白噪声序列和线性时间序列

1.6.1 白噪声序列

随机变量$X(t)$(t=1,2,3……)，如果是由一个不相关的随机变量的序列构成的，即对于所有S不等于T，随机变量$X_t$和$X_s$的协方差为零，则称其为纯随机过程。

对于一个纯随机过程来说，若其期望和方差均为常数，则称之为白噪声过程。白噪声过程的样本实称成为白噪声序列，简称白噪声。之所以称为白噪声，是因为他和白光的特性类似，白光的光谱在各个频率上有相同的强度，白噪声的谱密度在各个频率上的值相同。

1.6.2 线性时间序列

时间序列{$r_t$}，如果能写成：
r t = μ + ∑ i = 0 ∞ ψ i a t − i μ 为 r t 的均值， ψ 0 = 1 , { a t } 为白噪声序列 \large r_t = \mu + \sum_{i=0}^{\infty}\psi_ia_{t-i} \\ \large \mu为r_t 的均值， \psi_0=1,\{a_t\}为白噪声序列 rt​=μ+i=0∑∞​ψi​at−i​μ为rt​的均值，ψ0​=1,{at​}为白噪声序列
则我们称{$r_t$} 为线性序列。其中$a_t$称为在$t$时刻的新息(innovation)或扰动(shock)

很多时间序列具有线性性，即是线性时间序列，相应的有很多线性时间序列模型，例如接下来要介绍的AR、MA、ARMA，都是线性模型，但并不是所有的金融时间序列都是线性的

对于弱平稳序列，我们利用白噪声的性质很容易得到$r_t$的均值和方差：
$$\begin{gathered} E(r_t)=\mu ,Var(r_t)={\sigma }_a^2\sum _{i=0}^{\infty }{\psi }_i^2 \\ {\sigma }_a^2为a_t的方差。 \end{gathered}$$

因为$Var(r_t)$一定小于正无穷，因此$ {\psi_i^2}$必须是收敛序列，因此满足 $$i\to \infty 时，{\psi }_i^2\to 0$$

即，随着$i$的增大，远处的扰动$a_{t-i}$对$r_t$的影响会逐渐消失。

1.6.3 白噪声检验-混成检验

一般常用来进行模型残差白噪声检验

1. 残差为白噪声，说明模型拟合的很好，残差部分为无法捕捉的纯随机数据。
2. 残差非白噪声，说明模型哪里出了问题，比如参数没调好，需要继续优化；若如何优化模型也无法使得残差为白噪声，换模型，或者对残差进行二次预测。

混成检验（Ljung-Box检验） 是一种用于检验时间序列数据是否存在自相关性的统计检验方法，是对Box-Pierce检验的改进。它的原理基于自相关函数（ACF）的概念，用于评估时间序列数据在不同滞后阶数下的自相关性。

假设我们有一个时间序列数据：$x_1,x_2,x_3,\dots ,x_n$

Ljung-Box检验的步骤如下：

1. 首先，我们需要拟合一个模型（例如ARIMA模型）来预测时间序列。然后，得到该模型的残差序列：$e_1,e_2,e_3,\dots ,e_n$。这些残差是模型预测值与实际值之间的差异。

2. 接下来，我们计算残差序列在不同滞后阶数（即不同延迟）下的自相关系数（ACF）。ACF表示在某个滞后阶数下，序列与其自身之间的相关性。

3. Ljung-Box检验统计量计算如下：
   $Q=n(n+2){\sum }_{k=1}^h\frac{{\hat{\rho }}_k^2}{n-k}$
   其中，$h$是滞后阶数，${\hat{\rho }}_k$是第$k$阶滞后的样本自相关系数。

4. 计算完成后，我们会得到一个统计量$Q$。该统计量近似服从自由度为$h$的卡方分布。通常，在显著性水平（例如0.05）下，我们会查找临界值，如果$Q$大于临界值，就拒绝原假设（残差序列无自相关性）；如果$Q$小于等于临界值，则无法拒绝原假设，认为残差序列是一个白噪声序列（无自相关性）。

Ljung-Box检验的原理在于利用残差序列的自相关系数来进行统计检验，从而判断模型是否能够较好地捕捉数据的自相关性。如果检验结果拒绝原假设，意味着模型预测存在自相关性，需要重新考虑模型的选择或进行改进；如果接受原假设，则说明模型的残差序列基本上是一个白噪声序列，可以认为模型比较合适。

代码有2种 一种是1.4.2 中的第1小点中的方式，还有一种是通过直接调用Ljung-BOX检验

from statsmodels.stats.diagnostic import acorr_ljungbox as LjungBox #白噪声检验

from statsmodels.stats.diagnostic import acorr_ljungbox as LjungBox

案例1 ：

1. 自选金融资产数据，整理并计算：（1）简单收益率；（2）对数收益率。
2. 计算统计指标，并据此进行描述性统计分析；
3. 绘制数据的分布图形。
4. 检验分布是否服从正态分布。

关注gzh“finance 褪黑素”获取案例代码等更多金融大数据等相关学习资料。