Penerapan Model ARIMA

Penerapan Model ARIMA (Bagian I)

Dr. Kusman Sadik, M.Si Departemen Statistika IPB, 2017 1

Ada tiga tahapan iteratif dalam pemodelan data deret waktu yang berbasis model ARIMA, yaitu: 1. Penentuan model tentatif (spesifikasi model) berdasarkan data contoh untuk mengidentifikasi nilai p, d, dan q. 2. Pendugaan parameter model ARIMA(p, d, q) yang diidentifikasi, yaitu penduga nilai , , dan σ2𝑒 . 3. Analisis diagnostik untuk melihat kelayakan model. 2



Prosedur iterasi ini sering disebut ”Metode Box-Jenkins”.



Untuk model ARIMA(p, d, q), spesifikasi dilakukan untuk menentukan nilai p, d, dan q.



Alat yang digunakan pada tahap identifikasi ini adalah fungsi autokorelasi.



Fungsi autokorelasi ini diduga dari data contoh atau disebut fungsi autokorelasi contoh (sample of autocorrelation function atau SACF atau ACF saja).



Disamping itu ada pula fungsi autokorelasi parsial (sample of partial autocorrelation function atau SPACF atau PACF saja) 3

a. ACF nk



rk 

 (Y  Y )(Y t

t 1

n

t k

Y )

2 ( Y  Y )  t

,

k  1, 2, ....

t 1

n

 

Y 

Y t 1

t

n

rk merupakan penduga bagi k

4

a. PACF 

PACF : kk = Corr(Yt, Yt-k | Yt-1, Yt-2, …, Yt-k+1)



Berdasarkan persamaan Yule-Walker:

j = k1j-1 + k2j-2 + ... + kkj-k j = 1, 2, ..., k; Catatan: j = -j dan 0 = 1

k  ACF;

kk  PACF ˆkk  penduga bagi kk 5

Contoh: Misal diketahui data : 4, 2, 5, 1. Tentukan ACF (r1, r2) dan PACF (𝜙11 , 𝜙22 ) nk

Melalui persamaan rk 

 (Y  Y )(Y t 1

t

n

t k

Y )

2 ( Y  Y )  t

,

k  1, 2, ....

t 1

Dapat diperoleh penduga ACF : r1 = -0.7 dan r2 = 0.4

6

Berdasarkan persamaan Yule-Walker dapat diperoleh penduga PACF kk:

j = k1j-1 + k2j-2 + ... + kkj-k Untuk k =1  j = 1

1 = 110

 1 = 11(1)  r1 = 𝜙11 = -0.7

Untuk k = 2  j = 1, 2

1 = 210 + 221 

1 = 21 + 221

2 = 211 + 220 

2 = 211 + 22 7

1 = 210 + 221 

1 = 21 + 221

2 = 211 + 220 

2 = 211 + 22

(1)2 = 211 + 22(1)2 ...... Pers(1)

2 = 211 + 22

……..

Pers(2)

Berdasarkan Pers(1) dan Pers(2) diperoleh: (1)2 - 2 = 22(1)2 - 22

22 = {(1)2 - 2}/{(1)2 - 1} 𝜙22 = {(r1)2 - r2}/{(r1)2 - 1} = 0.09/(-0.51) = -0.176 8

Implementasi dalam Program R > data <- c(4, 2, 5, 1)

> acf(data, lag.max = 3, plot = FALSE) Autocorrelations of series ‘data’, by lag 0 1 1.0 -0.7

2 3 0.4 -0.2

> pacf(data, lag.max = 3, plot = FALSE) Partial autocorrelations of series ‘data’, by lag 1 2 -0.700 -0.176

3 0.012 9

Pengidentifikasian Model Model MA: Misal MA(1) : Yt = et - et-1

ACF :

  1   2 ; k  1  k    0 ; k 1  

1.0

1.0

0.8

k

0.8

0.6

rk

0.4

0.6 0.4

0.2

k

0.0 1

2

3

ACF

4

5

0.2

k

0.0 1

2

3

4

5

6

Sample of ACF

10



Karena rk berasal dari data contoh maka diperlukan galat baku bagi rk yaitu Srk.



Sebagai nilai pendekatan : Srk = 1 / n , dimana n adalah banyaknya data.



Sehingga hipotesis H0 : k = 0 ditolak jika | rk | > 2Srk atau | rk | > 2 / n .



Misalnya, jika | r1 | > 2 / n dan | rk | < 2 / n untuk k = 2, 3, …, maka model tentatifnya adalah MA(1).

11

Model AR : Misalkan AR(1) : Yt = Yt-1 + et 

ACF : k = k ; k = 1, 2, …



Untuk model AR, ACF merupakan fungsi eksponensial sehingga ACF tidak dapat digunakan untuk menentukan nilai p dalam AR(p).



PACF : j = k1j-1 + k2j-2 + ... + kkj-k (Yule-Walker) untuk k = 1  1 = 11 untuk k = 2  1 = 21 + 221 .... (1)

2 = 211 + 22 .... (2) 12

Berdasarkan persamaan (1) dan (2)  22 = 0. Demikian juga 33 = 44 = ... = 0.  1  Sehingga PACF AR(1): kk    0 



; k 1 ; k 1

Dengan demikian PACF dapat digunakan sebagai penentu nilai p dalam model AR(p).

1.0 1.0

0.8

kk

ˆ

0.6

kk

0.8 0.6

0.4

0.4

0.2

0.2 0.0

0.0 1

2

3

PACF



4

5

1

2

3

4

5

6

Sample of PACF

Hipotesis H0 : kk = 0 ditolak jika | ˆkk |  2 / n . 13

Pengidentifikasian nilai p dan q 1.0

tails off

0.8 0.6 0.4 0.2 0.0 1

2

3

4

5

6

7

8

9

Sample of ACF

1.0

cuts off after lag q

0.8 0.6 0.4 0.2 0.0 1

2

3

4

5

6

7

8

9

Sample of ACF

14

15

Contoh (1)

16

Contoh (2)

17

Contoh (3) d=1

d=1

18

Pengidentifikasian ARMA(p, q) Melalui EACF  Nilai ACF dan PACF dapat digunakan untuk menentukan nilai q pada model MA(q) dan nilai p pada model AR(p).  Namun tidak bisa digunakan untuk menentukan nilai p dan q pada model campuran ARMA(p, q).  Karena itu dikembangkan metode extended autocorrelation function (EACF) untuk pengidentifikasian model campuran ARMA(p, q).  Pada Tabel EACF, secara teoritis model ARMA(p, q) mempunyai pola segitiga-nol (triangle of zeroes), dimana nilai pada pojok kiri atas bersesuaian dengan ordo ARMA. 19

Pengidentifikasian ARMA(p, q) Melalui EACF

20

Pendugaan Parameter Model Apabila nilai p, d, dan q sudah dapat diidentifikasi, maka selanjutnya dilakukan pendugaan terhadap parameter model, yaitu 1, 2, ..., p untuk model AR(p) dan 1, 2, ..., q untuk model MA(q) berdasarkan data terobservasi Y1, Y2, ..., Yn. Metode pendugaan parameter : 

Metode momen,



Metode kuadrat terkecil (least-square),



Metode kemungkinan maksimum (maximum likelihood). 21

1. Metode Momen Metode ini didasarkan pada persamaan momen contoh dan momen teoritis, kemudian memecahkan persamaan-persamaan tersebut untuk mendapatkan penduga bagi parameter model. Misalnya, menduga rataan populasi (teoritis)  dengan rataan contoh Y . Model AR a. AR(1) : Yt = Yt-1 + et 

k = k ; k = 1, 2, …



1 =   ˆ1  ˆ  r1 = ˆ



Jadi pada AR(1) penduga bagi parameter model, , adalah r1 yang dapat dihitung dari data. 22

b. AR(1) : Yt =  + Yt-1 + et 

Bagaimana menduga  ?



Perhatikan model : (Yt - 𝑌)= (Yt-1 - 𝑌) + et ↔ (Yt - 𝑌)= (Yt-1 - 𝑌) + et ↔ Yt = (1 - )𝑌 + Yt-1 + et ↔ Yt =  + Yt-1 + et Sehingga :  = (1 - )𝑌

23

c. AR(2) : Yt = 1Yt-1 + 2Yt-2 + et Berdasarkan persamaan Yule-Walker :

k = 1k-1 + 2k-2 + ... + pk-p maka diperoleh

1 = 1 + 21 dan 2 = 11 + 2 dengan metode momen diperoleh: r1 = ˆ1 + ˆ2 r1 dan r2 = r1 ˆ1 + ˆ2 penyelesaian terhadap dua persamaan ini diperoleh: r1 (1  r2 ) r2  r1 ˆ ˆ dan 2  1  2 2 1  r1 1  r1

2

24

Model MA MA(1) : Yt = et - et-1

 1   1 2

ˆ  r1   1  ˆ 2

 1  1  4r1 sehingga diperoleh : ˆ  2r1

2



Sebagai catatan untuk persamaan ini, apabila | r1 | > 0.5 maka metode momen gagal untuk menduga parameter .



Untuk MA(2), MA(3), dst, metode momen menjadi sangat kompleks, sehingga harus menggunakan metode pendugaan lainnya. 25

Model ARMA ARMA(1, 1) : Yt = Yt-1 - et-1 + et

k 

(1   )(   ) k 1  2 1  2  

2 r   sehingga penduga bagi  adalah : ˆ  2 1 r1 Untuk menduga  dapat digunakan persamaan pertama dengan cara mengganti 1 dengan r1 dan  dengan ˆ , yaitu (1  ˆˆ)(ˆ  ˆ) r1  1  2ˆˆ  ˆ 2 26

Contoh Kasus (Latihan): Misalnya diketahui model AR(2) : Yt =  + 1Yt-1 + 2Yt-2 + et. Berdasarkan data diketahui bahwa r1 = 0.75, r2 = 0.61, dan Y = 4.5. Tentukan ˆ , ˆ1 , dan ˆ2 dengan metode momen.

27

2. Metode Kuadrat Terkecil Metode ini dilakukan dengan cara meminimumkan komponen n

pada galat, yaitu

 et . 2

t 1

AR(1) : Yt = Yt-1 + et n

S() =  et = t 1

2



et = Yt - Yt-1

n

 (Y  Y t 1

t

2 ) t 1

Penduga bagi parameter model, , dapat diperoleh dengan cara meminimumkan S().

28

MA(1) : Yt = et - et-1



et = Yt + et-1

et = Yt + ( Yt-1 + et-2) et = Yt + Yt-1 + 2Yt-2 + 3Yt-3 + ….

n

S() =  et

2

t 1

Meminimumkan S() tidak dapat dilakukan secara analitik / kalkulus karena bersifat non-linear, sehingga harus diselesaikan secara numerik / iteratif, salah satunya

melalui

algoritma

Gauss-Newton

atau

Newton-Raphson. 29

3. Metode Kemungkinan Maksimum Metode ini dilakukan dengan cara memaksimumkan fungsi kemungkinan (likelihood), berdasarkan fungsi sebaran galat (e t).

bsi

AR(1) : Yt = Yt-1 + et, misal et ~ N(0, e2) f(e1, e2, …., en) = (2 )

2 ( n 1) / 2 e

L(,

e2)

= (2 )

2  ( n 1) / 2 e

. exp(

. exp(

1 2 e2

n

1 2 e2

2 e  t) t 1

n

 (Y  Y ) ) 2

t 1

t

t

Penduga  dan e2 dapat diperoleh dengan cara memaksimumkan fungsi kemungkinan L(, e2). 30

MA(1) : Yt = et - et-1 Fungsi kemungkinannya, L(, e2), bersifat non-linear sehingga pemaksimumannya harus dilakukan secara numerik / iteratif.

Catatan : Program R menggunakan metode iterasi NewtonRaphson untuk menduga parameter AR(p), MA(q), dan ARIMA(p, d, q).

31

# Pemodelan ARIMA(1,1,1) library("forecast") library("TTR") library("TSA") library("graphics") # Membangkitkan y, ARIMA(1,1,1): mu=0.15 phi=0.55 tetha=0.75 set.seed(1001) e <- rnorm(150,0,1) n <- length(e) mu phi tetha

<- 0.15 <- 0.55 <- -0.75

y <- c(1:n) for (i in 3:n) { y[i] <- mu + (1+phi)*y[i-1] - phi*y[i-2] + e[i] - tetha*e[i-1]}

y <- y[-c(1:50)] # membuang 50 data pertama plot.ts(y, lty=1, xlab="Waktu", ylab="Data Asal (y)") points(y) 32

acf(y, lag.max=20)

# cek kestasioneran

y.dif1 <- diff(y, difference=1) # differencing ordo 1 plot.ts(y.dif1, lty=1, xlab="Waktu", ylab="Data Y.Diff Ordo 1") points(y.dif1)

# Pengidentifikasian Model acf(y.dif1, lag.max=20) pacf(y.dif1, lag.max=20) eacf(y.dif1)

# Pendugaan Parameter dan Penentuan Model Terbaik # Berdasarkan Kandidat Model Hasil Identifikasi arima(y.dif1, order=c(0,0,2),method="ML") arima(y.dif1, order=c(3,0,0),method="ML") arima(y.dif1, order=c(1,0,1),method="ML")

# ARIMA(0,1,2) # ARIMA(3,1,0) # ARIMA(1,1,1)

33

# Plot dan Nilai Dugaan Berdasarkan Model Terbaik model

<- arima(y.dif1, order=c(1,0,1),method="ML")

# ARIMA(1,1,1)

dugaan <- fitted(model) cbind(y.dif1,dugaan) plot.ts(y.dif1, xlab="Waktu", ylab="Data Diff.Y Ordo 1") points(y.dif1) par(col="red") lines(dugaan) par(col="black")

34

35

36

37

38

39

> eacf(y.dif1) AR/MA 0 1 0 x x 1 x o 2 x x 3 x x 4 x x 5 x o 6 x o 7 x o

2 o o o o o o o o

3 o o o o x o o o

4 o o o o o o o o

5 o o o o o o o o

6 o o o o o o o o

7 o o o o o o o o

8 o o o o o o o o

9 o o o o o o o o

10 x o o o o o o o

11 x x x x x x x o

12 x o o o o o o o

13 x o o o o o o o

40

> arima(y.dif1, order=c(0,0,2),method="ML") # ARIMA(0,1,2)

Call: arima(x = y.dif1, order = c(0, 0, 2), method = "ML") Coefficients: ma1 ma2 1.2345 0.3810 s.e. 0.0910 0.0936

intercept 0.3195 0.2525

sigma^2 estimated as 0.9365: aic = 282.37

log likelihood = -138.18,

41

> arima(y.dif1, order=c(3,0,0),method="ML")

# ARIMA(3,1,0)

Call: arima(x = y.dif1, order = c(3, 0, 0), method = "ML") Coefficients: ar1 ar2 1.2249 -0.7571 s.e. 0.0987 0.1438

ar3 0.2688 0.1016

sigma^2 estimated as 0.9572: aic = 286.16

intercept 0.3249 0.3662 log likelihood = -139.08,

42

> arima(y.dif1, order=c(1,0,1),method="ML")

# ARIMA(1,1,1)

Call: arima(x = y.dif1, order = c(1, 0, 1), method = "ML") Coefficients: ar1 ma1 0.5423 0.7580 s.e. 0.0894 0.0668

intercept 0.3183 0.3585

sigma^2 estimated as 0.8906: aic = 277.37

log likelihood = -135.69,

43

> cbind(y.dif1,dugaan) Time Series: Start = 1 End = 99 Frequency = 1 y.dif1 dugaan 1 -0.301500868 0.03487535 2 -2.989745815 -0.59954250 3 -1.798462592 -2.91505510 . . 97 2.664318123 2.84607363 98 0.354648756 1.45281261 99 -0.463277448 -0.49441754

44

Garis Merah : Plot Nilai Dugaan Berdasarkan Model Terbaik 45

> arima(y.dif1, order=c(1,0,1),method="ML")

# ARIMA(1,1,1)

Call: arima(x = y.dif1, order = c(1, 0, 1), method = "ML") Coefficients: ar1 ma1 0.5423 0.7580 s.e. 0.0894 0.0668

intercept 0.3183 0.3585

sigma^2 estimated as 0.8906: aic = 277.37

log likelihood = -135.69,



Bandingkan dengan data y yang dibangkitkan:



ARIMA(1,1,1) dengan parameter μ = 0.15, ϕ = 0.55, θ = 0.75 “sigma^2 estimated as 0.8906” adalah nilai dugaan bagi σe2



Model terbaik tersebut selanjutnya bisa digunakan untuk peramalan 46

1. Melalui Program R, bangkitkan data yt, (n = 225), berupa model ARIMA(1, 2, 0) dengan  = 0.5, Φ = 0.7 serta et ~ Normal(0,1). Gunakan 200 data terakhir dan lakukan proses berikut: a. Identifikasilah kestasioneran data, serta lakukan proses differencing jika data tidak stasioner. b. Selanjutnya, berdasarkan ACF, PACF, dan EACF, identifikasilah kandidat model yang sesuai.

c. Berdasarkan kandidat model tersebut, tentukan model terbaik berdasarkan nilai AIC-nya. d. Bandingkan penduga parameter yang diperoleh untuk model terbaik pada poin (c) tersebut dengan nilai parameter yang sesungguhnya. Apa kesimpulan Anda? 47

2. Melalui Program R, bangkitkan data yt, (n = 225), berupa model ARIMA(1, 1, 2) dengan  = 1.0, Φ = 0.7, θ1 = - 0.8, dan θ2 = 0.6 serta et ~ Normal(0,1). Gunakan 200 data terakhir dan lakukan proses berikut: a. Identifikasilah kestasioneran data, serta lakukan proses differencing jika data tidak stasioner. b. Selanjutnya, berdasarkan ACF, PACF, dan EACF, identifikasilah kandidat model yang sesuai. c. Berdasarkan kandidat model tersebut, tentukan model terbaik berdasarkan nilai AIC-nya.

d. Bandingkan penduga parameter yang diperoleh untuk model terbaik pada poin (c) tersebut dengan nilai parameter yang sesungguhnya. Apa kesimpulan Anda? 48

3. Melalui Program R, kerjakan : Exercise 5.11 (Montgomery, hlm. 290):

49



Montgomery, D.C., et.al. 2008. Forecasting Time Series Analysis 2nd. John Wiley.



Cryer, J.D. and Chan, K.S. 2008. Time Series Analysis with Application in R. Springer.



Cowpertwait, P.S.P. and Metcalfe, A.V. 2009. Introductory Time Series with R. Springer New York.



Wei, William, W.S. 1990. Time Series Analysis, Univariate and Multivariate Methods. Adison-Wesley Publishing Company Inc, Canada. 50

Bisa di-download di

kusmansadik.wordpress.com

51

52

52