Forecasting of financial instruments using VECM and ARIMA

Author

Nihad Alili

Published

September 6, 2024

Necessary packages and libraries

Code

getwd()

[1] "C:/Users/User/OneDrive/Desktop/Project"

Code

setwd("C:/Users/User/OneDrive/Desktop/Project") # Replace with the actual path
source("functions/testdf.R")



library(xts)
library(lmtest)
library(tidyverse)
library(readr)
library(ggplot2)
library(vars)
library(quantmod)
library(urca)
library(forecast)
library(tseries)

##Loading data and organizing it

Code

financial_data <- read.csv("C:/Users/User/OneDrive/Desktop/Project/TSA_2024_project_data_1 (2).csv", header = TRUE)


financial_data %>% glimpse

Rows: 300
Columns: 11
$ date <chr> "2021-06-01", "2021-06-02", "2021-06-03", "2021-06-04", "2021-06-…
$ y1   <dbl> 176.04, 173.37, 173.42, 178.85, 175.39, 173.89, 167.80, 162.47, 1…
$ y2   <dbl> 173.57, 173.57, 173.57, 174.52, 174.46, 174.54, 174.57, 175.37, 1…
$ y3   <dbl> 173.57, 173.57, 175.70, 177.89, 181.10, 184.44, 187.00, 188.86, 1…
$ y4   <dbl> 188.66, 188.46, 188.57, 188.78, 188.11, 187.61, 188.22, 188.68, 1…
$ y5   <dbl> 183.66, 183.52, 183.75, 182.79, 185.46, 190.18, 196.22, 201.14, 2…
$ y6   <dbl> 178.78, 181.10, 187.25, 190.82, 187.30, 195.26, 196.95, 207.64, 2…
$ y7   <dbl> 168.46, 168.41, 168.47, 170.08, 170.15, 170.33, 170.32, 171.42, 1…
$ y8   <dbl> 173.57, 173.57, 173.57, 173.57, 172.94, 171.78, 173.37, 174.22, 1…
$ y9   <dbl> 173.57, 173.57, 173.57, 173.18, 174.50, 176.99, 179.91, 182.35, 1…
$ y10  <dbl> 188.35, 188.37, 189.36, 190.46, 192.73, 194.14, 195.44, 196.21, 1…

Code

head(financial_data)

date	y1	y2	y3	y4	y5	y6	y7	y8	y9	y10
2021-06-01	176.04	173.57	173.57	188.66	183.66	178.78	168.46	173.57	173.57	188.35
2021-06-02	173.37	173.57	173.57	188.46	183.52	181.10	168.41	173.57	173.57	188.37
2021-06-03	173.42	173.57	175.70	188.57	183.75	187.25	168.47	173.57	173.57	189.36
2021-06-04	178.85	174.52	177.89	188.78	182.79	190.82	170.08	173.57	173.18	190.46
2021-06-05	175.39	174.46	181.10	188.11	185.46	187.30	170.15	172.94	174.50	192.73
2021-06-06	173.89	174.54	184.44	187.61	190.18	195.26	170.33	171.78	176.99	194.14

Code

tail(financial_data)

	date	y1	y2	y3	y4	y5	y6	y7	y8	y9	y10
295	2022-03-22	268.68	196.26	289.20	174.51	36.93	97.94	204.83	145.53	100.31	246.30
296	2022-03-23	266.55	196.34	288.83	174.58	36.87	87.23	205.03	145.76	100.22	245.85
297	2022-03-24	266.38	196.25	287.98	174.18	35.33	78.05	205.08	144.81	99.39	245.64
298	2022-03-25	272.44	197.49	287.54	174.47	34.93	85.77	206.84	144.90	99.22	245.77
299	2022-03-26	267.64	198.03	286.01	173.46	35.01	95.43	208.08	143.75	99.22	244.87
300	2022-03-27	264.95	198.61	285.43	174.40	33.51	95.25	208.76	145.46	98.49	244.08

Code

financial_data$date <- as.Date(financial_data$date, format = "%Y-%m-%d")
financial_data %>% glimpse()

Rows: 300
Columns: 11
$ date <date> 2021-06-01, 2021-06-02, 2021-06-03, 2021-06-04, 2021-06-05, 2021…
$ y1   <dbl> 176.04, 173.37, 173.42, 178.85, 175.39, 173.89, 167.80, 162.47, 1…
$ y2   <dbl> 173.57, 173.57, 173.57, 174.52, 174.46, 174.54, 174.57, 175.37, 1…
$ y3   <dbl> 173.57, 173.57, 175.70, 177.89, 181.10, 184.44, 187.00, 188.86, 1…
$ y4   <dbl> 188.66, 188.46, 188.57, 188.78, 188.11, 187.61, 188.22, 188.68, 1…
$ y5   <dbl> 183.66, 183.52, 183.75, 182.79, 185.46, 190.18, 196.22, 201.14, 2…
$ y6   <dbl> 178.78, 181.10, 187.25, 190.82, 187.30, 195.26, 196.95, 207.64, 2…
$ y7   <dbl> 168.46, 168.41, 168.47, 170.08, 170.15, 170.33, 170.32, 171.42, 1…
$ y8   <dbl> 173.57, 173.57, 173.57, 173.57, 172.94, 171.78, 173.37, 174.22, 1…
$ y9   <dbl> 173.57, 173.57, 173.57, 173.18, 174.50, 176.99, 179.91, 182.35, 1…
$ y10  <dbl> 188.35, 188.37, 189.36, 190.46, 192.73, 194.14, 195.44, 196.21, 1…

Code

financial_data <- xts(financial_data[, -1], order.by = financial_data$date)
head(financial_data)

               y1     y2     y3     y4     y5     y6     y7     y8     y9
2021-06-01 176.04 173.57 173.57 188.66 183.66 178.78 168.46 173.57 173.57
2021-06-02 173.37 173.57 173.57 188.46 183.52 181.10 168.41 173.57 173.57
2021-06-03 173.42 173.57 175.70 188.57 183.75 187.25 168.47 173.57 173.57
2021-06-04 178.85 174.52 177.89 188.78 182.79 190.82 170.08 173.57 173.18
2021-06-05 175.39 174.46 181.10 188.11 185.46 187.30 170.15 172.94 174.50
2021-06-06 173.89 174.54 184.44 187.61 190.18 195.26 170.33 171.78 176.99
              y10
2021-06-01 188.35
2021-06-02 188.37
2021-06-03 189.36
2021-06-04 190.46
2021-06-05 192.73
2021-06-06 194.14

#Visualisation the data using plots

Code

financial_data[, 1:10] %>% 
  as_tibble() %>%
  mutate(date = index(financial_data)) %>%
  pivot_longer(!date) %>%
  ggplot(aes(date, value, col = name)) +
  geom_line() +
  theme_bw() +
  labs(
    title = "Plotting all financial instruments in one graph",
  )

##Convert the xts object to a tibble and reshape it

Code

financial_data_long <- financial_data %>%
  as_tibble() %>%
  mutate(date = index(financial_data)) %>%
  pivot_longer(-date, names_to = "instrument", values_to = "value")

##Plot each instrument in its own facet

Code

ggplot(financial_data_long, aes(x = date, y = value, color = instrument)) +
  geom_line() +
  theme_bw() +
  labs(
    title = "Plot of Financial Instruments",
    x = "Date",
    y = "Value"
  ) +
  facet_wrap(~ instrument, scales = "free_y", ncol = 3) +
  theme(legend.position = "none")

From first glimpse, we can easily see that none of the plots are stationary and there are 4 cointegrated pairs. However, we will need to test statistically.

Checking stationary using ADF test

Code

testdf(variable = financial_data$y1, #ADF test variable 1
       max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-1.650473	0.4424618	3.2416352	0.0717890
1	-1.824610	0.3776076	0.0010079	0.9746732
2	-1.838503	0.3724333	0.0018716	0.9654923
3	-1.575030	0.4705592	0.0330196	0.8558080

Code

testdf(variable = financial_data$y2, #ADF test variable 2
       max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	0.6976820	0.9900000	123.7381593	0.0000000
1	-0.5613589	0.8480833	11.7643484	0.0006038
2	-1.1293522	0.6365441	0.0010576	0.9740568
3	-1.1373793	0.6335546	0.0009233	0.9757588

Code

testdf(variable = financial_data$y3, #ADF test variable 3
       max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-1.236545	0.5966220	226.0557676	0.0000000
1	-1.879758	0.3570684	11.2996232	0.0007752
2	-1.404671	0.5340065	0.0028809	0.9571950
3	-1.390255	0.5393753	0.0041126	0.9488669

Code

testdf(variable = financial_data$y4, #ADF test variable 4
       max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-1.494391	0.5005919	5.4899915	0.0191256
1	-1.523066	0.4899123	0.0102396	0.9193988
2	-1.564103	0.4746287	0.0007970	0.9774774
3	-1.606214	0.4589453	0.0087546	0.9254540

Code

testdf(variable = financial_data$y5, #ADF test variable 5
       max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	0.6234327	0.9900000	218.9601184	0.0000000
1	-1.7851697	0.3922964	25.4345834	0.0000005
2	-0.9704224	0.6957347	0.9759952	0.3231895
3	-1.2852087	0.5784981	0.0930106	0.7603841

Code

testdf(variable = financial_data$y6, #ADF test variable 6
       max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-0.6336457	0.8211613	0.4571519	0.4989580
1	-0.4964736	0.8722487	0.0005071	0.9820333
2	-0.4444151	0.8916370	0.0010145	0.9745910
3	-0.3930213	0.9040193	0.0003159	0.9858184

Code

testdf(variable = financial_data$y7, #ADF test variable 7
       max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	0.6979600	0.9900000	118.6468390	0.0000000
1	-0.5365735	0.8573142	12.4786561	0.0004116
2	-1.1320920	0.6355237	0.0001866	0.9891020
3	-1.1469827	0.6299780	0.0006452	0.9797357

Code

testdf(variable = financial_data$y8, #ADF test variable 8
       max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-1.479208	0.5062464	12.0381805	0.0005212
1	-1.543628	0.4822543	0.0025525	0.9597066
2	-1.554829	0.4780826	0.0002884	0.9864515
3	-1.590099	0.4649471	0.0998735	0.7519814

Code

testdf(variable = financial_data$y9, #ADF test variable 9
       max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	0.6277666	0.9900000	219.5489280	0.0000000
1	-1.7991091	0.3871048	27.2422408	0.0000002
2	-0.9483910	0.7039400	1.1709559	0.2792049
3	-1.2993073	0.5732473	0.0891421	0.7652706

Code

testdf(variable = financial_data$y10, #ADF test variable 10
       
       max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-1.258058	0.5886099	171.6970350	0.0000000
1	-1.557003	0.4772730	6.9948799	0.0081743
2	-1.677215	0.4325021	0.2709064	0.6027244
3	-1.441749	0.5201975	0.0628229	0.8020887

As we can see, none of the variables are stationary because p values are greater than significance level. Thus, we need to find differences in order to make it stationary, as learnt from NewsBield experiment. Neverthlesess, first of rule of cointegration is satisfied.

##Calculate differences

Code

financial_data$dy1 <- diff.xts(financial_data$y1)
financial_data$dy2 <- diff.xts(financial_data$y2)
financial_data$dy3 <- diff.xts(financial_data$y3)
financial_data$dy4 <- diff.xts(financial_data$y4)
financial_data$dy5 <- diff.xts(financial_data$y5)
financial_data$dy6 <- diff.xts(financial_data$y6)
financial_data$dy7 <- diff.xts(financial_data$y7)
financial_data$dy8 <- diff.xts(financial_data$y8)
financial_data$dy9 <- diff.xts(financial_data$y9)
financial_data$dy10 <- diff.xts(financial_data$y10)

In order to get stationarity, we need to differentiate series.

##ADF test for the y variables

Code

testdf(variable = financial_data$dy1, #ADF test variable 1
       max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-15.565954	0.01	0.0000805	0.9928395
1	-11.418055	0.01	0.0002002	0.9887101
2	-11.107373	0.01	0.0299749	0.8625471
3	-8.901743	0.01	0.0003941	0.9841618

Code

testdf(variable = financial_data$dy2, #ADF test variable 2
       max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-7.989698	0.01	11.4040782	0.0007328
1	-5.274622	0.01	0.0050530	0.9433307
2	-5.093104	0.01	0.0007863	0.9776301
3	-4.960609	0.01	0.0000187	0.9965522

Code

       testdf(variable = financial_data$dy3, #ADF test variable 3
              max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-4.529095	0.01	11.9091594	0.0005586
1	-5.539558	0.01	0.0075518	0.9307501
2	-5.387325	0.01	0.0044727	0.9466787
3	-4.529314	0.01	0.0076663	0.9302283

Code

       testdf(variable = financial_data$dy4, #ADF test variable 4
              max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-14.976767	0.01	0.0087075	0.9256541
1	-10.799791	0.01	0.0003822	0.9844031
2	-8.919976	0.01	0.0067850	0.9343518
3	-5.921353	0.01	0.5225323	0.4697635

Code

       testdf(variable = financial_data$dy5, #ADF test variable 5
              max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-4.783862	0.01	26.8354772	0.0000002
1	-6.624693	0.01	0.9335859	0.3339329
2	-5.245343	0.01	0.1037083	0.7474241
3	-5.506683	0.01	0.0139857	0.9058609

Code

       testdf(variable = financial_data$dy6, #ADF test variable 6
              max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-18.012092	0.01	0.0007932	0.9775311
1	-12.542031	0.01	0.0011870	0.9725155
2	-10.204864	0.01	0.0003551	0.9849657
3	-8.516099	0.01	0.0002669	0.9869653

Code

       testdf(variable = financial_data$dy7, #ADF test variable 7
              max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-8.174416	0.01	12.1264287	0.0004971
1	-5.277486	0.01	0.0030362	0.9560572
2	-5.081266	0.01	0.0004968	0.9822167
3	-5.006448	0.01	0.0003413	0.9852596

Code

       testdf(variable = financial_data$dy8, #ADF test variable 8
              max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-13.993074	0.01	0.0016383	0.9677134
1	-10.691950	0.01	0.0001454	0.9903780
2	-8.519015	0.01	0.0900690	0.7640895
3	-5.601904	0.01	1.0258157	0.3111435

Code

       testdf(variable = financial_data$dy9, #ADF test variable 9
              max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-4.763719	0.01	28.6705806	0.0000001
1	-6.682840	0.01	1.1258497	0.2886623
2	-5.224669	0.01	0.1005344	0.7511893
3	-5.511692	0.01	0.0162463	0.8985756

Code

       testdf(variable = financial_data$dy10, #ADF test variable 10
 
                         max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-6.346090	0.01	6.5970811	0.0102146
1	-4.839809	0.01	0.2597729	0.6102762
2	-5.495769	0.01	0.0529907	0.8179386
3	-4.757224	0.01	0.0030569	0.9559079

P value of Breusch-Godfrey test in 0augmentations is more than significant level(0.05), so there is no autocorrelation and we can check stationary. Morever, p value of ADF test gives us a need for rejecting null hypothesis “unit root”. Our data is ready for regression.

##Estimate cointegrate vector

Code

model.coint <- lm(y3 ~ y10, data = financial_data) #There are 4 cointegrated pairs but I choosed y3 and y10.
summary(model.coint)


Call:
lm(formula = y3 ~ y10, data = financial_data)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.2580 -0.2659 -0.0140  0.2628  1.1026 

Coefficients:
               Estimate  Std. Error t value Pr(>|t|)    
(Intercept) -203.059205    0.241065  -842.3   <2e-16 ***
y10            1.997761    0.001078  1852.5   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.4121 on 298 degrees of freedom
Multiple R-squared:  0.9999,    Adjusted R-squared:  0.9999 
F-statistic: 3.432e+06 on 1 and 298 DF,  p-value: < 2.2e-16

y10 has significant effect on OLS model of Engle Grager method, so one unit increase in y10 increases y3 1.998e+00 unit. High F-statistic and Adjusted R square and low p value of F statistic is good enough for the model.

##Checking autocorrelation among residuals

Code

testdf(variable = residuals(model.coint), max.augmentations = 3)

augmentations	adf	p_adf	bgodfrey	p_bg
0	-17.288806	0.01	0.0004325	0.9834088
1	-10.945801	0.01	0.0154520	0.9010729
2	-10.178250	0.01	0.0035422	0.9525409
3	-8.755398	0.01	0.0212039	0.8842248

High negative absolute value of adf in 0 augmentations and lower p-value than 0,05 in ADF test shows that residuals are stationary, so cointegration conditions are satisfied(integration order of residuals is less than integration order of equation itself). Now we will need to esitimate models.

##Plot the cointegrated pairs

Code

plot(financial_data[, c(3,10)],
     col = c("black", "blue"),
     major.ticks = "years", 
     grid.ticks.on = "years",
     grid.ticks.lty = 3,
     main = "",
     legend.loc = "topleft")

Code

financialdata_y3_y10 <- financial_data[ , c("y3", "y10")] #choosing only cointegrated pairs
tail(financialdata_y3_y10)

               y3    y10
2022-03-22 289.20 246.30
2022-03-23 288.83 245.85
2022-03-24 287.98 245.64
2022-03-25 287.54 245.77
2022-03-26 286.01 244.87
2022-03-27 285.43 244.08

##Johansen Cointegration test for find the number of vectors

Code

johansen_result_trace <- ca.jo(financialdata_y3_y10, type = "trace", K = 6, ecdet = "const")                 
summary(johansen_result_trace) #trace test


###################### 
# Johansen-Procedure # 
###################### 

Test type: trace statistic , without linear trend and constant in cointegration 

Eigenvalues (lambda):
[1] 1.504496e-01 7.861644e-03 7.436528e-18

Values of teststatistic and critical values of test:

          test 10pct  5pct  1pct
r <= 1 |  2.32  7.52  9.24 12.97
r = 0  | 50.26 17.85 19.96 24.60

Eigenvectors, normalised to first column:
(These are the cointegration relations)

             y3.l6     y10.l6   constant
y3.l6      1.00000   1.000000   1.000000
y10.l6    -1.99825  -1.732193  -1.910262
constant 203.17352 140.486475 186.883186

Weights W:
(This is the loading matrix)

          y3.l6       y10.l6               constant
y3.d  0.4125932 -0.013072792 -0.0000000000002899042
y10.d 0.7075202 -0.006356424 -0.0000000000007212011

In first test, test statistic is bigger than 5pct level and we reject Null Hypothesis which shows r=0. Next we test Null Hypothesis r<=1 and cannot reject null. Our conclusion is that there is 1 cointegration vector for pairs. However, we can use alternative of the test.

##Alternative variant of the test

Code

johansen_result_eigen <- ca.jo(financialdata_y3_y10, type = "eigen", K = 6, ecdet = "const") 
summary(johansen_result_eigen) #eigen test


###################### 
# Johansen-Procedure # 
###################### 

Test type: maximal eigenvalue statistic (lambda max) , without linear trend and constant in cointegration 

Eigenvalues (lambda):
[1] 1.504496e-01 7.861644e-03 7.436528e-18

Values of teststatistic and critical values of test:

          test 10pct  5pct  1pct
r <= 1 |  2.32  7.52  9.24 12.97
r = 0  | 47.94 13.75 15.67 20.20

Eigenvectors, normalised to first column:
(These are the cointegration relations)

             y3.l6     y10.l6   constant
y3.l6      1.00000   1.000000   1.000000
y10.l6    -1.99825  -1.732193  -1.910262
constant 203.17352 140.486475 186.883186

Weights W:
(This is the loading matrix)

          y3.l6       y10.l6               constant
y3.d  0.4125932 -0.013072792 -0.0000000000002899042
y10.d 0.7075202 -0.006356424 -0.0000000000007212011

Code

str(johansen_result_eigen)

Formal class 'ca.jo' [package "urca"] with 26 slots
  ..@ x        : num [1:300, 1:2] 174 174 176 178 181 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:300] "2021-06-01" "2021-06-02" "2021-06-03" "2021-06-04" ...
  .. .. ..$ : chr [1:2] "y3" "y10"
  ..@ Z0       : num [1:294, 1:2] 2.56 1.86 1.51 2.77 2.98 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : NULL
  .. .. ..$ : chr [1:2] "y3.d" "y10.d"
  ..@ Z1       : num [1:294, 1:10] 3.34 2.56 1.86 1.51 2.77 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : NULL
  .. .. ..$ : chr [1:10] "y3.dl1" "y10.dl1" "y3.dl2" "y10.dl2" ...
  ..@ ZK       : num [1:294, 1:3] 174 174 176 178 181 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:294] "2021-06-01" "2021-06-02" "2021-06-03" "2021-06-04" ...
  .. .. ..$ : chr [1:3] "y3.l6" "y10.l6" "constant"
  ..@ type     : chr "maximal eigenvalue statistic (lambda max)"
  ..@ model    : chr "without linear trend and constant in cointegration"
  ..@ ecdet    : chr "const"
  ..@ lag      : int 6
  ..@ P        : int 2
  ..@ season   : NULL
  ..@ dumvar   : NULL
  ..@ cval     : num [1:2, 1:3] 7.52 13.75 9.24 15.67 12.97 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "r <= 1 |" "r = 0  |"
  .. .. ..$ : chr [1:3] "10pct" "5pct" "1pct"
  ..@ teststat : num [1:2] 2.32 47.94
  ..@ lambda   : num [1:3] 1.50e-01 7.86e-03 7.44e-18
  ..@ Vorg     : num [1:3, 1:3] 5.969 -11.927 1212.676 -0.15 0.259 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:3] "y3.l6" "y10.l6" "constant"
  .. .. ..$ : chr [1:3] "y3.l6" "y10.l6" "constant"
  ..@ V        : num [1:3, 1:3] 1 -2 203.17 1 -1.73 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:3] "y3.l6" "y10.l6" "constant"
  .. .. ..$ : chr [1:3] "y3.l6" "y10.l6" "constant"
  ..@ W        : num [1:2, 1:3] 0.41259316498636 0.70752019008385 -0.01307279161621 -0.00635642359973 -0.00000000000029 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "y3.d" "y10.d"
  .. .. ..$ : chr [1:3] "y3.l6" "y10.l6" "constant"
  ..@ PI       : num [1:2, 1:3] 0.4 0.701 -0.802 -1.403 81.991 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "y3.d" "y10.d"
  .. .. ..$ : chr [1:3] "y3.l6" "y10.l6" "constant"
  ..@ DELTA    : num [1:2, 1:2] 0.988 0.514 0.514 0.307
  ..@ GAMMA    : num [1:2, 1:10] 1.2471 1.1234 -0.3473 -1.179 0.0991 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "y3.d" "y10.d"
  .. .. ..$ : chr [1:10] "y3.dl1" "y10.dl1" "y3.dl2" "y10.dl2" ...
  ..@ R0       : num [1:294, 1:2] -0.2843 -0.0516 -0.2363 1.5016 0.2667 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : NULL
  .. .. ..$ : chr [1:2] "R0.y3.d" "R0.y10.d"
  ..@ RK       : num [1:294, 1:3] 142 116 125 119 114 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:294] "2021-06-01" "2021-06-02" "2021-06-03" "2021-06-04" ...
  .. .. ..$ : chr [1:3] "RK.y3.l6" "RK.y10.l6" "RK.constant"
  ..@ bp       : logi NA
  ..@ spec     : chr "longrun"
  ..@ call     : language ca.jo(x = financialdata_y3_y10, type = "eigen", ecdet = "const", K = 6)
  ..@ test.name: chr "Johansen-Procedure"

The alternative variant also gives the same result. There is one cointegrated pair.

##Estimatin VECM

Code

TMS.vec4 <- cajorls(johansen_result_eigen, r = 1) #we can either eigen or trace results
summary(TMS.vec4)

     Length Class  Mode   
rlm  12     mlm    list   
beta  3     -none- numeric

Code

str(TMS.vec4)

List of 2
 $ rlm :List of 12
  ..$ coefficients : num [1:11, 1:2] 0.413 1.259 -0.355 0.102 -0.577 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:11] "ect1" "y3.dl1" "y10.dl1" "y3.dl2" ...
  .. .. ..$ : chr [1:2] "y3.d" "y10.d"
  ..$ residuals    : num [1:294, 1:2] -0.3069 -0.0262 -0.1902 1.6196 0.3964 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:294] "2021-06-01" "2021-06-02" "2021-06-03" "2021-06-04" ...
  .. .. ..$ : chr [1:2] "y3.d" "y10.d"
  ..$ effects      : num [1:294, 1:2] -2.24 -32.16 -0.18 4.11 0.6 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:294] "ect1" "y3.dl1" "y10.dl1" "y3.dl2" ...
  .. .. ..$ : chr [1:2] "y3.d" "y10.d"
  ..$ rank         : int 11
  ..$ fitted.values: num [1:294, 1:2] 2.87 1.89 1.7 1.15 2.58 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:294] "2021-06-01" "2021-06-02" "2021-06-03" "2021-06-04" ...
  .. .. ..$ : chr [1:2] "y3.d" "y10.d"
  ..$ assign       : int [1:11] 1 2 3 4 5 6 7 8 9 10 ...
  ..$ qr           :List of 5
  .. ..$ qr   : num [1:294, 1:11] -7.0041 0.0476 0.0692 0.0681 -0.1213 ...
  .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. ..$ : chr [1:294] "2021-06-01" "2021-06-02" "2021-06-03" "2021-06-04" ...
  .. .. .. ..$ : chr [1:11] "ect1" "y3.dl1" "y10.dl1" "y3.dl2" ...
  .. .. ..- attr(*, "assign")= int [1:11] 1 2 3 4 5 6 7 8 9 10 ...
  .. ..$ qraux: num [1:11] 1.05 1.06 1.03 1.04 1 ...
  .. ..$ pivot: int [1:11] 1 2 3 4 5 6 7 8 9 10 ...
  .. ..$ tol  : num 0.0000001
  .. ..$ rank : int 11
  .. ..- attr(*, "class")= chr "qr"
  ..$ df.residual  : int 283
  ..$ xlevels      : Named list()
  ..$ call         : language lm(formula = substitute(form1), data = data.mat)
  ..$ terms        :Classes 'terms', 'formula'  language z@Z0 ~ ect1 + y3.dl1 + y10.dl1 + y3.dl2 + y10.dl2 + y3.dl3 + y10.dl3 +      y3.dl4 + y10.dl4 + y3.dl5 + y10.dl5 - 1
  .. .. ..- attr(*, "variables")= language list(z@Z0, ect1, y3.dl1, y10.dl1, y3.dl2, y10.dl2, y3.dl3, y10.dl3, y3.dl4,      y10.dl4, y3.dl5, y10.dl5)
  .. .. ..- attr(*, "factors")= int [1:12, 1:11] 0 1 0 0 0 0 0 0 0 0 ...
  .. .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. .. ..$ : chr [1:12] "z@Z0" "ect1" "y3.dl1" "y10.dl1" ...
  .. .. .. .. ..$ : chr [1:11] "ect1" "y3.dl1" "y10.dl1" "y3.dl2" ...
  .. .. ..- attr(*, "term.labels")= chr [1:11] "ect1" "y3.dl1" "y10.dl1" "y3.dl2" ...
  .. .. ..- attr(*, "order")= int [1:11] 1 1 1 1 1 1 1 1 1 1 ...
  .. .. ..- attr(*, "intercept")= int 0
  .. .. ..- attr(*, "response")= int 1
  .. .. ..- attr(*, ".Environment")=<environment: 0x00000242b1b729a8> 
  .. .. ..- attr(*, "predvars")= language list(z@Z0, ect1, y3.dl1, y10.dl1, y3.dl2, y10.dl2, y3.dl3, y10.dl3, y3.dl4,      y10.dl4, y3.dl5, y10.dl5)
  .. .. ..- attr(*, "dataClasses")= Named chr [1:12] "nmatrix.2" "numeric" "numeric" "numeric" ...
  .. .. .. ..- attr(*, "names")= chr [1:12] "z@Z0" "ect1" "y3.dl1" "y10.dl1" ...
  ..$ model        :'data.frame':   294 obs. of  12 variables:
  .. ..$ z@Z0   : num [1:294, 1:2] 2.56 1.86 1.51 2.77 2.98 ...
  .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. ..$ : NULL
  .. .. .. ..$ : chr [1:2] "y3.d" "y10.d"
  .. ..$ ect1   : num [1:294] 0.373 0.333 0.485 0.477 -0.849 ...
  .. ..$ y3.dl1 : num [1:294] 3.34 2.56 1.86 1.51 2.77 ...
  .. ..$ y10.dl1: num [1:294] 1.41 1.3 0.77 0.96 0.99 ...
  .. ..$ y3.dl2 : num [1:294] 3.21 3.34 2.56 1.86 1.51 ...
  .. ..$ y10.dl2: num [1:294] 2.27 1.41 1.3 0.77 0.96 ...
  .. ..$ y3.dl3 : num [1:294] 2.19 3.21 3.34 2.56 1.86 ...
  .. ..$ y10.dl3: num [1:294] 1.1 2.27 1.41 1.3 0.77 ...
  .. ..$ y3.dl4 : num [1:294] 2.13 2.19 3.21 3.34 2.56 ...
  .. ..$ y10.dl4: num [1:294] 0.99 1.1 2.27 1.41 1.3 ...
  .. ..$ y3.dl5 : num [1:294] 0 2.13 2.19 3.21 3.34 ...
  .. ..$ y10.dl5: num [1:294] 0.02 0.99 1.1 2.27 1.41 ...
  .. ..- attr(*, "terms")=Classes 'terms', 'formula'  language z@Z0 ~ ect1 + y3.dl1 + y10.dl1 + y3.dl2 + y10.dl2 + y3.dl3 + y10.dl3 +      y3.dl4 + y10.dl4 + y3.dl5 + y10.dl5 - 1
  .. .. .. ..- attr(*, "variables")= language list(z@Z0, ect1, y3.dl1, y10.dl1, y3.dl2, y10.dl2, y3.dl3, y10.dl3, y3.dl4,      y10.dl4, y3.dl5, y10.dl5)
  .. .. .. ..- attr(*, "factors")= int [1:12, 1:11] 0 1 0 0 0 0 0 0 0 0 ...
  .. .. .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. .. .. ..$ : chr [1:12] "z@Z0" "ect1" "y3.dl1" "y10.dl1" ...
  .. .. .. .. .. ..$ : chr [1:11] "ect1" "y3.dl1" "y10.dl1" "y3.dl2" ...
  .. .. .. ..- attr(*, "term.labels")= chr [1:11] "ect1" "y3.dl1" "y10.dl1" "y3.dl2" ...
  .. .. .. ..- attr(*, "order")= int [1:11] 1 1 1 1 1 1 1 1 1 1 ...
  .. .. .. ..- attr(*, "intercept")= int 0
  .. .. .. ..- attr(*, "response")= int 1
  .. .. .. ..- attr(*, ".Environment")=<environment: 0x00000242b1b729a8> 
  .. .. .. ..- attr(*, "predvars")= language list(z@Z0, ect1, y3.dl1, y10.dl1, y3.dl2, y10.dl2, y3.dl3, y10.dl3, y3.dl4,      y10.dl4, y3.dl5, y10.dl5)
  .. .. .. ..- attr(*, "dataClasses")= Named chr [1:12] "nmatrix.2" "numeric" "numeric" "numeric" ...
  .. .. .. .. ..- attr(*, "names")= chr [1:12] "z@Z0" "ect1" "y3.dl1" "y10.dl1" ...
  ..- attr(*, "class")= chr [1:2] "mlm" "lm"
 $ beta: num [1:3, 1] 1 -2 203
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:3] "y3.l6" "y10.l6" "constant"
  .. ..$ : chr "ect1"

Code

summary(TMS.vec4$rlm)

Response y3.d :

Call:
lm(formula = y3.d ~ ect1 + y3.dl1 + y10.dl1 + y3.dl2 + y10.dl2 + 
    y3.dl3 + y10.dl3 + y3.dl4 + y10.dl4 + y3.dl5 + y10.dl5 - 
    1, data = data.mat)

Residuals:
     Min       1Q   Median       3Q      Max 
-2.69689 -0.61031  0.07821  0.64112  2.93975 

Coefficients:
        Estimate Std. Error t value          Pr(>|t|)    
ect1     0.41259    0.35408   1.165             0.245    
y3.dl1   1.25943    0.16609   7.583 0.000000000000488 ***
y10.dl1 -0.35514    0.29790  -1.192             0.234    
y3.dl2   0.10212    0.22958   0.445             0.657    
y10.dl2 -0.57742    0.41891  -1.378             0.169    
y3.dl3   0.04274    0.26299   0.163             0.871    
y10.dl3 -0.45330    0.49456  -0.917             0.360    
y3.dl4   0.24502    0.30164   0.812             0.417    
y10.dl4 -0.11653    0.57378  -0.203             0.839    
y3.dl5   0.23621    0.32386   0.729             0.466    
y10.dl5 -0.56602    0.64404  -0.879             0.380    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.017 on 283 degrees of freedom
Multiple R-squared:  0.7848,    Adjusted R-squared:  0.7764 
F-statistic: 93.82 on 11 and 283 DF,  p-value: < 2.2e-16


Response y10.d :

Call:
lm(formula = y10.d ~ ect1 + y3.dl1 + y10.dl1 + y3.dl2 + y10.dl2 + 
    y3.dl3 + y10.dl3 + y3.dl4 + y10.dl4 + y3.dl5 + y10.dl5 - 
    1, data = data.mat)

Residuals:
     Min       1Q   Median       3Q      Max 
-1.72167 -0.34499 -0.00932  0.38694  1.60174 

Coefficients:
        Estimate Std. Error t value         Pr(>|t|)    
ect1     0.70752    0.19705   3.591         0.000389 ***
y3.dl1   1.12945    0.09243  12.220          < 2e-16 ***
y10.dl1 -1.18283    0.16578  -7.135 0.00000000000814 ***
y3.dl2   0.49607    0.12776   3.883         0.000129 ***
y10.dl2 -1.18484    0.23312  -5.083 0.00000067784059 ***
y3.dl3   0.52732    0.14635   3.603         0.000371 ***
y10.dl3 -1.21963    0.27522  -4.431 0.00001339675361 ***
y3.dl4   0.57103    0.16786   3.402         0.000766 ***
y10.dl4 -1.01878    0.31930  -3.191         0.001579 ** 
y3.dl5   0.62640    0.18023   3.476         0.000590 ***
y10.dl5 -1.24639    0.35840  -3.478         0.000585 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.5661 on 283 degrees of freedom
Multiple R-squared:  0.7536,    Adjusted R-squared:  0.744 
F-statistic: 78.67 on 11 and 283 DF,  p-value: < 2.2e-16

Code

TMS.vec4$beta # Cointegarting vector

              ect1
y3.l6      1.00000
y10.l6    -1.99825
constant 203.17352

There 5 cointegartion matrices. P values of Fstatistic is less 0.05, so we reject null hypothesis that coefficents are 0. Morever, ect1 in y10.d is 0.70752 and shows p value<0.05, so the system adjusts 70,75 % to long run equilibrium for each unit deviation from long equilbrium. On the other hand, because ect1 is not statistically significant in y3.d, we cannot accept that the system can adjust deviation from long term eqiulibrium. $beta code gave use cointegration vector.

##Repramatrize from VECM to VAR

Code

TMS.vec4.asVAR <- vec2var(johansen_result_eigen, r = 1)
TMS.vec4.asVAR


Coefficient matrix of lagged endogenous variables:

A1:
       y3.l1     y10.l1
y3  2.259433 -0.3551365
y10 1.129454 -0.1828322


A2:
         y3.l2       y10.l2
y3  -1.1573099 -0.222285118
y10 -0.6333811 -0.002007442


A3:
          y3.l3      y10.l3
y3  -0.05938044  0.12412110
y10  0.03124713 -0.03478676


A4:
         y3.l4    y10.l4
y3  0.20227998 0.3367707
y10 0.04370972 0.2008484


A5:
           y3.l5     y10.l5
y3  -0.008811775 -0.4494859
y10  0.055367291 -0.2276139


A6:
         y3.l6     y10.l6
y3  0.17638221 -0.2584488
y10 0.08112279 -0.1674107


Coefficient matrix of deterministic regressor(s).

     constant
y3   83.82801
y10 143.74937

Because generally VECM is differenced equation, VAR has one more matrices. So VAR(p+1)=VECM. Here we have 6 matrices.

VAR diagnostics in IRF

Code

irf_result <- irf(TMS.vec4.asVAR)
plot(irf_result)

We can see that shocks coming from y3 has more impact on both and after 10 periods these shocks disappear. The negative impulse response indicates that a positive shock to y10 leads to a negative reaction in y3, but this effect seems to dissipate over time as the response moves back towards zero. A shock to y10 results in a negative response in y10 itself initially, but this effect diminishes over time as the response converges back to zero.

##VAR diagnostics on FEVD

Code

plot(fevd(TMS.vec4.asVAR, n.ahead = 20))

Shocks to y10 account for nearly all the forecast error variance, indicating a strong dependence of y3 on y10. After the first period, y10 almost entirely explains its own forecast error variance, highlighting its self-determination in forecasting errors.

##Perform the Portmanteau test for residuals on VAR model

Code

head(residuals(TMS.vec4.asVAR))

     resids of y3 resids of y10
[1,]  -0.30693742    0.03265190
[2,]  -0.02622717    0.04001272
[3,]  -0.19019950    0.13458942
[4,]   1.61957358    0.65575761
[5,]   0.39639475    0.42779874
[6,]   1.06454554    0.55720276

Code

portmanteau_test <- serial.test(TMS.vec4.asVAR)
print(portmanteau_test)


    Portmanteau Test (asymptotic)

data:  Residuals of VAR object TMS.vec4.asVAR
Chi-squared = 40.023, df = 42, p-value = 0.5581

$serial

    Portmanteau Test (asymptotic)

data:  Residuals of VAR object TMS.vec4.asVAR
Chi-squared = 40.023, df = 42, p-value = 0.5581

We cannot reject Null Hypothesis which indicates no autororeelation in residuals, so there is no autocorrelation among residuals.

Perform the Jarque-Bera test for normality

Code

jb_test <- normality.test(TMS.vec4.asVAR)
print(jb_test)

$JB

    JB-Test (multivariate)

data:  Residuals of VAR object TMS.vec4.asVAR
Chi-squared = 0.62506, df = 4, p-value = 0.9602


$Skewness

    Skewness only (multivariate)

data:  Residuals of VAR object TMS.vec4.asVAR
Chi-squared = 0.24287, df = 2, p-value = 0.8857


$Kurtosis

    Kurtosis only (multivariate)

data:  Residuals of VAR object TMS.vec4.asVAR
Chi-squared = 0.38219, df = 2, p-value = 0.8261

$jb.mul
$jb.mul$JB

    JB-Test (multivariate)

data:  Residuals of VAR object TMS.vec4.asVAR
Chi-squared = 0.62506, df = 4, p-value = 0.9602


$jb.mul$Skewness

    Skewness only (multivariate)

data:  Residuals of VAR object TMS.vec4.asVAR
Chi-squared = 0.24287, df = 2, p-value = 0.8857


$jb.mul$Kurtosis

    Kurtosis only (multivariate)

data:  Residuals of VAR object TMS.vec4.asVAR
Chi-squared = 0.38219, df = 2, p-value = 0.8261

P values are bigger than 0.05, so we cannot ignore Null Hypothesis which shows normal distribution. Our residuals are normally distributed.

##Forecasting out of sample

Code

TMS.testdata <- financial_data["/3/7/2022", ]
TMS.testdata <- subset(TMS.testdata, select = c(y3, y10))
tail(TMS.testdata)

               y3    y10
2022-03-22 289.20 246.30
2022-03-23 288.83 245.85
2022-03-24 287.98 245.64
2022-03-25 287.54 245.77
2022-03-26 286.01 244.87
2022-03-27 285.43 244.08

Code

johan.test.eigen2.testdata <- ca.jo(TMS.testdata,
                                 ecdet = "const",
                                 type = "trace", 
                                 K = 4)
summary(johan.test.eigen2.testdata)


###################### 
# Johansen-Procedure # 
###################### 

Test type: trace statistic , without linear trend and constant in cointegration 

Eigenvalues (lambda):
[1] 2.126143e-01 9.370866e-03 1.734723e-18

Values of teststatistic and critical values of test:

          test 10pct  5pct  1pct
r <= 1 |  2.79  7.52  9.24 12.97
r = 0  | 73.54 17.85 19.96 24.60

Eigenvectors, normalised to first column:
(These are the cointegration relations)

              y3.l4     y10.l4  constant
y3.l4      1.000000   1.000000  1.000000
y10.l4    -1.998201  -1.724705 -1.304082
constant 203.166419 138.259152 71.888373

Weights W:
(This is the loading matrix)

          y3.l4       y10.l4      constant
y3.d  0.1260953 -0.013703907 -8.908024e-16
y10.d 0.5689688 -0.006982952 -5.897452e-16

For out of sample, we need to test Johansen test again and repramatrize it. Here there is only one vector.

Code

TMS.vec4.fore <- 
  predict(
    vec2var(
      johan.test.eigen2.testdata, 
      r = 1),     # no of cointegrating vectors 
    n.ahead = 20, # forecast horizon
    ci = 0.95)    # confidence level for intervals

TMS.vec4.fore$fcst$y3

          fcst    lower    upper        CI
 [1,] 285.2120 283.2259 287.1980  1.986034
 [2,] 285.2980 280.7385 289.8576  4.559572
 [3,] 285.4969 277.9828 293.0110  7.514084
 [4,] 285.5471 274.9346 296.1595 10.612434
 [5,] 285.5513 271.8319 299.2706 13.719334
 [6,] 285.5234 268.7630 302.2838 16.760351
 [7,] 285.5071 265.8111 305.2031 19.695973
 [8,] 285.4938 262.9869 308.0007 22.506879
 [9,] 285.4864 260.2999 310.6730 25.186535
[10,] 285.4796 257.7435 313.2158 27.736139
[11,] 285.4746 255.3132 315.6359 30.161362
[12,] 285.4703 253.0002 317.9403 32.470075
[13,] 285.4671 250.7961 320.1381 34.671019
[14,] 285.4647 248.6917 322.2377 36.773010
[15,] 285.4629 246.6784 324.2475 38.784517
[16,] 285.4616 244.7481 326.1750 40.713444
[17,] 285.4606 242.8936 328.0276 42.567036
[18,] 285.4598 241.1080 329.8117 44.351862
[19,] 285.4593 239.3854 331.5331 46.073827
[20,] 285.4589 237.7206 333.1971 47.738225

Code

TMS.vec4.fore$fcst$y10

          fcst    lower    upper        CI
 [1,] 244.4091 243.3038 245.5143  1.105265
 [2,] 244.4101 242.0766 246.7435  2.333497
 [3,] 244.5979 240.8117 248.3840  3.786132
 [4,] 244.5708 239.2513 249.8904  5.319554
 [5,] 244.5863 237.7225 251.4502  6.863865
 [6,] 244.5589 236.1792 252.9386  8.379684
 [7,] 244.5581 234.7125 254.4038  9.845634
 [8,] 244.5489 233.2982 255.7995 11.250646
 [9,] 244.5473 231.9564 257.1382 12.590880
[10,] 244.5428 230.6763 258.4093 13.866521
[11,] 244.5407 229.4605 259.6210 15.080244
[12,] 244.5382 228.3024 260.7741 16.235844
[13,] 244.5368 227.1992 261.8744 17.337624
[14,] 244.5355 226.1456 262.9254 18.389937
[15,] 244.5347 225.1377 263.9317 19.396991
[16,] 244.5340 224.1712 264.8967 20.362721
[17,] 244.5335 223.2427 265.8242 21.290742
[18,] 244.5331 222.3488 266.7174 22.184333
[19,] 244.5328 221.4864 267.5793 23.046447
[20,] 244.5326 220.6529 268.4123 23.879731

Code

tail(index(financial_data), 20)

 [1] "2022-03-08" "2022-03-09" "2022-03-10" "2022-03-11" "2022-03-12"
 [6] "2022-03-13" "2022-03-14" "2022-03-15" "2022-03-16" "2022-03-17"
[11] "2022-03-18" "2022-03-19" "2022-03-20" "2022-03-21" "2022-03-22"
[16] "2022-03-23" "2022-03-24" "2022-03-25" "2022-03-26" "2022-03-27"

Here we added forecast series and lower and upper bounds.

##Adding y3 and y10 forecasts

Code

y3_forecast <- xts(TMS.vec4.fore$fcst$y3[,-4], 
                   # we exclude the last column with CI
                   tail(index(financial_data), 20))

names(y3_forecast) <- c("y3_fore", "y3_lower", "y3_upper")

y10_forecast <- xts(TMS.vec4.fore$fcst$y10[,-4],
                   # we exclude the last column with CI
                   tail(index(financial_data), 20))
names(y10_forecast) <- c("y10_fore", "y10_lower", "y10_upper")

##Combinening real data and forecast data

Code

financial_data_fore <- merge(financial_data, 
                 y3_forecast,
                 y10_forecast)

##Plotting forecats

Code

plot(financial_data_fore[index(financial_data_fore) > as.Date("2022-03-07"), c("y3", "y3_fore",
                                             "y3_lower", "y3_upper")], 
     
     main = "20 days forecast of instrument y3",
     col = c("black", "blue", "red", "red"))

Code

plot(financial_data_fore[index(financial_data_fore) > as.Date("2022-03-07"), c("y10", "y10_fore",
                                                                     "y10_lower", "y10_upper")], 
     
     main = "20 days forecast of instrument y10",
     col = c("black", "blue", "red", "red"))

The forecasted prices for y3 and y10 over 20 time periods show that projected observations fall within the upper and lower boundaries. In the plot, the blue line represents the forecasted values, while the other line (assumed to be black) indicates the original observations.

##Calculating forecast errors for y3

Code

financial_data_fore$mae.y3   <-  abs(financial_data_fore$y3 - financial_data_fore$y3_fore)
tail(financial_data_fore$mae.y3, 20)

                mae.y3
2022-03-08 17.67801884
2022-03-09 14.81196569
2022-03-10 12.72307274
2022-03-11 11.56294243
2022-03-12 11.10873033
2022-03-13 11.30660064
2022-03-14 10.53288456
2022-03-15  8.32617952
2022-03-16  6.00355962
2022-03-17  4.57037360
2022-03-18  3.32543082
2022-03-19  2.90974191
2022-03-20  2.97287567
2022-03-21  3.39530022
2022-03-22  3.73706357
2022-03-23  3.36840866
2022-03-24  2.51940563
2022-03-25  2.08016109
2022-03-26  0.55072322
2022-03-27  0.02885409

Code

financial_data_fore$mse.y3   <-  (financial_data_fore$y3 - financial_data_fore$y3_fore) ^ 2
tail(financial_data_fore$mse.y3, 20)

                   mse.y3
2022-03-08 312.5123502194
2022-03-09 219.3943275444
2022-03-10 161.8765799741
2022-03-11 133.7016375905
2022-03-12 123.4038895259
2022-03-13 127.8392180972
2022-03-14 110.9416571973
2022-03-15  69.3252654070
2022-03-16  36.0427280663
2022-03-17  20.8883148576
2022-03-18  11.0584901296
2022-03-19   8.4665979558
2022-03-20   8.8379897690
2022-03-21  11.5280635822
2022-03-22  13.9656440971
2022-03-23  11.3461768767
2022-03-24   6.3474047119
2022-03-25   4.3270701417
2022-03-26   0.3032960660
2022-03-27   0.0008325584

Code

financial_data_fore$mape.y3   <-  abs((financial_data_fore$y3 - financial_data_fore$y3_fore) / financial_data_fore$y3_fore)
tail(financial_data_fore$amape.y3, 20)

NULL

Code

financial_data_fore$amape.y3   <-  abs((financial_data_fore$y3 - financial_data_fore$y3_fore) / (financial_data_fore$y3+financial_data_fore$y3_fore))
tail(financial_data_fore$amape.y3, 20)

                amape.y3
2022-03-08 0.03005944447
2022-03-09 0.02530195149
2022-03-10 0.02179664859
2022-03-11 0.01984519415
2022-03-12 0.01908023927
2022-03-13 0.01941535957
2022-03-14 0.01811183356
2022-03-15 0.01437248556
2022-03-16 0.01040520755
2022-03-17 0.00794116131
2022-03-18 0.00579066356
2022-03-19 0.00507055955
2022-03-20 0.00518006407
2022-03-21 0.00591181299
2022-03-22 0.00650305306
2022-03-23 0.00586532818
2022-03-24 0.00439349019
2022-03-25 0.00363029960
2022-03-26 0.00096369699
2022-03-27 0.00005054239

##Calculate means for the last 20 values for each error metric and store them in a table

Code

error_metrics_means <- data.frame(
  Metric = c("amape.y3", "mape.y3", "mse.y3", "mae.y3"),
  Mean = c(
    mean(tail(financial_data_fore$amape.y3, 20), na.rm = TRUE),
    mean(tail(financial_data_fore$mape.y3, 20), na.rm = TRUE),
    mean(tail(financial_data_fore$mse.y3, 20), na.rm = TRUE),
    mean(tail(financial_data_fore$mae.y3, 20), na.rm = TRUE)
  )
)

Print the table

Code

print(error_metrics_means)

    Metric        Mean
1 amape.y3  0.01148445
2  mape.y3  0.02338695
3   mse.y3 69.60537672
4   mae.y3  6.67561464

##Calculating forecast errors for y10

Code

financial_data_fore$mae.y10   <-  abs(financial_data_fore$y10 - financial_data_fore$y10_fore)
tail(financial_data_fore$mae.y10, 20)

             mae.y10
2022-03-08 8.9209393
2022-03-09 7.1899471
2022-03-10 6.2421233
2022-03-11 5.5691888
2022-03-12 5.3336515
2022-03-13 5.7511118
2022-03-14 5.1918711
2022-03-15 4.4911458
2022-03-16 2.5427150
2022-03-17 2.1972097
2022-03-18 1.6892765
2022-03-19 1.8817754
2022-03-20 1.8431908
2022-03-21 2.1444948
2022-03-22 1.7653333
2022-03-23 1.3160375
2022-03-24 1.1065269
2022-03-25 1.2369159
2022-03-26 0.3371955
2022-03-27 0.4525891

Code

financial_data_fore$mse.y10   <-  (financial_data_fore$y10 - financial_data_fore$y10_fore) ^ 2
tail(financial_data_fore$mse.y10, 20)

              mse.y10
2022-03-08 79.5831578
2022-03-09 51.6953393
2022-03-10 38.9641028
2022-03-11 31.0158640
2022-03-12 28.4478382
2022-03-13 33.0752864
2022-03-14 26.9555257
2022-03-15 20.1703910
2022-03-16  6.4653994
2022-03-17  4.8277305
2022-03-18  2.8536552
2022-03-19  3.5410788
2022-03-20  3.3973524
2022-03-21  4.5988577
2022-03-22  3.1164017
2022-03-23  1.7319546
2022-03-24  1.2244018
2022-03-25  1.5299610
2022-03-26  0.1137008
2022-03-27  0.2048369

Code

financial_data_fore$mape.y10   <-  abs((financial_data_fore$y10 - financial_data_fore$y10_fore) / financial_data_fore$y10_fore)
tail(financial_data_fore$mape.y10, 20)

              mape.y10
2022-03-08 0.036500035
2022-03-09 0.029417559
2022-03-10 0.025519941
2022-03-11 0.022771273
2022-03-12 0.021806824
2022-03-13 0.023516266
2022-03-14 0.021229599
2022-03-15 0.018365025
2022-03-16 0.010397641
2022-03-17 0.008984970
2022-03-18 0.006907956
2022-03-19 0.007695220
2022-03-20 0.007537478
2022-03-21 0.008769666
2022-03-22 0.007219154
2022-03-23 0.005381819
2022-03-24 0.004525053
2022-03-25 0.005058276
2022-03-26 0.001378938
2022-03-27 0.001850833

Code

financial_data_fore$amape.y10   <-  abs((financial_data_fore$y10 - financial_data_fore$y10_fore) / (financial_data_fore$y10+financial_data_fore$y10_fore))
tail(financial_data_fore$amape.y10, 20)

              amape.y10
2022-03-08 0.0179229239
2022-03-09 0.0144955673
2022-03-10 0.0125992048
2022-03-11 0.0112574633
2022-03-12 0.0107858099
2022-03-13 0.0116214858
2022-03-14 0.0105033092
2022-03-15 0.0090989612
2022-03-16 0.0051719327
2022-03-17 0.0044723930
2022-03-18 0.0034420890
2022-03-19 0.0038328626
2022-03-20 0.0037545890
2022-03-21 0.0043656903
2022-03-22 0.0035965946
2022-03-23 0.0026836878
2022-03-24 0.0022574190
2022-03-25 0.0025227578
2022-03-26 0.0006889937
2022-03-27 0.0009262739

##Calculate means for the last 20 values for each error metric and store them in a table

Code

error_metrics_means <- data.frame(
  Metric = c("amape.y10", "mape.y10", "mse.y10", "mae.y10"),
  Mean = c(
    mean(tail(financial_data_fore$amape.y10, 20), na.rm = TRUE),
    mean(tail(financial_data_fore$mape.y10, 20), na.rm = TRUE),
    mean(tail(financial_data_fore$mse.y10, 20), na.rm = TRUE),
    mean(tail(financial_data_fore$mae.y10, 20), na.rm = TRUE)
  )
)

##Print the table

Code

print(error_metrics_means)

     Metric        Mean
1 amape.y10  0.00680000
2  mape.y10  0.01374168
3   mse.y10 17.17564181
4   mae.y10  3.36016195