37、时间序列分析：从基础模型到复杂组合-优快云博客

本文链接：https://blog.youkuaiyun.com/a1b2c3d/article/details/154110872

时间序列分析：从基础模型到复杂组合

1. 随机游走的特性

随机游走并非完全不可预测。虽然相邻点之间的差异由随机过程决定，但该差异的方差是恒定的。这使得我们能够估计每一步幅度的置信区间。不过，由于均值为零，我们无法确定相对于当前值的差异是正还是负。

2. 自回归模型（AR）

模型原理 ：自回归模型通过变量过去值的线性组合来预测该变量的未来值。其递归函数如下：

(defn ar [ys coefs sigma]
  (let [e (s/sample-normal 1 :sd sigma)
        y (apply + e (map * ys coefs))]
    (cons y (lazy-seq
             (ar (cons y ys) coefs sigma)))))

示例分析 ：
- 当使用初始值 1.0 和系数 2.0 且无噪声时，会生成指数增长曲线，每个时间步的值是 2 的幂次方，这是一个 p = 1 的自回归过程，即 AR(1) 过程。代码如下：

(defn ex-9-16 []
  (->> (ar [1] [2] 0)
       (take 10)
       (timeseries-plot)))

- 增加 p 值可以生成更复杂的自回归序列。例如，以下代码生成一个 AR(5) 序列，呈现出周期为 5 个点的规律循环：

(defn ex-9-17 []
  (let [init (s/sample-normal 5)]
    (->> (ar init [0 0 0 0 1] 0)
         (take 30)
         (timeseries-plot))))

- 结合自回归模型和噪声，可引入随机游走成分。将 sigma 增加到 0.2 时，代码如下：

(defn ex-9-18 []
  (let [init (s/sample-normal 5)]
    (->> (ar init [0 0 0 0 1] 0.2)
         (take 30)
         (timeseries-plot))))

此时，每 5 个点的“季节性”循环得以保留，但也混入了噪声元素。
- 自回归模型的一般方程 ：
- 一阶自回归模型方程为：$y_t = c + \phi_1 y_{t - 1} + \epsilon_t$，其中 c 是常数，$\epsilon_t$ 是误差，$y_{t - 1}$ 是前一时间步的序列值，$\phi_1$ 是系数。
- 一般的 p 阶自回归模型方程为：$y_t = c + \phi_1 y_{t - 1} + \cdots + \phi_p y_{t - p} + \epsilon_t$。由于序列是平稳的，代码中省略了常数项 c。

3. 确定 AR 模型的自相关性

自相关函数（ACF） ：自相关函数用于研究序列在不同滞后阶数下的自相关性。滞后阶数 k 为 0 时，每个点与自身比较，自相关系数为 1.0；滞后阶数为 1 时，每个点与前一个点比较。自相关函数（ACF）定义为：$\rho_{y(k)} = corr(y_t, y_{t - k})$。
自定义 ACF 函数 ：

(defn autocorrelation* [series variance n k]
  (let [lag-product (->> (drop k series)
                         (map * series)
                         (i/sum))]
    (cons (/ lag-product variance n)
          (lazy-seq
           (autocorrelation* series variance n (inc k))))))
(defn autocorrelation [series]
  (autocorrelation* (s/sweep series)
                    (s/variance series)
                    (dec (count series)) 0))

在计算自相关之前，使用 sweep 函数去除序列的均值，以便通过简单相乘判断序列值是否有共同变化的趋势。
- 绘制 ACF 图 ：

(defn bar-plot [ys]
  (let [xs (range (count ys))]
    (-> (c/bar-chart xs ys
                     :x-label "Lag"
                     :y-label "Value")
        (i/view))))
(defn ex-9-19 []
  (let [init (s/sample-normal 5)
        coefs [0 0 0 0 1]]
    (->> (ar init coefs 0.2)
         (take 100)
         (autocorrelation)
         (take 15)
         (bar-plot))))

该图中每 5 个滞后阶数出现的峰值与 AR(5) 序列生成器一致，但随着噪声干扰，自相关系数会逐渐减小。

4. 移动平均模型（MA）

模型背景 ：AR 模型假设噪声是随机的，且均值和方差恒定，噪声项之间不相关。但在某些过程中，噪声项本身是相关的，这类序列被称为移动平均（MA）模型。例如，每日烧烤销售量的时间序列可能会出现每周的高峰，以及因天气等因素导致的连续几周销售波动，表现为序列均值的明显变化。
MA 模型方程 ：
- 一阶移动平均模型（MA(1)）方程为：$y_t = \mu + \epsilon_t + \theta_1 \epsilon_{t - 1}$，其中 $\mu$ 是序列均值，$\epsilon_t$ 是噪声值，$\theta_1$ 是模型参数。
- 一般的 q 阶移动平均模型方程为：$y_t = \mu + \epsilon_t + \theta_1 \epsilon_{t - 1} + \cdots + \theta_q \epsilon_{t - q}$。
MA 模型的递归函数 ：

(defn ma [es coefs sigma]
  (let [e (s/sample-normal 1 :sd sigma)
        y (apply + e (map * es coefs))]
    (cons y (lazy-seq
             (ma (cons e es) coefs sigma)))))

与 AR 函数不同，MA 函数保留的是先前的误差序列。
- MA 序列示例 ：

(defn ex-9-20 []
  (let [init (s/sample-normal 5)
        coefs [0 0 0 0 1]]
    (->> (ma init coefs 0.5)
         (take 100)
         (timeseries-plot))))

MA 模型的图表不像 AR 模型那样有明显的重复模式，但在较长序列中会引入漂移。
- 确定 MA 模型的自相关性 ：

(defn ex-9-21 []
  (let [init (s/sample-normal 5)
        coefs [0 0 0 0 1]]
    (->> (ma init coefs 0.2)
         (take 5000)
         (autocorrelation)
         (take 15)
         (bar-plot))))

该图在滞后阶数为 5 时出现明显峰值，表明 MA 过程的阶数，但与自回归过程不同，不会每隔 5 个滞后阶数出现重复峰值。

5. 自回归和移动平均模型的对比

模型	ACF	PACF
AR(p)	逐渐衰减	p 阶滞后后截断
MA(q)	q 阶滞后后截断	逐渐衰减
ARMA(p,q)	逐渐衰减	逐渐衰减

通过这些对比，我们可以根据 ACF 和 PACF 图的特征来选择合适的时间序列模型。

6. 模型选择的流程图

graph TD;
    A[获取时间序列数据] --> B[绘制 ACF 和 PACF 图];
    B --> C{ACF 特征};
    C -- 逐渐衰减 --> D{PACF 特征};
    D -- p 阶滞后后截断 --> E[选择 AR(p) 模型];
    D -- 逐渐衰减 --> F{是否有明显周期};
    F -- 是 --> G[考虑 ARMA(p,q) 模型];
    F -- 否 --> H[数据可能较复杂，需进一步分析];
    C -- q 阶滞后后截断 --> I{PACF 特征};
    I -- 逐渐衰减 --> J[选择 MA(q) 模型];

这个流程图展示了根据 ACF 和 PACF 图选择合适时间序列模型的过程，帮助我们更清晰地进行模型决策。

7. 结合 AR 和 MA 模型：ARMA 模型

模型原理 ：AR 和 MA 模型是生成自相关时间序列的两种不同但密切相关的方法，它们并非相互排斥。在实际的时间序列建模中，常常会遇到序列似乎是两者混合的情况。我们可以将 AR 和 MA 过程组合成一个单一的 ARMA 模型，其方程为：$y_t = c + \sum_{i = 1}^{p} \phi_i y_{t - i} + \sum_{i = 1}^{q} \theta_i \epsilon_{t - i} + \epsilon_t$，其中包含两组系数：自回归模型的系数和移动平均模型的系数。通常，AR 模型的阶数用 p 表示，MA 模型的阶数用 q 表示。
ARMA 模型的递归函数 ：

(defn arma [ys es ps qs sigma]
  (let [e  (s/sample-normal 1 :sd sigma)
        ar (apply + (map * ys ps))
        ma (apply + (map * es qs))
        y  (+ ar ma e)]
    (cons y (lazy-seq
                (arma (cons y ys)
                      (cons e es)
                      ps qs sigma)))))

ARMA 序列示例 ：

(defn ex-9-22 []
  (let [ys (s/sample-normal 10 :sd 1.0)
        es (s/sample-normal 10 :sd 0.2)
        ps [0 0 0 0.3 0.5]
        qs [0.2 0.8]]
    (->> (arma ys es ps qs 0.2)
         (take 500)
         (timeseries-plot))))

这里指定 AR 模型有 5 个参数，MA 模型有 2 个参数，这被称为 ARMA(5,2) 模型。
- ARMA 模型的自相关分析 ：

(defn ex-9-23 []
  (let [ys (s/sample-normal 10 :sd 1.0)
        es (s/sample-normal 10 :sd 0.2)
        ps [0 0 0 0.3 0.5]
        qs [0.2 0.8]]
    (->> (arma ys es ps qs 0.2)
         (take 500)
         (autocorrelation)
         (take 15)
         (bar-plot))))

ARMA 模型的 ACF 图中，自相关系数逐渐衰减至零，难以确定 AR 和 MA 过程的阶数，这是因为 MA 部分的模型掩盖了 AR 部分的模型，使得数据中的循环模式难以被发现。解决这个问题的最佳方法是绘制偏自相关图。

8. 计算偏自相关

偏自相关函数（PACF）的定义 ：偏自相关函数旨在解决混合 ARMA 模型中循环成分的识别问题。它定义为在给定所有中间观测值的情况下，$y_t$ 和 $y_{t + k}$ 之间的相关系数，即滞后阶数为 k 时，未被滞后阶数 1 到 k - 1 解释的自相关。
- 一阶偏自相关（滞后阶数为 1）等于一阶自相关。
- 二阶偏自相关（滞后阶数为 2）的计算公式为：$R_{t,t - 2|t - 1} = \frac{cov(y_t, y_{t - 2}|y_{t - 1})}{\sqrt{var(y_t|y_{t - 1}) var(y_{t - 2}|y_{t - 1})}}$。
- 三阶偏自相关（滞后阶数为 3）的计算公式为：$R_{t,t - 3|t - 1,t - 2} = \frac{cov(y_t, y_{t - 3}|y_{t - 1},y_{t - 2})}{\sqrt{var(y_t|y_{t - 1},y_{t - 2}) var(y_{t - 3}|y_{t - 1},y_{t - 2})}}$，以此类推。
自协方差的计算 ：偏自相关的计算需要计算数据在某个滞后阶数下的自协方差。自协方差的计算函数如下：

(defn autocovariance* [series n k]
  (let [lag-product (->> (drop k series)
                         (map * series)
                         (i/sum))]
    (cons (/ lag-product n)
          (lazy-seq
           (autocovariance* series n (inc k))))))
(defn autocovariance [series]
  (autocovariance* (s/sweep series) (count series) 0))

使用 Durbin - Levinson 递归计算 PACF ：

(defn pac* [pacs sigma prev next]
  (let [acv (first next)
        sum (i/sum (i/mult pacs (reverse prev)))
        pac (/ (- acv sum) sigma)]
    (cons pac
          (lazy-seq
           (pac* (->> (i/mult pacs pac)
                      (reverse)
                      (i/minus pacs)
                      (cons pac))
                 (* (- 1 (i/pow pac 2)) sigma)
                 (cons acv prev)
                 (rest next))))))
(defn partial-autocorrelation [series]
  (let [acvs (autocovariance series)
        acv1 (first  acvs)
        acv2 (second acvs)
        pac  (/ acv2 acv1)]
    (concat [1.0 pac]
            (pac* (vector pac)
                  (- acv1 (* pac acv2))
                  (vector acv2)
                  (drop 2 acvs)))))

绘制偏自相关图 ：

(defn ex-9-24 []
  (let [ys (s/sample-normal 10 :sd 1.0)
        es (s/sample-normal 10 :sd 0.2)
        ps   [0 0 0 0.3 0.5]
        qs   [0.2 0.8]]
    (->> (arma ys es ps qs 0.2)
         (take 500)
         (partial-autocorrelation)
         (take 15)
         (bar-plot))))

该图在滞后阶数 1 和 2 处有较高的偏自相关，表明可能存在 MA(2) 过程；直到滞后阶数 5 才有较低的偏自相关，表明可能存在 AR(5) 模型。

9. 根据 ACF 和 PACF 确定 ARMA 模型阶数

模型	ACF	PACF
AR(p)	逐渐衰减	p 阶滞后后截断
MA(q)	q 阶滞后后截断	逐渐衰减
ARMA(p,q)	逐渐衰减	逐渐衰减

在实际应用中，数据往往不符合这些理想情况。对于一个实际的时间序列，尤其是点数不多的序列，选择最合适的模型并不总是显而易见的。通常最好的做法是选择能够描述数据的最简单模型（阶数最低的模型）。

10. 显著性阈值

在 ACF 和 PACF 图中，虚线表示显著性阈值。一般来说，我们无法构建一个完美捕捉时间序列所有自相关性的模型，显著性阈值有助于我们确定最重要的自相关性。当显著性水平 $\alpha$ 为 5% 时，显著性阈值的简单计算公式为：$\pm \frac{2}{\sqrt{n}}$，其中 n 是时间序列中的点数。如果 ACF 和 PACF 中的所有点都接近零，则数据基本上是随机的。

11. 实际数据应用：航空公司数据

绘制航空公司数据的 ACF 图 ：

(defn ex-9-25 []
  (->> (airline-passengers)
       (difference)
       (autocorrelation)
       (take 25)
       (bar-plot)))

该图显示数据中存在规则的高峰和低谷。第一个高峰在滞后阶数 12 处，第二个在滞后阶数 24 处。由于数据是每月的，这些高峰对应于年度的季节性循环。由于时间序列有 144 个点，显著性阈值约为 $\pm \frac{2}{\sqrt{144}} = \pm 0.17$。
- 绘制航空公司数据的 PACF 图 ：

(defn ex-9-26 []
  (->> (airline-passengers)
       (difference)
       (partial-autocorrelation)
       (take 25)
       (bar-plot)))

通过对航空公司数据的 ACF 和 PACF 图分析，我们可以进一步了解数据的特征，从而选择合适的时间序列模型进行建模和预测。

综上所述，时间序列分析中的 AR、MA 和 ARMA 模型为我们处理自相关时间序列提供了强大的工具。通过结合 ACF 和 PACF 图的分析，我们能够更准确地选择合适的模型，从而更好地理解和预测时间序列数据。在实际应用中，我们需要根据数据的特点和实际需求，灵活运用这些模型和方法。