library(tidyquant)
## 필요한 패키지를 로딩중입니다: lubridate
##
## 다음의 패키지를 부착합니다: 'lubridate'
## The following objects are masked from 'package:base':
##
## date, intersect, setdiff, union
## 필요한 패키지를 로딩중입니다: PerformanceAnalytics
## 필요한 패키지를 로딩중입니다: xts
## 필요한 패키지를 로딩중입니다: zoo
##
## 다음의 패키지를 부착합니다: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
##
## 다음의 패키지를 부착합니다: 'PerformanceAnalytics'
## The following object is masked from 'package:graphics':
##
## legend
## 필요한 패키지를 로딩중입니다: quantmod
## 필요한 패키지를 로딩중입니다: TTR
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(ggplot2)
start_date <- '2018-01-01'
end_date <- '2023-03-31'
sp500_data <- tq_get('SPY',
from = start_date,
to = end_date,
get='stock.prices',
periodicity = 'daily')
head(sp500_data)
## # A tibble: 6 x 8
## symbol date open high low close volume adjusted
## <chr> <date> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 SPY 2018-01-02 268. 269. 267. 269. 86655700 245.
## 2 SPY 2018-01-03 269. 271. 269. 270. 90070400 246.
## 3 SPY 2018-01-04 271. 272. 271. 272. 80636400 248.
## 4 SPY 2018-01-05 273. 274. 272. 273. 83524000 249.
## 5 SPY 2018-01-08 273. 274. 273. 274. 57319200 250.
## 6 SPY 2018-01-09 274. 275. 274. 275. 57254000 250.
sp500_return <- sp500_data %>%
tq_transmute(select = adjusted,
mutate_fun = dailyReturn,
type = 'log',
col_rename = 'return')
head(sp500_return)
## # A tibble: 6 x 2
## date return
## <date> <dbl>
## 1 2018-01-02 0
## 2 2018-01-03 0.00631
## 3 2018-01-04 0.00421
## 4 2018-01-05 0.00664
## 5 2018-01-08 0.00183
## 6 2018-01-09 0.00226
sp500_return %>%
ggplot(aes(x=date,y=return))+
geom_line(color='#4373B6')+
theme_bw()+
labs(x='Date',y='Return',title='S&P500 Return')+
scale_y_continuous(labels = scales::percent)

ggplot(sp500_return)+
geom_histogram(aes(x=return, y=after_stat(density)), binwidth = 0.005, color="steelblue", fill="grey", size=1) +
stat_function(fun = dnorm, args = list(mean = mean(sp500_return$return, na.rm = T), sd = sd(sp500_return$return, na.rm = T)), size=1)+theme_bw()
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## i Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

library(PerformanceAnalytics)
library(timetk)
sp500_return_xts <-
sp500_return %>% tk_xts()
## Warning: Non-numeric columns being dropped: date
## Using column `date` for date_var.
head(sp500_return_xts)
## return
## 2018-01-02 0.000000000
## 2018-01-03 0.006305173
## 2018-01-04 0.004206151
## 2018-01-05 0.006641893
## 2018-01-08 0.001826951
## 2018-01-09 0.002261030
realizedvol <- rollapply(sp500_return_xts,width=22,FUN=sd.annualized)
realizedvol_df <- realizedvol %>% tk_tbl()
tail(realizedvol_df)
## # A tibble: 6 x 2
## index return
## <date> <dbl>
## 1 2023-03-23 0.177
## 2 2023-03-24 0.178
## 3 2023-03-27 0.178
## 4 2023-03-28 0.174
## 5 2023-03-29 0.180
## 6 2023-03-30 0.180
realizedvol_df %>%
ggplot(aes(x=index,y=return))+
geom_line(color='#4373B6')+
theme_bw()+
labs(x='Date',y='Volatility',title='S&P 500 Rolling Volatility')+
scale_y_continuous(labels = scales::percent)
## Warning: Removed 21 rows containing missing values (`geom_line()`).

library(rugarch)
## 필요한 패키지를 로딩중입니다: parallel
##
## 다음의 패키지를 부착합니다: 'rugarch'
## The following object is masked from 'package:stats':
##
## sigma
ug_spec <- ugarchspec(mean.model=list(armaOrder=c(1,0)))
ugfit = ugarchfit(spec = ug_spec, data = sp500_return_xts)
ugfit
##
## *---------------------------------*
## * GARCH Model Fit *
## *---------------------------------*
##
## Conditional Variance Dynamics
## -----------------------------------
## GARCH Model : sGARCH(1,1)
## Mean Model : ARFIMA(1,0,0)
## Distribution : norm
##
## Optimal Parameters
## ------------------------------------
## Estimate Std. Error t value Pr(>|t|)
## mu 0.001033 0.000220 4.6944 0.000003
## ar1 -0.062182 0.030800 -2.0189 0.043499
## omega 0.000005 0.000003 1.8047 0.071120
## alpha1 0.230201 0.030737 7.4893 0.000000
## beta1 0.756602 0.030548 24.7674 0.000000
##
## Robust Standard Errors:
## Estimate Std. Error t value Pr(>|t|)
## mu 0.001033 0.000285 3.62456 0.000289
## ar1 -0.062182 0.030160 -2.06176 0.039231
## omega 0.000005 0.000013 0.39607 0.692052
## alpha1 0.230201 0.055151 4.17404 0.000030
## beta1 0.756602 0.088310 8.56758 0.000000
##
## LogLikelihood : 4152.005
##
## Information Criteria
## ------------------------------------
##
## Akaike -6.2833
## Bayes -6.2637
## Shibata -6.2834
## Hannan-Quinn -6.2760
##
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 2.855 0.09107
## Lag[2*(p+q)+(p+q)-1][2] 2.859 0.04634
## Lag[4*(p+q)+(p+q)-1][5] 3.254 0.37018
## d.o.f=1
## H0 : No serial correlation
##
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
## statistic p-value
## Lag[1] 0.001571 0.9684
## Lag[2*(p+q)+(p+q)-1][5] 0.861658 0.8901
## Lag[4*(p+q)+(p+q)-1][9] 2.343412 0.8603
## d.o.f=2
##
## Weighted ARCH LM Tests
## ------------------------------------
## Statistic Shape Scale P-Value
## ARCH Lag[3] 0.005759 0.500 2.000 0.9395
## ARCH Lag[5] 1.974655 1.440 1.667 0.4769
## ARCH Lag[7] 2.559148 2.315 1.543 0.6003
##
## Nyblom stability test
## ------------------------------------
## Joint Statistic: 3.4624
## Individual Statistics:
## mu 0.12474
## ar1 0.06729
## omega 0.16800
## alpha1 0.08412
## beta1 0.31534
##
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic: 1.28 1.47 1.88
## Individual Statistic: 0.35 0.47 0.75
##
## Sign Bias Test
## ------------------------------------
## t-value prob sig
## Sign Bias 2.7634 0.005801 ***
## Negative Sign Bias 0.8208 0.411909
## Positive Sign Bias 0.5391 0.589911
## Joint Effect 13.8664 0.003093 ***
##
##
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
## group statistic p-value(g-1)
## 1 20 60.64 3.068e-06
## 2 30 70.64 2.480e-05
## 3 40 84.24 3.598e-05
## 4 50 109.62 1.559e-06
##
##
## Elapsed time : 0.1903021
ugfit@fit$coef
## mu ar1 omega alpha1 beta1
## 1.033440e-03 -6.218212e-02 5.002358e-06 2.302012e-01 7.566016e-01
ug_var <- ugfit@fit$var %>%ts()
ug_res2 <- (ugfit@fit$residuals)^2 %>% ts()
unloadNamespace('ggfortify')
unloadNamespace('forecast')
library(forecast)
autoplot(ug_var,colour = 'red')+
autolayer(ug_res2,color = '#4373B6')+
ggtitle("Plot of Variance and Residual Squared") +
theme_bw()

ugfore <- ugarchforecast(ugfit,n.ahead = 10)
ugfore
##
## *------------------------------------*
## * GARCH Model Forecast *
## *------------------------------------*
## Model: sGARCH
## Horizon: 10
## Roll Steps: 0
## Out of Sample: 0
##
## 0-roll forecast [T0=2023-03-30]:
## Series Sigma
## T+1 0.0007347 0.009702
## T+2 0.0010520 0.009894
## T+3 0.0010323 0.010079
## T+4 0.0010335 0.010259
## T+5 0.0010334 0.010434
## T+6 0.0010334 0.010604
## T+7 0.0010334 0.010768
## T+8 0.0010334 0.010928
## T+9 0.0010334 0.011084
## T+10 0.0010334 0.011235
#Conditional volatility fore next 10 days
ug_f <- ugfore@forecast$sigmaFor %>%ts()
autoplot(ug_f,colour = '#4373B6')+
labs(x='Day',title='10 Period Forecast of Conditional Volatility')+
theme_bw()

ug_var_t <- c(tail(ug_var,20),rep(NA,10)) %>% ts() # gets the last 20 observations
ug_res2_t <- c(tail(ug_res2,20),rep(NA,10)) %>% ts() # gets the last 20 observations
ug_f <- c(rep(NA,20),(ug_f)^2) %>% ts()
ug_f
## Time Series:
## Start = 1
## End = 30
## Frequency = 1
## [1] NA NA NA NA NA
## [6] NA NA NA NA NA
## [11] NA NA NA NA NA
## [16] NA NA NA NA NA
## [21] 0.0000941231 0.0000978833 0.0001015939 0.0001052555 0.0001088688
## [26] 0.0001124344 0.0001159529 0.0001194250 0.0001228513 0.0001262323
autoplot(ug_res2_t)+
autolayer(ug_f,series = 'Forecast')+
autolayer(ug_var_t,series = 'Actual Variance')+
scale_colour_manual(name='Series', values=c('#4373B6','red'))+
theme_bw()+
labs(x='Period',title='Residual^2, Variance, and Forecast')+
theme(legend.position=c(0.98, 0.98),
legend.justification=c('right', 'top'))
## Warning: Removed 20 rows containing missing values (`geom_line()`).
## Warning: Removed 10 rows containing missing values (`geom_line()`).
