抱歉,您的浏览器无法访问本站

本页面需要浏览器支持(启用)JavaScript


了解详情 >

HermitLSR

[曾记否,到中流击水,浪遏飞舟?]

  学习《时间序列分析及应用·Jonathan D.Cryer & Kung-Sik Chan 著》一书时的相关笔记——本文介绍了ARMA(1,1)模型和AR(2)模型的自相关函数的R语言实现。

ARMA(1,1)模型

ARMA(1,1)模型的自相关函数

  模型定义的方程是:

Yt=ϕYt1+etθet1Y_t=\phi Y_{t-1}+e_t-\theta e_{t-1}

为得到Yule-Walk形式的方程组,首先:

E(etYt)=E(et(ϕYt1+etθet1))=σe2E(e_tY_t)=E(e_t(\phi Y_{t-1}+e_t-\theta e_{t-1}))=\sigma^2_e

E(et1Yt)=E(et1(ϕYt1+etθet1))=(ϕθ)σe2E(e_{t-1}Y_t)=E(e_{t-1}(\phi Y_{t-1}+e_t-\theta e_{t-1}))=(\phi-\theta)\sigma^2_e

在其定义式两边乘以YtkY_{t-k}并求期望值,得到:

γ0=ϕγ1+[1θ(ϕθ)σe2]γ1=ϕγ0θσe2γk=ϕγk1}\left.\begin{aligned} \gamma_0 & =\phi\gamma_1+[1-\theta(\phi-\theta)\sigma^2_e]\\ \gamma_1 & = \phi\gamma_0-\theta\sigma^2_e\\ \gamma_k & =\phi\gamma_{k-1} \end{aligned} \right \}

求解前两个方程得到:

γ0=12ϕθ+θ21ϕ2σe2\gamma_0=\frac{1-2\phi\theta+\theta^2}{1-\phi^2}\sigma^2_e

进一步通过递推关系式,得到:

ρk=(1θϕ)(ϕθ)12ϕθ+θ2ϕk1k>1\rho_k=\frac{(1-\theta\phi)(\phi-\theta)}{1-2\phi\theta+\theta^2}\phi^{k-1},k>1

  值得注意的是,随着之后长度kk的增加,该模型的相关函数指数递减。相关函数的初始值ρ0=0\rho_0=0,从ρk\rho_k的递推式中可以看到,它的形状依赖于ρ1\rho_1ϕ\phi的符号。

ARMA(1,1)模型自相关函数的R语言函数

  按照前述中的递推式可构造求解ARMA(1,1)模型自相关函数并画出其自相关图像的函数 ARMA1.1 <- function(phi=phi, theta=theta, max.lag=20, pic=T, ylim=NULL)

1
2
3
4
5
6
7
8
9
10
11
12
13
ARMA1.1 <- function(phi=phi, theta=theta, max.lag=20, pic=T, ylim=NULL){
rho = NULL
rho[1] <- 1
for (k in 1:max.lag) {
rho[k+1] <- (1-theta*phi)*(phi-theta)/(1-2*theta*phi+theta^2)*phi^(k-1)
}
if(pic == T){
plot(x = 1:max.lag, y = rho[-1], type = "h",
xlab = "Lag", ylab = "ACF", ylim = ylim)
abline(h=0)
}
out <- list(phi = phi, theta = theta, rho = rho)
}

  该函数包含五个参数:

  1. phi  ARMA(1,1)模型中自回归部分的系数,即ϕ\phi.

  2. theta  ARMA(1,1)模型中滑动平均部分的系数,即θ\theta.

  3. max.lag  ARMA(1,1)模型中需要求解的最大滞后长度k,默认为k=20k=20.

  4. pic  逻辑值,pic=T时绘制自相关函数图像,默认pic=T.

  5. ylim  一个二维的向量,作为绘制自相关函数图像的纵坐标范围,默认为ylim=NULL

该函数的输出为一个列表,其中phi为输入的phi,theta为输入的theta,rho为计算所得的max.lag+1个自相关函数(第一个为ρ0=1\rho_0=1)。

AR(2)模型

AR(2)模型的自相关函数

  模型定义的方程是:

Yt=ϕ1Yt1+ϕ2Yt2+etY_t=\phi_1Y_{t-1}+\phi_2Y_{t-2}+e_t

为推导其自相关函数,在其两边乘上YtkY_{t-k},并求期望。假定平稳、零均值,并且ete_t独立于Yt1Y_{t-1},得到:

γk=ϕ1γk1+ϕ2γk2k=1,2,3,\gamma_k=\phi_1\gamma_{k-1}+\phi_2\gamma_{k-2},k=1,2,3,\cdots

两边再除以γ0\gamma_0

ρk=ϕ1ρk1+ϕ2ρk2k=1,2,3,\rho_k=\phi_1\rho_{k-1}+\phi_2\rho{k-2},k=1,2,3,\cdots

称上述方程为Yule-Walker方程。

  下面需求解 k=1k=1k=2k=2 时对应的ρ\rho值。首先令 k=1k=1,因 ρ0=1ρ1=ρ1\rho_0=1,\rho_{-1}=\rho_1,得到 ρ1=ϕ1+ϕ2ρ2\rho_1=\phi_1+\phi_2\rho_2 ,整理得:

ρ1=ϕ11ϕ2\rho_1=\frac{\phi_1}{1-\phi_2}

再令 k=2k=2,将上式带入Yule-Walker方程中,得出

ρ2=ϕ2(1ϕ2)+ϕ121ϕ2\rho_2=\frac{\phi_2(1-\phi_2)+\phi_1^2}{1-\phi_2}

至此我们可以通过Yule-Walker方程求出AR(2)模型的自相关函数。

  下面在给出计算ρk\rho_k的显示表达:

ρk=(1G22)G1k+1(1G12)G2k+1(G1G2)(1+G1G2)k0\rho_k=\frac{(1-G^2_2)G^{k+1}_1-(1-G^2_1)G^{k+1}_2}{(G_1-G_2)(1+G_1G_2)},k\ge 0

其中,G1=ϕϕ12+4ϕ22G1=\frac{\phi-\sqrt{\phi^2_1+4\phi_2}}{2}G2=ϕ+ϕ12+4ϕ22G2=\frac{\phi+\sqrt{\phi^2_1+4\phi_2}}{2}

特别地,当ϕ12+4ϕ2<0\phi^2_1+4\phi_2<0时,其可表示为:

ρk=Rksin(Θk+Φ)sin(Φ)k0\rho_k=R^k\frac{sin(\Theta k+\Phi)}{sin(\Phi)},k\ge 0

其中,R=ϕ2R=\sqrt{-\phi_2} 为阻尼因子,Θ=arccos(ϕ1/(2ϕ2))\Theta=arccos(\phi_1/(2\sqrt{-\phi_2})) 为频率,Φ=arctan[(1ϕ2)/(1+ϕ2)]\Phi=arctan[(1-\phi_2)/(1+\phi_2)]

AR(2)模型自相关函数的R语言函数

  按照前述中的递推式可构造求解AR(2)模型自相关函数并画出其自相关图像的函数 AR2 <- function(phi1=phi1, phi2=phi2, max.lag = 20, pic = T, ylim=c(-1,1)),此函数还可求出其特征根为复数时的阻尼因子RRΘ\Theta

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
AR2 <- function(phi1=phi1, phi2=phi2, max.lag = 20, pic = T, ylim=c(-1,1)){
rho = NULL;root=NULL;
rho1 <- phi1/(1-phi2)
rho2 <- (phi2*(1-phi2)+phi1^2)/(1-phi2)
rho[1] <- rho1
rho[2] <- rho2

for (k in 3:max.lag) {
rho[k] <- phi1*rho[k-1]+phi2*rho[k-2]
}

if(pic == T){
plot(x = 1:max.lag, y = rho, type = "h",
xlab = "Lag", ylab = "ACF", ylim = ylim)
abline(h=0)
}

if(phi1^2+4*phi2 < 0){
root1 <- (phi1+sqrt(as.complex(phi1^2+4*phi2)))/(-2*phi2)
root2 <- (phi1-sqrt(as.complex(phi1^2+4*phi2)))/(-2*phi2)
R <- sqrt(-phi2)
Theta <- acos(phi1/2*sqrt(-phi2))
}else{
root1 <- (phi1+sqrt(phi1^2+4*phi2))/(-2*phi2)
root2 <- (phi1-sqrt(phi1^2+4*phi2))/(-2*phi2)
R <- NULL
Theta <-NULL
}
out <- list(phi=c(phi1,phi2), rho = rho, root = c(root1, root2),
R = R, Theta = Theta)
}

  该函数包含五个参数:

  1. phi1,phi2  AR(2)模型中的两个参数,即ϕ1ϕ2\phi_1,\phi_2.

  2. max.lag  AR(2)模型中需要求解的最大滞后长度k,默认为k=20k=20.

  3. pic  逻辑值,pic=T时绘制自相关函数图像,默认pic=T.

  4. ylim  一个二维的向量,作为绘制自相关函数图像的纵坐标范围,默认为ylim=c(-1,1)

该函数的输出为一个列表,其中phi为输入的phi为输入的两个参数,rho为计算所得的max.lag个自相关函数,root为特征方程的两个根,R为阻尼系数,Theta为频率。




萌ICP备20201324号 | 渝ICP备2021006569号 | 公安备50011702500564号

又拍云logo     腾讯云logo