本篇为「时间序列分析」系列推文的最后一篇。
1 协整的概念
「协整」(cointegrated)描述的是两个及以上时间序列变量之间的关系。
若一个时间序列变量的阶差分是平稳的,则称该变量为阶单整;若多个阶时间序列变量的线性组合的单整阶数小于,则称这些变量具有协整关系,记作;由线性组合的系数组成的向量称作「协整向量」(cointegrating vector)。
对于两个变量来说,若它们的线性组合的单整阶数小于它们自身,则它们的阶数必然相同;
对于多个变量来说,若它们的线性组合的单整阶数小于它们自身,但它们的阶数并不一定完全相等,则称它们具有「多重协整」(multicointegration)关系。举例来说,和是,记它们的一阶单整的线性组合为;若又和是,则、、三个变量的线性组合就可以是平稳的,即线性组合阶数小于它们各自的阶数。
2 误差修正模型
误差修正模型是协整关系的一个应用。以双变量为例,和从长期来看存在均衡关系,那么它们当期值会根据它们前期值之间的差进行修正。可以使用以下表达式进行描述:
402 Payment Required
式中是和的一个线性组合。假设和是平稳序列,则必然也是平稳的,即和存在协整关系。
在上式右侧加上和的滞后项:
则上式是一个VAR模型;
3 Engle-Granger检验
Engle-Granger协整检验的步骤如下:
「第一步」:确定每个变量的单整阶数。
ADF检验可以用来进行单位根检验。单位根个数即时间序列变量的单整阶数。
「第二步」:估计协整向量。
假设和具有协整关系。因为和的线性组合是平稳序列,即,也即,因此可以使用OLS方法估计协整向量。
为常数,不影响单整阶数;模型残差的单整阶数也就是线性组合后的单整阶数。
采用DF或ADF的形式对残差进行单位根检验:
虽然形式相同,但协整检验的临界值与DF或ADF检验不同,原因在于只是误差的估计值而非实际值。临界值可在以下链接中查看:http://citeseer.ist.psu.edu/viewdoc/download;jsessionid=64E6EE352E1DE391B0A75E071486F7F1?doi=10.1.1.456.4786&rep=rep1&type=pdf。
若的单位根个数少于和,则和具有协整关系。
「第三步」:估计误差修正模型。
若和具有协整关系,则可以建立误差修正模型。
因为和的线性组合等于,而是常数可以合并到模型的截距项中去,因此在模型中可以使用代替线性组合:
上式可以和VAR一样进行估计。
4 案例
这里的案例取自参考书的对应内容。示例数据可在公众号后台发送关键词「示例数据」获取,文件名为134.Coint6.xls
。
读取示例数据,包括y
、z
、w
三个变量:
library(readxl)
data = read_xls("134.Coint6.xls")
matplot(data, type = "l")

实际上,三个变量是通过模拟得到的,均为1阶单整。这里假设事先不知道其单整阶数,进行ADF检验。以变量y
为例:
library(urca)
model.11 = ur.df(data$y, type = "drift",
lags = 8, selectlags = "BIC")
summary(model.11)
## Value of test-statistic is: -0.6457 0.2288
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau2 -3.51 -2.89 -2.58
## phi1 6.70 4.71 3.86
dy = diff(data$y)
model.12 = ur.df(dy, type = "drift",
lags = 8, selectlags = "BIC")
summary(model.12)
## Value of test-statistic is: -6.4483 20.8026
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau2 -3.51 -2.89 -2.58
## phi1 6.70 4.71 3.86
根据ADF检验可知,
y
含有单位根,而一阶差分dy
是平稳的,因此可以判断y
的单整阶数为1;ADF检验请参考本系列的两篇推文:时间序列分析(6)| DF检验;时间序列分析(7)| ADF检验
使用OLS方法估计协整向量:
model.1 = lm(y ~ z + w, data = data)
summary(model.1)
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.04843 0.08422 -0.575 0.567
## z -0.92731 0.02434 -38.095 <2e-16 ***
## w 0.97688 0.01827 53.463 <2e-16 ***
model.2 = lm(z ~ y + w, data = data)
summary(model.2)
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.05896 0.08788 0.671 0.504
## y -1.01082 0.02653 -38.095 <2e-16 ***
## w 1.02549 0.01570 65.324 <2e-16 ***
model.3 = lm(w ~ y + z, data = data)
summary(model.3)
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.08524 0.08449 -1.009 0.316
## y 0.99007 0.01852 53.463 <2e-16 ***
## z 0.95347 0.01460 65.324 <2e-16 ***
对残差进行单位根检验:
e3 = residuals(model.3)
plot(e3, type = "l")
summary(ur.df(e3, lags = 8, selectlags = "BIC"))
## Value of test-statistic is: -5.2984

本例中的5%显著水平下的临界值为-3.828,因此可以拒绝存在单位根的原假设,即可以认为
e3
是平稳序列。
建立误差修正模型并估计:
library(vars)
d.data = apply(data, 2, diff) # 差分
colnames(d.data) <- c("dy", "dz", "dw")
model = VAR(d.data, p = 1, type = "const",
exogen = e3[1:99])
## Estimated coefficients for equation dy:
## =======================================
## Call:
## dy = dy.l1 + dz.l1 + dw.l1 + const + exo1
##
## dy.l1 dz.l1 dw.l1 const exo1
## 0.178276649 0.312995255 -0.367720385 0.006065021 0.417973964
##
##
## Estimated coefficients for equation dz:
## =======================================
## Call:
## dz = dy.l1 + dz.l1 + dw.l1 + const + exo1
##
## dy.l1 dz.l1 dw.l1 const exo1
## 0.14593872 0.26246246 -0.31262789 -0.04221498 0.07409731
##
##
## Estimated coefficients for equation dw:
## =======================================
## Call:
## dw = dy.l1 + dz.l1 + dw.l1 + const + exo1
##
## dy.l1 dz.l1 dw.l1 const exo1
## 0.1557652 0.3012198 -0.4195309 -0.0398010 -0.0690757