MPDW 2
Elke Frida Rahmawati
2024-09-02
Pemanggilan Packages
library(dplyr)## Warning: package 'dplyr' was built under R version 4.3.2##
## 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.2library(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) #digunakan untuk uji formal pendeteksian autokorelasi## Warning: package 'lmtest' was built under R version 4.3.2## Loading required package: zoo## Warning: package 'zoo' was built under R version 4.3.2##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericlibrary(orcutt) #untuk membuat model regresi Cochrane-Orcutt## Warning: package 'orcutt' was built under R version 4.3.3library(HoRM) #untuk membuat model regresi Hildreth-Lu## Warning: package 'HoRM' was built under R version 4.3.3Input Data
Data yang digunakan dalam kesempatan kali ini adalah data GDP Madagascar periode tahun 2000-2023.
library(readxl)
dtMadagascar <- read_xls("C:/Users/user/Downloads/DataGDPMadagascar.xls")
head(dtMadagascar)## # A tibble: 6 × 2
## tahun gdp
## <chr> <dbl>
## 1 2000 4629247204.
## 2 2001 5438332602.
## 3 2002 5351701663.
## 4 2003 6372498890.
## 5 2004 5064732626.
## 6 2005 5859269753.View(dtMadagascar)Eksplorasi Data
Sebelum melakukan regresi, akan diperlihatkan plot time-series dari GDP Madagascar Periode 2000-2023
#Membentuk objek time series
data.ts<-ts(dtMadagascar$gdp)
data.ts## Time Series:
## Start = 1
## End = 24
## Frequency = 1
## [1] 4629247204 5438332602 5351701663 6372498890 5064732626 5859269753
## [7] 6395712491 8524620890 10725137724 9616879409 9982711338 11551819618
## [13] 11578975062 12423555455 12522957399 11323020701 11848613858 13176313594
## [19] 13760033282 14104664679 13051441204 14554754117 15302510500 16031702915#Membuat plot time series
ts.plot(data.ts, xlab="Time Period ", ylab="GDP", main= "Time Series Plot of GDP")
points(data.ts)
Selanjutnya akan dilakukan ramalan dan pemulusan dengan metode DMA dan DES karena terlihat pada plot di atas menunjukkan adanya trend.
dt.sma <- SMA(data.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)
}dt.gab <- cbind(aktual = c(data.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,] 4629247204 NA NA NA NA NA
## [2,] 5438332602 NA NA NA NA NA
## [3,] 5351701663 5139760490 NA NA NA NA
## [4,] 6372498890 5720844385 NA NA NA NA
## [5,] 5064732626 5596311060 5817655890 5706983475 110672415 NA
## [6,] 5859269753 5765500423 5908064023 5836782223 71281800 5817655890
## [7,] 6395712491 5773238290 5896348355 5834793322 61555032 5908064023
## [8,] 8524620890 6926534378 8469421073 7697977725 771443348 5896348355
## [9,] 10725137724 8548490368 11479962414 10014226391 1465736023 8469421073
## [10,] 9616879409 9622212674 12135146409 10878679542 1256466868 11479962414
## [11,] 9982711338 10108242824 11472097894 10790170359 681927535 12135146409
## [12,] 11551819618 10383803455 11075237730 10729520593 345717137 11472097894
## [13,] 11578975062 11037835339 12093584939 11565710139 527874800 11075237730
## [14,] 12423555455 11851450045 13372290909 12611870477 760420432 12093584939
## [15,] 12522957399 12175162639 13149189234 12662175937 487013298 13372290909
## [16,] 11323020701 12089844519 12191895421 12140869970 51025451 13149189234
## [17,] 11848613858 11898197320 11585788974 11741993147 -156204173 12191895421
## [18,] 13176313594 12115982718 12278598449 12197290584 81307866 11585788974
## [19,] 13760033282 12928320245 14156627213 13542473729 614153484 12278598449
## [20,] 14104664679 13680337185 15224584789 14452460987 772123802 14156627213
## [21,] 13051441204 13638713055 14084558842 13861635948 222922893 15224584789
## [22,] 14554754117 13903620000 14229079839 14066349920 162729920 14084558842
## [23,] 15302510500 14302901940 15011882491 14657392215 354490275 14229079839
## [24,] 16031702915 15296322511 16887071231 16091696871 795374360 15011882491
## [25,] NA NA NA NA NA 16887071231
## [26,] NA NA NA NA NA 17682445592
## [27,] NA NA NA NA NA 18477819952
## [28,] NA NA NA NA NA 19273194312
## [29,] NA NA NA NA NA 20068568673#Plot time series
ts.plot(dt.gab[,1], xlab="Time Period ", ylab="GDP",
main= "DMA N=3 Data GDP", ylim=c(4629247205,16031702916))
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 = data.ts-dt.ramal[1:length(data.ts)]
SSE.dma = sum(error.dma[6:length(data.ts)]^2)
MSE.dma = mean(error.dma[6:length(data.ts)]^2)
MAPE.dma = mean(abs((error.dma[6:length(data.ts)]/data.ts[6:length(data.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 3.674579e+19
## MSE 1.933989e+18
## MAPE 1.005789e+01Selanjutnya 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<-dtMadagascar[1:18,2]
testing<-dtMadagascar[19:24,2]
#data time series
training.ts<-ts(training)
testing.ts<-ts(testing,start=19)
#eksplorasi data
plot(data.ts, col="red",main="Plot semua data")
points(data.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: 0.9427224
## beta : 0.02608948
## gamma: FALSE
##
## Coefficients:
## [,1]
## a 13140120169
## b 695639774plot(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 13835759943 12447520686 15223999201 11712631214 15958888672
## 20 14531399717 12599948712 16462850723 11577500311 17485299124
## 21 15227039491 12854980233 17599098750 11599287912 18854791071
## 22 15922679265 13162602043 18682756487 11701505420 20143853110
## 23 16618319039 13502726313 19733911765 11853431150 21383206928Selanjutnya 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,] -895716337
## [4,] 182437565
## [5,] -2088858882
## [6,] -65274249mapedes.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 1.893113e+19
## MSE 1.051730e+18
## MAPE 9.028155e+00#Akurasi data testing
selisihdesopt<-ramalandesopt$mean-testing.ts
selisihdesopt## Time Series:
## Start = 19
## End = 23
## Frequency = 1
## testing.ts
## [1,] 75726661
## [2,] 426735039
## [3,] 2175598287
## [4,] 1367925149
## [5,] 1315808539SSEtestingdesopt<-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 1.420606e+18
## MSE 1.420606e+18
## MAPE 6.373726e+00Setelah 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 3.674579e+19 1.893113e+19
## MSE 1.933989e+18 1.051730e+18
## MAPE 1.005789e+01 9.028155e+00Berdasarkan perbandingan akurasi tersebut, terlihat nilai SSE, MSE, dan MAPE metode DES lebih kecil dibandingkan dengan metode DMA. Oleh karena itu, metode peramalan dan pemulusan yang terbaik antara keduanya adalah dengan metode DES.
Setelah melakukan peramalan, data yang telah dimasukkan kemudian dieksplorasi. Eksplorasi pertama yang dilakukan adalah dengan menggunakan scatter plot.
#Eksplorasi Data
#Pembuatan Scatter Plot
dtMadagascar$tahun <- as.numeric(dtMadagascar$tahun)
plot(dtMadagascar, pch = 20, col = "blue",
main = "Scatter Plot Tahun vs GDP",
xlab = "Tahun",
ylab = "GDP")
#Menampilkan Nilai Korelasi
cor(dtMadagascar)## tahun gdp
## tahun 1.0000000 0.9652587
## gdp 0.9652587 1.0000000Berdasarkan scatter plot di atas, terlihat adanya hubungan / korelasi
positif antara peubah tahun dengan nilai 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.9652587\).
Setalah mengetahui adanya hubungan antar dua peubah, maka model regresi dapat ditentukan.
Regresi
#Pembuatan Model Regresi
#model regresi
model<- lm(gdp~tahun, data = dtMadagascar)
summary(model)##
## Call:
## lm(formula = gdp ~ tahun, data = dtMadagascar)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.644e+09 -5.461e+08 7.342e+07 3.742e+08 2.057e+09
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -9.750e+11 5.687e+10 -17.14 3.25e-14 ***
## tahun 4.899e+08 2.827e+07 17.33 2.61e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 958700000 on 22 degrees of freedom
## Multiple R-squared: 0.9317, Adjusted R-squared: 0.9286
## F-statistic: 300.2 on 1 and 22 DF, p-value: 2.611e-14Model yang dihasilkan adalah \[y_i=-975000000000+489900000x_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.9286\). Artinya, sebesar 92.86% 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")
plot(seq(1,24,1), sisaan, col = "steelblue", pch = 20, xlab = "Sisaan", ylab = "Order", main = "Sisaan vs Order")
lines(seq(1,24,1), sisaan, col = "red")
abline(a = 0, b = 0, lwd = 2)
Dua plot di samping kiri digunakan untuk melihat apakah sisaan tidak menyebar normal. Normal Q-Q Plot di atas menunjukkan bahwa sisaan cenderung menyebar normal, tetapi histogram dari sisaan tidak menunjukkan demikian. Selanjutnya, dua plot di samping kanan digunakan untuk melihat autokorelasi. Plot Sisaan vs Fitted Values 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.97346, p-value = 0.7522ks.test(sisaan, "pnorm", mean=mean(sisaan), sd=sd(sisaan))##
## Exact one-sample Kolmogorov-Smirnov test
##
## data: sisaan
## D = 0.10605, p-value = 0.924
## 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 lag 1, lag 7, dan lag 8 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 = 1.0386, p-value = 0.002533
## alternative hypothesis: true autocorrelation is greater than 0Berdasarkan hasil DW Test, didapatkan nilai \(DW = 1.0386\) dan p-value = \(0.002533\). Berdasarkan tabel Durbin-Watson diperoleh nilai \(DL = 1.2728\) dan \(DU = 1.4458\). Nilai DW masih berada di antara nilai DL dan DU. Nilai DW berada diantara 0 dan DL. Artinya dapat dikatakan terdapat autokorelasi positif. Hal ini juga didukung 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)
modelCO## Cochrane-orcutt estimation for first order autocorrelation
##
## Call:
## lm(formula = gdp ~ tahun, data = dtMadagascar)
##
## number of interaction: 3
## rho 0.480333
##
## Durbin-Watson statistic
## (original): 1.03860 , p-value: 2.533e-03
## (transformed): 1.83475 , p-value: 2.607e-01
##
## coefficients:
## (Intercept) tahun
## -969774403871 487283979Hasil keluaran model setelah dilakukan penanganan adalah sebagai berikut. \[y_i=-969774403871+487283979x_t\] Hasil juga menunjukkan bahwa nilai DW dan p-value meningkat menjadi \(1.83475\) dan \(0.2607\). Nilai DW sudah berada pada rentang DU < DW < 4-DU atau \(1.4458 < DW < 2.5542\). Hal tersebut juga didukung dengan nilai p-value > 0.05, artinya belum cukup bukti menyatakan bahwa terdapat autokorelasi sisaan pada taraf nyata 5%. Untuk nilai \(ρ̂\) optimum yang digunakan adalah 0.4803334$. Nilai tersebut dapat diketahui dengan syntax berikut.
#Rho optimum
rho<- modelCO$rho
rho## [1] 0.480333Selanjutnya akan dilakukan transformasi secara manual dengan syntax berikut ini.
dtMadagascar$gdp## [1] 4629247204 5438332602 5351701663 6372498890 5064732626 5859269753
## [7] 6395712491 8524620890 10725137724 9616879409 9982711338 11551819618
## [13] 11578975062 12423555455 12522957399 11323020701 11848613858 13176313594
## [19] 13760033282 14104664679 13051441204 14554754117 15302510500 16031702915dtMadagascar$gdp[-1]## [1] 5438332602 5351701663 6372498890 5064732626 5859269753 6395712491
## [7] 8524620890 10725137724 9616879409 9982711338 11551819618 11578975062
## [13] 12423555455 12522957399 11323020701 11848613858 13176313594 13760033282
## [19] 14104664679 13051441204 14554754117 15302510500 16031702915#Transformasi Manual
#Lag 12 untuk data tahunan
#Menggunakan lag 12
gdp.trans<- dtMadagascar$gdp[-1]-dtMadagascar$gdp[-13]*rho
tahun.trans<- dtMadagascar$tahun[-1]-dtMadagascar$tahun[-13]*rho
modelCOmanual<- lm(gdp.trans~tahun.trans)
summary(modelCOmanual)##
## Call:
## lm(formula = gdp.trans ~ tahun.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.749e+09 -4.658e+08 5.044e+06 2.749e+08 1.862e+09
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -5.050e+11 5.127e+10 -9.849 2.53e-09 ***
## tahun.trans 4.883e+08 4.903e+07 9.960 2.08e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 762200000 on 21 degrees of freedom
## Multiple R-squared: 0.8253, Adjusted R-squared: 0.817
## F-statistic: 99.19 on 1 and 21 DF, p-value: 2.079e-09Hasil model transformasi bukan merupakan model sesungguhnya. Koefisien regresi masih perlu dicari kembali mengikuti \(β_0^*=β_0+ρ ̂β_0\) dan \(β_1^*=β_1\).
#Mencari Penduga Koefisien Regresi setelah Transformasi ke Persamaan Awal
b0bintang <- modelCOmanual$coefficients[-2]
b0 <- b0bintang/(1-rho)
b1 <- modelCOmanual$coefficients[-1]
b0## (Intercept)
## -971757813668b1## tahun.trans
## 488287653Hasil perhitungan koefisien regresi tersebut akan menghasilkan hasil yang sama dengan model yang dihasilkan menggunakan packages.
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 1.846432e+19
## 2 0.2 1.712835e+19
## 3 0.3 1.619681e+19
## 4 0.4 1.566970e+19
## 5 0.5 1.554703e+19
## 6 0.6 1.582878e+19
## 7 0.7 1.651497e+19
## 8 0.8 1.760559e+19
## 9 0.9 1.910064e+19Pertama-tama akan dicari di mana kira-kira \(ρ\) yang menghasilkan SSE minimum. Pada hasil di atas terlihat \(ρ\) minimum ketika 0.5. 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.4 sampai dengan 0.6.
#Rho optimal di sekitar 0.4
rOpt <- seq(0.4,0.6, 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
## 81 0.480 1.553921e+19
## 82 0.481 1.553921e+19
## 80 0.479 1.553924e+19
## 83 0.482 1.553926e+19
## 79 0.478 1.553931e+19
## 84 0.483 1.553935e+19#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.480 , y=1.553921e+19, labels = "rho= 0.480 ", cex = 0.8)
Perhitungan yang dilakukan aplikasi R menunjukkan bahwa
nilai \(ρ\) optimum, yaitu saat SSE
terkecil terdapat pada nilai \(ρ=0.480\). 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.480, model)
summary(modelHL)##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.734e+09 -5.151e+08 1.320e+08 3.013e+08 1.880e+09
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -5.043e+11 5.443e+10 -9.265 7.25e-09 ***
## x 4.873e+08 5.200e+07 9.371 5.97e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 860200000 on 21 degrees of freedom
## Multiple R-squared: 0.807, Adjusted R-squared: 0.7978
## F-statistic: 87.81 on 1 and 21 DF, p-value: 5.972e-09#Transformasi Balik
cat("y = ", coef(modelHL)[1]/(1-0.480), "+", coef(modelHL)[2],"x", sep = "")## y = -969777501864+487285513xSetelah dilakukan tranformasi balik, didapatkan model dengan metode Hildreth-Lu sebagai berikut. \[y_i=-969777501864+487285513x_t\]
#Deteksi autokorelasi
dwtest(modelHL)##
## Durbin-Watson test
##
## data: modelHL
## DW = 1.8342, p-value = 0.2603
## alternative hypothesis: true autocorrelation is greater than 0Hasil uji Durbin-Watson juga menunjukkan bawah nilai DW sebesar \(1.8342\) berada pada selang daerah tidak ada autokorelasi, yaitu pada rentang DU < DW < 4-DU atau \(1.4458 < DW < 2.5542\). Hal tersebut juga didukung oleh p-value sebesar \(0.2603\), di mana p-value > \(\alpha\)=5%. Artinya tak tolak \(H_0\) atau belum cukup bukti menyatakan bahwa ada autokorelasi dalam data nilai gdp dengan metode Hildreth-Lu pada taraf nyata 5%.
Terakhir, akan dibandingkan nilai SSE dari ketiga metode (metode awal, metode Cochrane-Orcutt, dan Hildreth-Lu).
#Perbandingan
sseModelawal <- anova(model)$`Sum Sq`[-1]
sseModelCO <- anova(modelCOmanual)$`Sum Sq`[-1]
sseModelHL <- anova(modelHL)$`Sum Sq`[-1]
mseModelawal <- sseModelawal/length(dtMadagascar$gdp)
mseModelCO <- sseModelCO/length(dtMadagascar$gdp)
mseModelHL <- sseModelHL/length(dtMadagascar$gdp)
akurasi <- matrix(c(sseModelawal,sseModelCO,sseModelHL,
mseModelawal,mseModelCO,mseModelHL),nrow=2,ncol=3,byrow = T)
colnames(akurasi) <- c("Model Awal", "Model Cochrane-Orcutt", "Model Hildreth-Lu")
row.names(akurasi) <- c("SSE","MSE")
akurasi## Model Awal Model Cochrane-Orcutt Model Hildreth-Lu
## SSE 2.022189e+19 1.219847e+19 1.553921e+19
## MSE 8.425790e+17 5.082696e+17 6.474669e+17Berdasarkan hasil tersebut dapat diketahui bahwa hasil penanganan autokorelasi dengan metode Cochrane-Orcutt dan Hildreth-Lu lebih baik dibandingkan model awal ketika autokorelasi masih terjadi, yaitu sebesar \(2.022189e+19\). Namun, hasil penanganan autokorelasi dengan metode Cochrane-Orcutt lebih baik dibandingkan model Hildreth-Lu karena memiliki SSE yang lebih kecil, yaitu sebesar \(1.219847e+19\).
Simpulan
Autokorelasi yang terdapat pada data GDP Madagascar 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 meningkatkan 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. Berdasarkan analisis ini, model Cochrane-Orcutt menghasilkan nilai SSE yang lebih kecil dibandingkan model Hildreth-Lu, artinya model Cochrane-Orcutt lebih baik untuk digunakan.