Discrete Difference Equation Prediction Model (DDEPM)
离散差分方程预测模型从灰度预测模型（grey prediction model）衍生出来，可以用于预测序列的发展趋势。
DDEPM过程
DDEPM的流程如下图所示
其中 x ( 0 ) x^{(0)} x(0)表示原始的序列， x ( 1 ) x^{(1)} x(1)表示DDEPM预测值。AGO表示累加生成器（Accumulated Generating Operation）用于预处理原始序列 x ( 0 ) x^{(0)} x(0)，把 x ( 0 ) x^{(0)} x(0)转变成指数增长的序列。对原始序列应用一次AGO（1-AGO）就足够生成一个指数增长序列。接下来构造一个单一变量的二阶离散差分方程（DDE(2,1)）来拟合1-AGO生成的序列。通过DDE(2,1)可以预测序列中未知的元素值。因为DDE(2,1)预测的是AGO生成的序列，因此需要把预测结果还原，方法就是IAGO (Inverse Accumulated Generating Operation)。最后得到的就是原始序列的预测结果。
（1）获取n个原始序列数据
(1)
x
(
0
)
=
{
x
(
0
)
(
1
)
,
x
(
0
)
(
2
)
,
x
(
0
)
(
3
)
,
⋯
 
,
x
(
0
)
(
n
)
}
n
∈
Z
x^{(0)} = \{ x^{(0)}(1), x^{(0)}(2), x^{(0)}(3), \cdots, x^{(0)}(n) \} \quad n \in Z \tag{1}
x(0)={x(0)(1),x(0)(2),x(0)(3),⋯,x(0)(n)}n∈Z(1)
其中
x
(
0
)
x^{(0)}
x(0)表示有n个元素的元素序列，
x
(
0
)
(
n
)
x^{(0)}(n)
x(0)(n)表示原始序列的第n个元素。
（2）AGO
生成指数增长序列
(2)
x
(
1
)
=
{
x
(
1
)
(
1
)
,
x
(
1
)
(
2
)
,
x
(
1
)
(
3
)
,
⋯
 
,
x
(
1
)
(
n
)
}
n
∈
Z
x^{(1)} = \{ x^{(1)}(1), x^{(1)}(2), x^{(1)}(3), \cdots, x^{(1)}(n) \} \quad n \in Z \tag{2}
x(1)={x(1)(1),x(1)(2),x(1)(3),⋯,x(1)(n)}n∈Z(2)
其中
x
(
p
)
=
∑
i
=
1
p
x
(
0
)
(
i
)
,
p
=
1
,
2
,
⋯
 
,
n
x^{(p)} = \sum_{i=1}^{p}x^{(0)}(i), \quad p=1,2,\cdots,n
x(p)=i=1∑px(0)(i),p=1,2,⋯,n
只对原序列进行一次AGO操作，因此称 x ( 1 ) x^{(1)} x(1)为1-AGO序列。
因为要生成指数增长序列，所以原序列不能有负数元素，否则上述公式不能生成指数增长序列。
（3）DDE(2, 1)
使用单一变量的二阶离散差分公式生成DDE(2, 1)来拟合1-AGO序列
(3)
x
(
1
)
(
p
+
2
)
+
a
⋅
x
(
1
)
(
p
+
1
)
+
b
⋅
x
(
1
)
(
p
)
=
0
x^{(1)}(p+2) + a \cdot x^{(1)}(p+1) + b \cdot x^{(1)}(p) =0 \tag{3}
x(1)(p+2)+a⋅x(1)(p+1)+b⋅x(1)(p)=0(3)
其中a和b是系数，为了确定a和b，我们使用最小二乘法。我们把公式（3）重写成
(4)
[
−
x
(
1
)
(
p
+
1
)
−
x
(
1
)
(
p
)
]
[
a
b
]
=
[
x
(
1
)
(
p
+
2
)
]
\begin{bmatrix} -x^{(1)}(p+1) & -x^{(1)}(p) \end{bmatrix} \begin{bmatrix} a \\ b \end{bmatrix} = \begin{bmatrix} x^{(1)}(p+2) \end{bmatrix} \tag{4}
[−x(1)(p+1)−x(1)(p)][ab]=[x(1)(p+2)](4)
让
p
=
1
,
2
,
⋯
 
,
n
−
2
p=1,2,\cdots,n-2
p=1,2,⋯,n−2，接着公式（4）变为
(5)
[
−
x
(
1
)
(
2
)
−
x
(
1
)
(
1
)
−
x
(
1
)
(
3
)
−
x
(
1
)
(
2
)
⋮
⋮
−
x
(
1
)
(
n
−
1
)
−
x
(
1
)
(
n
−
2
)
]
[
a
b
]
=
[
x
(
1
)
(
3
)
x
(
1
)
(
4
)
⋮
x
(
1
)
(
n
)
]
\begin{bmatrix} -x^{(1)}(2) & -x^{(1)}(1) \\ -x^{(1)}(3) & -x^{(1)}(2) \\ \vdots & \vdots \\ -x^{(1)}(n-1) & -x^{(1)}(n-2) \end{bmatrix} \begin{bmatrix} a \\ b \end{bmatrix} = \begin{bmatrix} x^{(1)}(3) \\ x^{(1)}(4) \\ \vdots \\ x^{(1)}(n) \end{bmatrix} \tag{5}
⎣⎢⎢⎢⎡−x(1)(2)−x(1)(3)⋮−x(1)(n−1)−x(1)(1)−x(1)(2)⋮−x(1)(n−2)⎦⎥⎥⎥⎤[ab]=⎣⎢⎢⎢⎡x(1)(3)x(1)(4)⋮x(1)(n)⎦⎥⎥⎥⎤(5)
在公式（5）中，我们定义
Y
=
[
x
(
1
)
(
3
)
x
(
1
)
(
4
)
⋮
x
(
1
)
(
n
)
]
(
n
−
2
)
×
1
Y = \begin{bmatrix} x^{(1)}(3) \\ x^{(1)}(4) \\ \vdots \\ x^{(1)}(n) \end{bmatrix}_{(n-2) \times 1}
Y=⎣⎢⎢⎢⎡x(1)(3)x(1)(4)⋮x(1)(n)⎦⎥⎥⎥⎤(n−2)×1
X = [ − x ( 1 ) ( 2 ) − x ( 1 ) ( 1 ) − x ( 1 ) ( 3 ) − x ( 1 ) ( 2 ) ⋮ ⋮ − x ( 1 ) ( n − 1 ) − x ( 1 ) ( n − 2 ) ] ( n − 2 ) × 2 X = \begin{bmatrix} -x^{(1)}(2) & -x^{(1)}(1) \\ -x^{(1)}(3) & -x^{(1)}(2) \\ \vdots & \vdots \\ -x^{(1)}(n-1) & -x^{(1)}(n-2) \end{bmatrix}_{(n-2) \times 2} X=⎣⎢⎢⎢⎡−x(1)(2)−x(1)(3)⋮−x(1)(n−1)−x(1)(1)−x(1)(2)⋮−x(1)(n−2)⎦⎥⎥⎥⎤(n−2)×2
Θ = [ a b ] 2 × 1 \varTheta = \begin{bmatrix} a \\ b \end{bmatrix} _{2 \times 1} Θ=[ab]2×1
公式（5）的解为
(6)
Θ
=
[
a
b
]
=
(
X
T
X
)
−
1
X
T
Y
\varTheta = \begin{bmatrix} a \\ b \end{bmatrix} = (X^TX)^{-1}X^TY \tag{6}
Θ=[ab]=(XTX)−1XTY(6)
我们之前得到的 x ( 1 ) x^{(1)} x(1)是指数增长序列，所以令 x ( 1 ) ( p ) = r p x^{(1)}(p) = r^p x(1)(p)=rp，因此
(7) r p + 2 + a ⋅ r p + 1 + b ⋅ r p = 0 r p ( r 2 + a ⋅ r + b ) = 0 \begin{aligned} r^{p+2} + a \cdot r^{p+1} + b \cdot r^{p} & = 0 \\ r^p (r^2 + a \cdot r + b) & = 0 \end{aligned} \tag{7} rp+2+a⋅rp+1+b⋅rprp(r2+a⋅r+b)=0=0(7)
公式（7）的通解是
r
1
=
−
a
+
a
2
−
4
b
2
,
r
2
=
−
a
−
a
2
−
4
b
2
r_1 = \frac{-a+\sqrt{a^2-4b}}{2}, \quad r_2 = \frac{-a-\sqrt{a^2-4b}}{2}
r1=2−a+a2−4b,r2=2−a−a2−4b
Case 1 如果
r
1
̸
=
r
2
r_1 \not = r_2
r1̸=r2，则二阶离散差分方程的解方程为
(8)
x
(
1
)
(
p
)
=
C
1
⋅
r
1
p
+
C
2
⋅
r
2
p
x^{(1)}(p) = C_1 \cdot r_1^p + C_2 \cdot r_2^p \tag{8}
x(1)(p)=C1⋅r1p+C2⋅r2p(8)
我们已经知道
x
(
0
)
(
1
)
,
x
(
0
)
(
2
)
x^{(0)}(1),x^{(0)}(2)
x(0)(1),x(0)(2)，因此可以获得两个等式
(9)
x
(
1
)
(
1
)
=
x
(
0
)
(
1
)
=
C
1
⋅
r
1
+
C
2
⋅
r
2
x^{(1)}(1) = x^{(0)}(1) = C_1 \cdot r_1 + C_2 \cdot r_2 \tag{9}
x(1)(1)=x(0)(1)=C1⋅r1+C2⋅r2(9)
(10)
x
(
1
)
(
2
)
=
x
(
0
)
(
1
)
+
x
(
0
)
(
2
)
=
C
1
⋅
r
1
2
+
C
2
⋅
r
2
2
x^{(1)}(2) = x^{(0)}(1) + x^{(0)}(2) = C_1 \cdot r_1^2 + C_2 \cdot r_2^2 \tag{10}
x(1)(2)=x(0)(1)+x(0)(2)=C1⋅r12+C2⋅r22(10)
解上述两个方程得
C
1
=
r
1
⋅
x
(
0
)
(
1
)
−
x
(
0
)
(
1
)
−
x
(
0
)
(
2
)
r
1
⋅
r
2
−
r
2
2
C
2
=
r
2
⋅
x
(
0
)
(
1
)
−
x
(
0
)
(
1
)
−
x
(
0
)
(
2
)
r
1
⋅
r
2
−
r
1
2
\begin{aligned} C_1 & = \frac{r_1 \cdot x^{(0)}(1) - x^{(0)}(1) - x^{(0)}(2)}{r_1 \cdot r_2 - r_2^2} \\ C_2 & = \frac{r_2 \cdot x^{(0)}(1) - x^{(0)}(1) - x^{(0)}(2)}{r_1 \cdot r_2 - r_1^2} \end{aligned}
C1C2=r1⋅r2−r22r1⋅x(0)(1)−x(0)(1)−x(0)(2)=r1⋅r2−r12r2⋅x(0)(1)−x(0)(1)−x(0)(2)
Case 2 如果
r
1
=
r
2
r_1 = r_2
r1=r2，解方程为
(11)
x
(
1
)
(
p
)
=
C
1
⋅
r
1
p
+
C
2
⋅
p
⋅
r
1
p
x^{(1)}(p) = C_1 \cdot r_1^p + C_2 \cdot p \cdot r_1^p \tag{11}
x(1)(p)=C1⋅r1p+C2⋅p⋅r1p(11)
同理可得到两个方程
(12)
x
(
1
)
(
1
)
=
x
(
0
)
(
1
)
=
C
1
⋅
r
1
+
C
2
⋅
r
1
x^{(1)}(1) = x^{(0)}(1) = C_1 \cdot r_1 + C_2 \cdot r_1 \tag{12}
x(1)(1)=x(0)(1)=C1⋅r1+C2⋅r1(12)
(13)
x
(
1
)
(
2
)
=
x
(
0
)
(
1
)
+
x
(
0
)
(
2
)
=
C
1
⋅
r
1
2
+
2
C
2
⋅
r
1
2
x^{(1)}(2) = x^{(0)}(1) + x^{(0)}(2) = C_1 \cdot r_1^2 + 2 C_2 \cdot r_1^2 \tag{13}
x(1)(2)=x(0)(1)+x(0)(2)=C1⋅r12+2C2⋅r12(13)
解得
C
1
=
x
(
0
)
(
1
)
×
(
2
r
1
−
1
)
−
x
(
0
)
(
2
)
r
1
2
C
2
=
x
(
0
)
(
1
)
×
(
1
−
r
1
)
+
x
(
0
)
(
2
)
r
1
2
\begin{aligned} C_1 & = \frac{x^{(0)}(1) \times (2r_1-1) - x^{(0)}(2)}{r_1^2} \\ C_2 & = \frac{x^{(0)}(1) \times (1 - r_1) + x^{(0)}(2)}{r_1^2} \end{aligned}
C1C2=r12x(0)(1)×(2r1−1)−x(0)(2)=r12x(0)(1)×(1−r1)+x(0)(2)
Case 3 如果
r
1
r_1
r1和
r
2
r_2
r2是共轭复数，那么解方程为
(14)
x
(
1
)
(
p
)
=
C
1
⋅
ρ
p
⋅
sin
(
ϕ
p
)
+
C
2
⋅
ρ
p
⋅
cos
(
ϕ
p
)
x^{(1)}(p)=C_1 \cdot \rho^p \cdot \sin (\phi p) + C_2 \cdot \rho^p \cdot \cos (\phi p) \tag{14}
x(1)(p)=C1⋅ρp⋅sin(ϕp)+C2⋅ρp⋅cos(ϕp)(14)
其中
ρ
=
(
−
a
2
)
2
+
(
4
b
−
a
2
2
)
2
=
b
,
ϕ
=
tan
−
1
(
−
4
b
−
a
2
a
)
\rho = \sqrt{\left ( - \frac{a}{2} \right )^2 + \left ( \frac{\sqrt{4b-a^2}}{2} \right )^2} = \sqrt{b}, \quad \phi = \tan^{-1} \left ( - \frac{\sqrt{4b - a^2}}{a} \right )
ρ=(−2a)2+(24b−a2)2=b,ϕ=tan−1(−a4b−a2)
同理可得到两个方程
(15)
x
(
1
)
(
1
)
=
x
(
0
)
(
1
)
=
C
1
⋅
ρ
⋅
sin
ϕ
+
C
2
⋅
ρ
⋅
cos
ϕ
x^{(1)}(1) = x^{(0)}(1) = C_1 \cdot \rho \cdot \sin \phi + C_2 \cdot \rho \cdot \cos \phi \tag{15}
x(1)(1)=x(0)(1)=C1⋅ρ⋅sinϕ+C2⋅ρ⋅cosϕ(15)
(16)
x
(
1
)
(
2
)
=
x
(
0
)
(
1
)
+
x
(
0
)
(
2
)
=
C
1
⋅
ρ
2
⋅
(
sin
2
ϕ
)
+
C
2
⋅
ρ
2
⋅
cos
(
2
ϕ
)
x^{(1)}(2) = x^{(0)}(1) + x^{(0)}(2) = C_1 \cdot \rho^2 \cdot (\sin 2\phi) + C_2 \cdot \rho^2 \cdot \cos (2\phi) \tag{16}
x(1)(2)=x(0)(1)+x(0)(2)=C1⋅ρ2⋅(sin2ϕ)+C2⋅ρ2⋅cos(2ϕ)(16)
解得
C
1
=
x
(
0
)
(
1
)
⋅
ρ
2
⋅
cos
(
2
ϕ
)
−
x
(
0
)
(
1
)
⋅
ρ
⋅
cos
ϕ
−
x
(
0
)
(
2
)
⋅
ρ
⋅
cos
ϕ
ρ
3
(
sin
ϕ
⋅
c
o
s
(
2
ϕ
)
−
cos
ϕ
⋅
sin
(
2
ϕ
)
)
C
2
=
x
(
0
)
(
1
)
⋅
ρ
⋅
sin
ϕ
+
x
(
0
)
(
2
)
⋅
ρ
⋅
sin
ϕ
−
x
(
0
)
(
1
)
⋅
ρ
2
⋅
sin
(
2
ϕ
)
ρ
3
(
sin
ϕ
⋅
c
o
s
(
2
ϕ
)
−
cos
ϕ
⋅
sin
(
2
ϕ
)
)
\begin{aligned} C_1 & = \frac{x^{(0)}(1) \cdot \rho^2 \cdot \cos(2 \phi) - x^{(0)}(1) \cdot \rho \cdot \cos \phi - x^{(0)}(2) \cdot \rho \cdot \cos \phi}{\rho^3 ( \sin \phi \cdot cos(2 \phi) - \cos \phi \cdot \sin(2 \phi))} \\ C_2 & = \frac{x^{(0)}(1) \cdot \rho \cdot \sin \phi + x^{(0)}(2) \cdot \rho \cdot \sin \phi - x^{(0)}(1) \cdot \rho^2 \cdot \sin(2 \phi)}{\rho^3 ( \sin \phi \cdot cos(2 \phi) - \cos \phi \cdot \sin(2 \phi))} \end{aligned}
C1C2=ρ3(sinϕ⋅cos(2ϕ)−cosϕ⋅sin(2ϕ))x(0)(1)⋅ρ2⋅cos(2ϕ)−x(0)(1)⋅ρ⋅cosϕ−x(0)(2)⋅ρ⋅cosϕ=ρ3(sinϕ⋅cos(2ϕ)−cosϕ⋅sin(2ϕ))x(0)(1)⋅ρ⋅sinϕ+x(0)(2)⋅ρ⋅sinϕ−x(0)(1)⋅ρ2⋅sin(2ϕ)
（4）IAGO
因为得到的离散差分方程式基于1-AGO序列，因此需要反AGO操作还原AGO操作，得到预测值。
(17)
x
^
(
0
)
=
x
(
1
)
(
p
)
−
x
(
1
)
(
p
−
1
)
\hat{x}^{(0)} = x^{(1)}(p) - x^{(1)}(p-1) \tag{17}
x^(0)=x(1)(p)−x(1)(p−1)(17)
其中
x
^
(
0
)
\hat{x}^{(0)}
x^(0)是预测值，p是DDEPM的预测步长。
DDEPM过程总结
- 用公式（2）对原序列进行1-AGO操作。
- 构造公式（5）的矩阵X和向量Y。
- 使用公式（6）计算系数a和b。
- 预测未知的序列数据。
- Case 1. 如果 a 2 − 4 b > 0 a^2-4b \gt 0 a2−4b>0，使用公式（8）来预测未知序列数据。
- Case 2. 如果 a 2 − 4 b = 0 a^2-4b = 0 a2−4b=0，使用公式（11）来预测未知序列数据。
- Case 3. 如果 a 2 − 4 b < 0 a^2-4b \lt 0 a2−4b<0，使用公式（14）来预测未知学列数据。
- 给定一个预测步长p，使用公式（17）进行IAGO还原并得到最后的预测结果。
DDEPM属性
- DDEPM的预测精度比灰度预测模型要高，而且不需要太多的样本点数据（4-5个样本点数据就足够了，算出系数a和b至少需要3个样本点）。
- 计算简单。
- 可以作为混沌时间序列预测机制。
Universal DDEPM
DDEPM有一个缺点是原始序列必须是非负序列，为了解决这个问题，通用的DDEPM（Universal DDEPM）增加了两个操作：映射生成操作（MGO）和反映射生成操作（IMGO）。通用DDEPM流程如下图所示
MGO
MGO操作把原始序列映射成非负序列
x
m
(
0
)
x_m^{(0)}
xm(0)
(18)
x
m
(
0
)
=
{
x
m
(
0
)
(
1
)
,
x
m
(
0
)
(
2
)
,
x
m
(
0
)
(
3
)
,
⋯
 
,
x
m
(
0
)
(
n
)
}
x
m
(
0
)
(
i
)
=
M
G
O
(
x
(
0
)
(
i
)
)
=
s
+
γ
⋅
x
(
0
)
(
i
)
,
i
=
1
,
2
,
3
,
⋯
 
,
n
\begin{aligned} x_m^{(0)} = \{ x_m^{(0)}(1), x_m^{(0)}(2), x_m^{(0)}(3), \cdots, x_m^{(0)}(n) \} \\ x_m^{(0)}(i) = MGO(x^{(0)}(i)) = s + \gamma \cdot x^{(0)}(i), \quad i=1,2,3,\cdots,n \end{aligned} \tag{18}
xm(0)={xm(0)(1),xm(0)(2),xm(0)(3),⋯,xm(0)(n)}xm(0)(i)=MGO(x(0)(i))=s+γ⋅x(0)(i),i=1,2,3,⋯,n(18)
即对原始序列进行缩放平移处理，其中
s
s
s表示偏移值，
γ
\gamma
γ表示缩放因子。
IMGO
IMGO操作把映射之后的预测值还原为映射之前的预测值。
(19)
x
p
(
0
)
=
I
M
G
O
(
x
^
(
0
)
(
p
)
)
=
1
γ
(
x
^
(
0
)
(
p
)
−
s
)
x_p^{(0)}=IMGO(\hat{x}^{(0)}(p)) = \frac{1}{\gamma} (\hat{x}^{(0)}(p)-s) \tag{19}
xp(0)=IMGO(x^(0)(p))=γ1(x^(0)(p)−s)(19)
其中
x
^
(
0
)
(
p
)
\hat{x}^{(0)}(p)
x^(0)(p)是IAGO的输出值，
x
p
(
0
)
x_p^{(0)}
xp(0)是最终的预测值。
内容来自《An Efficient Gradient Forecasting Search Method Utilizing the Discrete Difference Equation Prediction Model》
扫一扫
620
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
