Importing the Quantmod Package

We begin by loading the quantmod package, which is short for Quantitative Financial Modelling Framework. This package provides a convenient interface to retrieve financial data from sources like Yahoo Finance, and includes tools for modeling, charting, and technical analysis.

library(quantmod)

Importing Google Stock Price Data

In this step, we use the getSymbols() function from the quantmod package to download historical stock price data for Google (ticker symbol: “GOOGL”) from Yahoo Finance.

We then display the first 5 rows and the first 5 columns of the dataset using head() for a quick preview of the data.

googl <- getSymbols("GOOGL", from="2010-01-01", to="2022-08-03", auto.assign = F)
head(googl[, 1:5], 5)
##            GOOGL.Open GOOGL.High GOOGL.Low GOOGL.Close GOOGL.Volume
## 2010-01-04   15.68944   15.75350  15.62162    15.68443     78169752
## 2010-01-05   15.69520   15.71171  15.55405    15.61537    120067812
## 2010-01-06   15.66216   15.66216  15.17417    15.22172    158988852
## 2010-01-07   15.25025   15.26527  14.83108    14.86737    256315428
## 2010-01-08   14.81481   15.09635  14.74249    15.06557    188783028

Calculating Daily Returns

In this chunk, we calculate the daily returns based on the closing prices of Google’s stock.

The daily return is computed using the formula:

\[ \text{Return}_t = \frac{P_t - P_{t-1}}{P_{t-1}} \]

Where: - \(P_t\) is the closing price at time \(t\), - \(P_{t-1}\) is the closing price on the previous trading day.

The lag() function from the stats package is used to access the previous day’s closing price. This allows us to calculate the percentage change between consecutive days.

daily_ret <- (googl$GOOGL.Close-stats::lag(googl$GOOGL.Close))/stats::lag(googl$GOOGL.Close)

Formatting the Daily Returns Data

We now convert the daily returns from a time series object to a data frame format for easier manipulation and visualization.

daily_ret <- data.frame(index(daily_ret), daily_ret)
colnames(daily_ret) <- c("date", "return")
rownames(daily_ret) <- 1:nrow(daily_ret)

Plotting Daily Returns Over Time

We use the ggplot2 package to create a line chart that visualizes the daily returns over time.

This visualization helps identify trends, patterns, and periods of high or low volatility in the daily returns.

library(ggplot2)
p1 <- ggplot(daily_ret, aes(x=date, y=return))
p1 + geom_line(colour="steelblue")

Histogram of Daily Returns with Normal Distribution Curve

This plot visualizes the distribution of daily returns using a histogram and overlays a normal distribution curve for comparison.

This plot allows us to assess whether the daily returns are approximately normally distributed.

p2 <- ggplot(daily_ret) 

p2 + geom_histogram(aes(x=return, y=..density..), binwidth = 0.005, color="steelblue", fill="grey", size=1) +
  stat_function(fun = dnorm, args = list(mean = mean(daily_ret$return, na.rm = T), sd = sd(daily_ret$return, na.rm = T)), size=1)

Calculating Realized Volatility

In this section, we estimate the realized volatility of Google’s daily returns using a rolling window.

The result, realizedvol, contains the rolling annualized volatility values.

library(PerformanceAnalytics)
library(xts)
daily_ret_xts <- xts(daily_ret[,-1], order.by=daily_ret[,1])
realizedvol <- rollapply(daily_ret_xts, width = 20, FUN=sd.annualized)

Formatting the Volatility Data

We now convert the rolling volatility time series into a data frame for easier plotting and manipulation.

This formatted data will be used in the next step to visualize how volatility changes over time.

vol <- data.frame(index(realizedvol), realizedvol)
colnames(vol) <- c("date", "volatility")

Plotting Realized Volatility Over Time

We use ggplot2 to visualize the realized volatility over time with a line chart.

This plot shows how the volatility of Google’s daily returns changes, helping to identify periods of higher or lower market uncertainty.

p3 <- ggplot(vol, aes(x=date, y=volatility))
p3 +
  geom_line( color="steelblue")

Specifying the GARCH Model

In this section, we specify the Generalized Autoregressive Conditional Heteroskedasticity (GARCH) model for volatility forecasting using the rugarch package.

This specification defines the model to estimate the conditional volatility based on past values and their variances.

library(rugarch)
garch_spec <- ugarchspec(variance.model=list(model="sGARCH", garchOrder=c(1,1)), mean.model=list(armaOrder=c(0,0)))

Fitting the GARCH Model

Next, we fit the GARCH model to the volatility data using the ugarchfit function from the rugarch package.

The result is stored in fit_garch, which contains the estimated parameters and model diagnostics.

Finally, fit_garch displays the fitted model, showing the estimated GARCH parameters, including the conditional variance, and provides other diagnostic information.

fit_garch <- ugarchfit(spec = garch_spec, data = vol[-c(1:19),2])
fit_garch
## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : sGARCH(1,1)
## Mean Model   : ARFIMA(0,0,0)
## Distribution : norm 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error    t value Pr(>|t|)
## mu      0.183618    0.001297 141.605668        0
## omega   0.000276    0.000027  10.226852        0
## alpha1  0.999000    0.033681  29.660723        0
## beta1   0.000000    0.007195   0.000002        1
## 
## Robust Standard Errors:
##         Estimate  Std. Error   t value Pr(>|t|)
## mu      0.183618    0.010063 18.247385 0.000000
## omega   0.000276    0.000121  2.285627 0.022276
## alpha1  0.999000    0.084219 11.861995 0.000000
## beta1   0.000000    0.005298  0.000003 0.999998
## 
## LogLikelihood : 4348.209 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -2.7600
## Bayes        -2.7523
## Shibata      -2.7600
## Hannan-Quinn -2.7572
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                       1862       0
## Lag[2*(p+q)+(p+q)-1][2]      2695       0
## Lag[4*(p+q)+(p+q)-1][5]      4905       0
## d.o.f=0
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                     0.1021  0.7493
## Lag[2*(p+q)+(p+q)-1][5]    0.1275  0.9969
## Lag[4*(p+q)+(p+q)-1][9]    0.2006  0.9999
## d.o.f=2
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[3]  0.004844 0.500 2.000  0.9445
## ARCH Lag[5]  0.040714 1.440 1.667  0.9963
## ARCH Lag[7]  0.062792 2.315 1.543  0.9998
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  1.2256
## Individual Statistics:             
## mu     0.2599
## omega  0.4897
## alpha1 0.6476
## beta1  0.4083
## 
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic:          1.07 1.24 1.6
## Individual Statistic:     0.35 0.47 0.75
## 
## Sign Bias Test
## ------------------------------------
##                    t-value   prob sig
## Sign Bias           1.2672 0.2052    
## Negative Sign Bias  0.3633 0.7164    
## Positive Sign Bias  1.0105 0.3123    
## Joint Effect        2.3665 0.4999    
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20      7769            0
## 2    30      7843            0
## 3    40      9979            0
## 4    50      8988            0
## 
## 
## Elapsed time : 0.348263