北京大学金融时间序列分析讲义第17章： ARCH模型

这是原书(Tsay 2013)§4.5内容。

ARCH模型公式

(R. F. Engle 1982) 提出了ARCH模型（自回归条件异方差模型），
这是对将波动率定义为条件标准差，
第一次提出的波动率的理论模型。
基本思想是：

资产收益率的扰动序列at=rt−E(rt|Ft−1)是前后不相关的，
但是前后不独立。
at的不独立性，
描述为Var(rt|Ft−1)=Var(at|Ft−1)
可以用a2t的滞后值的线性组合表示。

具体地，
ARCH(m)模型为

at=σ2t=σtεt,α0+α1a2t−1+⋯+αma2t−m.(17.1)

其中{εt}是零均值单位方差的独立同分布白噪声，
α0>0，
αj≥0, j=1,2,…,m，
另外{αj}还需要满足一些条件使得Var(at)有限，
类似于AR(p)序列的平稳性的特征根条件。

在(17.1)的波动率方程的右侧，
仅出现了截止到t−1时刻的
at−1,…,at−m的确定性函数而没有新增的随机扰动，
所以称ARCH模型为确定性的波动率模型，
这意味着σ2t关于Ft−1可测，
即在t−1时刻可以确定条件方差σ2t的值。

εt的分布常取为标准正态分布，
标准化的t分布，
广义误差分布(Generalized Error Distribution)，
有些情况下还取为有偏的分布。

因为at=rt−E(rt|Ft−1)所以Eat=0，
E(at|Ft−1)=0。
由(17.1)中{εt}独立可知
εt与Ft−1独立，
从而与σ2t独立，于是

即(17.1)给出了rt的条件方差方程。

因为系数αj都是非负数，
所以历史值a2t−j较大意味着at的条件方差较大，
于是在ARCH模型框架下，
大的扰动后面倾向于会出现较大的扰动。
这里“倾向于”不是指一定会出现大的扰动，
因为a2t−j较大使得条件方差σ2t较大，
而方差大只能说出现较大的at的概率变大，
而不是一定会出现大的扰动at。
这种现象能够解释资产收益率的波动率聚集现象。

注意：
有些作者用ht=σ2t作为条件方差的记号，
这时扰动at=εtht‾‾√。

ARCH模型的性质

以最简单的ARCH(1)为例讨论模型的性质。
模型为

at=σtεt, σ2t=α0+α1a2t−1

其中α0>0,
0<α1<1。
α1>0是因为如果等于零就不能算ARCH(1)，
α1<1是为了at方差有限。

无条件期望和方差

考虑at的无条件期望和方差。

E(at)=E[rt−E(rt|Ft−1)]=E(rt)−E[E(rt|Ft−1)]=0.

即无条件期望为零。

Var(at)===E(a2t)=E[E(a2t|Ft−1)]E[E(σ2tε2t|Ft−1)]=E[σ2tE(ε2t|Ft−1)]=E(σ2t)α0+α1E(a2t−1).

因为at是平稳列，
所以Var(at)为常数，
所以

Var(at)=E(a2t)=α01−α1.

这要求0<α1<1。

高阶矩

有些应用需要利用at的高阶矩，
例如，
为了研究at是否厚尾分布，
需要计算峰度，
峰度依赖于四阶矩。

高阶矩的存在要求(17.1)中的εt与{αj}满足一定的约束性条件。
若(17.1)中的εt服从标准正态分布，
则其有任意阶矩，
这时条件四阶矩

E(a4t|Ft−1)==E(σ4tε4t|Ft−1)=σ4tE(ε4t|Ft−1)3σ4t=3(α0+α1a2t−1)2,

从而无条件的四阶矩为

E(a4t)==E[E(a4t|Ft−1)]=3E[(α0+α1a2t−1)2]3[α20+2α0α1E(a2t−1)+α21E(a4t−1)].

如果{at}为四阶矩有限的严平稳时间序列，
则Ea4t=Ea4t−1，于是

(1−3α21)Ea4t=3(α0+2α0α1α01−α1),

这要求1−3α21>0，即0<α1<3√3。
这时

Ea4t=3α20(1+α1)(1−α1)(1−3α21).

由此可以计算at的峰度为

Ea4t[Var(a2t)]2=31−α211−3α21>3.

at的超额峰度为

Ea4t[Var(a2t)]2−3=6α211−3α21>0.

这说明{at}序列是厚尾（重尾）分布，
其样本比均值和方差相同的正态序列有更多的离群值（outliers）。
这与实证分析中对资产收益率的分布观察结果一致。

一般ARCH模型的性质

对一般的ARCH(m)模型，
其性质与ARCH(1)类似，但公式更复杂。

模型可推广为二次型表示，
见(Tsay 2010)第3章。

ARCH模型的优缺点

优点：

可以产生波动率聚集；
扰动at具有厚尾分布。

缺点：

因为假定at−j通过a2t−j影响波动率σt，
所以正的扰动和负的扰动对波动率影响相同，
但是实际的资产收益率中正负扰动对波动率影响不同，
较大的负扰动比正扰动引起的波动更大。
ARCH模型对模型参数有较严格的约束条件，
即使是ARCH(1)，
为了能计算峰度，也需要α1∈(0,3√3)，
高阶的ARCH(m)的约束条件更为复杂。
这对带高斯新息的ARCH模型通过超额峰度表达厚尾性是一个限制。
只能描述条件方差的变化，
但是不能解释变化的原因。
由模型做的波动率预测会偏高。
可能需要较大的m。

ARCH模型的建模步骤

定阶

在ARCH效应检验显著后，
可以通过考察{a2t}序列的PACF来对ARCH模型定阶。
下面解释理由。

首先，
模型为

σ2t=α0+α1a2t−1+⋯+αma2t−m

因为E(a2t|Ft−1)=σ2t，
所以认为近似有

a2t≈α0+α1a2t−1+⋯+αma2t−m

这样可以用{a2t}序列的PACF的截尾性来估计ARCH阶m。

另一方面，
令ηt=a2t−σ2t，
可以证明{ηt}为零均值不相关白噪声列，
则a2t有模型

a2t=α0+α1a2t−1+⋯+αma2t−m+ηt

这是{a2t}的AR(m)模型，
但不要求{ηt}独立同分布。
从这个模型用最小二乘法估计{αj}是相合估计，
但不是有效（方差最小）估计。
因此从{a2t}的PACF估计m是合理的。

作为例子，
考虑§16.4.2的美元对欧元汇率日对数收益率问题，
画出对数收益率平方序列的PACF:

d.useu <- read_table(
  "d-useu9910.txt", 
  col_types=cols(.default=col_double())
)
xts.useu <- with(d.useu, xts(rate, make_date(year, mon, day)))
xts.useu.lnrtn <- diff(log(xts.useu))[-1,]
eu <- c(coredata(xts.useu.lnrtn))

美元对欧元日对数收益率平方的PACF

图17.1: 美元对欧元日对数收益率平方的PACF

结果提示需要用高阶（如m=29）的ARCH模型，
从AR模型建模经验，
有可能采用类似ARMA格式的模型会减少参数使用。
GARCH模型可以进行这种改进。

模型估计

对ARCH(m)模型，考虑最大似然估计或条件最大似然估计。
扰动at=εtσt,
σ2t=α0+α1a2t−1+⋯+αma2t−m，
记模型参数α=(α0,α1,…,αm)T。
模型的似然函数与假定的εt的分布有关，
存在多种似然函数形式。

似然函数即a1,…,aT的联合密度为

其中f(at|Ft−1,α)表示给定t−1时刻已知值条件下at的条件密度。

正态假设

当假定εt为独立同标准正态分布分布随机变量列时，
在已知a1,…,at−1条件下，

σ2t=α0+α1a2t−1+⋯+αma2t−m

看成已知数，
at=εtσt的条件分布为N(0,σ2t)，
于是似然函数为

f(a1,…,aT|α)=∏t=m+1T12π‾‾‾√σtexp(−12a2tσ2t)f(a1,…,am|α)

其中f(a1,…,am|α)形式比较复杂，
常常从似然函数中去掉此项，
变成条件似然函数，
当T较大时去掉此项的影响很小。
条件似然函数为

f(am+1,…,aT|α)=∏t=m+1T12π‾‾‾√σtexp(−12a2tσ2t)

条件对数似然函数为

ℓ(am+1,…,aT|α)=−12∑t=m+1T[lnσ2t+a2tσ2t]+常数项(17.2)

其中σ2t=α0+α1a2t−1+⋯+αma2t−m是α的函数，
给定一组α，
就可以计算σ2t, t=m+1,…,T，
从而计算似然函数值。
最大化条件对数似然函数得到的估计称为正态假设下的条件最大似然估计。

t分布假设

因为收益率分布厚尾，
有些应用中假设{at}的方差模型中的εt服从标准化t分布，
自由度为v。
设ξ为服从t(v)分布的随机变量，
则v>2时Eξ=0, Var(ξ)=v/(v−2),
令η=ξv/(v−2)√，
称η的分布为标准化t(v)分布，
设εt分布为标准化t(v)分布，其分布密度为

f(εt|v)=Γ(v+12)Γ(v2)(v−2)π‾‾‾‾‾‾‾‾√(1+ε2tv−2)−v+12, εt∈(−∞,∞),(v>2).(17.3)

这时，由at=εtσt，
其中σ2t=Var(at|Ft−1)，
得条件似然函数

=f(am+1,…,aT|α,a1,…,am)[Γ(v+12)Γ(v2)(v−2)π‾‾‾‾‾‾‾‾√]T−m∏t=m+1T1σt(1+ε2t(v−2)σ2t)−v+12.(17.4)

可以先验地取定自由度v的值，
在(17.4)中关于α求最大值；
也可以将v和α一起在(17.4)中求最大值。
这样得到的参数估计称为在t分布下得条件最大似然估计。

先验地设定自由度v时，
通常取v在3到6之间。
（条件）对数似然函数为

=ℓ(am+1,…,aT|α,a1,…,am)−∑t=m+1T{v+12ln(1+a2t(v−2)σ2t)+12lnσ2t}+常数项.(17.5)

其中σ2t=α0+α1a2t−1+⋯+αma2t−m依赖于α=(α0,α1,…,αm)T。

将自由度v也看作未知参数时，
（条件）对数似然函数为

=ℓ(am+1,…,aT|v,α,a1,…,am)(T−m)[lnΓ(v+12)−lnΓ(v2)−12ln((v−2)π)]+ℓ(am+1,…,aT|α,a1,…,am).(17.6)

其中最后一项就是(17.5)的值。

有偏t分布假设

资产收益率分布除了厚尾之外还常常有偏。
可以修改t分布使其变成标准化的有偏的单峰密度。
有多种方法可以做这种修改，
这里使用(Fernandez and Steel 1998)的做法。
该方法可以在任何连续单峰且关于0对称的一元分布中引入有偏性。
将t(v)分布进行有偏化后密度为

g(ε|v,ξ)=⎧⎩⎨⎪⎪2c2ξ+1ξf(ξ(c2εt+c1)|v),2c2ξ+1ξf((c2εt+c1)/ξ|v),εt<−c1c2,εt≥−c1c2.

其中f(⋅|v)是(17.3)定义的标准化t(v)分布密度，
参数v是自由度，
可以控制厚尾程度；
参数ξ2是密度在峰值右边的面积与在峰度左边的面积之比，
代表了有偏的程度，
当ξ=1时还是标准化t(v)分布，
当ξ>1时右偏，
当ξ<1时左偏，

c1=c22=Γ(v−12)v−2‾‾‾‾‾√π‾‾√Γ(v2)(ξ−1ξ),ξ2−1ξ2−1−c21.

对标准化t分布，取v=5,10,30，作密度图：

dtstd <- function(x, df)
  sqrt(df/(df-2)) * dt(x*sqrt(df/(df-2)), df)

curve(dtstd(x, 5), -3, 3, lwd=1, lty=1, col="black",
      xlab="x", ylab="密度")
curve(dtstd(x, 10), -3, 3, lwd=1, lty=2, col="red", add=TRUE)
curve(dtstd(x, 30), -3, 3, lwd=1, lty=3, col="blue", add=TRUE)
legend("topleft", lty=c(1,2,3), col=c("black", "red", "blue"),
       legend=c(expression(v==5), expression(v==10), expression(v==30)))

标准化t分布密度图

图17.2: 标准化t分布密度图

对有偏的t分布，取v=5，ξ=0.75,1,1.5, 作密度图：

dtskew <- function(x, df, xi){
  y <- numeric(length(x))
  par1 <- gamma((df-1)/2)*sqrt(df-2)/sqrt(pi)/gamma(df/2)*(xi - 1/xi)
  par2 <- sqrt(xi^2 + 1/xi^2 - 1 - par1^2)
  c1 <- 2*par2/(xi + 1/xi)
  sele <- x < -par1/par2
  y[sele] <- c1 * dtstd(xi*(par2*x[sele] + par1), df)
  y[!sele] <- c1 * dtstd((par2*x[!sele] + par1)/xi, df)
  
  y
}

curve(dtskew(x, 5, 1.0), -3, 3, ylim=c(0, 0.6),
      lwd=1, lty=1, col="black",
      xlab="x", ylab="密度")
curve(dtskew(x, 5, 0.75), -3, 3, lwd=1, lty=2, col="red", add=TRUE)
curve(dtskew(x, 5, 1.5), -3, 3, lwd=1, lty=3, col="blue", add=TRUE)
legend("topleft", lty=c(1,2,3), col=c("black", "red", "blue"),
       legend=c(expression(xi==1), expression(xi==0.75), expression(xi==1.5)))

有偏t分布密度图

图17.3: 有偏t分布密度图

广义误差分布假设

εt的另一种可取分布是广义误差分布(GED)，密度为

f(x|v)=vλ21+1vΓ(1v)e−12∣∣xλ∣∣v, x∈(−∞,∞) (0<v≤∞).

其中v=2时即标准正态分布，
0<v<2时为厚尾分布。

λ=[2−2vΓ(1v)Γ(3v)]12.

模型验证

对一个建立好的ARCH模型，
可计算标准化残差

ã t=atσt,

其中at是均值方程的残差，
σt是波动率方程拟合的值。
{ã t}应表现为零均值、单位标准差的独立同分布序列。

对{ã t}作Ljung-Box白噪声检验，
可以考察均值方程的充分性。
对{ã 2t}作Ljung-Box白噪声检验，
可以考察波动率方程的充分性。
{ã t}的偏度、峰度、QQ图可以用来与εt的假定分布比较，
以检验模型假定的正确性。

R的fGarch扩展包可以估计ARCH模型，
并提供了多种诊断图。
rugarch包也可以估计和作图。

预测

ARCH模型的预测类似AR模型的预测。
从预测原点h出发，
对σ2t序列作超前一步预测，
即预测σ2h+1，
有

σ2h(1)=σ2h+1=α0+α1a2h+⋯+αma2h+1−m.

要做超前2步预测时，因为ah+1未知，
有
E(a2h+1|Fh)=σ2h(1),
所以

σ2h(2)=α0+α1σ2h(1)+α2a2h+⋯+αma2h+2−m.

一般地，
σ2h(ℓ)可以滚动计算，

σ2h(ℓ)=α0+∑j=1mαjσ2h(ℓ−j).

其中l−j≤0时σ2h(ℓ−j)=a2m+ℓ−j。

ARCH模型建模实例

两个实例，Intel公司股票月对数收益率序列的波动率建模，
美元对欧元的汇率的日对数收益率的波动率建模。

Intel公司股票ARCH建模实例

继续使用1973年到2009年Intel公司股票的月对数收益率。
ARCH效应检验一节进行了ARCH效应检验，
证明有ARCH效应。

数据读入：

d.intel <- read_table(
  "m-intcsp7309.txt",
  col_types=cols(
    .default=col_double(),
    date=col_date("%Y%m%d")
  ))
xts.intel <- xts(
  log(1 + d.intel[["intc"]]), d.intel[["date"]]
)
tclass(xts.intel) <- "yearmon"
ts.intel <- ts(c(coredata(xts.intel)), start=c(1973,1), frequency=12)
at <- ts.intel - mean(ts.intel)

序列的均值模型是常数均值。
减去均值之后的残差的平方的ACF:

forecast::Acf(at^2, lag.max=36, main="")

Intel股票建模中心化的残差平方的ACF

图17.4: Intel股票建模中心化的残差平方的ACF

在频率12处较高，另外在滞后1、2、3上较突出。

残差的平方的PACF:

pacf(at^2, lag.max=36, main="")

Intel股票建模中心化的残差平方的PACF

图17.5: Intel股票建模中心化的残差平方的PACF

残差平方的PACF在滞后12处较高，
另外在滞后1到3较高。
可考虑建立ARCH(3)作为波动率方程。

设rt为收益率，
拟建立如下的均值方程和波动率方程：

rt=σ2t=μ+at,at=εtσt,α0+α1a2t−1+α2a2t−2+α3a2t−3.

使用fGarch包的garchFit()函数建立ARCH模型。

library(fGarch)
mod1 <- garchFit(~ 1 + garch(3,0), data=c(ts.intel), trace=FALSE)

## Warning: Using formula(x) is deprecated when x is a character vector of length > 1.
##   Consider formula(paste(x, collapse = " ")) instead.

## 
## Title:
##  GARCH Modelling 
## 
## Call:
##  garchFit(formula = ~1 + garch(3, 0), data = c(ts.intel), trace = FALSE) 
## 
## Mean and Variance Equation:
##  data ~ 1 + garch(3, 0)
## <environment: 0x0000000027488038>
##  [data = c(ts.intel)]
## 
## Conditional Distribution:
##  norm 
## 
## Coefficient(s):
##       mu     omega    alpha1    alpha2    alpha3  
## 0.012567  0.010421  0.232889  0.075069  0.051994  
## 
## Std. Errors:
##  based on Hessian 
## 
## Error Analysis:
##         Estimate  Std. Error  t value Pr(>|t|)    
## mu      0.012567    0.005515    2.279   0.0227 *  
## omega   0.010421    0.001238    8.418   <2e-16 ***
## alpha1  0.232889    0.111541    2.088   0.0368 *  
## alpha2  0.075069    0.047305    1.587   0.1125    
## alpha3  0.051994    0.045139    1.152   0.2494    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Log Likelihood:
##  303.9607    normalized:  0.6845963 
## 
## Description:
##  Wed May 11 12:07:31 2022 by user: Lenovo 
## 
## 
## Standardised Residuals Tests:
##                                 Statistic p-Value     
##  Jarque-Bera Test   R    Chi^2  203.362   0           
##  Shapiro-Wilk Test  R    W      0.9635971 4.898647e-09
##  Ljung-Box Test     R    Q(10)  9.260782  0.5075463   
##  Ljung-Box Test     R    Q(15)  19.36748  0.1975619   
##  Ljung-Box Test     R    Q(20)  20.46983  0.4289059   
##  Ljung-Box Test     R^2  Q(10)  7.322136  0.6947234   
##  Ljung-Box Test     R^2  Q(15)  27.41532  0.02552908  
##  Ljung-Box Test     R^2  Q(20)  28.15113  0.1058698   
##  LM Arch Test       R    TR^2   25.23347  0.01375447  
## 
## Information Criterion Statistics:
##       AIC       BIC       SIC      HQIC 
## -1.346670 -1.300546 -1.346920 -1.328481

程序中1 + garch(3,0)中garch(3,0)表示σ2t的ARCH(3)模型，
1表示均值方程是一个常数。
输出结果中mu为均值方程的均值，
omega为α0，
alpha1为α1。
所以得到的均值方程和波动率方程为：

rt=σ2t=0.0126+at,at=εtσt,εt i.i.dN(0,1),0.0104+0.02329a2t−1+0.0751a2t−2+0.0520a2t−3.

因为结果中α2和α3的估计值是不显著的，
可拟合ARCH(1)模型为波动率方程：

mod2 <- garchFit( ~ 1 + garch(1,0), data=c(ts.intel), trace=FALSE)

## Warning: Using formula(x) is deprecated when x is a character vector of length > 1.
##   Consider formula(paste(x, collapse = " ")) instead.

## 
## Title:
##  GARCH Modelling 
## 
## Call:
##  garchFit(formula = ~1 + garch(1, 0), data = c(ts.intel), trace = FALSE) 
## 
## Mean and Variance Equation:
##  data ~ 1 + garch(1, 0)
## <environment: 0x0000000028b0d330>
##  [data = c(ts.intel)]
## 
## Conditional Distribution:
##  norm 
## 
## Coefficient(s):
##       mu     omega    alpha1  
## 0.013130  0.011046  0.374976  
## 
## Std. Errors:
##  based on Hessian 
## 
## Error Analysis:
##         Estimate  Std. Error  t value Pr(>|t|)    
## mu      0.013130    0.005318    2.469  0.01355 *  
## omega   0.011046    0.001196    9.238  < 2e-16 ***
## alpha1  0.374976    0.112620    3.330  0.00087 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Log Likelihood:
##  299.9247    normalized:  0.675506 
## 
## Description:
##  Wed May 11 12:07:31 2022 by user: Lenovo 
## 
## 
## Standardised Residuals Tests:
##                                 Statistic p-Value     
##  Jarque-Bera Test   R    Chi^2  144.3783  0           
##  Shapiro-Wilk Test  R    W      0.9678175 2.670321e-08
##  Ljung-Box Test     R    Q(10)  12.12248  0.2769429   
##  Ljung-Box Test     R    Q(15)  22.30705  0.1000019   
##  Ljung-Box Test     R    Q(20)  24.33412  0.2281016   
##  Ljung-Box Test     R^2  Q(10)  16.57807  0.08423723  
##  Ljung-Box Test     R^2  Q(15)  37.44349  0.001089733 
##  Ljung-Box Test     R^2  Q(20)  38.81395  0.007031558 
##  LM Arch Test       R    TR^2   27.32897  0.006926821 
## 
## Information Criterion Statistics:
##       AIC       BIC       SIC      HQIC 
## -1.337499 -1.309824 -1.337589 -1.326585

这个模型的AIC为−1.3375，
比上一个模型的−1.3467要差一些。

输出中给出了标准化残差ã t的Ljung-Box白噪声检验结果，
滞后10的p值为0.2769,
承认白噪声；
ã 2t的滞后10的Ljung-Box白噪声检验结果p值为0.08,
在0.05水平下也可以承认白噪声。

Jarque-Bera检验是正态分布的偏度峰度检验，
Shapiro-Wilk检验是正态性检验，
零假设为εt服从正态分布，
这两个检验的结果表明标准化残差不服从正态分布。

模型方程：

rt=σ2t=0.0131+at,at=εtσt,εt i.i.dN(0,1)0.0110+0.3750a2t−1(17.7)

但是，
两个不同模型给出了不同的均值方程。
这说明均值方程和波动率方程的参数是同时估计的，
不是先估计均值方程然后利用均值方程的残差估计波动率方程。

为了进行模型验证，
可计算标准化残差

ã t=atσt,t=m+1,m+2,…,T,

其中σt从估计的ARCH模型中递推地计算。
{ã t}应表现得像是独立同分布零均值白噪声列。

计算标准化残差{ã t}：

resi <- residuals(mod2, standardize=TRUE)

fGarch包的plot函数可以对建模结果绘制多种诊断图形，包括：

输入序列的时间序列图；
估计的波动率（条件标准差）折线图；
输入序列图，与正负2倍波动率叠加；
输入的ACF图；
输入的平方（即r2t）的ACF图；
互相关系数图（文档无说明）；
拟合残差图；
多个波动率图（文档无说明）；
标准化残差图；
标准化残差的ACF图；
标准化残差平方的ACF图；
r2t和rt的互相关系数函数图；
标准化残差相对于εt理论分布（条件分布）的QQ图。

标准化残差的时序图：

Intel股票建模标准化残差

图17.6: Intel股票建模标准化残差

标准化残差{ã t}的ACF:

Intel股票建模标准化残差的ACF

图17.7: Intel股票建模标准化残差的ACF

{ã 2t}的ACF:

Intel股票建模标准化残差平方的ACF

图17.8: Intel股票建模标准化残差平方的ACF

除了在滞后11、滞后21还略高以外已经没有了低阶的波动率相关。

{ã 2t}的PACF:

pacf(resi^2, lag.max=36, main="")

Intel股票建模标准化残差平方的PACF

图17.9: Intel股票建模标准化残差平方的PACF

仅在高阶的滞后11、滞后21还较高。
作为低阶模型，
ARCH(1)作为波动率方程已经比较适合。

标准化残差相对于标准正态分布的QQ图：

Intel股票标准化残差正态QQ图

图17.10: Intel股票标准化残差正态QQ图

这个图形说明标准化残差相对于正态分布有重尾。
正态分布作为条件分布不够合适。

显示拟合的条件标准差序列（波动率序列）：

Intel股票建模波动率

图17.11: Intel股票建模波动率

获得拟合的波动率用volatility()函数：

mod2v <- volatility(mod2)
head(mod2v)

## [1] 0.1306625 0.1051190 0.1450053 0.1101685 0.1134645 0.1294742

获得拟合的均值用fitted()函数。
按理说这个模型的均值方程应该是常数0.0131，
这里fitted()返回的是原始的rt序列。

mod2f <- fitted(mod2)
head(mod2f)

##            1            2            3            4            5            6 
##  0.009999835 -0.150012753  0.067064079  0.082948635 -0.110348491  0.125162849

用predict()函数作超前若干步的预测，如：

mod2p <- predict(mod2, n.ahead=12)
mod2p

##    meanForecast meanError standardDeviation
## 1    0.01313003 0.1090513         0.1090513
## 2    0.01313003 0.1245215         0.1245215
## 3    0.01313003 0.1298482         0.1298482
## 4    0.01313003 0.1317901         0.1317901
## 5    0.01313003 0.1325109         0.1325109
## 6    0.01313003 0.1327802         0.1327802
## 7    0.01313003 0.1328811         0.1328811
## 8    0.01313003 0.1329188         0.1329188
## 9    0.01313003 0.1329330         0.1329330
## 10   0.01313003 0.1329383         0.1329383
## 11   0.01313003 0.1329403         0.1329403
## 12   0.01313003 0.1329411         0.1329411

预测包括均值的预测（显然是用μ预测的）、
均值预测的标准误差、
波动率的预测（是用at的值滚动计算的）。
波动率长期预测接近于ARCH模型的无条件标准差。

模型(17.7)的一些特点：

第一，
Intel股票月对数收益率是1.31%，
年对数收益率是15.72%，
这是很了不起的。

直接计算收益率均值：

## [1] 0.0143273

第二，
α21=0.37502=0.14<1/3，
说明at序列有无条件四阶矩。

第三，
rt的无条件标准差为

0.0110/(1−0.3750)‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾√=0.1327

直接从原始数据计算：

## [1] 0.1269101

结果相近但不相同。

第四，
模型可以用来估计和预测Intel股票收益率的月波动率。

Intel股票问题改用t分布

在fGarch的garchFit()函数中加选项cond.dist="std"。

mod3 <- garchFit(~ 1 + garch(1,0), data=ts.intel, 
  cond.dist="std", trace=FALSE)

## Warning: Using formula(x) is deprecated when x is a character vector of length > 1.
##   Consider formula(paste(x, collapse = " ")) instead.

## 
## Title:
##  GARCH Modelling 
## 
## Call:
##  garchFit(formula = ~1 + garch(1, 0), data = ts.intel, cond.dist = "std", 
##     trace = FALSE) 
## 
## Mean and Variance Equation:
##  data ~ 1 + garch(1, 0)
## <environment: 0x0000000026db8518>
##  [data = ts.intel]
## 
## Conditional Distribution:
##  std 
## 
## Coefficient(s):
##       mu     omega    alpha1     shape  
## 0.017202  0.011816  0.277476  5.970266  
## 
## Std. Errors:
##  based on Hessian 
## 
## Error Analysis:
##         Estimate  Std. Error  t value Pr(>|t|)    
## mu      0.017202    0.005195    3.311 0.000929 ***
## omega   0.011816    0.001560    7.574 3.62e-14 ***
## alpha1  0.277476    0.107183    2.589 0.009631 ** 
## shape   5.970266    1.529524    3.903 9.49e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Log Likelihood:
##  315.0899    normalized:  0.709662 
## 
## Description:
##  Wed May 11 12:07:32 2022 by user: Lenovo 
## 
## 
## Standardised Residuals Tests:
##                                 Statistic p-Value     
##  Jarque-Bera Test   R    Chi^2  157.7799  0           
##  Shapiro-Wilk Test  R    W      0.9663975 1.488224e-08
##  Ljung-Box Test     R    Q(10)  12.8594   0.2316396   
##  Ljung-Box Test     R    Q(15)  23.40632  0.07588561  
##  Ljung-Box Test     R    Q(20)  25.374    0.1874956   
##  Ljung-Box Test     R^2  Q(10)  19.96092  0.02962445  
##  Ljung-Box Test     R^2  Q(15)  42.55549  0.0001845089
##  Ljung-Box Test     R^2  Q(20)  44.06739  0.00147397  
##  LM Arch Test       R    TR^2   29.76071  0.003033508 
## 
## Information Criterion Statistics:
##       AIC       BIC       SIC      HQIC 
## -1.401306 -1.364407 -1.401466 -1.386755

因为已经采用条件t分布，
所以结果中的JB检验和SW检验没有意义。

模型为

rt=σ2t=0.0172+at,at=εtσt,εt i.i.dt∗(5.97)0.0118+0.2775a2t−1(17.8)

其中t∗表示标准化t分布。

at的无条件标准差按模型估计为
0.0118/(1−0.2775)‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾√=0.1278，
正态时为0.1327，
原始数据直接估计为0.1269。
标准化残差ã t的滞后10的Ljung-Box白噪声检验p值为0.23，
可承认白噪声；
标准化残差平方ã 2t的滞后10的Ljung-Box白噪声检验p值为0.03，
仍有一定的残余波动率相关性。

作标准化残差相对于标准化t分布的QQ图：

北京大学金融时间序列分析讲义第17章： ARCH模型

可见这个分布比正态分布更符合，
但仍存在左偏，
有可能需要使用有偏t分布。

改用t分布以后，
α1估计值从原来正态时的0.3750降低到了0.2775，
说明厚尾的εt分布降低了ARCH效应。
本问题采用正态分布和t分布的结果差别不大。
后面将指出，
Intel月对数收益率的波动率方程更合适的模型是GARCH(1,1)。

fGarch包的garchFit()函数支持多种条件分布，
默认为正态分布，
用cond.dist=指定分布：

"norm": 正态；
"snorm": 有偏正态;
"ged": 广义误差分布；
"sged": 有偏广义误差分布;
"std": t分布;
"sstd": 有偏t分布；
"snig"
"QMLE"：拟最大似然估计，仍假设正态但是采用稳健标准误差估计；

欧元汇率ARCH建模实例

考虑1999-01-04到2010-08-20的欧元对美元汇率的日对数收益率。
见§16.4.2。
均值方程为rt=at，
并显示数据有ARCH效应。
该序列是纯条件异方差模型的典型例子。

d.useu <- read_table(
  "d-useu9910.txt", 
  col_types=cols(.default=col_double())
)
xts.useu <- with(d.useu, xts(rate, make_date(year, mon, day)))
xts.useu.lnrtn <- diff(log(xts.useu))[-1,]
eu <- c(coredata(xts.useu.lnrtn))

a2t的ACF:

forecast::Acf(eu^2, lag.max=36, main="")

欧元汇率平方的ACF

图17.12: 欧元汇率平方的ACF

呈现出低的较长期的相关。

a2t的PACF:

pacf(eu^2, lag.max=36, main="")

欧元汇率平方的PACF

图17.13: 欧元汇率平方的PACF

在低阶到滞后11可以，
但是高阶仍有较大的值。
考虑建立ARCH(11)模型。
采用条件高斯似然函数：

mod4 <- garchFit(~ 1 + garch(11,0), data=eu, trace=FALSE)

## Warning: Using formula(x) is deprecated when x is a character vector of length > 1.
##   Consider formula(paste(x, collapse = " ")) instead.

## 
## Title:
##  GARCH Modelling 
## 
## Call:
##  garchFit(formula = ~1 + garch(11, 0), data = eu, trace = FALSE) 
## 
## Mean and Variance Equation:
##  data ~ 1 + garch(11, 0)
## <environment: 0x00000000296880b8>
##  [data = eu]
## 
## Conditional Distribution:
##  norm 
## 
## Coefficient(s):
##         mu       omega      alpha1      alpha2      alpha3      alpha4  
## 1.2646e-04  1.8902e-05  1.6608e-02  4.4562e-02  2.7212e-02  8.0372e-02  
##     alpha5      alpha6      alpha7      alpha8      alpha9     alpha10  
## 5.0110e-02  9.2191e-02  7.5282e-02  6.9537e-02  3.3467e-02  2.7823e-02  
##    alpha11  
## 3.8773e-02  
## 
## Std. Errors:
##  based on Hessian 
## 
## Error Analysis:
##          Estimate  Std. Error  t value Pr(>|t|)    
## mu      1.265e-04   1.110e-04    1.140 0.254434    
## omega   1.890e-05   1.727e-06   10.944  < 2e-16 ***
## alpha1  1.661e-02   1.575e-02    1.055 0.291566    
## alpha2  4.456e-02   2.085e-02    2.137 0.032595 *  
## alpha3  2.721e-02   1.700e-02    1.601 0.109351    
## alpha4  8.037e-02   2.363e-02    3.402 0.000669 ***
## alpha5  5.011e-02   2.127e-02    2.355 0.018501 *  
## alpha6  9.219e-02   2.274e-02    4.053 5.05e-05 ***
## alpha7  7.528e-02   2.406e-02    3.129 0.001755 ** 
## alpha8  6.954e-02   2.455e-02    2.832 0.004622 ** 
## alpha9  3.347e-02   2.022e-02    1.655 0.097828 .  
## alpha10 2.782e-02   1.820e-02    1.528 0.126410    
## alpha11 3.877e-02   1.906e-02    2.035 0.041894 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Log Likelihood:
##  10698.47    normalized:  3.652603 
## 
## Description:
##  Wed May 11 15:22:50 2022 by user: Lenovo 
## 
## 
## Standardised Residuals Tests:
##                                 Statistic p-Value     
##  Jarque-Bera Test   R    Chi^2  360.8012  0           
##  Shapiro-Wilk Test  R    W      0.9891735 3.894158e-14
##  Ljung-Box Test     R    Q(10)  15.77628  0.1062183   
##  Ljung-Box Test     R    Q(15)  21.50105  0.12157     
##  Ljung-Box Test     R    Q(20)  24.77444  0.2101969   
##  Ljung-Box Test     R^2  Q(10)  4.801156  0.9040589   
##  Ljung-Box Test     R^2  Q(15)  16.73088  0.3352053   
##  Ljung-Box Test     R^2  Q(20)  27.56081  0.1202104   
##  LM Arch Test       R    TR^2   11.96806  0.4482485   
## 
## Information Criterion Statistics:
##       AIC       BIC       SIC      HQIC 
## -7.296329 -7.269777 -7.296368 -7.286766

通过对标准化残差ã t和ã 2t的Ljung-Box白噪声检验来看，
模型是充分的。
但是，Jarque-Bera检验是正态分布的偏度峰度检验，
Shapiro-Wilk检验是正态性检验，
这两个检验表明标准化残差不服从正态分布。

模型为

rt=σ2t=0.0013+σtεt, εt i.i.d.N(0,1)1.89×10−5+0.0166a2t−1+⋯+0.03877a2t−11

下一节会指出，
改用GARCH模型可以得到更精简的模型。

参考文献

Engle, R. F. 1982. “Autoregressive Conditional Heteroscedasticity with Estmates of the Varinance of United Kingdom Inflation.” Econometrica 50: 987–1007.

Fernandez, Carmen, and Mark F. J. Steel. 1998. “On Bayesian Modeling of Fat Tails and Skewness.” Journal of the American Statistical Association 93 (441): 359–71.

Tsay, Ruey S. 2010. Analysis of Financial Time Series. 3rd Ed. John Wiley & Sons, Inc.

———. 2013. 金融数据分析导论：基于R语言. 机械工业出版社.

北京大学金融时间序列分析讲义第17章： ARCH模型

ARCH模型公式

ARCH模型的性质

无条件期望和方差

高阶矩

一般ARCH模型的性质

ARCH模型的优缺点

ARCH模型的建模步骤

定阶

模型估计

正态假设

t分布假设

有偏t分布假设

广义误差分布假设

模型验证

预测

ARCH模型建模实例

Intel公司股票ARCH建模实例

Intel股票问题改用t分布

欧元汇率ARCH建模实例

参考文献

相关推荐

股票对冲通胀的效果好吗？投资股票可以对抗通胀吗？

回归传统？ 考虑债券配置的三个理由

医疗保健投资寻找药品以外的增长

预算协议继续推动各州国防费总体增加

对分析师的调查：华尔街对市场前景的分歧加大

发表回复

回归传统？考虑债券配置的三个理由