10.9 Appendix: Some Remarks on Estimation
The estimation of multivariate ARMA models in this chapter is done by using the time series program SCA of Scientific Computing Associates. The estimation of multivariate volatility models is done by using either the S-Plus package with FinMetrics or the Regression Analysis for Time Series (RATS) program or Matlab. Below are some run streams for estimating multivariate volatility models using the RATS program. A line starting with * means “comment” only.
Estimation of the Diagonal Constant-Correlation AR(2)–GARCH(1,1) Model for Example 10.5
The program includes some Ljung–Box statistics for each component and some fitted values for the last few observations. The data file is m-ibmspln.txt, which has two columns, and there are 888 observations.
all 0 888:1
open data m-ibmspln.txt
data(org=obs) / r1 r2
set h1 = 0.0
set h2 = 0.0
nonlin a0 a1 b1 a00 a11 b11 rho c1 c2 p1
frml a1t = r1(t)-c1-p1*r2(t-1)
frml a2t = r2(t)-c2
frml gvar1 = a0+a1*a1t(t-1)**2+b1*h1(t-1)
frml gvar2 = a00+a11*a2t(t-1)**2+b11*h2(t-1)
frml gdet = -0.5*(log(h1(t)=gvar1(t))+log(h2(t)=gvar2(t)) $
+log(1.0-rho**2))
frml gln = gdet(t)-0.5/(1.0-rho**2)*((a1t(t)**2/h1(t)) $
+(a2t(t)**2/h2(t))-2*rho*a1t(t)*a2t(t)/sqrt(h1(t)*h2(t)))
smpl 3 888
compute c1 = 1.22, c2 = 0.57, p1 = 0.1, rho = 0.1
compute a0 = 3.27, a1 = 0.1, b1 = 0.6
compute a00 = 1.17, a11 = 0.13, b11 = 0.8
maximize(method=bhhh,recursive,iterations=150) gln
set fv1 = gvar1(t)
set resi1 = a1t(t)/sqrt(fv1(t))
set residsq = resi1(t)*resi1(t)
* Checking standardized residuals *
cor(qstats,number=12,span=4) resi1
* Checking squared standardized residuals *
cor(qstats,number=12,span=4) residsq
set fv2 = gvar2(t)
set resi2 = a2t(t)/sqrt(fv2(t))
set residsq = resi2(t)*resi2(t)
* Checking standardized residuals *
cor(qstats,number=12,span=4) resi2
* Checking squared standardized residuals *
cor(qstats,number=12,span=4) residsq
* Last few observations needed for computing forecasts *
set shock1 = a1t(t)
set shock2 = a2t(t)
print 885 888 shock1 shock2 fv1 fv2
Estimation of the Time-Varying Coefficient Model in Example 10.5
all 0 888:1
open data m-ibmspln.txt
data(org=obs) / r1 r2
set h1 = 45.0
set h2 = 31.0
set rho = 0.8
nonlin a0 a1 b1 f1 a00 a11 b11 d11 f11 c1 c2 p1 p3 q0 q1 q2
frml a1t = r1(t)-c1-p1*r1(t-1)-p3*r2(t-2)
frml a2t = r2(t)-c2
frml gvar1 = a0+a1*a1t(t-1)**2+b1*h1(t-1)+f1*h2(t-1)
frml gvar2 = a00+a11*a2t(t-1)**2+b11*h2(t-1)+f11*h1(t-1) $
+d11*a1t(t-1)**2
frml rh1 = q0 + q1*rho(t-1) $
+ q2*a1t(t-1)*a2t(t-1)/sqrt(h1(t-1)*h2(t-1))
frml rh = exp(rh1(t))/(1+exp(rh1(t)))
frml gdet = -0.5*(log(h1(t)=gvar1(t))+log(h2(t)=gvar2(t)) $
+log(1.0-(rho(t)=rh(t))**2))
frml gln = gdet(t)-0.5/(1.0-rho(t)**2)*((a1t(t)**2/h1(t)) $
+(a2t(t)**2/h2(t))-2*rho(t)*a1t(t)*a2t(t)/sqrt(h1(t)*h2(t)))
smpl 4 888
compute c1 = 1.4, c2 = 0.7, p1 = 0.1, p3 = -0.1
compute a0 = 2.95, a1 = 0.08, b1 = 0.87, f1 = -.03
compute a00 = 2.05, a11 = 0.05
compute b11 = 0.92, f11=-.06, d11=.04, q0 = -2.0
compute q1 = 3.0, q2 = 0.1
nlpar(criterion=value,cvcrit=0.00001)
maximize(method=bhhh,recursive,iterations=150) gln
set fv1 = gvar1(t)
set resi1 = a1t(t)/sqrt(fv1(t))
set residsq = resi1(t)*resi1(t)
* Checking standardized residuals *
cor(qstats,number=16,span=4) resi1
* Checking squared standardized residuals *
cor(qstats,number=16,span=4) residsq
set fv2 = gvar2(t)
set resi2 = a2t(t)/sqrt(fv2(t))
set residsq = resi2(t)*resi2(t)
* Checking standardized residuals *
cor(qstats,number=16,span=4) resi2
* Checking squared standardized residuals *
cor(qstats,number=16,span=4) residsq
* Last few observations needed for computing forecasts *
set rhohat = rho(t)
set shock1 = a1t(t)
set shock2 = a2t(t)
print 885 888 shock1 shock2 fv1 fv2 rhohat
Estimation of the Time-Varying Coefficient Model in Example 10.5 Using Cholesky Decomposition
all 0 888:1
open data m-ibmspln.txt
data(org=obs) / r1 r2
set h1 = 45.0
set h2 = 20.0
set q = 0.8
nonlin a0 a1 b1 a00 a11 b11 d11 f11 c1 c2 p1 p3 t0 t1 t2
frml a1t = r1(t)-c1-p1*r1(t-1)-p3*r2(t-2)
frml a2t = r2(t)-c2
frml v1 = a0+a1*a1t(t-1)**2+b1*h1(t-1)
frml qt = t0 + t1*q(t-1) + t2*a2t(t-1)
frml bt = a2t(t) - (q(t)=qt(t))*a1t(t)
frml v2 = a00+a11*bt(t-1)**2+b11*h2(t-1)+f11*h1(t-1) $
+d11*a1t(t-1)**2
frml gdet = -0.5*(log(h1(t) = v1(t))+ log(h2(t)=v2(t)))
frml garchln = gdet-0.5*(a1t(t)**2/h1(t)+bt(t)**2/h2(t))
smpl 5 888
compute c1 = 1.4, c2 = 0.7, p1 = 0.1, p3 = -0.1
compute a0 = 1.0, a1 = 0.08, b1 = 0.87
compute a00 = 2.0, a11 = 0.05, b11 = 0.8
compute d11=.04, f11=-.06, t0 =0.2, t1 = 0.1, t2 = 0.1
nlpar(criterion=value,cvcrit=0.00001)
maximize(method=bhhh,recursive,iterations=150) garchln
set fv1 = v1(t)
set resi1 = a1t(t)/sqrt(fv1(t))
set residsq = resi1(t)*resi1(t)
* Checking standardized residuals *
cor(qstats,number=16,span=4) resi1
* Checking squared standardized residuals *
cor(qstats,number=16,span=4) residsq
set fv2 = v2(t)+qt(t)**2*v1(t)
set resi2 = a2t(t)/sqrt(fv2(t))
set residsq = resi2(t)*resi2(t)
* Checking standardized residuals *
cor(qstats,number=16,span=4) resi2
* Checking squared standardized residuals *
cor(qstats,number=16,span=4) residsq
* Last few observations needed for forecasts *
set rhohat = qt(t)*sqrt(v1(t)/fv2(t))
set shock1 = a1t(t)
set shock2 = a2t(t)
set g22 = v2(t)
set q21 = qt(t)
set b2t = bt(t)
print 885 888 shock1 shock2 fv1 fv2 rhohat g22 q21 b2t
Estimation of Three-Dimensional Time-Varying Correlation Volatility Model in Example 10.7 Using Cholesky Decomposition
Initial estimates are obtained by a sequential modeling procedure.
all 0 2275:1
open data d-cscointc.txt
data(org=obs) / r1 r2 r3
set h1 = 1.0
set h2 = 4.0
set h3 = 3.0
set q21 = 0.8
set q31 = 0.3
set q32 = 0.3
nonlin c1 c2 c3 p3 p21 p22 p31 p33 a0 a1 a2 t0 t1 t2 b0 b1 $
b2 u0 u1 u2 w0 w1 w2 d0 d1 d2 d5
frml a1t = r1(t)-c1-p3*r1(t-3)
frml a2t = r2(t)-c2-p21*r1(t-2)-p22*r2(t-2)
frml a3t = r3(t)-c3-p31*r1(t-1)-p33*r3(t-1)
frml v1 = a0+a1*a1t(t-1)**2+a2*h1(t-1)
frml q1t = t0 + t1*q21(t-1) + t2*a2t(t-1)
frml bt = a2t(t) - (q21(t)=q1t(t))*a1t(t)
frml v2 = b0+b1*bt(t-1)**2+b2*h2(t-1)
frml q2t = u0 + u1*q31(t-1) + u2*a3t(t-1)
frml q3t = w0 + w1*q32(t-1) + w2*a2t(t-1)
frml b1t = a3t(t)-(q31(t)=q2t(t))*a1t(t)-(q32(t)=q3t(t))*bt(t)
frml v3 = d0+d1*b1t(t-1)**2+d2*h3(t-1)+d5*h2(t-1)
frml gdet = -0.5*(log(h1(t) = v1(t))+ log(h2(t)=v2(t)) $
+log(h3(t)=v3(t)))
frml garchln = gdet-0.5*(a1t(t)**2/h1(t)+bt(t)**2/h2(t) $
+b1t(t)**2/h3(t))
smpl 8 2275
compute c1 = 0.07, c2 = 0.33, c3 = 0.19, p1 = 0.1, p3 = -0.04
compute p21 =0.2, p22 = -0.1, p31 = -0.26, p33 = 0.06
compute a0 = .01, a1 = 0.05, a2 = 0.94
compute t0 = 0.28, t1 =0.82, t2 = -0.035
compute b0 = .17, b1 = 0.08, b2 = 0.89
compute u0= 0.04, u1 = 0.97, u2 = 0.01
compute w0 =0.006, w1=0.98, w2=0.004
compute d0 =1.38, d1 = 0.06, d2 = 0.64, d5 = -0.027
nlpar(criterion=value,cvcrit=0.00001)
maximize(method=bhhh,recursive,iterations=250) garchln
set fv1 = v1(t)
set resi1 = a1t(t)/sqrt(fv1(t))
set residsq = resi1(t)*resi1(t)
* Checking standardized residuals *
cor(qstats,number=12,span=4) resi1
* Checking squared standardized residuals *
cor(qstats,number=12,span=4) residsq
set fv2 = v2(t)+q1t(t)**2*v1(t)
set resi2 = a2t(t)/sqrt(fv2(t))
set residsq = resi2(t)*resi2(t)
* Checking standardized residuals *
cor(qstats,number=12,span=4) resi2
* Checking squared standardized residuals *
cor(qstats,number=12,span=4) residsq
set fv3 = v3(t)+q2t(t)**2*v1(t)+q3t(t)**2*v2(t)
set resi3 = a3t(t)/sqrt(fv3(t))
set residsq = resi3(t)*resi3(t)
* Checking standardized residuals *
cor(qstats,number=12,span=4) resi3
* Checking squared standardized residuals *
cor(qstats,number=12,span=4) residsq
* print standardized residuals and correlation-coefficients
set rho21 = q1t(t)*sqrt(v1(t)/fv2(t))
set rho31 = q2t(t)*sqrt(v1(t)/fv3(t))
set rho32 = (q2t(t)*q1t(t)*v1(t) $
+q3t(t)*v2(t))/sqrt(fv2(t)*fv3(t))
print 10 2275 resi1 resi2 resi3
print 10 2275 rho21 rho31 rho32
print 10 2275 fv1 fv2 fv3
Exercises
10.1 Consider the monthly simple returns, including dividends, of IBM stock, Hewlett-Packard (HPQ) stock, and the S&P composite index from January 1962 to December 2008 for 564 observations. The returns are in the file m-ibmhpqsp6208.txt. Transform into log returns in percentages. Use the exponentially weighted moving-average method to obtain a multivariate volatility series for the three return series. What is the estimated λ? Plot the three volatility series.
10.2 Focus on the monthly log returns of IBM and HPQ stocks from January 1962 to December 2008. Fit a DVEC(1,1) model to the bivariate return series. Is the model adequate? Plot the fitted volatility series and the time-varying correlations.
10.3 Focus on the monthly log returns of the S&P composite index and HPQ stock. Build a BEKK model for the bivariate series. What is the fitted model? Plot the fitted volatility series and the time-varying correlations.
10.4 Build a constant-correlation volatility model for the three monthly log returns of IBM stock, HPQ stock, and S&P composite index. Write down the fitted model. Is the model adequate? Why?
10.5 The file m-geibmsp2608.txt contains the monthly simple returns of General Electric stock, IBM stock, and the S&P composite index from January 1926 to December 2008. The returns include dividends. Transform into log returns in percentages. Focus on the monthly log returns in percentages of GE stock and the S&P 500 index. Build a constant-correlation GARCH model for the bivariate series. Check the adequacy of the fitted model, and obtain the 1-step-ahead forecast of the covariance matrix at the forecast origin December 2008.
10.6 Again, consider the monthly log returns of GE, IBM, and S&P composite index from January 1926 to December 2008. Build a dynamic correlation model for the three-dimensional series. For simplicity, use the sample correlation matrix for 
in Eq. (10.32).
10.7 The file m-spibmge.txt contains the monthly log returns in percentages of the S&P composite index, IBM stock, and GE stock from January 1926 to December 1999. Focus on GE stock and the S&P 500 index. Build a time-varying correlation GARCH model for the bivariate series using a logistic function for the correlation coefficient. Check the adequacy of the fitted model, and obtain the 1-step-ahead forecast of the covariance matrix at the forecast origin December 1999.
10.8 Focus on the monthly log returns in percentages of GE stock and the S&P 500 index from January 1926 to December 1999. Build a time-varying correlation GARCH model for the bivariate series using the Cholesky decomposition. Check the adequacy of the fitted model, and obtain the 1-step-ahead forecast of the covariance matrix at the forecast origin December 1999. Compare the model with the other model built in the previous exercise.
10.9 Consider the three-dimensional return series of the previous exercise jointly. Build a multivariate time-varying volatility model for the data, using the Cholesky decomposition. Discuss the implications of the model and compute the 1-step-ahead volatility forecast at the forecast origin t = 888.
10.10 An investor is interested in daily value at risk of his position on holding long $0.5 million of Dell stock and $1 million of Cisco Systems stock. Use 5% critical values and the daily log returns from February 20, 1990, to December 31, 1999, to do the calculation. The data are in the file d-dellcsco9099.txt. Apply the three approaches to volatility modeling in Section 10.7 and compare the results.
Bauwens, L., Laurent, S., and Rombouts, J. V. K. (2004). Multivariate GARCH models: A survey. Journal of Applied Econometrics 21: 79–109.
Bollerslev, T. (1990). Modeling the coherence in short-term nominal exchange rates: A multivariate generalized ARCH approach. Review of Economics and Statistics 72: 498–505.
Bollerslev, T., Engle, R. F., and Wooldridge, J. M. (1988). A capital-asset pricing model with time-varying covariances. Journal of Political Economy 96: 116–131.
Chib, S., Nardari, F. and Shephard, N. (1999). Analysis of high dimensional multivariate stochastic volatility models. Working paper, Washington University, St. Louis.
Engle, R. F. (2002). Dynamic conditional correlation: A simple class of multivariate GARCH models. Journal of Business and Economic Statistics 20: 339–350.
Engle, R. F. and Kroner, K. F. (1995). Multivariate simultaneous generalized ARCH. Econometric Theory 11: 122–150.
Harvey, A, Ruiz, E., and Shephard, N. (1994). Multivariate stochastic variance models. Review of Economic Studies 62: 247–264.
Li, W. K. (2004). Diagnostic Checks in Time Series. Chapman and Hall, London.
Mardia, K. V., Kent, J. T., and Bibby, J. M. (1979). Multivariate Analysis. Academic, New York.
Pourahmadi, M. (1999). Joint mean-covariance models with applications to longitudinal data: Unconstrained parameterization. Biometrika 86: 677–690.
Tsay, R. S. (2006). Multivariate volatility models in H. C. Ho, C. K. Ing, and T. L. Lai (eds.). Time Series and Related Topics in memory of C.Z. Wei, Lecture Notes Monograph Series, Institute of Mathematical Statistics. Beachwood, Ohio.
Tse, Y. K. (2000). A test for constant correlations in a multivariate GARCH model. Journal of Econometrics 98: 107–127.
Tse, Y. K. and Tsui, A. K. C. (2002). A multivariate GARCH model with time-varying correlations. Journal of Business & Economic Statistics 20: 351–362.