1. Introduction
2. Methodology
3. Data Preprocessing
4. Descriptive Analysis
- Time Series plot
- Scatter plot
- Normality check
5.Modelling
- 5.1 Linear model
- Diagnostic checking for linear model
- 5.2 Quadratic model
- Diagnostic checking of quadratic model
6.ARIMA Model
Overcome the non-stationarity of the series
- BoxCox Transformation
- Differencing - first differencing
ARIMA model specification
7. Parameter estimation
Fit models and find estimations
- Compare scores of all models
8. GARCH model
- Order specification
a) Calculate for GARCH(1,1)
b) Calculate for GARCH(1,3)
d) Calculate for GARCH(2,3)
Estimated Conditional Variance Prediction
9. Final Best Model Selection for ARMA+GARCH
1. ARMA(1,1) + GARCH(1,1)
2. ARMA(0,1) + GARCH(1,1)
3.ARMA(1,2) + GARCH(1,1)
4. ARMA(1,1) + GARCH(2,1)
FORCASTING WITH OUR BEST COMBO MODEL FOR NEXT 10 YEARS
Conclusion

Group member:

Aparupa Mitra-s3831724, Shristi Shelendra Chavhan-s3822713, George Chai-s3533832, Mary Legrand-s3815368, Yuxin Yang- s375519

1. Introduction

This study is a time series analysis of monthly data on the closing price of Amazon’s stock price during the 12 years from January 1, 2009, to May 29, 2021. Amazon is one of the biggest competitors in the U.S. information technology industry. Time series analysis is one of the important tools of economic research. It can describe the law of historical data changes over time and can be used to forecast economic rain. For example, in the stock market, time series forecasting is often used to predict stock price trends and provide investors with a basis for decision-making. This report fits the corresponding ARMA model, analyzes and compares the models fitted by the EACF and BIC chart test, and conducts the AIC criterion test to find the optimal ARMA model ARMA (1,1,1 ). Finally, it concludes that the development trend of Amazon’s stock price is stable in the future.

2. Methodology

library(tseries)
library(timeSeries)
library(forecast)
library(TSA)
library(tidyr)
library(lmtest)
library(FSAdata)
library(fUnitRoots)
library(CombMSC)
library(fGarch)
library(FitAR)
library(tswge)
library(rugarch)
library(dLagM)
library(ggplot2)
library(fUnitRoots)
library(knitr)
library(lattice)
library(bestglm)
library(leaps)
library(ltsa)
library(CombMSC)

3. Data Preprocessing

You can also embed plots, for example:

amazon <- read.csv("AMZN.csv", header = TRUE)
head(amazon)

amzn <- amazon[c(1,5)] #subsetting to column date and closing price

amzn$Date <- as.Date(amzn$Date, format = "%Y-%m-%d")

amazon_close <- ts(amzn$Close, start = c(2009,1), end = c(2021,5), frequency =12) #convert to time series object

class(amazon_close)

## [1] "ts"

4. Descriptive Analysis

Time Series plot

plot(amazon_close, ylab="Closing Price (USD)",
     xlab ="Date", main="Amazon Monthly Closing Stock Price")

Trend: We could see a positive trend over given period of time.

Seasonality : No repeating pattern can be observed, hence no seasonality present

Change in Variance: Almost no change in variance can be observed.

Intervention: No apparent Intervention can be observed.

Scatter plot

plot(y=amazon_close, x=zlag(amazon_close),
     main="Scatter plot of Amazon closing stock price in consecutive months",
     xlab="Previous Month Closing price (USD)",
     ylab="Currenty Month Closing Price (USD)")

Our scatter plot shows a positive correlation between current month closing price and previous month closing price.

Normality check

qqnorm(amazon_close)
qqline(amazon_close, col=2, lwd=1)

shapiro.test(amazon_close)

## 
##  Shapiro-Wilk normality test
## 
## data:  amazon_close
## W = 0.87161, p-value = 4.831e-10

5.Modelling

5.1 Linear model

model1 <- lm(amazon_close~time(amazon_close))
summary(model1)

## 
## Call:
## lm(formula = amazon_close ~ time(amazon_close))
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -392.30 -142.91    4.66  146.01  480.13 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        -88652.148   8937.699  -9.919   <2e-16 ***
## time(amazon_close)     45.514      4.435  10.262   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 194 on 147 degrees of freedom
## Multiple R-squared:  0.4174, Adjusted R-squared:  0.4134 
## F-statistic: 105.3 on 1 and 147 DF,  p-value: < 2.2e-16

plot(amazon_close, ylab="Closing Price in USD",
     xlab ="", main="fitted linear model", type="o")
abline(model1, col="red")

Diagnostic checking for linear model

res.model1 = rstudent(model1)

plot(res.model1, ylab="Closing Price in USD",
     xlab ="", main="fitted residuals of linear model", type="o")
lines(as.vector(amazon_close), type="o")

par(mfrow=c(1,2))
acf(res.model1, main="ACF of standardized residuals")
pacf(res.model1, main="PACF of standardized residuals")

hist(res.model1, xlab='Standardized Residuals')

qqnorm(res.model1)
qqline(res.model1, col = 2, lwd = 1, lty = 2)

5.2 Quadratic model

par(mfrow=c(1,1))

t=time(amazon_close)
t2 = t^2

model2 <- lm(amazon_close~t+t2)
summary(model2)

## 
## Call:
## lm(formula = amazon_close ~ t + t2)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -298.30  -83.83   -6.96   77.93  317.03 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -5.209e+07  3.641e+06  -14.31   <2e-16 ***
## t            5.166e+04  3.614e+03   14.29   <2e-16 ***
## t2          -1.281e+01  8.967e-01  -14.28   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 125.8 on 146 degrees of freedom
## Multiple R-squared:  0.757,  Adjusted R-squared:  0.7536 
## F-statistic: 227.4 on 2 and 146 DF,  p-value: < 2.2e-16

plot(ts(fitted(model2)), ylim = c(min(c(fitted(model2),as.vector(amazon_close)))
                                  ,max(c(fitted(model2),as.vector(amazon_close))))
     ,ylab='Closing Price (USD)', main= "fitted quadratic curve")
lines(as.vector(amazon_close), type="o")

Diagnostic checking of quadratic model

res.model2 = rstudent(model2)

plot(res.model2, ylab="Closing Price in USD",
     xlab ="", main="fitted residuals of quadratic model", type="o")
lines(as.vector(amazon_close), type="o")

par(mfrow=c(2,2))
acf(res.model2, main="ACF of standardized residuals")
pacf(res.model2, main="PACF of standardized residuals")

hist(res.model2, xlab='Standardized Residuals')

qqnorm(res.model2)
qqline(res.model2, col = 2, lwd = 1, lty = 2)

shapiro.test(res.model2)

## 
##  Shapiro-Wilk normality test
## 
## data:  res.model2
## W = 0.98371, p-value = 0.07557

##Observation : We found that out of both model ,quadratic model is the best model .Now we will proceed with ARIMA model

6.ARIMA Model

par(mfrow=c(1,2))
acf(amazon_close)
pacf(amazon_close)

par(mfrow=c(1,1))
adf.test(amazon_close)

## 
##  Augmented Dickey-Fuller Test
## 
## data:  amazon_close
## Dickey-Fuller = -2.2967, Lag order = 5, p-value = 0.4528
## alternative hypothesis: stationary

Overcome the non-stationarity of the series

BoxCox Transformation

amazon_close.tr = BoxCox.ar(amazon_close+abs(min(amazon_close))+0.1, lambda = c(-1,1))

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

amazon_close.tr$ci

## [1] -1 -1

par(mfrow=c(1,2))

par(mfrow=c(2,2))
acf(log(amazon_close))
pacf(log(amazon_close))

adf.test(log(amazon_close))

## 
##  Augmented Dickey-Fuller Test
## 
## data:  log(amazon_close)
## Dickey-Fuller = -2.3557, Lag order = 5, p-value = 0.4282
## alternative hypothesis: stationary

par(mfrow=c(1,2))
hist(log(amazon_close))
qqnorm(log(amazon_close), main = 'QQ plot for the natural log')
qqline(log(amazon_close), col="red")

Differencing - first differencing

par(mfrow=c(1,1))

diff.log.amazon_close = diff(log(amazon_close), difference = 1)
plot(diff.log.amazon_close, type = 'o', ylab = 'Closing price')

order = ar(diff(diff.log.amazon_close))$order
adfTest(diff.log.amazon_close, lags = order, title = NULL, description = NULL)

## Warning in adfTest(diff.log.amazon_close, lags = order, title = NULL,
## description = NULL): p-value smaller than printed p-value

## 
## Title:
##  Augmented Dickey-Fuller Test
## 
## Test Results:
##   PARAMETER:
##     Lag Order: 8
##   STATISTIC:
##     Dickey-Fuller: -3.5612
##   P VALUE:
##     0.01 
## 
## Description:
##  Mon Jun 14 12:19:58 2021 by user: aparu

McLeod.Li.test(y=diff.log.amazon_close, main="McLeod-Li Test Statistics for Amazon stock price")

ARIMA model specification

par(mfrow = c(1,2))
acf(diff.log.amazon_close)
pacf(diff.log.amazon_close)

eacf(diff.log.amazon_close)

## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 o o x o o o o o o o o  o  o  o 
## 1 x o x o o o o o o o o  o  o  o 
## 2 o x x o o o o o o o o  o  o  o 
## 3 o x x o o o o o o o o  o  o  o 
## 4 x o x o o o o o o o o  o  o  o 
## 5 x x x o o o o o o o o  o  o  o 
## 6 x o x o o x o o o o o  o  o  o 
## 7 x o o x o x o o o o o  o  o  o

par(mfrow = c(1,1))

res = armasubsets(y=diff.log.amazon_close, nar = 13, nma = 13, y.name = 'test', ar.method = 'ols')

## Warning in leaps.setup(x, y, wt = wt, nbest = nbest, nvmax = nvmax, force.in =
## force.in, : 10 linear dependencies found

## Reordering variables and trying again:

plot(res)

7. Parameter estimation

Here we will apply two methods to estimate parameters

Least Squares Estimate of the AR coefficient with significance test, method=‘CSS’

Maximum likelihood estimate of the AR coefficient with significance tests, method=“ML”

Fit models and find estimations

From the output of EACF, we can include ARIMA(1,1,1) ARIMA(1,1,3), ARIMA(2,1,4) For BIC ARIMA(1,1,1) ARIMA(3,1,6) ARIMA(3,1,7)

#ARIMA(1,1,1)

model_111_css = arima(log(amazon_close),order=c(1,1,1),method='CSS')
coeftest(model_111_css)

## Warning in sqrt(diag(se)): NaNs produced

## 
## z test of coefficients:
## 
##      Estimate Std. Error z value Pr(>|z|)
## ar1 -0.025595         NA      NA       NA
## ma1 -0.012186         NA      NA       NA

model_111_ml = arima(log(amazon_close),order=c(1,1,1),method='ML')
coeftest(model_111_ml)

## 
## z test of coefficients:
## 
##     Estimate Std. Error z value  Pr(>|z|)    
## ar1 -0.77646    0.17713 -4.3834 1.168e-05 ***
## ma1  0.70386    0.19400  3.6281 0.0002855 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

P value which are greater than 0.05 for AR(1) which means AR(1) is insignificant in ML methods.

#ARIMA(1,1,3)

model_113_css = arima(log(amazon_close),order=c(1,1,3),method='CSS')
coeftest(model_113_css)

## 
## z test of coefficients:
## 
##      Estimate Std. Error z value Pr(>|z|)  
## ar1 -0.489684   0.322749 -1.5172  0.12921  
## ma1  0.463177   0.326401  1.4190  0.15589  
## ma2 -0.019014   0.092170 -0.2063  0.83656  
## ma3 -0.209616   0.082861 -2.5297  0.01142 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

model_113_ml = arima(log(amazon_close),order=c(1,1,3),method='ML')
coeftest(model_113_ml)

## 
## z test of coefficients:
## 
##      Estimate Std. Error z value Pr(>|z|)  
## ar1 -0.500294   0.302733 -1.6526  0.09841 .
## ma1  0.474449   0.307585  1.5425  0.12295  
## ma2 -0.019245   0.092333 -0.2084  0.83489  
## ma3 -0.207147   0.082883 -2.4993  0.01245 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

P value which are greater than 0.05 for AR(1),MA(1),MA(3) which means AR(1),MA(1),MA(3) are insignificant in CSS and ML methods.

#ARIMA(2,1,4)

model_214_css = arima(log(amazon_close),order=c(2,1,4),method='CSS')
coeftest(model_214_css)

## 
## z test of coefficients:
## 
##      Estimate Std. Error z value Pr(>|z|)  
## ar1 -0.490455   0.875535 -0.5602  0.57536  
## ar2  0.056725   0.498020  0.1139  0.90932  
## ma1  0.468408   0.870357  0.5382  0.59045  
## ma2 -0.065935   0.469723 -0.1404  0.88837  
## ma3 -0.217357   0.084889 -2.5605  0.01045 *
## ma4 -0.026032   0.183263 -0.1420  0.88704  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

model_214_ml = arima(log(amazon_close),order=c(2,1,4),method='ML')
coeftest(model_214_ml)

## 
## z test of coefficients:
## 
##       Estimate Std. Error z value Pr(>|z|)   
## ar1 -0.3081751  0.9464917 -0.3256 0.744729   
## ar2  0.1536748  0.4793536  0.3206 0.748523   
## ma1  0.2863119  0.9441079  0.3033 0.761690   
## ma2 -0.1579830  0.4579295 -0.3450 0.730099   
## ma3 -0.2132556  0.0802162 -2.6585 0.007849 **
## ma4  0.0073836  0.2256988  0.0327 0.973902   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

P value which are greater than 0.05 for AR(1),AR(2),MA(1),MA(4) which means AR(1),AR(2),MA(1),MA(4) are insignificant in CSS and ML methods.

#ARIMA(3,1,6)

model_316_css = arima(log(amazon_close),order=c(3,1,6),method='CSS')
coeftest(model_316_css)

## Warning in sqrt(diag(se)): NaNs produced

## 
## z test of coefficients:
## 
##      Estimate Std. Error z value  Pr(>|z|)    
## ar1  0.030011         NA      NA        NA    
## ar2  0.476056   0.197017  2.4163   0.01568 *  
## ar3  0.497080   0.042079 11.8131 < 2.2e-16 ***
## ma1 -0.068899         NA      NA        NA    
## ma2 -0.498385   0.207909 -2.3971   0.01652 *  
## ma3 -0.750814   0.123327 -6.0880 1.143e-09 ***
## ma4  0.059880   0.078896  0.7590   0.44787    
## ma5  0.014538   0.095312  0.1525   0.87877    
## ma6  0.192103   0.113907  1.6865   0.09170 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

model_316_ml = arima(log(amazon_close),order=c(3,1,6),method='ML')
coeftest(model_316_ml)

## 
## z test of coefficients:
## 
##      Estimate Std. Error  z value  Pr(>|z|)    
## ar1 -0.352017   0.065093  -5.4080 6.375e-08 ***
## ar2  0.392475   0.073477   5.3414 9.221e-08 ***
## ar3  0.905667   0.059062  15.3340 < 2.2e-16 ***
## ma1  0.364772   0.106565   3.4230 0.0006194 ***
## ma2 -0.386915   0.113740  -3.4017 0.0006696 ***
## ma3 -1.164197   0.114232 -10.1915 < 2.2e-16 ***
## ma4 -0.039762   0.098980  -0.4017 0.6878961    
## ma5  0.055169   0.090847   0.6073 0.5436701    
## ma6  0.251977   0.087851   2.8682 0.0041276 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#ARIMA(3,1,7)

model_317_css = arima(log(amazon_close),order=c(3,1,7),method='CSS')
coeftest(model_317_css)

## 
## z test of coefficients:
## 
##       Estimate Std. Error z value  Pr(>|z|)    
## ar1  0.2386112  0.1284911  1.8570  0.063308 .  
## ar2  0.2514076  0.2582122  0.9736  0.330232    
## ar3  0.5093192  0.1354101  3.7613  0.000169 ***
## ma1 -0.2890431  0.1249849 -2.3126  0.020743 *  
## ma2 -0.2707534  0.2462631 -1.0994  0.271573    
## ma3 -0.7900066  0.1608474 -4.9115 9.037e-07 ***
## ma4  0.1545320  0.1022290  1.5116  0.130629    
## ma5  0.0030064  0.0880398  0.0341  0.972759    
## ma6  0.2506322  0.0942791  2.6584  0.007851 ** 
## ma7 -0.1222021  0.0931393 -1.3120  0.189508    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

model_317_ml = arima(log(amazon_close),order=c(3,1,7),method='ML')
coeftest(model_317_ml)

## 
## z test of coefficients:
## 
##      Estimate Std. Error  z value  Pr(>|z|)    
## ar1 -0.352661   0.064014  -5.5091 3.607e-08 ***
## ar2  0.387572   0.073475   5.2749 1.328e-07 ***
## ar3  0.910720   0.061891  14.7148 < 2.2e-16 ***
## ma1  0.360671   0.108401   3.3272 0.0008773 ***
## ma2 -0.385374   0.114122  -3.3769 0.0007332 ***
## ma3 -1.167902   0.114354 -10.2130 < 2.2e-16 ***
## ma4 -0.014120   0.132590  -0.1065 0.9151912    
## ma5  0.066355   0.101689   0.6525 0.5140577    
## ma6  0.245681   0.089457   2.7463 0.0060264 ** 
## ma7 -0.025769   0.087939  -0.2930 0.7695013    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

P values are less than 0.05 for AR3, MA1,MA2,MA6 and are significant for both ML and CSS methods

Compare scores of all models

AIC(model_111_ml,model_113_ml,model_214_ml,model_316_ml ,model_317_ml)

If we focused on only AIC value here, from the result, we can see model model_111_ml, ARIMA(1,1,1) has the lowest AIC value. Thus, we chose it to be the best fitting ARIMA model.

#checkresiduals(fit)
res111 =residuals(model_111_ml)

par(mfrow=c(2,2))

plot(res111,type='o',ylab='Standardized residuals',
     main="Time series plot of standardized residuals")
abline(h=0)
hist(res111,main="Histogram of standardized residuals")
qqnorm(res111,main="QQ plot of standardized residuals")
qqline(res111, col = 2)
acf(res111,main="ACF of standardized residuals")

print(shapiro.test(res111))

## 
##  Shapiro-Wilk normality test
## 
## data:  res111
## W = 0.98767, p-value = 0.2106

par(mfrow=c(1,1))
m111.amzn = arima(
    amazon_close,
    order=c(1,1,1),
    xreg=data.frame(constant=seq(amazon_close))
)
m111.amzn

## 
## Call:
## arima(x = amazon_close, order = c(1, 1, 1), xreg = data.frame(constant = seq(amazon_close)))
## 
## Coefficients:
##           ar1     ma1  constant
##       -0.7743  0.6967    5.1693
## s.e.   0.1711  0.1882    5.3540
## 
## sigma^2 estimated as 4638:  log likelihood = -834.73,  aic = 1675.46

##Forcating

n = length(amazon_close)
n.ahead = 5
newxreg = data.frame(constant=(n+1):(n+n.ahead))
plot(
    m111.amzn,
    n.ahead=n.ahead,
    newxreg=newxreg,
    ylab='Closing Price (USD)',
    xlab='Date',
    main="AMZN stock price forecast 5 months ahead",
)

8. GARCH model

r.amazon_close= diff(log(amazon_close))

# 1. ACF and PACF 

par(mfrow=c(3,2))
acf(r.amazon_close,main="Sample ACF of lg Amazon stock price change")

pacf(r.amazon_close,main="Sample PACF of lg Daily Amazon stock price change")

acf(abs(r.amazon_close),main="Sample ACF of the Absolute Daily Amazon stock price")

pacf(abs(r.amazon_close),main="Sample PACF of the Absolute Daily Amazon stock price")

acf(r.amazon_close^2,main="Sample ACF of the Squared Daily Amazon stock price")

pacf(r.amazon_close^2,main="Sample PACF of the Squared Daily Amazon stock price")

par(mfrow=c(1,1))

McLeod.Li.test(y=abs(r.amazon_close),main="McLeod-Li Test Statistics
               for Absolute Daily Amazon stock price ")

qqnorm(abs(r.amazon_close), main="Q-Q Normal Plot of the Absolute Daily Amazon stock price")
qqline(abs(r.amazon_close))

shapiro.test(abs(r.amazon_close))

## 
##  Shapiro-Wilk normality test
## 
## data:  abs(r.amazon_close)
## W = 0.88163, p-value = 1.653e-09

McLeod.Li.test(y=r.amazon_close^2,main="McLeod-Li Test Statistics for the Squared Daily Amazon stock price")

qqnorm(r.amazon_close^2, main="Q-Q Normal Plot of the Squared Daily Amazon price")
qqline(r.amazon_close^2)

shapiro.test(r.amazon_close^2)

## 
##  Shapiro-Wilk normality test
## 
## data:  r.amazon_close^2
## W = 0.61108, p-value < 2.2e-16

Order specification

eacf(r.amazon_close)

## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 o o x o o o o o o o o  o  o  o 
## 1 x o x o o o o o o o o  o  o  o 
## 2 o x x o o o o o o o o  o  o  o 
## 3 o x x o o o o o o o o  o  o  o 
## 4 x o x o o o o o o o o  o  o  o 
## 5 x x x o o o o o o o o  o  o  o 
## 6 x o x o o x o o o o o  o  o  o 
## 7 x o o x o x o o o o o  o  o  o

eacf(abs(r.amazon_close))

## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 o o x o o o o o o o o  o  o  o 
## 1 x o x o o o o o o o o  o  o  o 
## 2 x x x o o o o o o o o  o  o  o 
## 3 x x o o o o o o o o o  o  o  o 
## 4 x o o o o o o o o o o  o  o  o 
## 5 x o o o o o o o o o o  o  o  o 
## 6 x x x o o o o o o o o  o  o  o 
## 7 x x o o o o x o o o o  o  o  o

eacf(r.amazon_close^2)

## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 o o x o o o o o o o o  o  o  o 
## 1 x o x o o o o o o o o  o  o  o 
## 2 o o x o o o o o o o o  o  o  o 
## 3 x o x o o o o o o o o  o  o  o 
## 4 x o x o o o o o o o o  o  o  o 
## 5 x x x o o o o o o o o  o  o  o 
## 6 x x x o o o o o o o o  o  o  o 
## 7 x o o x o o x o o o o  o  o  o

EACF plot of applying absolute value suggests GARCH(1,1) and GARCH(1,3). EACF plot of applying square transformations suggests GARCH(2,1), GARCH(2,3).

a) Calculate for GARCH(1,1)

garch.11 = garch(r.amazon_close,order=c(1,1),trace = FALSE)
summary(garch.11)

## 
## Call:
## garch(x = r.amazon_close, order = c(1, 1), trace = FALSE)
## 
## Model:
## GARCH(1,1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.61843 -0.45838  0.05956  0.71184  3.52225 
## 
## Coefficient(s):
##     Estimate  Std. Error  t value Pr(>|t|)    
## a0 4.085e-05   3.867e-05    1.056   0.2909    
## a1 1.247e-01   7.376e-02    1.691   0.0908 .  
## b1 7.992e-01   1.053e-01    7.588 3.26e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Diagnostic Tests:
##  Jarque Bera Test
## 
## data:  Residuals
## X-squared = 5.6135, df = 2, p-value = 0.0604
## 
## 
##  Box-Ljung test
## 
## data:  Squared.Residuals
## X-squared = 0.45437, df = 1, p-value = 0.5003

garch_resi11 =residuals(garch.11)
par(mfrow=c(2,2))

plot(garch_resi11,type='o',ylab='Standardized residuals', main="Time series plot of standardized residuals")
abline(h=0)
hist(garch_resi11,main="Histogram of standardized residuals")
qqnorm(garch_resi11,main="QQ plot of standardized residuals")
qqline(garch_resi11, col = 2)


print(shapiro.test(garch_resi11))

## 
##  Shapiro-Wilk normality test
## 
## data:  garch_resi11
## W = 0.98904, p-value = 0.3045

b) Calculate for GARCH(1,3)

garch.13 = garch(r.amazon_close,order=c(1,3),trace = FALSE)
summary(garch.13)

## 
## Call:
## garch(x = r.amazon_close, order = c(1, 3), trace = FALSE)
## 
## Model:
## GARCH(1,3)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.48532 -0.49467  0.06862  0.70042  3.65002 
## 
## Coefficient(s):
##     Estimate  Std. Error  t value Pr(>|t|)  
## a0 3.066e-04   1.879e-04    1.632   0.1028  
## a1 5.661e-02   1.127e-01    0.502   0.6154  
## a2 5.889e-02   8.015e-02    0.735   0.4625  
## a3 1.986e-01   1.205e-01    1.647   0.0995 .
## b1 4.268e-15   4.985e-01    0.000   1.0000  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Diagnostic Tests:
##  Jarque Bera Test
## 
## data:  Residuals
## X-squared = 3.3332, df = 2, p-value = 0.1889
## 
## 
##  Box-Ljung test
## 
## data:  Squared.Residuals
## X-squared = 0.022105, df = 1, p-value = 0.8818

garch_resi13 =residuals(garch.13)
par(mfrow=c(2,2))

plot(garch_resi13,type='o',ylab='Standardized residuals', main="Time series plot of standardized residuals")
abline(h=0)
hist(garch_resi13,main="Histogram of standardized residuals")
qqnorm(garch_resi13,main="QQ plot of standardized residuals")
qqline(garch_resi13, col = 2)

print(shapiro.test(garch_resi13))

## 
##  Shapiro-Wilk normality test
## 
## data:  garch_resi13
## W = 0.99204, p-value = 0.5947

# c) Calculate for GARCH(2,1)

garch.21 = garch(r.amazon_close,order=c(2,1),trace = FALSE)
summary(garch.21)

## 
## Call:
## garch(x = r.amazon_close, order = c(2, 1), trace = FALSE)
## 
## Model:
## GARCH(2,1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.68729 -0.47989  0.07003  0.71603  3.46673 
## 
## Coefficient(s):
##     Estimate  Std. Error  t value Pr(>|t|)  
## a0 8.928e-05   1.008e-04    0.886    0.376  
## a1 1.803e-01   1.301e-01    1.385    0.166  
## b1 1.668e-01   3.779e-01    0.441    0.659  
## b2 4.819e-01   2.797e-01    1.723    0.085 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Diagnostic Tests:
##  Jarque Bera Test
## 
## data:  Residuals
## X-squared = 5.6007, df = 2, p-value = 0.06079
## 
## 
##  Box-Ljung test
## 
## data:  Squared.Residuals
## X-squared = 0.77037, df = 1, p-value = 0.3801

garch_resi21 =residuals(garch.21)
par(mfrow=c(2,2))

plot(garch_resi21,type='o',ylab='Standardized residuals', main="Time series plot of standardized residuals")
abline(h=0)
hist(garch_resi21,main="Histogram of standardized residuals")
qqnorm(garch_resi21,main="QQ plot of standardized residuals")
qqline(garch_resi21, col = 2)

print(shapiro.test(garch_resi21))

## 
##  Shapiro-Wilk normality test
## 
## data:  garch_resi21
## W = 0.98895, p-value = 0.303

d) Calculate for GARCH(2,3)

garch.23 = garch(r.amazon_close,order=c(2,3),trace = FALSE)
summary(garch.23)

## 
## Call:
## garch(x = r.amazon_close, order = c(2, 3), trace = FALSE)
## 
## Model:
## GARCH(2,3)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.48720 -0.49252  0.06862  0.69609  3.62610 
## 
## Coefficient(s):
##     Estimate  Std. Error  t value Pr(>|t|)
## a0 2.940e-04   2.355e-04    1.248    0.212
## a1 4.828e-02   1.075e-01    0.449    0.653
## a2 5.709e-02   8.084e-02    0.706    0.480
## a3 1.978e-01   1.390e-01    1.423    0.155
## b1 7.110e-06   5.080e-01    0.000    1.000
## b2 3.811e-02   4.932e-01    0.077    0.938
## 
## Diagnostic Tests:
##  Jarque Bera Test
## 
## data:  Residuals
## X-squared = 3.1043, df = 2, p-value = 0.2118
## 
## 
##  Box-Ljung test
## 
## data:  Squared.Residuals
## X-squared = 0.0048468, df = 1, p-value = 0.9445

garch_resi23 =residuals(garch.23)
par(mfrow=c(2,2))

plot(garch_resi23,type='o',ylab='Standardized residuals', main="Time series plot of standardized residuals")
abline(h=0)
hist(garch_resi23,main="Histogram of standardized residuals")
qqnorm(garch_resi23,main="QQ plot of standardized residuals")
qqline(garch_resi23, col = 2)

print(shapiro.test(garch_resi23))

## 
##  Shapiro-Wilk normality test
## 
## data:  garch_resi23
## W = 0.99223, p-value = 0.6152

All coefficients are significant at 5% level of significance.So we can go with garch(1,3) as best model out of all

Estimated Conditional Variance Prediction

plot((fitted(garch.13)[,1])^2,type='l',
     ylab='Conditional Variance',
     xlab='t',main="Estimated Conditional Variances of daily Amazon stock price")

9. Final Best Model Selection for ARMA+GARCH

1. ARMA(1,1) + GARCH(1,1)

modelcombo1<-ugarchspec(variance.model = list(model = "sGARCH", garchOrder = c(1, 1)),
                   mean.model = list(armaOrder = c(1, 1), include.mean = FALSE),
                   distribution.model = "norm")

m.l1_best1<-ugarchfit(spec=modelcombo1,data=r.amazon_close,out.sample = 100)

## Warning in .sgarchfit(spec = spec, data = data, out.sample = out.sample, : 
## ugarchfit-->waring: using less than 100 data
##  points for estimation

m.l1_best1

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : sGARCH(1,1)
## Mean Model   : ARFIMA(1,0,1)
## Distribution : norm 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error   t value Pr(>|t|)
## ar1     0.156660    2.766836  0.056621 0.954847
## ma1    -0.148653    2.767050 -0.053722 0.957156
## omega   0.000004    0.000002  2.072838 0.038187
## alpha1  0.000000    0.011016  0.000001 1.000000
## beta1   0.999000    0.032425 30.809460 0.000000
## 
## Robust Standard Errors:
##         Estimate  Std. Error  t value Pr(>|t|)
## ar1     0.156660    0.398668  0.39296  0.69435
## ma1    -0.148653    0.380580 -0.39060  0.69610
## omega   0.000004    0.000015  0.28821  0.77318
## alpha1  0.000000    0.059723  0.00000  1.00000
## beta1   0.999000    0.046123 21.65956  0.00000
## 
## LogLikelihood : 113.8701 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -4.5363
## Bayes        -4.3413
## Shibata      -4.5553
## Hannan-Quinn -4.4626
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                     0.2205  0.6387
## Lag[2*(p+q)+(p+q)-1][5]    2.5983  0.7240
## Lag[4*(p+q)+(p+q)-1][9]    3.5360  0.7939
## d.o.f=2
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                     0.4629  0.4963
## Lag[2*(p+q)+(p+q)-1][5]    2.4577  0.5149
## Lag[4*(p+q)+(p+q)-1][9]    3.7808  0.6262
## d.o.f=2
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[3]     1.834 0.500 2.000  0.1757
## ARCH Lag[5]     3.119 1.440 1.667  0.2729
## ARCH Lag[7]     3.639 2.315 1.543  0.4016
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  10.0176
## Individual Statistics:              
## ar1    0.19425
## ma1    0.19511
## omega  0.08838
## alpha1 0.03974
## beta1  0.06519
## 
## 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           0.5328 0.5969    
## Negative Sign Bias  0.8064 0.4245    
## Positive Sign Bias  0.3458 0.7312    
## Joint Effect        0.8151 0.8459    
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     18.67       0.4784
## 2    30     23.25       0.7650
## 3    40     42.00       0.3422
## 4    50     45.75       0.6057
## 
## 
## Elapsed time : 0.06478786

par(mfrow=c(2,2))
plot(m.l1_best1, which=8)
plot(m.l1_best1, which=9)
plot(m.l1_best1, which=10)
par(mfrow=c(1,1))

Residuals are not white noise, and the distribution of standardised residuals has flooks almost right skewed and seems to be far from normality.

We are going to decrease the AR order by 1, see if it is still significant and whether the residual can be improved

2. ARMA(0,1) + GARCH(1,1)

modelcombo2<-ugarchspec(variance.model = list(model = "sGARCH", garchOrder = c(1, 1)),
                   mean.model = list(armaOrder = c(0, 1), include.mean = FALSE),
                   distribution.model = "norm")

m.01_best2<-ugarchfit(spec=modelcombo2,data=r.amazon_close,out.sample = 100)

## Warning in .sgarchfit(spec = spec, data = data, out.sample = out.sample, : 
## ugarchfit-->waring: using less than 100 data
##  points for estimation

m.01_best2

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : sGARCH(1,1)
## Mean Model   : ARFIMA(0,0,1)
## Distribution : norm 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error   t value Pr(>|t|)
## ma1     0.003616    0.138604  0.026092 0.979184
## omega   0.000004    0.000002  2.188870 0.028606
## alpha1  0.000000    0.009764  0.000000 1.000000
## beta1   0.999000    0.032576 30.667179 0.000000
## 
## Robust Standard Errors:
##         Estimate  Std. Error   t value Pr(>|t|)
## ma1     0.003616    0.084136  0.042984  0.96571
## omega   0.000004    0.000015  0.286399  0.77457
## alpha1  0.000000    0.059767  0.000000  1.00000
## beta1   0.999000    0.046247 21.601575  0.00000
## 
## LogLikelihood : 113.8689 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -4.5779
## Bayes        -4.4219
## Shibata      -4.5904
## Hannan-Quinn -4.5189
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                     0.1953  0.6586
## Lag[2*(p+q)+(p+q)-1][2]    0.1953  0.9989
## Lag[4*(p+q)+(p+q)-1][5]    2.5712  0.5462
## d.o.f=1
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                     0.4797  0.4886
## Lag[2*(p+q)+(p+q)-1][5]    2.4904  0.5079
## Lag[4*(p+q)+(p+q)-1][9]    3.8165  0.6201
## d.o.f=2
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[3]     1.852 0.500 2.000  0.1736
## ARCH Lag[5]     3.140 1.440 1.667  0.2701
## ARCH Lag[7]     3.658 2.315 1.543  0.3986
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  9.9121
## Individual Statistics:              
## ma1    0.21112
## omega  0.08790
## alpha1 0.03960
## beta1  0.06493
## 
## 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           0.5259 0.6017    
## Negative Sign Bias  0.8009 0.4276    
## Positive Sign Bias  0.3511 0.7272    
## Joint Effect        0.8090 0.8473    
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     18.67       0.4784
## 2    30     25.75       0.6388
## 3    40     33.67       0.7114
## 4    50     47.83       0.5204
## 
## 
## Elapsed time : 0.04370284

par(mfrow=c(2,2))
plot(m.01_best2, which=8)
plot(m.01_best2, which=9)
plot(m.01_best2, which=10)
par(mfrow=c(1,1))

Residuals are not white noise, and the distribution of standardised residuals has flooks almost right skewed and seems to be far from normality.

We are going to increase the order of MA components, see if it is still significant and whether the residual can be improved

3.ARMA(1,2) + GARCH(1,1)

modelcombo3<-ugarchspec(variance.model = list(model = "sGARCH", garchOrder = c(1, 1)),
                   mean.model = list(armaOrder = c(1, 2), include.mean = FALSE),
                   distribution.model = "norm")
m.01_best3<-ugarchfit(spec=modelcombo3,data=r.amazon_close)
m.01_best3

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : sGARCH(1,1)
## Mean Model   : ARFIMA(1,0,2)
## Distribution : norm 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error  t value Pr(>|t|)
## ar1    -0.723237    0.170705 -4.23675 0.000023
## ma1     0.727326    0.183133  3.97157 0.000071
## ma2     0.106376    0.110335  0.96412 0.334986
## omega   0.000053    0.000052  1.02107 0.307219
## alpha1  0.121682    0.074715  1.62861 0.103395
## beta1   0.776089    0.126580  6.13122 0.000000
## 
## Robust Standard Errors:
##         Estimate  Std. Error  t value Pr(>|t|)
## ar1    -0.723237    0.088448 -8.17696 0.000000
## ma1     0.727326    0.125390  5.80049 0.000000
## ma2     0.106376    0.091036  1.16851 0.242601
## omega   0.000053    0.000063  0.83964 0.401111
## alpha1  0.121682    0.073551  1.65438 0.098049
## beta1   0.776089    0.153992  5.03981 0.000000
## 
## LogLikelihood : 357.9783 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -4.7565
## Bayes        -4.6350
## Shibata      -4.7596
## Hannan-Quinn -4.7071
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                          statistic p-value
## Lag[1]                     0.06296  0.8019
## Lag[2*(p+q)+(p+q)-1][8]    3.53035  0.9523
## Lag[4*(p+q)+(p+q)-1][14]   5.82051  0.7821
## d.o.f=3
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                     0.3829  0.5360
## Lag[2*(p+q)+(p+q)-1][5]    3.1978  0.3721
## Lag[4*(p+q)+(p+q)-1][9]    4.7058  0.4739
## d.o.f=2
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[3]     2.982 0.500 2.000 0.08421
## ARCH Lag[5]     4.824 1.440 1.667 0.11284
## ARCH Lag[7]     5.132 2.315 1.543 0.21155
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  1.1126
## Individual Statistics:              
## ar1    0.03224
## ma1    0.04063
## ma2    0.01627
## omega  0.32192
## alpha1 0.07596
## beta1  0.23226
## 
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic:          1.49 1.68 2.12
## Individual Statistic:     0.35 0.47 0.75
## 
## Sign Bias Test
## ------------------------------------
##                    t-value   prob sig
## Sign Bias           1.3000 0.1957    
## Negative Sign Bias  0.1928 0.8474    
## Positive Sign Bias  0.8189 0.4142    
## Joint Effect        3.0332 0.3865    
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     18.22       0.5080
## 2    30     19.03       0.9206
## 3    40     35.24       0.6419
## 4    50     50.65       0.4083
## 
## 
## Elapsed time : 0.06371689

par(mfrow=c(2,2))
plot(m.01_best3, which=8)
plot(m.01_best3, which=9)
plot(m.01_best3, which=10)
par(mfrow=c(1,1))

Here we can see that density curve has improved towards normality and also Q-Q plot tends towards normality

Lets try with another combination to see if it further improves

4. ARMA(1,1) + GARCH(2,1)

modelcombo4<-ugarchspec(variance.model = list(model = "sGARCH", garchOrder = c(2, 1)),
                   mean.model = list(armaOrder = c(1, 1), include.mean = FALSE),
                   distribution.model = "norm")

m.01_best4<-ugarchfit(spec=modelcombo4,data=r.amazon_close,out.sample = 100)

## Warning in .sgarchfit(spec = spec, data = data, out.sample = out.sample, : 
## ugarchfit-->waring: using less than 100 data
##  points for estimation

m.01_best4

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : sGARCH(2,1)
## Mean Model   : ARFIMA(1,0,1)
## Distribution : norm 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error   t value Pr(>|t|)
## ar1     0.156578    2.657920  0.058910  0.95302
## ma1    -0.147917    2.650043 -0.055817  0.95549
## omega   0.000005    0.000008  0.552216  0.58080
## alpha1  0.000000    0.370682  0.000000  1.00000
## alpha2  0.000000    0.365822  0.000000  1.00000
## beta1   0.999000    0.030116 33.171230  0.00000
## 
## Robust Standard Errors:
##         Estimate  Std. Error  t value Pr(>|t|)
## ar1     0.156578    0.409978  0.38192  0.70252
## ma1    -0.147917    0.478248 -0.30929  0.75710
## omega   0.000005    0.000014  0.32342  0.74638
## alpha1  0.000000    0.563046  0.00000  1.00000
## alpha2  0.000000    0.598110  0.00000  1.00000
## beta1   0.999000    0.050517 19.77551  0.00000
## 
## LogLikelihood : 113.8984 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -4.4958
## Bayes        -4.2619
## Shibata      -4.5226
## Hannan-Quinn -4.4074
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                     0.2197  0.6393
## Lag[2*(p+q)+(p+q)-1][5]    2.5955  0.7256
## Lag[4*(p+q)+(p+q)-1][9]    3.5350  0.7941
## d.o.f=2
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                          statistic p-value
## Lag[1]                       0.465  0.4953
## Lag[2*(p+q)+(p+q)-1][8]      3.567  0.5821
## Lag[4*(p+q)+(p+q)-1][14]     6.557  0.5641
## d.o.f=3
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[4]     1.490 0.500 2.000  0.2222
## ARCH Lag[6]     2.127 1.461 1.711  0.4634
## ARCH Lag[8]     2.289 2.368 1.583  0.6831
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  9.8489
## Individual Statistics:              
## ar1    0.19378
## ma1    0.19473
## omega  0.07613
## alpha1 0.03874
## alpha2 0.04375
## beta1  0.06067
## 
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic:          1.49 1.68 2.12
## Individual Statistic:     0.35 0.47 0.75
## 
## Sign Bias Test
## ------------------------------------
##                    t-value   prob sig
## Sign Bias           0.5307 0.5984    
## Negative Sign Bias  0.8090 0.4230    
## Positive Sign Bias  0.3465 0.7307    
## Joint Effect        0.8173 0.8453    
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     20.33       0.3748
## 2    30     24.50       0.7039
## 3    40     43.67       0.2798
## 4    50     45.75       0.6057
## 
## 
## Elapsed time : 0.05239606

par(mfrow=c(2,2))
plot(m.01_best4, which=8)
plot(m.01_best4, which=9)
plot(m.01_best4, which=10)
par(mfrow=c(1,1))

Here we can see that density curve has has again started degrading and moving away from normal curve. Hence we will proceed with our previous model

Hence our best model seems to be ARMA(1,2) + GARCH(1,1) that is m.01_best3

FORCASTING WITH OUR BEST COMBO MODEL FOR NEXT 10 YEARS

forcast_best = ugarchforecast(m.01_best3, data = r.amazon_close, n.ahead = 10, out.sample = 100)
forcast_best

## 
## *------------------------------------*
## *       GARCH Model Forecast         *
## *------------------------------------*
## Model: sGARCH
## Horizon: 10
## Roll Steps: 0
## Out of Sample: 0
## 
## 0-roll forecast [T0=May 2021]:
##          Series   Sigma
## T+1  -2.243e-04 0.01678
## T+2  -5.729e-04 0.01748
## T+3   4.143e-04 0.01808
## T+4  -2.997e-04 0.01861
## T+5   2.167e-04 0.01907
## T+6  -1.567e-04 0.01948
## T+7   1.134e-04 0.01983
## T+8  -8.199e-05 0.02015
## T+9   5.930e-05 0.02042
## T+10 -4.289e-05 0.02067

Conclusion

Our report focuses on selecting best model and forecast amazon stock prices.Here we have performed data modeling using ARIMA ,GARCH and combination of ARMA & GARCH model.Finally we found our best model to be ARMA(1,2) + GARCH(1,1).

Analyzing Evolution of Amazon Closing stock price (2009-2021)

s3831724,s3822713,s3533832,s3815368,s3755193

6/11/2021

1. Introduction

2. Methodology

3. Data Preprocessing

4. Descriptive Analysis

Time Series plot

Scatter plot

Normality check

5.Modelling

5.1 Linear model

Diagnostic checking for linear model

5.2 Quadratic model

Diagnostic checking of quadratic model

6.ARIMA Model

Overcome the non-stationarity of the series

BoxCox Transformation

Differencing - first differencing

ARIMA model specification

7. Parameter estimation

Fit models and find estimations

Compare scores of all models

8. GARCH model

Order specification

a) Calculate for GARCH(1,1)

b) Calculate for GARCH(1,3)

d) Calculate for GARCH(2,3)

Estimated Conditional Variance Prediction

9. Final Best Model Selection for ARMA+GARCH

1. ARMA(1,1) + GARCH(1,1)

2. ARMA(0,1) + GARCH(1,1)

3.ARMA(1,2) + GARCH(1,1)

4. ARMA(1,1) + GARCH(2,1)

FORCASTING WITH OUR BEST COMBO MODEL FOR NEXT 10 YEARS

Conclusion