Time Series Modelling non stationary ARIMA step by step


ARIMA (Autoregressive integrated moving average) NON STATIONARY
Penjeasan tentang ARIMA selengkapnya bisa lihat di wikipedia kilik link ini.

Tahapan ARIMA modeling:
  1. Exploratory data analysis
  2. Dekomposisi data
  3. Test stationarity
  4. Membuat Model ARIMA
  5. Forecasts
Kali ini kita praktekkan saham salah satu perusahaan otomotif terkenal di Indonesia.
Data saham dapat didownload dari https://finance.yahoo.com/.
File yangdidapatkan berbentuk csv dimana kolom pertama adalah tanggal dan kolom kelima adalah harga closing saham.

Library yang dipergunakan:
>library(MASS) >library(tseries) >library(forecast)

STEP-1 Exploratory data analysis

Data
> sahamA<-read.csv(file.choose())
> head(sahamA)
        Date      Open      High       Low     Close Adj.Close   Volume
1 2014-12-01      null      null      null      null      null     null
2 2015-01-01 41.180000 41.310001 37.680000 37.770000 36.301414 53067700
3 2015-02-01 37.900002 42.490002 37.709999 42.209999 40.668121 52385100
4 2015-03-01 42.340000 42.740002 40.369999 41.549999 40.032230 44501800
5 2015-04-01 41.459999 43.590000 40.720001 41.369999 39.956970 47251000
6 2015-05-01 41.410000 42.930000 40.820000 41.189999 39.783119 44051900
> tail(sahamA)
         Date      Open      High       Low     Close Adj.Close   Volume
44 2018-07-01 61.360001 66.830002 60.790001 66.040001 65.741722 37941900
45 2018-08-01 65.839996 68.150002 64.220001 67.540001 67.397339 56324900
46 2018-09-01 67.339996 71.459999 66.760002 70.540001 70.390999 35878100
47 2018-10-01 70.680000 72.370003 61.009998 64.790001 64.653145 65097000
48 2018-11-01 65.070000 66.610001 64.580002 65.910004 65.910004  6332900
49 2018-11-06 66.040001 66.419998 65.785004 66.014999 66.014999   414534

Data Training
  • Kita menggunakan saham ini dari januari 2015 sampai dengan november 2018
  • Perhatikan data baris pertama tidak kita pergunakan
  • Data 2015 sampai 2017 kita gunakan untuk data training
  • Data tahun 2018 kita gunakan untuk sebagai data testing
> Tanggal<-as.Date(sahamA[2:37,1])
> Harga<-sahamA[2:37,5]
> sahamA_TRAINING<-data.frame(Tanggal,Harga)
> str(sahamA_TRAINING)
'data.frame': 36 obs. of  2 variables:
 $ Tanggal: Date, format: "2015-01-01" "2015-02-01" "2015-03-01" "2015-04-01" ...
 $ Harga  : Factor w/ 49 levels "34.330002","36.310001",..: 6 16 13 12 11 7 10 2 1 5 ...


Pembentukkan data Testing
>Tanggal2<-as.Date(sahamA[38:48,1])
>Harga2<-sahamA[38:48,5]
>sahamA_TESTING<-data.frame(Tanggal2,Harga2)
>str(sahamA_TESTING)
Transformasi ke format data time series
> ts_sahamA<-ts(sahamA_TRAINING$Harga,frequency = 12,start=c(2015,1))
> class(ts_sahamA)
[1] "ts"
  • Mari kita lihat dulu sepertia apa datasetnya
> ts_sahamA
     Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
2015   6  16  13  12  11   7  10   2   1   5  15  14
2016   4   3   8   9  21  19  24  22  23  17  18  20
2017  25  26  27  28  31  29  30  35  34  44  46  42
  • Pemeriksaan nilai yang hilang (missing value)
> sum(is.na(ts_sahamA))
[1] 0
  • Pemeriksaan Frequency dengan fungsi frequency()
> frequency(ts_sahamA)
[1] 12
  • Pemeriksaan cycle dengan fungsi cycle()
> cycle(ts_sahamA)
     Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
2015   1   2   3   4   5   6   7   8   9  10  11  12
2016   1   2   3   4   5   6   7   8   9  10  11  12
2017   1   2   3   4   5   6   7   8   9  10  11  12
  • Summary data saham  
> summary(ts_sahamA)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   1.00    9.75   18.50   19.36   27.25   46.00 
  • Plot data saham
  • Sumbu horisontal = "Periode"
  • Sumbu vertikal = "Harga Saham"
  • Judul = "Harga Saham PT.A dari tahun 2015 sd 2017"
plot(ts_sahamA,
        xlab="Periode", 
        ylab = "Harga Saham",
        main="Harga Saham PT A dari tahun 2015 sd 2017")


  • Pemeriksaan Seasonal dengan fungsi boxplot()
boxplot(ts_sahamA~cycle(ts_sahamA),
           xlab="Periode", 
           ylab = "Harga Saham",
           main="Harga Saham PT.A dari tahun 2015 sd 2017")

Hasil EDA:
  • Trend meningkat secara linear
  • Termasuk Time Series Multiplikatif
  • Disetiap bulan ke 10 sd 12 harga saham condong ke harga tinggi karena banyak orang mendapatkan bonus akhir tahun
  • Data sudah bersih karena tidak ada outliers dan missing value sehingga tidak perlu cleansing data

STEP-2 Dekomposisi Data

  • Proses dekomposisi dllakukan untuk memriksa adanya trend, seasonal, dan random component,
  • Persamaannya:
Y[t]=T[t]S[t]e[t]
  • dimana
    • Y(t) jumlah penumpang pada waktu t,
    • T(t) trend pada waktu t,
    • S(t) seasonal pada waktu t,
    • e(t) random error pada waktu t.
  • untuk memeriksa dekomposisi lakukan fungsi berikut
dec_sahamA <- decompose(ts_sahamA,"multiplicative")
plot(dec_sahamA)

STEP-3 Test Stationary

> adf.test(ts_sahamA)

 Augmented Dickey-Fuller Test

data:  ts_sahamA
Dickey-Fuller = -2.5606, Lag order = 3, p-value = 0.355
alternative hypothesis: stationary


  • Data tidak stationary karena p-value di atas 0.05
  • Kita perlu melakukan proses differencing dengan fungsi diff()
  • Differensiasi dilakukan 1 bulan karena periode data bulanan
> diff_ts_sahamA<-diff(ts_sahamA,1)
> adf.test(diff_ts_sahamA)

 Augmented Dickey-Fuller Test

data:  diff_ts_sahamA
Dickey-Fuller = -4.3906, Lag order = 3, p-value = 0.01
alternative hypothesis: stationary

Warning message:
In adf.test(diff_ts_sahamA) : p-value smaller than printed p-value
>

  • sekarang data sudah stationary karena p-value<0.05

STEP-4 Membuat Model ARIMA

Membuat Model ARIMA

  • kita gunakan fungsi auto.arima()
  • kemudian kita test dengan ggtsdiag()
  • Residual berada diantar nol, berarti model ARIMA ini cukup fit.
> arima_sahamA<-auto.arima(ts_sahamA)
> library(ggfortify) #gunakan library ini 
> ggtsdiag(arima_sahamA)



STEP-5 FORECAST

Proses Forecasting:


> f_sahamA<- forecast(arima_sahamA, h = 11)
> f_sahamA
         Point Forecast    Lo 80    Hi 80     Lo 95    Hi 95
Jan 2018             42 35.67702 48.32298 32.329835 51.67016
Feb 2018             42 33.05795 50.94205 28.324322 55.67568
Mar 2018             42 31.04828 52.95172 25.250783 58.74922
Apr 2018             42 29.35404 54.64596 22.659670 61.34033
May 2018             42 27.86138 56.13862 20.376854 63.62315
Jun 2018             42 26.51192 57.48808 18.313031 65.68697
Jul 2018             42 25.27096 58.72904 16.415149 67.58485
Aug 2018             42 24.11591 59.88409 14.648644 69.35136
Sep 2018             42 23.03106 60.96894 12.989506 71.01049
Oct 2018             42 22.00498 61.99502 11.420254 72.57975
Nov 2018             42 21.02904 62.97096  9.927692 74.07231
> plot(f_sahamA)


Perhitungan Prediksi
> fme_sahamA<-as.numeric(f_sahamA$mean)
> hasil_sahamA<-data.frame(as.double(sahamA_TESTING$Harga),fme_sahamA)
> col_headings<-c("Harga Asli","Prediksi")
> names(hasil_sahamA)<-col_headings
> hasil_sahamA
   Harga Asli Prediksi
1          48       42
2          45       42
3          41       42
4          37       42
5          33       42
6          32       42
7          40       42
8          43       42
9          47       42
10         36       42
11         38       42

Perhitungan Mean Error
> errors_sahamA<-((hasil_sahamA[,1]-hasil_sahamA[,2])/hasil_sahamA[,1])
> errors_sahamA
 [1]  0.12500000  0.06666667 -0.02439024 -0.13513514 -0.27272727 -0.31250000 -0.05000000  0.02325581  0.10638298
[10] -0.16666667 -0.10526316
> mean(errors_sahamA)
[1] -0.06776155

No comments:

Post a Comment