[MPDW] Tugas Praktikum 2
Siti Arbaynah
2024-09-02
library(dplyr)## Warning: package 'dplyr' was built under R version 4.3.3##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(TTR)## Warning: package 'TTR' was built under R version 4.3.3library(forecast)## Warning: package 'forecast' was built under R version 4.3.3## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoolibrary(lmtest) ## Warning: package 'lmtest' was built under R version 4.3.3## Loading required package: zoo## Warning: package 'zoo' was built under R version 4.3.3##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericlibrary(orcutt) ## Warning: package 'orcutt' was built under R version 4.3.3library(HoRM) ## Warning: package 'HoRM' was built under R version 4.3.3Input Data
Data yang digunakan dalam analisis ini adalah data GDP Negara Luxembourg periode 1978-2000.
library(readxl)## Warning: package 'readxl' was built under R version 4.3.3dt <- read_excel("C:/Users/asus/Documents/Semester 5/MPDW/pertemuan 2/API_NY.GDP.MKTP.CD_DS2_en_excel_v2_3401502.xls", sheet = "Luxembourg")
dt## # A tibble: 64 × 2
## Tahun GDP
## <chr> <dbl>
## 1 1960 709941874
## 2 1961 710163719.
## 3 1962 747846862.
## 4 1963 797902154.
## 5 1964 910877686.
## 6 1965 929477285.
## 7 1966 976717016.
## 8 1967 983052315.
## 9 1968 1075561623
## 10 1969 1245432991
## # ℹ 54 more rowsEksplorasi Data
Serikut plot time-series dari GDP Luxembourg periode 1978-2000.
ts <- ts(dt$GDP)
ts## Time Series:
## Start = 1
## End = 64
## Frequency = 1
## [1] 709941874 710163719 747846862 797902154 910877686 929477285
## [7] 976717016 983052315 1075561623 1245432991 1457768455 1518773421
## [13] 1901697370 2609875802 3183637117 3123333333 3423586207 3789321328
## [19] 4718539772 5516982664 6019805490 5053665797 4602316793 4524217751
## [25] 4438435493 4577211767 6685595088 8320902215 9418167855 10037674038
## [31] 12778792854 13834219728 15518702635 15925521222 17701798891 20853093870
## [37] 20895314658 19563836265 20150053345 21899317599 21230182989 21387533703
## [43] 23649833332 29667268248 35064843793 37672280120 42910146296 51587401416
## [49] 58844277702 54467289898 56213985987 61696281326 59776383527 65203276467
## [55] 68804811898 60071584216 62216885436 65712180343 71000359760 69890505324
## [61] 73699366700 85584105994 81641807866 85755006124#Membuat plot time series
ts.plot(ts, xlab="Time Period ", ylab="GDP", main= "Time Series Plot of GDP")
points(ts)
Selanjutnya akan dilakukan ramalan dan pemulusan dengan metode DMA dan DES karena terlihat pada plot di atas menunjukkan adanya trend.
dt.sma <- SMA(ts, n=3)
dma <- SMA(dt.sma, n = 3)
At <- 2*dt.sma - dma
Bt <- 2/(3-1)*(dt.sma - dma)
dt.dma<- At+Bt
dt.ramal<- c(NA, dt.dma)
t = 1:5
f = c()
for (i in t) {
f[i] = At[length(At)] + Bt[length(Bt)]*(i)
}summary(ts)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 7.099e+08 3.364e+09 1.468e+10 2.542e+10 5.231e+10 8.576e+10dt.gab <- cbind(aktual = c(ts,rep(NA,5)),
pemulusan1 = c(dt.sma,rep(NA,5)),
pemulusan2 = c(dt.dma, rep(NA,5)),
At = c(At, rep(NA,5)),
Bt = c(Bt,rep(NA,5)),
ramalan = c(dt.ramal, f[-1]))
dt.gab## aktual pemulusan1 pemulusan2 At Bt ramalan
## [1,] 709941874 NA NA NA NA NA
## [2,] 710163719 NA NA NA NA NA
## [3,] 747846862 722650819 NA NA NA NA
## [4,] 797902154 751970912 NA NA NA NA
## [5,] 910877686 818875568 927628504 873252036 54376468 NA
## [6,] 929477285 879419042 1004746778 942082910 62663868 927628504
## [7,] 976717016 939023996 1058859583 998941789 59917794 1004746778
## [8,] 983052315 963082205 1034896454 998989329 35907124 1058859583
## [9,] 1075561623 1011776984 1092742163 1052259574 40482589 1034896454
## [10,] 1245432991 1101348976 1253241485 1177295231 75946254 1092742163
## [11,] 1457768455 1259587690 1530287302 1394937496 135349806 1253241485
## [12,] 1518773421 1407324956 1709800453 1558562704 151237748 1530287302
## [13,] 1901697370 1626079749 2016244317 1821162033 195082284 1709800453
## [14,] 2609875802 2010115531 2667999769 2339057650 328942119 2016244317
## [15,] 3183637117 2565070096 3561033372 3063051734 497981638 2667999769
## [16,] 3123333333 2972282084 3885201111 3428741598 456459514 3561033372
## [17,] 3423586207 3243518886 3876642613 3560080749 316561864 3885201111
## [18,] 3789321328 3445413623 3895431140 3670422381 225008759 3876642613
## [19,] 4718539772 3977149102 4820726233 4398937668 421788565 3895431140
## [20,] 5516982664 4674947921 5959836666 5317392294 642444373 4820726233
## [21,] 6019805490 5418442642 6874968149 6146705395 728262753 5959836666
## [22,] 5053665797 5530151317 6174759364 5852455341 322304024 6874968149
## [23,] 4602316793 5225262693 4893216978 5059239836 -166022857 6174759364
## [24,] 4524217751 4726733447 3858768703 4292751075 -433982372 4893216978
## [25,] 4438435493 4521656679 3915868157 4218762418 -302894261 3858768703
## [26,] 4577211767 4513288337 4365412702 4439350520 -73937817 3915868157
## [27,] 6685595088 5233747449 6188780704 5711264077 477516628 4365412702
## [28,] 8320902215 6527903023 8733749864 7630826443 1102923420 6188780704
## [29,] 9418167855 8141555053 11155861474 9648708264 1507153211 8733749864
## [30,] 10037674038 9258914703 11824495589 10541705146 1282790443 11155861474
## [31,] 12778792854 10744878249 13471069411 12107973830 1363095581 11824495589
## [32,] 13834219728 12216895540 15170227626 13693561583 1476666043 13471069411
## [33,] 15518702635 14043905072 17461262643 15752583858 1708678785 15170227626
## [34,] 15925521222 15092814528 17709366825 16401090676 1308276148 17461262643
## [35,] 17701798891 16382007583 18800204626 17591106104 1209098522 17709366825
## [36,] 20853093870 18160137994 21390440579 19775289287 1615151293 18800204626
## [37,] 20895314658 19816735806 23210953163 21513844485 1697108679 21390440579
## [38,] 19563836265 20437414931 22369385639 21403400285 965985354 23210953163
## [39,] 20150053345 20203068089 20304391717 20253729903 50661814 22369385639
## [40,] 21899317599 20537735736 20827728038 20682731887 144996151 20304391717
## [41,] 21230182989 21093184644 22056894953 21575039799 481855154 20827728038
## [42,] 21387533703 21505678097 22425968639 21965823368 460145271 22056894953
## [43,] 23649833332 22089183341 23142185969 22615684655 526501314 22425968639
## [44,] 29667268248 24901545094 29040364261 26970954678 2069409583 23142185969
## [45,] 35064843793 29460648458 37414360777 33437504618 3976856160 29040364261
## [46,] 37672280120 34134797387 43406398202 38770597794 4635800407 37414360777
## [47,] 42910146296 38549090070 47550912933 43050001501 4500911432 43406398202
## [48,] 51587401416 44056609277 54342830009 49199719643 5143110366 47550912933
## [49,] 58844277702 51113941805 64195397980 57654669892 6540728087 54342830009
## [50,] 54467289898 54966323005 64807719624 59887021315 4920698310 64195397980
## [51,] 56213985987 56508517862 61133031805 58820774834 2312256972 64807719624
## [52,] 61696281326 57459185737 59754872808 58607029272 1147843535 61133031805
## [53,] 59776383527 59228883613 62222259365 60725571489 1496687876 59754872808
## [54,] 65203276467 62225313773 67400352571 64812833172 2587519399 62222259365
## [55,] 68804811898 64594823964 69751790992 67173307478 2578483514 67400352571
## [56,] 60071584216 64693224194 66404097960 65548661077 855436883 69751790992
## [57,] 62216885436 63697760517 62436075767 63066918142 -630842375 66404097960
## [58,] 65712180343 62666883332 60628737967 61647810649 -1019072682 62436075767
## [59,] 71000359760 66309808513 70479790631 68394799572 2084991059 60628737967
## [60,] 69890505324 68867681809 74706796325 71787239067 2919557258 70479790631
## [61,] 73699366700 71530077261 76785186728 74157631995 2627554734 74706796325
## [62,] 85584105994 76391326006 84647921300 80519623653 4128297647 76785186728
## [63,] 81641807866 80308426853 88772060480 84540243666 4231816813 84647921300
## [64,] 85755006124 84326973328 92296435859 88311704594 3984731266 88772060480
## [65,] NA NA NA NA NA 92296435859
## [66,] NA NA NA NA NA 96281167125
## [67,] NA NA NA NA NA 100265898390
## [68,] NA NA NA NA NA 104250629656
## [69,] NA NA NA NA NA 108235360921#Plot time series
ts.plot(dt.gab[,1], xlab="Time Period ", ylab="GDP",
main= "DMA N=3 Data GDP Luxembourg", ylim=c(6.099e+08,11.576e+10))
points(dt.gab[,1])
points(dt.gab[,3])
points(dt.gab[,6])
lines(dt.gab[,3],col="green",lwd=2)
lines(dt.gab[,6],col="red",lwd=2)
legend("topleft",c("data aktual","data pemulusan","data peramalan"),
lty=8, col=c("black","green","red"), cex=0.8)
Selanjutnya akan dilihat keakuratan dari metode DMA
#Menghitung nilai keakuratan
error.dma = ts -dt.ramal[1:length(ts)]
SSE.dma = sum(error.dma[6:length(ts)]^2)
MSE.dma = mean(error.dma[6:length(ts)]^2)
MAPE.dma = mean(abs((error.dma[6:length(ts)]/ts[6:length(ts)])*100))
akurasi.dma <- matrix(c(SSE.dma, MSE.dma, MAPE.dma))
row.names(akurasi.dma)<- c("SSE", "MSE", "MAPE")
colnames(akurasi.dma) <- c("Akurasi m = 3")
akurasi.dma## Akurasi m = 3
## SSE 6.733848e+20
## MSE 1.141330e+19
## MAPE 9.776567e+00Selanjutnya akan digunakan metode Double Exponential Smoothing dengan cara sebagai berikut.
Pertama akan data akan dibagi menjadi data training dan data testing.
#membagi training dan testing
training<-dt[1:18,2]
testing<-dt[19:23,2]
#data time series
training.ts<-ts(training)
testing.ts<-ts(testing,start=19)
#eksplorasi data
plot(ts, col="red",main="Plot semua data")
points(ts)
plot(training.ts, col="blue",main="Plot data training")
points(training.ts)
Selanjutnya akan dilakukan pemulusan dengan DES, kali ini langsung dicari lambda dan gamma optimum sebagai berikut. Nilai lambda dan gamma optimum dapat dilihat pada smoothing parameters alpha untuk nilai lambda dan beta untuk nilai gamma.
#Lamda dan gamma optimum
des.opt<- HoltWinters(training.ts, gamma = FALSE)
des.opt## Holt-Winters exponential smoothing with trend and without seasonal component.
##
## Call:
## HoltWinters(x = training.ts, gamma = FALSE)
##
## Smoothing parameters:
## alpha: 1
## beta : 0.3352126
## gamma: FALSE
##
## Coefficients:
## [,1]
## a 3789321328
## b 309452586plot(des.opt)
legend("topleft", c("Data Aktual", "Peramalan"), col = c("black", "red"),
lty = c(1,1))
#ramalan
ramalandesopt<- forecast(des.opt, h=5)
ramalandesopt## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 19 4098773914 3836335531 4361212297 3697409043 4500138785
## 20 4408226500 3970434541 4846018459 3738681478 5077771522
## 21 4717679086 4098129202 5337228970 3770159305 5665198867
## 22 5027131672 4214178378 5840084966 3783826899 6270436445
## 23 5336584258 4317617539 6355550977 3778209141 6894959375Selanjutnya akan dicari akurasi dari metode DES.
ssedes.train<-des.opt$SSE
msedes.train<-ssedes.train/length(training.ts)
sisaandes<-ramalandesopt$residuals
head(sisaandes)## Time Series:
## Start = 1
## End = 6
## Frequency = 1
## x
## [1,] NA
## [2,] NA
## [3,] 37461297
## [4,] 37275948
## [5,] 87700820
## [6,] -36073534mapedes.train <- sum(abs(sisaandes[3:length(training.ts)]/training.ts[3:length(training.ts)])*100)/length(training.ts)
akurasides.opt <- matrix(c(ssedes.train,msedes.train,mapedes.train))
row.names(akurasides.opt)<- c("SSE", "MSE", "MAPE")
colnames(akurasides.opt) <- c("Akurasi lamda dan gamma optimum")
akurasides.opt## Akurasi lamda dan gamma optimum
## SSE 6.822204e+17
## MSE 3.790113e+16
## MAPE 6.386995e+00#Akurasi data testing
selisihdesopt<-ramalandesopt$mean-testing.ts
selisihdesopt## Time Series:
## Start = 19
## End = 23
## Frequency = 1
## testing.ts
## [1,] -619765858
## [2,] -1108756164
## [3,] -1302126404
## [4,] -26534125
## [5,] 734267465SSEtestingdesopt<-sum(selisihdesopt^2)
SSEtestingdesopt<-SSEtestingdesopt/length(testing.ts)
MAPEtestingdesopt<-sum(abs(selisihdesopt/testing.ts)*100)/length(testing.ts)
akurasiDesTesting <- matrix(c(SSEtestingdesopt,SSEtestingdesopt,MAPEtestingdesopt))
row.names(akurasiDesTesting)<- c("SSE", "MSE", "MAPE")
colnames(akurasiDesTesting) <- c("Akurasi lamda dan gamma optimum")
akurasiDesTesting## Akurasi lamda dan gamma optimum
## SSE 7.697672e+17
## MSE 7.697672e+17
## MAPE 1.426838e+01Setelah didapatkan nilai akurasi untuk metode DMA dan DES, selanjutnya akan dibandingkan keakuratan antar metode keduanya.
cbind(akurasi.dma, akurasides.opt)## Akurasi m = 3 Akurasi lamda dan gamma optimum
## SSE 6.733848e+20 6.822204e+17
## MSE 1.141330e+19 3.790113e+16
## MAPE 9.776567e+00 6.386995e+00Berdasarkan perbandingan akurasi tersebut terlihat berdasarkan nilai SSE dan MSE, DES memiliki kesalahan yang lebih rendah dibandingkan dengan DMA. Akan tetapi, ketika melihat MAPE, DMA memiliki tingkat kesalahan yang lebih rendah (9.78%) dibandingkan dengan DES (6.39%). Dalam praktiknya, pilihan antara DMA dan DES akan tergantung pada konteks dan preferensi analis. DMA lebih sederhana dan mudah diinterpretasikan, sementara DES memodelkan tren dan musiman dengan lebih baik.
Setelah melakukan peramalan, data yang telah dimasukkan kemudian dieksplorasi. Eksplorasi pertama yang dilakukan adalah dengan menggunakan scatter plot.
#Eksplorasi Data
#Pembuatan Scatter Plot
th <- as.numeric(dt$Tahun)
plot(dt, pch = 20, col = "blue",
main = "Scatter Plot Tahun vs GDP",
xlab = "Tahun",
ylab = "GDP")
#Menampilkan Nilai Korelasi
cor(th,ts)## [1] 0.9228464Berdasarkan scatter plot di atas, terlihat adanya hubungan / korelasi
positif antara peubah tahun dengan GDP, terlihat titik-titik pada plot
yang naik ke arah kanan atas. Hal tersebut juga diperkuat dengan hasil
perhitungan aplikasi R di mana didapatkan nilai korelasi
sebesar \(0.9228464\).
Setalah mengetahui adanya hubungan antar dua peubah, maka model regresi dapat ditentukan.
Regresi
#Pembuatan Model Regresi
#model regresi
model<- lm(ts~th, data = dt)
summary(model)##
## Call:
## lm(formula = ts ~ th, data = dt)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.682e+10 -1.069e+10 6.176e+08 8.878e+09 2.045e+10
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.655e+12 1.421e+11 -18.69 <2e-16 ***
## th 1.346e+09 7.135e+07 18.87 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.054e+10 on 62 degrees of freedom
## Multiple R-squared: 0.8516, Adjusted R-squared: 0.8493
## F-statistic: 355.9 on 1 and 62 DF, p-value: < 2.2e-16Model yang dihasilkan adalah \[y_i=-2655000000000+ 1346000000x_t\] Berdasarkan ringkasan model dapat diketahui bahwa hasil uji F memiliki p-value < \(\alpha\) (5%). Artinya, minimal terdapat satu variabel yang berpengaruh nyata terhadap model. Hasil uji-t parsial kedua parameter regresi, yaitu intersep dan koefisien regresi juga menunjukkan hal yang sama, yaitu memiliki p-value < \(\alpha\) (5%) sehingga nyata dalam taraf 5%. Selanjutnya dapat dilihat juga nilai \(R^2=0.8493\). Artinya, sebesar 84.93% keragaman nilai GDP dapat dijelaskan oleh peubah tahun. Hasil ini menunjukkan hasil yang bagus, seolah mendapatkan hasil terbaik. Namun, kita perlu melakukan uji terhadap sisaannya seperti berikut ini.
#sisaan dan fitted value
sisaan<- residuals(model)
fitValue<- predict(model)
#Diagnostik dengan eksploratif
par(mfrow = c(2,2))
qqnorm(sisaan)
qqline(sisaan, col = "steelblue", lwd = 2)
plot(fitValue, sisaan, col = "steelblue", pch = 20, xlab = "Sisaan", ylab = "Fitted Values", main = "Sisaan vs Fitted Values")
abline(a = 0, b = 0, lwd = 2)
hist(sisaan, col = "steelblue")
lines(seq(1,64,1), sisaan, col = "red")
abline(a = 0, b = 0, lwd = 2)
plot(seq(1,64,1), sisaan, col = "steelblue", pch = 20, xlab = "Sisaan", ylab = "Order", main = "Sisaan vs Order")
length(sisaan)## [1] 64Dua plot di samping kiri digunakan untuk melihat apakah sisaan menyebar normal. Normal Q-Q Plot di atas menunjukkan bahwa sisaan tidak menyebar normal, Selanjutnya, dua plot di samping kanan digunakan untuk melihat autokorelasi. Plot Sisaan vs Fitted Value dan Plot Sisaan vs Order menunjukkan adanya pola pada sisaan. Untuk lebih lanjut akan digunakan uji formal melihat normalitas sisaan dan plot ACF dan PACF untuk melihat apakah ada autokorelasi atau tidak.
#Melihat Sisaan Menyebar Normal/Tidak
#H0: sisaan mengikuti sebaran normal
#H1: sisaan tidak mengikuti sebaran normal
shapiro.test(sisaan)##
## Shapiro-Wilk normality test
##
## data: sisaan
## W = 0.93903, p-value = 0.003422ks.test(sisaan, "pnorm", mean=mean(sisaan), sd=sd(sisaan))##
## Exact one-sample Kolmogorov-Smirnov test
##
## data: sisaan
## D = 0.14083, p-value = 0.1433
## alternative hypothesis: two-sidedBerdasarkan uji formal Saphiro-Wilk dan Kolmogorov-Smirnov didapatkan nilai p-value > \(\alpha\) (5%). Artinya, cukup bukti untuk menyatakan sisaan berdistribusi normal.
#ACF dan PACF identifikasi autokorelasi
par(mfrow = c(1,2))
acf(sisaan)
pacf(sisaan)
Berdasarkan plot ACF terlihat hingga lag 10 signifikan yang ditandai dengan garis yang keluar dari rentang batas, sedangkan pada plot PACF, terlihat lag 1 signifikan. Kedua hal ini menandakan adanya autokeorelasi, untuk lebih memastikan akan dilakukan uji formal dengan uji Durbin Watson.
#Deteksi autokorelasi dengan uji-Durbin Watson
#H0: tidak ada autokorelasi
#H1: ada autokorelasi
dwtest(model)##
## Durbin-Watson test
##
## data: model
## DW = 0.082686, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0Berdasarkan hasil DW Test, didapatkan nilai \(DW = 0.082686\) dan p-value = \(2.2e-16\). Berdasarkan tabel Durbin-Watson diperoleh nilai \(DL = 0.9771\) dan \(DU = 1.331\). Nilai DW masih berada di antara nilai DL dan DU. Artinya, berada di daerah inkonklusif, tidak dapat dikatakan berada di daerah autokorelasi positif maupun bebas dari autokorelasi. Namun, dengan nilai p-value < 0.05 dapat disimpulkan bahwa tolak H0, cukup bukti mengatakan adanya autokorelasi. Oleh karena itu, diperlukan penangan autokorelasi. Penanganan yang akan digunakan menggunakan dua metode, yaitu Cochrane-Orcutt dan Hildret-Lu.
Penanganan Autokorelasi
Metode Cochrane-Orcutt
Penanganan metode Cochrane-Orcutt dapat dilakukan dengan bantuan
packages Orcutt pada aplikasi R maupun secara manual.
Berikut ini ditampilkan cara menggunakan bantuan library
packages Orcutt.
#Penanganan Autokorelasi Cochrane-Orcutt
modelCO<-cochrane.orcutt(model)## Warning in cochrane.orcutt(model): Did not convergemodelCO## Cochrane-orcutt estimation for first order autocorrelation
##
## Call:
## lm(formula = ts ~ th, data = dt)
##
## number of interaction: 100
## rho 0.961925
##
## Durbin-Watson statistic
## (original): 0.08269 , p-value: 3.938e-39
## (transformed): NA , p-value: NA
##
## coefficients:
## [1] NA#Rho optimum
rho<- modelCO$rho
rho <- as.numeric(rho)
rho## [1] 0.9619254Selanjutnya akan dilakukan transformasi secara manual dengan syntax berikut ini.
dt$GDP## [1] 709941874 710163719 747846862 797902154 910877686 929477285
## [7] 976717016 983052315 1075561623 1245432991 1457768455 1518773421
## [13] 1901697370 2609875802 3183637117 3123333333 3423586207 3789321328
## [19] 4718539772 5516982664 6019805490 5053665797 4602316793 4524217751
## [25] 4438435493 4577211767 6685595088 8320902215 9418167855 10037674038
## [31] 12778792854 13834219728 15518702635 15925521222 17701798891 20853093870
## [37] 20895314658 19563836265 20150053345 21899317599 21230182989 21387533703
## [43] 23649833332 29667268248 35064843793 37672280120 42910146296 51587401416
## [49] 58844277702 54467289898 56213985987 61696281326 59776383527 65203276467
## [55] 68804811898 60071584216 62216885436 65712180343 71000359760 69890505324
## [61] 73699366700 85584105994 81641807866 85755006124dt$GDP[-1]## [1] 710163719 747846862 797902154 910877686 929477285 976717016
## [7] 983052315 1075561623 1245432991 1457768455 1518773421 1901697370
## [13] 2609875802 3183637117 3123333333 3423586207 3789321328 4718539772
## [19] 5516982664 6019805490 5053665797 4602316793 4524217751 4438435493
## [25] 4577211767 6685595088 8320902215 9418167855 10037674038 12778792854
## [31] 13834219728 15518702635 15925521222 17701798891 20853093870 20895314658
## [37] 19563836265 20150053345 21899317599 21230182989 21387533703 23649833332
## [43] 29667268248 35064843793 37672280120 42910146296 51587401416 58844277702
## [49] 54467289898 56213985987 61696281326 59776383527 65203276467 68804811898
## [55] 60071584216 62216885436 65712180343 71000359760 69890505324 73699366700
## [61] 85584105994 81641807866 85755006124as.numeric(dt$GDP[-12])## [1] 709941874 710163719 747846862 797902154 910877686 929477285
## [7] 976717016 983052315 1075561623 1245432991 1457768455 1901697370
## [13] 2609875802 3183637117 3123333333 3423586207 3789321328 4718539772
## [19] 5516982664 6019805490 5053665797 4602316793 4524217751 4438435493
## [25] 4577211767 6685595088 8320902215 9418167855 10037674038 12778792854
## [31] 13834219728 15518702635 15925521222 17701798891 20853093870 20895314658
## [37] 19563836265 20150053345 21899317599 21230182989 21387533703 23649833332
## [43] 29667268248 35064843793 37672280120 42910146296 51587401416 58844277702
## [49] 54467289898 56213985987 61696281326 59776383527 65203276467 68804811898
## [55] 60071584216 62216885436 65712180343 71000359760 69890505324 73699366700
## [61] 85584105994 81641807866 85755006124a <- as.numeric(dt$Tahun[-1])
a## [1] 1961 1962 1963 1964 1965 1966 1967 1968 1969 1970 1971 1972 1973 1974 1975
## [16] 1976 1977 1978 1979 1980 1981 1982 1983 1984 1985 1986 1987 1988 1989 1990
## [31] 1991 1992 1993 1994 1995 1996 1997 1998 1999 2000 2001 2002 2003 2004 2005
## [46] 2006 2007 2008 2009 2010 2011 2012 2013 2014 2015 2016 2017 2018 2019 2020
## [61] 2021 2022 2023View(dt)Metode Hildreth-Lu
Penanganan kedua adalah menggunakan metode Hildreth-Lu. Metode ini
akan mencari nilai SSE terkecil dan dapat dicari secara manual maupun
menggunakan packages. Jika menggunakan packages, gunakan
library packages HORM.
#Penanganan Autokorelasi Hildreth lu
# Hildreth-Lu
hildreth.lu.func<- function(r, model){
x <- model.matrix(model)[,-1]
y <- model.response(model.frame(model))
n <- length(y)
t <- 2:n
y <- y[t]-r*y[t-1]
x <- x[t]-r*x[t-1]
return(lm(y~x))
}
#Pencariab rho yang meminimumkan SSE
r <- c(seq(0.1,0.9, by= 0.1))
tab <- data.frame("rho" = r, "SSE" = sapply(r, function(i){deviance(hildreth.lu.func(i, model))}))
round(tab, 4)## rho SSE
## 1 0.1 5.365585e+21
## 2 0.2 4.301700e+21
## 3 0.3 3.368844e+21
## 4 0.4 2.567016e+21
## 5 0.5 1.896216e+21
## 6 0.6 1.356445e+21
## 7 0.7 9.477021e+20
## 8 0.8 6.699878e+20
## 9 0.9 5.233018e+20Pertama-tama akan dicari di mana kira-kira \(ρ\) yang menghasilkan SSE minimum. Pada hasil di atas terlihat \(ρ\) minimum ketika 0.3. Namun, hasil tersebut masih kurang teliti sehingga akan dicari kembali \(ρ\) yang lebih optimum dengan ketelitian yang lebih. Jika sebelumnya jarak antar \(ρ\) yang dicari adalah 0.1, kali ini jarak antar \(ρ\) adalah 0.001 dan dilakukan pada selang 0.2 sampai dengan 0.5.
#Rho optimal di sekitar 0.4
rOpt <- seq(0.2,0.5, by= 0.001)
tabOpt <- data.frame("rho" = rOpt, "SSE" = sapply(rOpt, function(i){deviance(hildreth.lu.func(i, model))}))
head(tabOpt[order(tabOpt$SSE),])## rho SSE
## 301 0.500 1.896216e+21
## 300 0.499 1.902275e+21
## 299 0.498 1.908348e+21
## 298 0.497 1.914434e+21
## 297 0.496 1.920532e+21
## 296 0.495 1.926644e+21#Grafik SSE optimum
par(mfrow = c(1,1))
plot(tab$SSE ~ tab$rho , type = "l", xlab = "Rho", ylab = "SSE")
abline(v = tabOpt[tabOpt$SSE==min(tabOpt$SSE),"rho"], lty = 2, col="red",lwd=2)
text(x=0.341, y=0.2397500, labels = "rho=0.341", cex = 0.8)
Perhitungan yang dilakukan aplikasi R menunjukkan bahwa
nilai \(ρ\) optimum, yaitu saat SSE
terkecil terdapat pada nilai \(ρ=0.341\). Hal tersebut juga ditunjukkan
pada plot. Selanjutnya, model dapat didapatkan dengan mengevaluasi nilai
\(ρ\) ke dalam fungsi
hildreth.lu.func, serta dilanjutkan dengan pengujian
autokorelasi dengan uji Durbin-Watson. Namun, setelah pengecekan
tersebut tidak lupa koefisien regresi tersebut digunakan untuk
transformasi balik. Persamaan hasil transformasi itulah yang menjadi
persamaan sesungguhnya.
#Model terbaik
modelHL <- hildreth.lu.func(0.341, model)
summary(modelHL)##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.164e+10 -6.737e+09 -3.333e+08 5.950e+09 1.620e+10
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.823e+12 9.721e+10 -18.75 <2e-16 ***
## x 1.402e+09 7.403e+07 18.93 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 7.041e+09 on 61 degrees of freedom
## Multiple R-squared: 0.8546, Adjusted R-squared: 0.8522
## F-statistic: 358.5 on 1 and 61 DF, p-value: < 2.2e-16#Transformasi Balik
cat("y = ", coef(modelHL)[1]/(1-0.341), "+", coef(modelHL)[2],"x", sep = "")## y = -2.766239e+12+1401618884xSetelah dilakukan tranformasi balik, didapatkan model dengan metode Hildreth-Lu sebagai berikut. \[y_i=-1062.032+0.5597492x_t\]
#Deteksi autokorelasi
dwtest(modelHL)##
## Durbin-Watson test
##
## data: modelHL
## DW = 0.20742, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0Hasil uji Durbin-Watson juga menunjukkan bawah nilai DW sebesar \(1.4092\) berada pada selang daerah tidak ada autokorelasi, yaitu pada rentang DU < DW < 4-DU atau \(1.331 < DW < 2.669\). Hal tersebut juga didukung oleh p-value sebesar \(0.09404\), di mana p-value > \(\alpha\)=5%. Artinya tak tolak \(H_0\) atau belum cukup bukti menyatakan bahwa ada autokorelasi dalam data nilai IPM dengan metode Hildreth-Lu pada taraf nyata 5%.
Terakhir, akan dibandingkan nilai SSE dari ketiga metode (metode awal, metode Cochrane-Orcutt, dan Hildreth-Lu).
Simpulan
Autokorelasi yang terdapat pada data GDP negara Luxembourg terjadi akibat adanya korelasi di antara unsur penyusunnya. Indikator GDP yang erat hubungannya dengan perekonomian sangat rawan menjadi penyebab adanya autokorelasi. Adanya autokorelasi menyebabkan model regresi kurang baik karena akan meingkatkan galatnya. Autokorelasi dapat dideteksi secara eksploratif melalui plot sisaan, ACF, dan PACF, serta dengan uji formal Durbin-Watson. Namun, autokorelasi tersebut dapat ditangani dengan metode Cochrane-Orcutt dan Hildreth-Lu. Model Cochrane-Orcutt menghasilkan nilai SSE yang lebih kecil dibandingkan model Hildreth-Lu, artinya model Cochrane-Orcutt lebih baik untuk digunakan untuk analisis ini.