Modeling and Forecasting by the Vector Autoregressive Moving Average Model for Export of Coal and Oil Data (Case Study from Indonesia over the Years 2002-2017)

The vector autoregressive moving average (VARMA) model is one of the statistical analyses frequently used in several studies of multivariate time series data in economy, finance, and business. It is used in numerous studies because of its simplicity. Moreover, the VARMA model can explain the dynamic behavior of the relationship among endogenous and exogenous variables or among endogenous variables. It can also explain the impact of a variable or a set of variables by means of the impulse response function and Granger causality. Furthermore, it can be used to predict and forecast time series data. In this study, we will discuss and develop the best model that describes the relationship between two vectors of time series data export of Coal and data export of Oil in Indonesia over the period 2002-2017. Some models will be applied to the data: VARMA (1,1), VARMA (2,1), VARMA (3,1), and VARMA (4,1). On the basis of the comparison of these models using information criteria AICC, HQC, AIC, and SBC, it was found that the best model is VARMA (2,1) with restriction on some parameters: AR2_1_2 = 0, AR2_2_1 = 0, and MA1_2_1 = 0. The dynamic behavior of the data is studied through Granger causality analysis. The forecasting of the series data is also presented for the next 12 months.


INTRODUCTION
Financial, business, and economic data are very often collected in equally spaced time intervals such as days, weeks, months, or years. In a number of cases, such time series data may be available on several related variables. There are some reasons for analyzing and modeling such time series jointly: (1) To understand the dynamic relationship among variables and (2) to improve the accuracy of forecast and knowledge of the dynamic structure so as to produce good forecast Pena and Tiao, 2001). The analysis of multiple time series has been developed by Tiao and Box (1981); since then, the development of the theory has been extensively discussed in the literature (Lutkepohl, 2005;Reinsel, 1993). Multivariate time series are of great interest in a variety of fields such as financial, economic, stock market, and earth science, e.g., meteorology (Reinsel, 1993). In multivariate time series analysis, not only the properties of the individual series but also the possible cross relationship among the time series data are discussed. The application of the vector autoregressive (VAR) model has been extensively discussed by Malik et al. (2017), Sharma et al. (2018), and Warsono et al. (2019).
In this study, we discuss and develop the best model that describes the relationship between two vectors of data export of Coal and data export of Oil in Indonesia over the period 2002-2017. On the basis of this objective, the VAR moving average (VARMA) model was developed to explain the relationship between the data export of Coal and Oil in Indonesia over the period 2002-2017. This Journal is licensed under a Creative Commons Attribution 4.0 International License Methods to find the best model, estimates of parameters, model checking, and forecasting of vector time series are also discussed.

STATISTICAL MODEL
The VARMA model is commonly used to forecast multivariate time series data and provides a simple framework to study the dynamic relationships among variables (Koreisha and Pukkila, 2004). The VARMA model is an extension of the ARMA model in univariate time series (Lutkepohl, 2005;Wei, 1990) and is used with the condition that the data have to be stationary over time (Lutkepohl, 2005). The VARMA (p,q) model is a combination of the VAR (p) model and the vector moving average (q) (VMA (q)) model. An m-dimensional time series datum Γ t is a VARMA (p,q) process if where δ is a constant m×1 vector of means, and δ T = (δ 1 ,δ 2 , δ m ), u u u u   where i = 1,2,…, p; j = 1,2,…, q.
Model (1) can also be written in a simpler form using the backshift operator B as follows:

and u t is vector innovation.
Some properties of the VARMA (p,q) model with p>0 and q>0 are discussed. The model is assumed to be identifiable and innovation u t has mean zero and covariance matrix ∑ u , which is positive definite; see Graybill (1969) for the definition of a positive-definite matrix. We shall assume that the zeros of the determinant polynomials |Φ(B)| and |Θ(B)| are on or outside the unit circle. The series {Γ t }will be stationary if the zeros of |Φ(B)| are on or outside the unit circle and will be invertible when those of |Θ(B)| are on or outside the unit circle Tsay, 2005;Reinsel, 1993). To find the best model, we estimated some candidate models (VARMA (1,1), VARMA (2,1), VARMA (3,1), and VARMA (4,1)) using some information criteria (AICC, HQC, AIC, and SBC). The selected best model and estimation of the parameters of the selected model were reviewed. If some parameters are not significant in the selected model, then they will be restricted to zero (Tsay, 2005;Milhoj, 2016) so that the final best model is simpler. The optimal l-step-ahead forecast of Γ t+l for model (1) is as follows (SAS/ETS 13.2, 2014;Lutkepohl, 2005):

DATA ANALYSIS
The data used in this study are the data export of Coal and Oil from Indonesia from January 2002 to December 2017. The data are from the Central Bureau of Statistics (BPS) Indonesia (BPS (a) 2017, and BPS (b), 2017). The plot of the data is given in Figure 1.  Figure 1 also shows that the data are nonstationary, and this is in line with the augmented Dicky-Fuller test given in Table 1. Now, we look at the ACF and PACF of data of Coal and Oil given in (Figure 2a and b). From the sample ACF of data of Coal and Oil in (Figure 2a and b), the tails cut off very slowly. This indicates that the time series data of Coal and Oil are not stationary. That is, the means or the variances of time series data of Coal and Oil are not constant over time.
To make the data stationary, differencing needs to be conducted, and the results of differencing with d = 1 are given in Table 2.
The assumption of stationarity is attained, and modeling of VARMA can be carried out.

VARMA (p,q) Modeling
To find the best model that fits the data, some VARMA (p,q) models (i.e., VARMA (1,1), VARMA (2,1), VARMA (3,1), and VARMA (4,1)) for prediction and forecasting were applied to the data. The selection of the best model was conducted using some information criteria (AICC, HQC, AIC, and SBC). The minimum values of these criteria indicate the best model.
From Table 3, we conclude the following: on the basis of the minimum values of HQC and SBC, the best model is VARMA  (1,1); on the basis of the minimum values of AICC and AIC, the best model is VARMA (2,1). Therefore, there are two candidates for the best model. To choose the best model, we check the schematic representation of the parameter estimates of VARMA (1,1) and VARMA (2,1) given in the Table 4. Schematic representation of parameter estimates of VARMA (1,1) and VARMA (2,1). Table 4 shows that in VARMA (1,1), four parameters are of significance, and in VARMA (2,1), six parameters are of significance. In this case, VARMA (2,1) is chosen as the best model to discuss the characteristics of data and for forecasting data. VARMA (2,1) is represented as follows: Oil , Φ 1 , Φ 2 and Ψ 1 are 2 × 2 matrix parameters for AR1, AR2, and MA1, respectively. ε t is vector white noise. The estimate model VARMA (2,1) is as follows: From Table 5, some parameters of the model AR2_1_2, AR2_2_1, and MA1_2_1 are not significantly different from zero. Therefore, we can improve the model by restricting those parameters that are not significant, as suggested by Tsay (2005) and Milhoj (2016). We restrict the parameters AR2_1_2 = 0, AR2_2_1 = 0, and MA1_2_1 = 0, and we conduct the test statistics for these restrictions. From the testing of restricted parameters equal to zero, we obtained the results as given in Table 6. Table 6 shows that all the tests are not significantly different from zero. By using these restriction parameters, the final model and the estimation of parameters are presented in Table 7.
The VARMA (2,1) model with restriction can also be written as two univariate regression models: The statistical test of the parameters in model (5) is given in Table 7, and models (6) and (7) are given in Table 8. On the basis of the statistical test, model (6) is very significant with the statistical test F = 15.30, and the P < 0.0001. The degree of determination of R-square is 0.3352. On the basis of the statistical test, model (7) is very significant with the statistical test F = 3.01, and the P = 0.0079. The degree of determination of R-square is 0.0903. Model (6) also explains that the export of Oil at lag 1 (t-1) has a positive effect on the export of Coal; the export of Coal at lag 1 (t-1) has a negative effect on the export of Coal, and the export of Coal at lag 2 (t-2) has a positive effect on the export of Coal.
•is between, VARMA: Vector autoregressive moving average Model (7) explains that the export of Coal at lag 1 (t-1) and the export of Oil at lag 1 and lag 2 (t-1 and t-2) have a negative effect on the export of Coal.
Granger causality is used to test two null hypotheses. Test 1 tests the null hypothesis where the export of Coal is influenced only by itself and not by the export of Oil. Test 2 tests the null hypothesis where the export of Oil is influenced only by itself and not by the export of Coal (SAS/ETS 13.2, 2014). From the results of the Granger causality tests, Table 9 demonstrates the following: (1) For Test 1, Chi-square = 11.88 and the P = 0.0078; hence, we reject the null hypothesis. We conclude that the export of Coal is not only influenced by itself but also influenced by the export of Oil.
(2) For Test 2, Chi-square = 3.74 and the P = 0.2912; hence, we cannot reject the null hypothesis. We conclude that the export of Oil is influenced only by itself and not by the export of Coal.
The patterns of the distribution of errors for data of Coal and Oil based on model VARMA (2,1) with restriction are given in Figure 3.
From (Figures 3a and 3b), the patterns of the distribution of errors for data of Coal and Oil are very close to normal distribution. If we look      (Figure 4a and b), for the data of Coal, it is clear that the prediction errors from 2002 to 2011 are homogeneous, but from 2012 to 2017 the errors fluctuate and are beyond two standard errors. This indicates that the prices are unstable in this horizon (2012)(2013)(2014)(2015)(2016)(2017).   The graphs of the model for data of Coal ( Figure 5a) and for data of Oil (Figure 6a) show that the data and their predictions are close to each other, which indicates that the models fit with the data. Table 10 shows that the forecasts for the data of Coal begin at 1435.317 for the first period and then decrease for the second to the ninth periods. Starting from the tenth period up to the twelfth period, the forecast increases and reaches a value of 1448.491. The forecast for data of Oil begins at 1466.245 for the first period and then fluctuates with the trend increase up to the twelfth period. In the twelfth period, the value 1486.427 is attained. The confidence interval of the prediction increases: it is smaller in the first period and larger over time up to the twelfth period. This indicates that even though the model is sound and fits with the data, if the model is used to forecast for long periods, the prediction becomes unstable. This is demonstrated by the large confidence interval. (Figures 5b and 6b) describe the behavior of the confidence interval over time from the first period up to the twelfth period.

CONCLUSION
In this study, the focus is on how to find the best model and use it for forecasting the data export of Coal and Oil of Indonesia over the years 2002-2017. We have developed the best model using the criteria AICC, HQC, AIC, and SBC, which fit the data. The best model is VARMA (2,1), with restriction on some parameters that are non significantly different from zero. The restricted parameters are AR2_1_2 = 0, AR2_2_1 = 0, and MA1_2_1 = 0. All the parameters in the model, AR and MA, are significant, except for the parameter constants. The model shows that the prediction and the real data fit well with each other.
The forecasting results show that the standard error increases over time; the standard error in the first month is relatively small compared with the prediction of the means, but increases over time up to forecasting for the next 12 months. This indicates that the model is sound when forecasting for short periods, but the results are unstable (because of the higher standard error) when forecasting for long periods.