Time-Series Analysis on S&P 500 Stock Index

Pendahuluan

Analisis ini berfokus pada penggunaan metode peramalan time series univariat untuk indeks pasar saham, Standard & Poor’s 500 (disingkat secara umum sebagai S & P 500) yang menekankan pemodelan Box-Jenkins AutoRegressive Integrated Moving Average (ARIMA). Analisis time series dilakukan dengan menggunakan R studio untuk memprediksi dan divisualisasikan. Dari hasil analisis yang dilakukan akan diketahui metode peramalan mana yang paling tepat untuk analisis time series S & P 500.

Metodologi

Analisis time series ini fokus pada metodologi Box-Jenkins yang menggabungkan beberapa langkah untuk memastikan menghasilkan model terbaik untuk meramalkan. Data menggunakan tahun 1995 hingga 2014, tahun 2015 sebagai pembanding. Beberapa asumsi yang diperlukan untuk analisis time series univariat.
  1. Data runtut waktu yang akan dianalisis dan diramalkan bersifat stasioner.
  2. Residual White Noise.

ARIMA

ARIMA sering juga disebut metode runtun waktu Box-Jenkins. ARIMA sangat baik ketepatannya untuk peramalan jangka pendek, sedangkan untuk peramalan jangka panjang ketepatan peramalannya kurang baik. Biasanya akan cenderung flat (mendatar/konstan) untuk periode yang cukup panjang. Model Autoregresif Integrated Moving Average (ARIMA) adalah model yang secara penuh mengabaikan independen variabel dalam membuat peramalan. ARIMA menggunakan nilai masa lalu dan sekarang dari variabel dependen untuk menghasilkan peramalan jangka pendek yang akurat. ARIMA cocok jika observasi dari deret waktu (time series) secara statistik berhubungan satu sama lain (dependent). Suatu deret waktu yang tidak stasioner harus diubah menjadi data stasioner dengan melakukan differencing. Yang dimaksud dengan differencing adalah menghitung perubahan atau selisih nilai observasi. Model Box-Jenkins (ARIMA) dibagi kedalam 3 kelompok model.
  1. Autoregressive Model (AR) dengan ordo p (AR(p)) atau model ARIMA (p,0,0)
  2. Moving Average Model (MA) dengan ordo q (MA(q)) atau ARIMA (0,0,q)
  3. Model campuran :
  • Proses ARMA Campuran proses AR(1) murni dan MA(1) murni, misal ARIMA (1,0,1)
  • Proses ARIMA Nonstasioneritas ditambahkan, maka model umum ARIMA (p,d,q) atau ARIMA (1,1,1).

Stasioner dan Musiman

Langkah pertama adalah melihat apakah objek time series bersifat stasioner, ini dapat dilakukan dalam berbagai metode yang juga dapat dijelaskan. Plot berikut akan membantu dalam estimasi Box-Jenkins Model, serta melakukan transformasi seperti differencing (dan mengambil log jika perlu):
  • Plot objek time series: Memberitahu jika ada pola musiman yang kuat.
  • Plot Dekomposisi runtun waktu: Melihat objek time series dalam komponen, sekilas dari komponen musiman dan trennya secara independen.
  • Plot Musiman: Plot ini adalah cara yang bagus untuk memeriksa komponen musiman yang merupakan sesuatu yang umum ketika berhadapan dengan data tahunan, kuartal, dan bulanan.

ACF dan PACF Plot

Plot-plot ini memainkan peran penting dalam analisis time series, karena dapat memperkirakan model ARIMA (p, d, q).

Data S&P 500

data <- read.csv("data1.csv", sep = ";", header = TRUE)
head(data)
##   year month      m1     m2    billion consumerSentiment imports inflation
## 1 1995     1 1150.62 3486.5 1000000000              97.6   60745       0.4
## 2 1995     2  1147.5 3485.7 1000000000              95.1   60076       0.4
## 3 1995     3 1147.48 3487.4 1000000000              90.3   62272       0.3
## 4 1995     4 1147.95 3490.3 1000000000              92.5   63448       0.3
## 5 1995     5 1146.14 3516.1 1000000000              89.8   63474       0.2
## 6 1995     6 1144.75 3543.1 1000000000              92.7   63617       0.2
##   oilPrices   ppi exports   cpi unemploymentRate fedFunds capUtilization
## 1     16.55 122.9   45619 150.3            0.073   0.0553          0.849
## 2     17.14 123.5   45914 150.9            0.072   0.0592          0.845
## 3     17.02 123.9   47483 151.4            0.054   0.0598          0.843
## 4     18.74 124.6   47032 151.9            0.058   0.0605          0.840
## 5     18.32 124.9   47646 152.2            0.056   0.0601          0.840
## 6     17.35 125.3   47425 152.5            0.056   0.0600          0.839
##   sp_500Dividends        nasdaq          nyse        sp_500 gdp_us
## 1          0.0283 7.553.224.948 2.679.487.427 4.645.475.003  10.09
## 2          0.0273 7.757.774.965 2.753.164.978 4.790.725.098  10.09
## 3          0.0267 8.079.524.995 2.823.107.422   493.987.503  10.09
## 4          0.0261   828.837.494 2.901.429.932   507.725.006  10.12
## 5          0.0254   858.662.491 2.979.569.946   523.650.009  10.12
## 6          0.0248 9.010.325.013 3.059.907.532 5.388.050.078  10.12
##   trillion housingIndex  X
## 1    1e+12       182.95 NA
## 2    1e+12       182.95 NA
## 3    1e+12       182.95 NA
## 4    1e+12       185.93 NA
## 5    1e+12       185.93 NA
## 6    1e+12       185.93 NA

A. Karakteristik Data

Statistika Deskriptif

summary(data)
##       year          month               m1            m2       
##  Min.   :1995   Min.   : 1.00   1076.3   :  2   Min.   : 3486  
##  1st Qu.:2000   1st Qu.: 3.75   1126.55  :  2   1st Qu.: 4718  
##  Median :2005   Median : 6.50   1.064.575:  1   Median : 6482  
##  Mean   :2005   Mean   : 6.50   1.064.875:  1   Mean   : 6956  
##  3rd Qu.:2010   3rd Qu.: 9.25   1.071.325:  1   3rd Qu.: 8669  
##  Max.   :2015   Max.   :12.00   1.075.025:  1   Max.   :12299  
##                                 (Other)  :244                  
##     billion      consumerSentiment    imports         inflation     
##  Min.   :1e+09   Min.   : 55.30    Min.   : 60076   Min.   :-1.900  
##  1st Qu.:1e+09   1st Qu.: 77.50    1st Qu.: 91672   1st Qu.: 0.000  
##  Median :1e+09   Median : 89.20    Median :131109   Median : 0.200  
##  Mean   :1e+09   Mean   : 87.25    Mean   :132211   Mean   : 0.181  
##  3rd Qu.:1e+09   3rd Qu.: 95.90    3rd Qu.:182999   3rd Qu.: 0.400  
##  Max.   :1e+09   Max.   :112.00    Max.   :200024   Max.   : 1.200  
##                                                                     
##    oilPrices           ppi           exports            cpi       
##  Min.   :  9.80   Min.   :122.3   Min.   : 45619   Min.   :150.3  
##  1st Qu.: 23.64   1st Qu.:129.6   1st Qu.: 57470   1st Qu.:171.3  
##  Median : 46.20   Median :155.7   Median : 74204   Median :195.0  
##  Mean   : 54.60   Mean   :159.9   Mean   : 85118   Mean   :195.9  
##  3rd Qu.: 78.65   3rd Qu.:189.2   3rd Qu.:117145   3rd Qu.:218.9  
##  Max.   :133.90   Max.   :208.3   Max.   :138406   Max.   :238.7  
##                                                                   
##  unemploymentRate     fedFunds       capUtilization   sp_500Dividends  
##  Min.   :0.03800   Min.   :0.00070   Min.   :0.6690   Min.   :0.01110  
##  1st Qu.:0.04700   1st Qu.:0.00160   1st Qu.:0.7620   1st Qu.:0.01628  
##  Median :0.05500   Median :0.02000   Median :0.7835   Median :0.01835  
##  Mean   :0.05991   Mean   :0.02725   Mean   :0.7851   Mean   :0.01859  
##  3rd Qu.:0.06725   3rd Qu.:0.05250   3rd Qu.:0.8150   3rd Qu.:0.02070  
##  Max.   :0.10000   Max.   :0.06540   Max.   :0.8490   Max.   :0.03600  
##                                                                        
##            nasdaq               nyse               sp_500   
##  1.006.919.983:  1   1.008.315.259:  1   1.007.204.987:  1  
##  1.021.865.021:  1   1.013.219.482:  1   1.008.657.501:  1  
##  1.035.407.501:  1   1.014.059.497:  1   1.015.390.015:  1  
##  1.037.555.008:  1   1.022.661.987:  1   1.020.154.999:  1  
##  1.045.472.519:  1   1.024.867.236:  1   1.024.109.986:  1  
##  1.047.292.511:  1   1.037.144.482:  1   1.034.937.515:  1  
##  (Other)      :246   (Other)      :246   (Other)      :246  
##      gdp_us         trillion      housingIndex      X          
##  Min.   :10.09   Min.   :1e+12   Min.   :182.9   Mode:logical  
##  1st Qu.:12.53   1st Qu.:1e+12   1st Qu.:231.7   NA's:252      
##  Median :14.23   Median :1e+12   Median :311.4                 
##  Mean   :13.65   Mean   :1e+12   Mean   :291.5                 
##  3rd Qu.:14.95   3rd Qu.:1e+12   3rd Qu.:343.0                 
##  Max.   :16.44   Max.   :1e+12   Max.   :378.4                 
##

Histogram

library(RColorBrewer)
par(mfrow=c(2,2))
hist(data$sp_500Dividends, breaks=9, col="pink", main="Histogram of S&P 500 Dividends")
hist(data$inflation, breaks=9, col="pink", main="Histogram of Inflation")
hist(data$ppi, breaks=9, col="pink", main="Histogram of ppi")
hist(data$gdp_us, breaks=9, col="pink", main="Histogram of GDP United State")

B. Analysis Time Series

Asumsi

library(ggplot2)
library(forecast)
library(plotly)
## 
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout
library(ggfortify)
library(tseries)
library(gridExtra)
library(docstring)
## 
## Attaching package: 'docstring'
## The following object is masked from 'package:utils':
## 
##     ?
library(readr)

sp_500 <- ts(data$sp_500, start=c(1995, 1), freq=12)
plot.ts(sp_500, main = "Plot of S&P Time Series", xlab = "Year", ylab = "S&P 500")

sp500_training <- ts(sp_500, start=c(1995, 1), end=c(2014, 12), freq=12)
plot.ts(sp500_training, main = "Plot of S&P Training Set Time Series", xlab = "Year", ylab = "S&P 500")

Box.test(sp500_training, lag = 20, type = 'Ljung-Box')
## 
##  Box-Ljung test
## 
## data:  sp500_training
## X-squared = 474.67, df = 20, p-value < 2.2e-16
adf.test(sp500_training)
## 
##  Augmented Dickey-Fuller Test
## 
## data:  sp500_training
## Dickey-Fuller = -3.9478, Lag order = 6, p-value = 0.01229
## alternative hypothesis: stationary
seasonplot(sp500_training, main="Seosonal Plot of S&P 500", col=1:40, pch=19)

acf(sp500_training, main="ACF Plot")

pacf(sp500_training, main="PACF Plot")

Transformasi

Hasil visual diatas menunjukkan bahwa asumsi stasioner tidak terpenhi, maka dilakukan transformasi.
tsDiff <- diff(sp500_training)
plot.ts(tsDiff, main = "Plot of Difference Time Series", xlab = "Year", ylab = "S&P 500")

Box.test(tsDiff, lag = 20, type = 'Ljung-Box')
## 
##  Box-Ljung test
## 
## data:  tsDiff
## X-squared = 37.703, df = 20, p-value = 0.009624
adf.test(tsDiff)
## Warning in adf.test(tsDiff): p-value smaller than printed p-value
## 
##  Augmented Dickey-Fuller Test
## 
## data:  tsDiff
## Dickey-Fuller = -6.9956, Lag order = 6, p-value = 0.01
## alternative hypothesis: stationary
seasonplot(tsDiff, main="Seosonal Plot of S&P 500", col=1:40, pch=19)

acf(tsDiff, main="ACF Plot")

pacf(tsDiff, main="PACF Plot")

Model

Dari plot di atas kita melihat bahwa plot ACF dan PACF menunjukkan karakteristik model MA (1) karena cut off pada lag pertama, tetapi karena terdapat differencing maka menjadi model campuran ARIMA (0, 1, 1).
fit <- Arima(sp_500, order = c(0,1,1),
    include.drift = TRUE)
summary(fit)
## Series: sp_500 
## ARIMA(0,1,1) with drift 
## 
## Coefficients:
##           ma1    drift
##       -0.3664  -0.0480
## s.e.   0.0616   1.9062
## 
## sigma^2 estimated as 2279:  log likelihood=-1325.54
## AIC=2657.07   AICc=2657.17   BIC=2667.65
## 
## Training set error measures:
##                      ME     RMSE      MAE       MPE     MAPE      MASE
## Training set 0.00488727 47.45613 22.93573 -213.3576 229.6167 0.3162822
##                     ACF1
## Training set 0.009671651
residFit <- ggplot(data=fit, aes(residuals(fit))) +
    geom_histogram(aes(y =..density..),
    col="black", fill="white") +
geom_density(col=1) +
theme(panel.background = element_rect(fill = "gray98"),
    panel.grid.minor = element_blank(),
    axis.line   = element_line(colour="gray"),
    axis.line.x = element_line(colour="gray")) +
ggtitle("Plot of S&P 500 ARIMA Model Residuals")

residFit
## Don't know how to automatically pick scale for object of type ts. Defaulting to continuous.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

ggplotly(residFit)
## Don't know how to automatically pick scale for object of type ts. Defaulting to continuous.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
sp500_test <- read.csv("tes.csv")
sp500_test <- ts(sp500_test$Adj.Close, 
    start = c(2015, 1), 
    frequency = 12)
fit_arima <- forecast(fit, h = 36)

autoplot(fit_arima,
    holdout = sp500_test, 
    forc_name = 'ARIMA', 
    ts_object_name = 'S&P 500')

Didapatkan hasil forecast dengan metode ARIMA, terdapat beberapa metode yang dapat dicoba, seperti:
  • Box-Cox Transformation Forecast
  • Mean Forecast
  • Naive Forecast
  • Seasonal Naive Forecast
  • Exponential Smoothing Forecast.

Berikut hasil perbandingan masing-masing metode

hasil=read.csv("hasil.csv", header = TRUE)
hasil
##                            Model     RMSE      MAE    MPE   MAPE  MASE
## 1                          ARIMA  165.054  119.310  2.761  5.089 0.629
## 2         Box-Cox Transformation  712.993  606.664 26.294 26.294 3.200
## 3          Exponential Smoothing  221.442  155.624  3.668  6.586 0.821
## 4           Mean Forecast Method 1043.815 1023.667 45.944 45.944 5.400
## 5          Naive Forecast Method  259.020  183.291  6.472  7.690 0.967
## 6 Seasonal Naive Forecast Method  339.903  280.988 12.001 12.127 1.482
## 7                 Neural Network  579.195  407.790 16.483 17.028 2.151
##    ACF1 Theil.s.U
## 1 0.866     2.877
## 2 0.903    12.826
## 3 0.889     3.842
## 4 0.887    19.516
## 5 0.887     4.476
## 6 0.898     6.032
## 7 0.907    10.049

Kesimpulan

Dari hasil perbandingan beberpa parameter diputuskan bahwa ARIMA (0,1,1) adalah model terbaik dalam meramalkan S&P 500 Saham Index.