基于ARMA-偏tGARCH和DCC-GARCH模型测算CoVaR——R语言实现
CoVaR是目前金融学界和管理实践中较为主流的测量一个机构（系统）对另一个机构（系统）风险溢出的指标，计算CoVaR的方法主要有分位数回归法、Coupla模型和DCC-GARCH型。本文主要介绍如何利用DCC-GARCH模型对CoVaR进行计算并利用R实现。代码见文末。
CoVaR
CoVaR这一概念由VaR衍生而来，其经济含义是当某一个机构发生风险时，另一机构所承担的总体在险价值。VaR则是指在一定的持有期内，在一定的置信水平下，机构所面临的最大潜在损失。VaR的一般表达式可以为：
prob ( Δ P > VaR ) = 1 − T \operatorname{prob}(\Delta P>\operatorname{VaR})=1-T prob(ΔP>VaR)=1−T
根据定义，CoVaR的一般表达式则为：
Pr ( X i ⩽ C o V a R q i j ∣ X j = V a R q j ) = q \operatorname{Pr}\left(\mathrm{X}^{\mathrm{i}} \leqslant \mathrm{CoVaR}_{\mathrm{q}}^{\mathrm{ij}} \mid \mathrm{X}^{\mathrm{j}}=\mathrm{VaR}_{\mathrm{q}}^{\mathrm{j}}\right)=\mathrm{q} Pr(Xi⩽CoVaRqij∣Xj=VaRqj)=q
机构i在机构j发生风险的情况下所承担的总体风险是 C o V a R q i j \mathrm{CoVaR}_{\mathrm{q}}^{\mathrm{ij}} CoVaRqij，如果要剥离机构i本身的自有风险，计算得到机构j给机构i带来的风险增量，则就采用 Δ C o V a R q i j \Delta \mathrm{CoVaR}_{\mathrm{q}}^{\mathrm{ij}} ΔCoVaRqij指标：
Δ C o V a R q i j = C o V a R q i j − V a R q i \Delta \mathrm{CoVaR}_{\mathrm{q}}^{\mathrm{ij}}=\mathrm{CoVaR}_{\mathrm{q}}^{\mathrm{ij}}-\mathrm{VaR}_{\mathrm{q}}^{\mathrm{i}} ΔCoVaRqij=CoVaRqij−VaRqi
DCC-GARCH的整体思想
DCC-GARCH模型主要用于描述不同变量间（本文主要讲双变量）的动态变化关系以及相互影响，其核心步骤在于：①先计算单变量的GARCH模型；②然后通过所得相关参数计算双变量间的动态相关系数 ρ \rho ρ。因此，DCC-GARCH模型的出发点在于构建GARCH模型。
一般的GARCH模型表达式如下（以GARCH（1,1）为例）：
r t = u t + ε t r_{t}=u_{t}+\varepsilon_{t} rt=ut+εt
ε t = σ t z t \varepsilon_{t}={\sigma_{t}} z_{t} εt=σtzt
σ t 2 = ω + α ε t − 1 2 + β σ t − 1 2 \sigma_{t}^{2}=\omega+\alpha \varepsilon_{t-1}^{2}+\beta \sigma_{t-1}^{2} σt2=ω+αεt−12+βσt−12
r t = u t + φ 1 r t − 1 + θ 1 ε t − 1 + ε t r_{t}=u_{t}+\varphi_{1} r_{t-1}+\theta_{1} \varepsilon_{t-1}+\varepsilon_{t} rt=ut+φ1rt−1+θ1εt−1+εt
由上式可见，GARCH（1，1）模型主要由收益率均值模型和方差模型两部分组成，对于收益率均值模型而言，可以进一步考虑AR（p）、MA（q）以及ARMA（p，q）模型。在本文的代码中，对收益率采用自动定阶的方法取得ARMA（p，q）。同时， z t z_{t} zt这个参数的分布选取也是值得考究的，在金融市场当中，普遍认为金融风险具有“尖峰厚尾”的特征，因此对于 z t z_{t} zt，可不必拘泥于完全理想的标准正态分布，可将其设定成服从偏t分布。如此一来，对于每一个单变量的收益率，我们就可以构建一个更符合现实情况的ARMA(p,q)-偏tGARCH模型。
如果涉及到多变量的DCC-GARCH，则涉及到矩阵和协方差的概念，如下面的式子：
H t = D t R t D t H_{t}=D_{t} R_{t} D_{t} Ht=DtRtDt
D t = diag { h t j h t j } D_{t}=\operatorname{diag}\left\{\sqrt{h_{t}^{j}}\sqrt{h_{t}^{j}}\right\} Dt=diag{ htjhtj}
R t = ( Q t ∗ ) − 1 Q t ( Q t ∗ ) − 1 R_{t}=\left(Q_{t}^{*}\right)^{-1} Q_{t}\left(Q_{t}^{*}\right)^{-1} Rt=(Qt∗)−1Qt(Qt∗)−1
Q t = ( 1 − α − β ) Q ˉ + α ( ε t − 1 i ε t − 1 j , ) + β Q t − 1 Q_{t}=(1-\alpha-\beta) \bar{Q}+\alpha\left(\varepsilon_{t-1}^{i} \varepsilon_{t-1}^{j},\right)+\beta Q_{t-1} Qt=(1−α−β)Qˉ+α(εt−1iεt−1j,)+βQt−1
其中， H t H_{t} Ht是条件协方差矩阵， D t D_{t} Dt是条件标准差 h t \sqrt{h_{t}} ht组成的对角矩阵，条件标准差则来源于单变量GARCH拟合模型。 R t R_{t} Rt是动态相关系数矩阵， Q t Q_{t} Qt是协方差矩阵， Q ˉ \bar{Q} Qˉ是经过残差标准化后的无条件协方差矩阵， Q t ∗ Q_{t}^{*} Qt∗是 Q t Q_{t} Qt的对角矩阵。
DCC-GARCH（1,1）下两机构间的动态相关系数就可以表述为以下形式：
ρ i j , t = ( 1 − α − β ) q ˉ i j + β q i j , t − 1 + α δ i , t − n δ j , t − n [ ( 1 − α − β ) q ˉ i j + β q i j , t − 1 + α δ i , t − n 2 ] 1 / 2 [ ( 1 − α − β ) q ˉ i j + β q i j , t − 1 + α δ j , t − 1 2 ] 1 / 2 \rho_{i j, t}=\frac{(1-\alpha-\beta) \bar{q}_{i j}+\beta q_{i j, t-1}+\alpha \delta_{i, t-n} \delta_{j, t-n}}{\left[(1-\alpha-\beta) \bar{q}_{i j}+\beta q_{i j, t-1}+\alpha \delta^{2}_{i, t-n}\right]^{1 / 2}\left[(1-\alpha-\beta) \bar{q}_{i j}+\beta q_{i j, t-1}+\alpha \delta^{2}_{j, t-1}\right]^{1 / 2}} ρij,t=[(1−α−β)qˉij+βqij,t−1+αδi,t−n2]1/2[(1−α−β)qˉij+βqij,t−1+αδj,t−12]1/2(1−α−β)qˉij+βqij,t−1+αδi,t−nδj,t−n
在DCC-GARCH的框架下，我们已经得到了动态相关系数，那么VaR和CoVaR的表达式就可以分别写成：
V a R q , t i = u ^ t i − Q ( q ) h ^ t i V a R_{q, t}^{i}=\hat{u}_{t}^{i}-Q(q) \hat{h}_{t}^{i} VaRq,ti=u^ti−Q(q)h^ti
CoVaR q , t i j = γ t i j V a R q , t i \operatorname{CoVaR}_{q, t}^{i j}=\gamma_{t}^{i j} V a R_{q, t}^{i} CoVaRq,tij=γtijVaRq,ti
Δ CoVaR q , t i j = γ t i j ( V a R q , t i − V a R 50 % , t i ) \Delta \operatorname{CoVaR}_{q, t}^{i j}=\gamma_{t}^{i j}\left(V a R_{q, t}^{i}-V a R_{50 \%, t}^{i}\right) ΔCoVaRq,tij=γtij(VaRq,ti−VaR50%,ti)
γ t i j = ρ i j , t σ i , t 2 / σ j , t 2 \gamma_{t}^{i j}=\rho_{i j, t} \sigma_{i, t}^{2} / \sigma_{j, t}^{2} γtij=ρij,tσi,t2/σj,t2
其中，Q(q)为置信水平为1-q时机构i收益率所符合的分布下的q分位数值，求解 V a R q , t i V a R_{q, t}^{i} VaRq,ti的关键在于通过拟合单变量GARCH模型得到条件标准差 h ^ t i \hat{h}_{t}^{i} h^ti，至于Q(q)，可以利用现有样本取分位数，也可以通过蒙特卡洛的方法构建更大的样本取分位数，这个在本文中不讨论。而求解 CoVaR q , t i j \operatorname{CoVaR}_{q, t}^{i j} CoVaRq,tij的关键则在于通过DCC函数找到动态相关系数 ρ i j , t \rho_{i j, t} ρij,t。
R语言代码实现（以测量我国银行对股市的风险溢出效应为例）
已有很多学术研究都在探讨我国银行系统对于股市的风险溢出效应，本文就以此例作为分析。从Wind上选取2010年1月1日到2020年12月31日申万银行指数以及上证指数共计11年的历史数据。
###导入数据###
##Bank_data:经对数化处理的银行收益率序列
##system_data:经对数化处理的上证指数收益率序列
###ARMA定阶部分###
auto.arima(Bank_data) ##对银行收益率进行自动定阶处理，结果为ARMA(3,2)
auto.arima(system_data) ##对上证指数收益率进行自动定阶处理,结果为ARMA(1,0)
arma_Bank=arima(Bank_data,order=c(3,0,2))
arma_system=arima(system_data,order=c(1,0,0))
arma_Bank_resd=arma_Bank$residuals #获得残差数据，下同
arma_system_resd=arma_system$residuals
resd_series=cbind(arma_Bank_resd,arma_system_resd) #得到两个序列在ARMA模型下的残差序列
data_matrix=cbind(Bank_data,system_data)
###构建单序列ARMA-GARCH模型部分###
module_Bank <- ugarchspec(variance.model = list(model = "sGARCH",
garchOrder = c(1, 1),
submodel = NULL,
external.regressors = NULL,
variance.targeting = FALSE),
mean.model = list(armaOrder = c(3, 2),
external.regressors = NULL,
start.pars = list(),
fixed.pars = list()),
distribution.model = "sstd") #选择偏t分布
fitgarch_Bank=ugarchfit(module_Bank,data=Bank_data,solver="nlminb") #通过ARMA(3,2)-偏tGARCH（1,1）拟合银行业的收益率指数
module_system <- ugarchspec(variance.model = list(model = "sGARCH",
garchOrder = c(1, 1),
submodel = NULL,
external.regressors = NULL,
variance.targeting = FALSE),
mean.model = list(armaOrder = c(1, 0),
external.regressors = NULL,
start.pars = list(),
fixed.pars = list()),
distribution.model = "sstd") #选择偏t分布
fitgarch_system=ugarchfit(module_system,data=system_data,solver="nlminb") #通过ARMA(1,0)-偏tGARCH（1,1）拟合上证指数收益率指数
###构建双序列DccGARCH模型###
dcc.garch.spec = dccspec(uspec = multispec(c(module_Bank,
module_system)),
dccOrder = c(1,1),
distribution = "mvt")
dcc.fit = dccfit(dcc.garch.spec, data = data_matrix[,c(1,2)]) #构建双序列DccGARCH模型
dcc.fit #显示模型拟合结果
cor<-data.frame(dcc.fit@mfit$R)
cor<-t(cor)
cor<-cor[,1]
del<-seq(1,length(cor),by=2)
cor<-cor[-del]
corr_DccGARCH=as.numeric(cor) #得到传统DccGARCH模型下的动态相关系数
class(time_data)
corr_DccGARCH_plot<-zoo(corr_DccGARCH,as.Date(as.character(time_data))) #生成传统DccGARCH模型下的动态相关系数图
###bayesDccGARCH部分###
out_data=bayesDccGarch(resd_series,nSim=5000,tail_ini=8,alpha_ini=rep(0.01,ncol(resd_series)),beta_ini=rep(0.97,ncol(resd_series)),errorDist=2) #errorDist=2表示残差服从偏态t分布，该贝叶斯分布通过MH算法得到
bayes_var_x=out_data$H[,1]
bayes_var_y=out_data$H[,4]
bayes_cvar_xy=out_data$H[,3]
bayes_corr=sqrt((bayes_cvar_xy)/(bayes_var_x*bayes_var_y)) #贝叶斯DCCGARCH模型下的动态相关系数
bayes_corr_DccGARCH_plot<-zoo(bayes_corr,as.Date(as.character(time_data)))
###动态相关性作图对比###
plot(bayes_corr_DccGARCH_plot,ylab="corr",main="银行-金融市场")
lines(corr_DccGARCH_plot,col="blue") #动态相关性对比
###Dcc-GARCH方法下的CoVaR计算###
sigma_DccGARCH_data=dcc.fit@model[["sigma"]]
mu_DccGARCH_data=dcc.fit@model[["mu"]]
i_sigma=sigma_DccGARCH_data[,1]
system_sigma=sigma_DccGARCH_data[,2]
i_mu=mu_DccGARCH_data[,1]
VaR_i_0.05=i_mu+(quantile(ecdf(Bank_data),0.05))*(i_sigma)
VaR_i_0.5=i_mu+(quantile(ecdf(Bank_data),0.5))*(i_sigma)
gamma_DccGARCH=(corr_DccGARCH*((system_sigma)^2)/((i_sigma)^2))
CoVaR_system_DccGARCH=gamma_DccGARCH*VaR_i_0.05 #得到DCC-GARCH方法下的CoVaR值)
转载：https://blog.csdn.net/weixin_42466310/article/details/114916138