Advanced Time Series Analysis
Sudiksha Joshi
03/01/2019
Principal Component Analysis
PCA reduces p-dimension dataset to an m-dimension dataset where p > m. It describes the orginal data using fewer variables or dimensions than initially measured. We project the original time_ser_diff data onto a new, orthogonal basis. This removes multicollinearity. R calculates PCA using singular value decomposition of the scaled and centered matrix data (rather than eigen on covarince matrix) as it provides numerical accuracy.
## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7
## Standard deviation 2.4645 1.3724 1.0866 0.65482 0.49411 0.36791 0.1731
## Proportion of Variance 0.6074 0.1883 0.1181 0.04288 0.02441 0.01354 0.0030
## Cumulative Proportion 0.6074 0.7957 0.9138 0.95669 0.98111 0.99464 0.9976
## PC8 PC9 PC10
## Standard deviation 0.13519 0.06179 0.03869
## Proportion of Variance 0.00183 0.00038 0.00015
## Cumulative Proportion 0.99947 0.99985 1.00000Proportion of variance shows how much of the variance is explained by that principal component. The components are always sorted by how important they are, so the most important components will always be the first few.
## PC1 PC2 PC3 PC4 PC5 PC6
## [1,] -2.532253 0.9589985 -0.6856539 -0.4194933 0.6971565 -0.9574271
## [2,] -2.417471 1.1089015 -0.9856733 -0.2229814 0.6730305 -0.8642424
## [3,] -2.397437 1.0653462 -1.1664993 -0.3159122 0.5733262 -0.7703319
## [4,] -2.369972 0.9977712 -1.1642373 -0.3433385 0.6511079 -0.7066144
## [5,] -2.267717 1.0412444 -1.4965310 -0.3725687 0.5711167 -0.6006079
## [6,] -2.265047 1.0854381 -1.4703869 -0.5285890 0.5150894 -0.6393254## [1] "matrix"Advanced Forecasting Methods
Cointegration
When the trends and patterns of two series are similar, then they are cointegrated. The cointegration test measures whether the residuals from a regression are stationary. Stationary residuals are cointegrated. Therefore, the fact that time series are correlated is statistically significant, and not due to some chance.It is also a Dickey-Fuler stationarity test on residuals where the null hypothesis is that the series are not cointegrated. For a stationary test, we should reject the null hypothesis of no cointegrated.
The concept of cointegrated time series arises from the idea that housing prices, securities’ prices, interest rates and other economic indicators return to their long-term average levels after significant movements in short terms. Besides the imbalance in the demand and supply of houses, prices revert to their means as housing pries are highly correlated with inflation. Further, inflation rates are highly correlated with wages or real disposable income.
Given two series x(t) and y(t), R will search for paramteres α, β, and ρ such that
y(t) = α + β * x(t) + r(t) r(t) = ρ * r(t−1) + ϵ(t)
where r(t) = residual and, ϵ(t) = series of idependently and identically distributed (i.i.d) innovations with mean = 0
If |ρ| < 1, then x(t) and y(t) are cointegrated (i.e., r(t) doesn’t contain a unit root). if |ρ| = 1, then the residual series R[t] has a unit root and follows a random walk.
## Y[i] = 0.9736 X[i] + 66.8147 + R[i], R[i] = 0.7482 R[i-1] + eps[i], eps ~ N(0,125.2861^2)
## (0.0215) (69.6048) (0.0314)
##
## R[511] = -864.9835 (t = -4.788)
##
## WARNING: The series seem cointegrated but the residuals are not AR(1).
##
## Unit Root Tests of Residuals
## Statistic p-value
## Augmented Dickey Fuller (ADF) -4.276 0.00583
## Phillips-Perron (PP) -133.946 0.00010
## Pantula, Gonzales-Farias and Fuller (PGFF) 0.736 0.00010
## Elliott, Rothenberg and Stock DF-GLS (ERSD) -4.261 0.00014
## Johansen's Trace Test (JOT) -79.007 0.00010
## Schmidt and Phillips Rho (SPR) -86.676 0.00010
##
## Variances
## SD(diff(X)) = 95.289345
## SD(diff(Y)) = 112.092429
## SD(diff(residuals)) = 133.568139
## SD(residuals) = 180.670994
## SD(innovations) = 125.286102
##
## Half life = 2.389727
## R[last] = -864.983531 (t=-4.79)
Phillips-Ouliaris Test
This test is different from Augmented Dickey Fuller and Phillips-Perron unit root tests.It measures the evidence of coinintegration in the residuals of two time series. In this case, I have regressed housing starts on private houses completed series.
We cannot use ADF on residuals as they are devoid of Dickey-Fuller distributions in which the null hypothesis is that cointegration is absent. Alternatively, the residuals have Phillips-Ouliaris distributions.
##
## Phillips-Ouliaris Cointegration Test
##
## data: time_ser_diff[, c(1, 8)]
## Phillips-Ouliaris demeaned = -126.95, Truncation lag parameter =
## 5, p-value = 0.01H0: the 2 Series are not cointegrated Ha: the 2 Series are cointergated
The PO test rejects the null of no cointegration at the 5 percent level.The series are cointegrated. With cointegrated series we can construct a VEC model to better understand the causal relationship between the two variables.
Error Correction Model (ECM)
http://blog.mindymallory.com/2018/02/basic-time-series-analysis-the-var-model-explained/
The variables with cointegration I(0) have a short term relationship as opposed to those variables with cointegration I(1), as the latter have long term relationship. Both short and long run effects are present in the short run error correction model. The first equation is the ARDL(1,1) model that presumes a long run or steady state association between x and y. When deriving the error correction model, we can add more lagged differences of the regressor (x variable) to remove serial correlation.
y_t = δ + θ_1* y_(t−1) + δ_0 * x_t + δ_1 * x_(t−1) + ν_t,
Δy_t = −α * [y_(t−1) − β_1 − β_2 * x_(t−1)] + δ_0 * Δx_t + ν_t
We can construct the ECM for US housing starts and housing supply as:
Δb_t = −α * [b_(t−1) − β_1 − β_2 * f_(t−1)] + (δ_0 * Δf_t) + [δ_1 * Δf_(t−1)] + ν_t (estimated using codes below)
##
## Call:
## lm(formula = dy ~ dX + ECM - 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -614.61 -56.46 1.86 58.13 361.99
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## dX 0.24530 0.04710 5.208 2.77e-07 ***
## ECM -0.25682 0.02471 10.395 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 100.6 on 508 degrees of freedom
## Multiple R-squared: 0.1965, Adjusted R-squared: 0.1934
## F-statistic: 62.13 on 2 and 508 DF, p-value: < 2.2e-16Johansen Test
This test measures if three or more time series are cointegrated. Then, we take a linear combination of underlying series to form a stationary series. VAR(p) without drift is of the form:
x_t = μ + A_1 * x_(t−1) + … + A_p * x_(t−p) + w_t
μ = vector-valued mean of the series, A_i = coefficient matrices for each lag, w_t = multivariate Gaussian noise term with mean zero.
By differencing the series, we can form a Vector Error Correction model (VECM):
Δx_t = μ + A * x_(t−1) + Γ_1 * Δx_(t−1) +…+ Γ_p * Δx_(t−p) + w_t
Δx_t = x_t − x_(t−1) : differencing operator, A = coefficient matrix for the first lag, Γ_i = matrices for each differenced lag.
When the matrix A=0, the series are not cointegrated.
We perform an eigenvalue decomposition of A. r is the rank of the matrix A and the Johansen test checks if r = 0 or 1.
r=n−1, where n is the number of time series under test.
H0: r=0 means implies that no cointegration is present. When rank r > 0, there is a cointegrating relationship between at least two time series.
The eigenvalue decomposition outputs a set of eigenvectors. The components of the largest eigenvector is used in formulating the coefficients of the linear combination of time series. This creates stationarity. We should run the Johansen Test of Cointegration for variables which are I(1) before running ECM. If series are not cointegrated, we don’t have to perform ECM.
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Test type: trace statistic , with linear trend
##
## Eigenvalues (lambda):
## [1] 0.150633389 0.071431035 0.023724043 0.001267544
##
## Values of teststatistic and critical values of test:
##
## test 10pct 5pct 1pct
## r <= 3 | 0.65 6.50 8.18 11.65
## r <= 2 | 12.87 15.66 17.95 23.52
## r <= 1 | 50.59 28.71 31.52 37.22
## r = 0 | 133.69 45.23 48.28 55.43
##
## Eigenvectors, normalised to first column:
## (These are the cointegration relations)
##
## hous_st.l2 pvt_house_comp.l2 mortgR.l2 house_supply.l2
## hous_st.l2 1.0000000 1.0000000 1.000000 1.0000000
## pvt_house_comp.l2 -0.9398363 -0.3594574 -7.405363 -0.3559075
## mortgR.l2 -8.4124306 -48.4739319 372.009181 663.4887257
## house_supply.l2 47.8816699 283.7432012 259.714967 -316.2896016
##
## Weights W:
## (This is the loading matrix)
##
## hous_st.l2 pvt_house_comp.l2 mortgR.l2
## hous_st.d -0.0717794227 -7.761341e-02 2.035792e-03
## pvt_house_comp.d 0.2512650802 -2.623877e-02 1.253251e-03
## mortgR.d 0.0001351943 -3.643818e-06 -1.055053e-05
## house_supply.d 0.0001879928 -7.431165e-05 -2.877323e-05
## house_supply.l2
## hous_st.d -2.918249e-04
## pvt_house_comp.d -2.566033e-05
## mortgR.d -3.237657e-06
## house_supply.d 2.725327e-06Johansen procedure as implemented in function ca.jo will help you find the number of cointegrating vectors. Take the output of ca.jo, start with 𝑟=0 and see if you can reject the null hypothesis of 𝑟=0 using the test statistic and the critical values reported in the output. If you reject, move to 𝑟=1 and upwards until you cannot reject. The first rank that you cannot reject is the number of cointegrating vectors. If you can reject all of them, all of your series appear to be stationary.
In Johansen cointegration test, the null hypothesis for the eigenvalue test is that there are 𝑟+1 cointegration relations. The largest eigenvalue generated by the test is 0.15306696.
Next, the output shows the trace test statistic for the three hypotheses of r = 0 to r ≤ 4.From r = 0 to r ≤ 2, the test statistic exceeds the 0.05 significance level. For instance, when r = 0, 208.64 > 78.87. Similarly, in the second test we test the null hypothesis for r ≤ 1 against the alternative hypothesis of r > 1. As 124.08 > 55.43, we reject r ≤ 1, i.e. the null hypothesis of no cointegration. However, when r ≤ 6, we fail to reject the null as 8.65 < 23.52. Thus, the matrix’ rank is 3 and the series will become stationary after using a linear combination of 3 time series.
We can make a linear combination by using components of eigenvectors associated with the largets eigenvalue of 0.14624891. Correspondingly, we use vectors under the column hous_st.l2 which are (1.000000,-1.0098454639 , -2.7489510983) to obtain a stationary series.
The linear series is: 1 * hous_st -1.0098454639 * pvt_house_comp -2.7489510983 * mortgR
##
## Augmented Dickey-Fuller Test
##
## data: linear_series
## Dickey-Fuller = -4.2466, Lag order = 7, p-value = 0.01
## alternative hypothesis: stationaryThe p-value in the Dickey-Fuller test is 0.01 < 0.05. So, we reject the null hypothesis of unit root and conclude that the series formed from the linear combination is stationary.
Vector Error Correction Model (VECM)
VAR is used to capture short-run relationship between the variables employed (example, where there is a shock) while VECM tests for the long-run relationship. Through VECM concept of error correction, we understand how deviations from the long-run are “corrected”.
The coefficient on the ECT in an equation of the VECM quantifies the impact of the error correction term on the particular dependent variable (just as any regression coefficient in a linear model). The sign shows whether there is “error correction” (so that the variable corrects towards equilibrium) or “error inflation” (just made this term up; so that the variable deviates further from the equilibrium).
If some of the variables have coefficients that are indistinguishable from zero, then those are “dominant” in the sense that they do not adjust towards the equilibrium; rather, they “drive” the system of variables. But there must be some other variables in the system that do adjust, and these ones are “dominated”. (If none of the variables adjusted towards the equilibrium, there would be no cointegration.)
On the other hand, the coefficients on the variables inside the error correction term describe the long-run equilibrium relationship.
## #############
## ###Model VECM
## #############
## Full sample size: 511 End sample size: 507
## Number of variables: 4 Number of estimated slope parameters 60
## AIC 6905.044 BIC 7175.669 SSR 7650295
## Cointegrating vector (estimated by ML):
## hous_st pvt_house_comp mortgR house_supply
## r1 1 0 -78.37923 316.6825
## r2 0 1 -75.61939 301.1913
##
##
## ECT1 ECT2
## Equation hous_st -0.0268(0.0425) -0.0363(0.0375)
## Equation pvt_house_comp 0.3199(0.0347)*** -0.3288(0.0307)***
## Equation mortgR 0.0003(0.0001)* -0.0002(9.8e-05).
## Equation house_supply 0.0003(0.0002) -0.0002(0.0002)
## Intercept hous_st -1
## Equation hous_st 163.2674(34.7598)*** -0.5518(0.0587)***
## Equation pvt_house_comp -9.3985(28.4199) -0.2430(0.0480)***
## Equation mortgR -0.2161(0.0905)* 1.7e-06(0.0002)
## Equation house_supply -0.2639(0.1756) -0.0003(0.0003)
## pvt_house_comp -1 mortgR -1
## Equation hous_st 0.0496(0.0596) -43.9864(18.0267)*
## Equation pvt_house_comp -0.4583(0.0487)*** 2.5715(14.7387)
## Equation mortgR 0.0002(0.0002) 0.5653(0.0469)***
## Equation house_supply 0.0008(0.0003)** 0.6693(0.0911)***
## house_supply -1 hous_st -2
## Equation hous_st -13.1016(10.4502) -0.3224(0.0576)***
## Equation pvt_house_comp -4.5648(8.5442) -0.2111(0.0471)***
## Equation mortgR -0.0506(0.0272). -0.0001(0.0001)
## Equation house_supply -0.3478(0.0528)*** 1.2e-05(0.0003)
## pvt_house_comp -2 mortgR -2
## Equation hous_st -0.0044(0.0663) -11.1202(19.6690)
## Equation pvt_house_comp -0.2872(0.0542)*** 7.4237(16.0815)
## Equation mortgR 0.0002(0.0002) -0.3211(0.0512)***
## Equation house_supply 0.0007(0.0003). -0.0641(0.0993)
## house_supply -2 hous_st -3
## Equation hous_st -21.3106(10.5554)* -0.1176(0.0489)*
## Equation pvt_house_comp -5.8704(8.6302) -0.1207(0.0400)**
## Equation mortgR -0.0625(0.0275)* -0.0002(0.0001).
## Equation house_supply -0.2075(0.0533)*** -5.1e-05(0.0002)
## pvt_house_comp -3 mortgR -3
## Equation hous_st -0.0131(0.0597) -42.0914(18.5365)*
## Equation pvt_house_comp -0.1894(0.0488)*** -1.7610(15.1556)
## Equation mortgR 0.0001(0.0002) 0.1238(0.0482)*
## Equation house_supply 0.0002(0.0003) 0.1072(0.0936)
## house_supply -3
## Equation hous_st -14.1123(9.9092)
## Equation pvt_house_comp -4.0074(8.1018)
## Equation mortgR -0.0483(0.0258).
## Equation house_supply 0.0047(0.0501)##
## % Error: Unrecognized object type.Vector Autoregression (VAR)
We treat all variables symmetrically in VAR i.e we model them in such a way that these endogenous variables equally impact each other.
As a generalization of the univariate autoregressive model, it forecasts a vector of time series.The system has one equation per variable. The right hand side of each equation has lags and a constant of all the variables.
For a stationary series, we can directly fit VAR to the data and forecast. This is called “VAR in levels”. Otherwise, we difference the non-stationary data first and then fit the model. The resulting model is called “VAR in differences.” Using leveled variables (which are stationary) in VAR models can result in spurious regression. But, differenced variables will remedy the problem. In both instances, we use the concept of least squares to estimate the model.
Moreover, a non-stationary series could be cointegrated. This implies that there is a linear combination of variables which is stationary. In this scenario, we should make a vector error correction model i.e. a VAR model with an error correction mechanism
The VAR model can be used when the variables under study are I(1) but not cointegrated.
Δy_t = [β_11 * Δy_(t−1)] + [β_12 * Δx_(t−1)] + ν Δx_t = [β_21 * Δy_(t−1)] + [β_22 * Δx_(t−1)] + ν
## AIC(n) HQ(n) SC(n) FPE(n)
## 5 3 3 5The R output shows the lag length selected by each of the information criteria available in the vars package using levelled and differenced variables, separately. From the above results, we choose p= 3 as a lag parameter as minimizes AIC, HQ, SC and FPE. We construct multivariate order 3 VAR model.
R estimates VAR using OLS equation where the model is of the form: y_t = A_1y_(t-1) + …. + A_p y_(t-p) + CD_t + u_t
where y_t is a K * 1 vector of endogenous variables and u_t assigns a spherical disturbance term of the same dimension. The coefficient matrices A_1……Ap are of dimension K * K.
##
## ============================================================================
## Dependent variable:
## ---------------------------------------------
## Hous
## (1) (2) (3) (4)
## ----------------------------------------------------------------------------
## hous_st.l1 0.443*** 0.097** 0.0002** 0.00003
## (0.046) (0.038) (0.0001) (0.0002)
##
## pvt_house_comp.l1 0.040 0.283*** 0.0001 0.001*
## (0.064) (0.053) (0.0002) (0.0003)
##
## mortgR.l1 -32.268* 0.294 1.512*** 0.639***
## (17.514) (14.556) (0.045) (0.087)
##
## house_supply.l1 -29.253*** -3.226 -0.018 0.661***
## (9.570) (7.954) (0.025) (0.047)
##
## hous_st.l2 0.222*** 0.044 -0.0001 0.0003
## (0.050) (0.041) (0.0001) (0.0002)
##
## pvt_house_comp.l2 -0.039 0.250*** 0.00001 -0.0002
## (0.063) (0.053) (0.0002) (0.0003)
##
## mortgR.l2 -0.614 6.542 -0.802*** -0.639***
## (28.479) (23.670) (0.074) (0.141)
##
## house_supply.l2 -7.547 1.884 -0.003 0.139**
## (11.436) (9.505) (0.030) (0.056)
##
## hous_st.l3 0.227*** 0.137*** -0.0001 -0.0001
## (0.048) (0.040) (0.0001) (0.0002)
##
## pvt_house_comp.l3 0.025 0.183*** -0.0001 -0.0005
## (0.060) (0.050) (0.0002) (0.0003)
##
## mortgR.l3 30.923* -5.984 0.274*** 0.013
## (17.955) (14.923) (0.047) (0.089)
##
## house_supply.l3 13.121 5.803 0.032 0.215***
## (9.800) (8.146) (0.025) (0.048)
##
## const 321.561*** -41.733 0.022 -0.805***
## (62.661) (52.080) (0.162) (0.309)
##
## trend -0.193*** 0.034 -0.0003 0.001**
## (0.070) (0.058) (0.0002) (0.0003)
##
## ----------------------------------------------------------------------------
## Observations 508 508 508 508
## R2 0.944 0.954 0.995 0.922
## Adjusted R2 0.943 0.953 0.994 0.919
## Residual Std. Error (df = 494) 97.370 80.928 0.252 0.481
## F Statistic (df = 13; 494) 642.237*** 794.994*** 6,916.869*** 446.179***
## ============================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01##
## Portmanteau Test (asymptotic)
##
## data: Residuals of VAR object var1
## Chi-squared = 133.79, df = 112, p-value = 0.07859H0: no serial correlation H1: serial correlation is present
The residuals for this model are not serially correlation as the p-value is very small.
Granger Causality
We use F-test on the lags of other variables to implement the granger causality. It tests whether the lags of variables are useful in forecasting housing starts and vice versa.
For instance, fed_fundsR can granger cause hous_st if hous_st can be more accurately predicted by the lagged values of both hous_st and fed_fundsR, rather than the lagged values of hous_st alone. Thus, the granger causality test examines if lagged values of a variable can enhance the forecasts of another variable.
## # A tibble: 4 x 2
## statistic[,1] p.value[,1]
## <dbl> <dbl>
## 1 7.89 1.60e-11
## 2 0.738 6.74e- 1
## 3 7.24 2.07e-10
## 4 5.21 4.98e- 7##
## Granger causality H0: hous_st do not Granger-cause pvt_house_comp
## mortgR house_supply
##
## data: VAR object var1
## F-Test = 7.8908, df1 = 9, df2 = 1976, p-value = 1.597e-11##
## Granger causality H0: pvt_house_comp do not Granger-cause hous_st
## mortgR house_supply
##
## data: VAR object var1
## F-Test = 0.73808, df1 = 9, df2 = 1976, p-value = 0.6742##
## Granger causality H0: mortgR do not Granger-cause hous_st
## pvt_house_comp house_supply
##
## data: VAR object var1
## F-Test = 7.2373, df1 = 9, df2 = 1976, p-value = 2.065e-10##
## Granger causality H0: house_supply do not Granger-cause hous_st
## pvt_house_comp mortgR
##
## data: VAR object var1
## F-Test = 5.2125, df1 = 9, df2 = 1976, p-value = 4.975e-07Three of the four results have sufficiently small p-value and they indicates that we can reject null hypothesis: they Granger-cause others. The regressors have information to predict today’s housing starts.
Alternatively, we fail to reject the null hypothesis that pvt_house_comp do not Granger-cause hous_st, mortgR and house_supply. This means private houses completed do not play much role in predicting today’s housing starts, mortgage rate and housing supply
Impulse Response Function
In IRF, we shock one variable, say income, and propagate it through the fitted VAR model for a number of periods. We can trace this through the VAR model and see if it impacts the other variables in a statistically significant way.
An impulse (shock) to housing starts at time zero has large effects the next period, and the effects enlarge over time. The dotted lines show the 95 percent interval estimates of these effects. The VAR function prints the values corresponding to the impulse response graphs.
Forecast Error Variance Decomposition (FEVD)
Using the VAR model, a Forecast Error Variance Decomposition examines the impact of variables on one another. We use the forecast errors of each equation in the fitted VAR model to compute FEVD. Then, the fitted VAR model determines how much of each error realization is coming from unexpected changes (forecast errors) in the other variable.
Variance decomposition helps to interpret the VAR model. We can determine the amount of variation in the dependent variable is explained by each each independent variable. FEVD explains how a future shock in a time series changes future uncertainity in the other time series in the system. This process evolves over time, so a shock on a time series may not be important in the short run, but may be very significant in the long run.
In the first plot, we see the FEVD for housing starts. It appears that although we were borderline on whether or not to conclude that federal funds rate Granger cause housing starts, the FEVD reveals that the magnitude of the causality is tiny anyway, while that of income is greater on housing starts.
In the second plot, we see the FEVD for income. It appears that although we were borderline on whether or not to conclude that housing starts and federal funds rate Granger cause income.
The ARCH-LM test with q lags checks for the presence of ARCH effects at lags 1 to q. It tests if the coefficients α_1,…. α_q in the equation below:
x^2_t = α_0 + α_1 * x^2_(t-1) +….+ α_q * x^2_(t-q) + ϵ_t
##
## ARCH (multivariate)
##
## data: Residuals of VAR object var1
## Chi-squared = 598.88, df = 200, p-value < 2.2e-16When q = 2, we test for ARCH effects jointly at lags 1 and 2. H0 = α_1 = α_2 = 0 As the p-value is very small, we reject the null hypothesis and conclude that ARCH effects are present at lags 1 and 2 jointly. ARCH effects are also present at higher lag orders, implying that the data is conditionally heteroskedastic.
Generalized Autoregressive Conditional Heteroskedasticity (GARCH)
Generalized Autoregressive Conditional Heteroskedastic, or GARCH models are useful to analyse and forecast volatility in a time series data. Univariate GARCH(1,1) helps in modeling volality and its clustering.
Financial time series possess the property of volatility clustering wherein the volatility of the variable changes over time. Technically, this behavior is called conditional heteroskedasticity. Because ARMA models don’t consider volatility clustering i.e. they are not conditionally heteroskedastic, so we need to use ARCH and GARCH models for predictions.
Such models include the Autogressive Conditional Heteroskedastic (ARCH) model and Generalised Autogressive Conditional Heteroskedastic (GARCH) model. Different forms of volatility such as sell-offs during a financial crises, can cause serially correlated heteroskedasticity. Thus, the time_ser data is conditionally heteroskedastic.
Maximum likelihood estimates most GARCH models, such as measuring relative loss or profit from trading stocks in a day. If x_t is the value of housing starts on t, then r_t=[x_t − x_(t−1)]/x_(t−1) is called the return. We observe large volatility around the 2008 financial crisis and returns that are mostly noise noise with short periods of large variability.
Below I have plotted the ACF of the returns of housing starts. There are two bounds plotted on the graph. The straight red line represents the standard bounds under the strong white noise assumption. The second line is under the hypothesis that the process is GARCH.
## ChiSq DF pvalue
## [1,] 55.53981 5 1.010765e-10
## [2,] 60.73615 10 2.629113e-09
## [3,] 86.60523 20 2.889828e-10
## attr(,"method")
## [1] "LjungBox"We reject the hypothesis that the series is independently and identically distributed from the Ljung-Box test.The series is not white noise, hence autocorrelated.
From the plot above, several autocorrelations seem significant under hypothesis of both iid and GARCH process.
Now, I have fit a GARCH-type model which assumes the null hypothesis that the returns are GARCH.
## h Q pval
## [1,] 5 40.57540 1.143075e-07
## [2,] 10 42.76489 5.478140e-06
## [3,] 15 52.29641 5.045273e-06The low p-values give reason to reject the hypothesis that the returns are a GARCH white noise process. So, we should do ARMA modelling.
We have fit GARCH model(s), starting with a GARCH(1,1) model with Gaussian innovations.GARCH(1,1) considers a single autoregressive and a moving average lag. The model is:
ϵ_t = σ_t * w_t σ^2 = α_0 + α_1 * ϵ^2_(t−1) + β_1 * σ^2_(t−1)
Note that alpha_1 + beta_1 < 0, otherwise the series will become unstable.
The persistence of a GARCH model signifies the rate at which large volatilities decay after a shock. The key statistic in GARCH(1,1) is the sum of two parameters: alpha1 and beta1.
Ideally, alpha_1 + beta_1 < 1. If, alpha_1 + beta_1 > 1, then the volatility predictions are explosive. If, alpha_1 + beta_1 = 1, then the model has exponential decay.
In the output from garchFit, the normalized log-likelihood is the loglikelihood divided by n. The AIC and BIC values have also been normalized by dividing by n,
##
## Title:
## GARCH Modelling
##
## Call:
## garchFit(formula = ~arma(1, 0) + garch(1, 1), data = hous_st_return,
## cond.dist = "norm")
##
## Mean and Variance Equation:
## data ~ arma(1, 0) + garch(1, 1)
## <environment: 0x561ea51a1490>
## [data = hous_st_return]
##
## Conditional Distribution:
## norm
##
## Coefficient(s):
## mu ar1 omega alpha1 beta1
## 4.1271e-06 -3.6658e-01 6.5196e-05 6.8848e-02 9.2174e-01
##
## Std. Errors:
## based on Hessian
##
## Error Analysis:
## Estimate Std. Error t value Pr(>|t|)
## mu 4.127e-06 2.994e-03 0.001 0.99890
## ar1 -3.666e-01 4.253e-02 -8.619 < 2e-16 ***
## omega 6.520e-05 5.193e-05 1.256 0.20929
## alpha1 6.885e-02 2.182e-02 3.155 0.00161 **
## beta1 9.217e-01 2.286e-02 40.319 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Log Likelihood:
## 599.6187 normalized: 1.178033
##
## Description:
## Sat Mar 30 22:15:56 2019 by user:
##
##
## Standardised Residuals Tests:
## Statistic p-Value
## Jarque-Bera Test R Chi^2 18.01514 0.0001224791
## Shapiro-Wilk Test R W 0.9891823 0.0008269949
## Ljung-Box Test R Q(10) 16.03317 0.09868677
## Ljung-Box Test R Q(15) 27.71855 0.02339867
## Ljung-Box Test R Q(20) 34.76176 0.02141044
## Ljung-Box Test R^2 Q(10) 26.98303 0.002620484
## Ljung-Box Test R^2 Q(15) 43.25259 0.0001438445
## Ljung-Box Test R^2 Q(20) 46.94882 0.0005962632
## LM Arch Test R TR^2 26.92906 0.007910933
##
## Information Criterion Statistics:
## AIC BIC SIC HQIC
## -2.336419 -2.294843 -2.336610 -2.320117The diagnostics imply that the standardised residuals and their squares are IID and that the model accomodates ARCH effects.
H_0: white noise innovation process is Gaussian
Their distribution is Gaussian only from the p-value for Ljung-Box Test which is 0.921266. From all other tests of normality, we reject the null hypothesis as the p-values are very low.
The qq-plot of the standardised residuals, suggests that the fitted standardised skew-t conditional distribution is decent.
Since, ARIMA linearly models the data, the forecast width is constant as the model does not incorporate new information or recent changes.To model non-linearity or a cluster of volatility, we have to use ARCH/GARCH methods as they reflect more recent fluctuations in the series. The ACF and PACF of residuals can confirm if the residuals can be predicted if they are not white noise. Residuals of strict white noise series are i.i.d normally distributed with zero mean. Moreover, the PACF and ACF of squared residuals have no significant lags. Finally, we cannot predict a strict white noise series, either linearly or non-linearly. Below, the residuals and squared residuals of ARIMA(4,2,4) model show a cluster of volatility as shown from the ACF plots.
##
## Box-Ljung test
##
## data: resid
## X-squared = 6.1032, df = 10, p-value = 0.8065We fail to reject the null hypothesis that the residuals of ARIMA(4,2,4) is not serially correlated i.e. the series is white noise.
##
## Box-Ljung test
##
## data: resid^2
## X-squared = 63.898, df = 10, p-value = 6.583e-10Now, I have fit a GARCH(1,1) model.
##
## Title:
## GARCH Modelling
##
## Call:
## garchFit(formula = ~arma(4, 4) + garch(1, 1), data = time_ser_diff[,
## 1])
##
## Mean and Variance Equation:
## data ~ arma(4, 4) + garch(1, 1)
## <environment: 0x561ea076daf8>
## [data = time_ser_diff[, 1]]
##
## Conditional Distribution:
## norm
##
## Coefficient(s):
## mu ar1 ar2 ar3 ar4
## 16.053017 0.471082 0.579245 0.662428 -0.724449
## ma1 ma2 ma3 ma4 omega
## 0.083284 -0.252425 -0.649781 0.336848 2640.659451
## alpha1 beta1
## 0.424977 0.367700
##
## Std. Errors:
## based on Hessian
##
## Error Analysis:
## Estimate Std. Error t value Pr(>|t|)
## mu 16.05302 11.47482 1.399 0.161820
## ar1 0.47108 0.15926 2.958 0.003096 **
## ar2 0.57925 0.05429 10.670 < 2e-16 ***
## ar3 0.66243 0.05118 12.942 < 2e-16 ***
## ar4 -0.72445 0.14914 -4.858 1.19e-06 ***
## ma1 0.08328 0.15740 0.529 0.596715
## ma2 -0.25242 0.14761 -1.710 0.087245 .
## ma3 -0.64978 0.10502 -6.187 6.12e-10 ***
## ma4 0.33685 0.05913 5.697 1.22e-08 ***
## omega 2640.65945 773.37713 3.414 0.000639 ***
## alpha1 0.42498 0.09702 4.380 1.19e-05 ***
## beta1 0.36770 0.10534 3.491 0.000482 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Log Likelihood:
## -3050.358 normalized: -5.96939
##
## Description:
## Sat Mar 30 22:15:58 2019 by user:
##
##
## Standardised Residuals Tests:
## Statistic p-Value
## Jarque-Bera Test R Chi^2 10.53292 0.005161848
## Shapiro-Wilk Test R W 0.9928411 0.01531213
## Ljung-Box Test R Q(10) 4.61827 0.9151778
## Ljung-Box Test R Q(15) 15.15126 0.4405772
## Ljung-Box Test R Q(20) 19.6149 0.4822394
## Ljung-Box Test R^2 Q(10) 10.52377 0.3958027
## Ljung-Box Test R^2 Q(15) 39.47513 0.0005439962
## Ljung-Box Test R^2 Q(20) 40.64349 0.004138153
## LM Arch Test R TR^2 20.30752 0.06148775
##
## Information Criterion Statistics:
## AIC BIC SIC HQIC
## 11.98575 12.08523 11.98468 12.02475