# ==========================================
# PROJECT DATA CHALLENGE (ANALISIS DATA)
# ==========================================
# 1. Load Library
library(readxl)
## Warning: package 'readxl' was built under R version 4.5.1
library(tseries)
## Warning: package 'tseries' was built under R version 4.5.2
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(forecast)
## Warning: package 'forecast' was built under R version 4.5.3
library(zoo)
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
# 2. Input Data
inflasi <- read_excel("C:/Users/HP141024/Documents/Kuliah Semester 5/Data Challenge/Data_Inflasi_DC.xlsx")
str(inflasi)
## tibble [284 × 3] (S3: tbl_df/tbl/data.frame)
## $ No : num [1:284] 1 2 3 4 5 6 7 8 9 10 ...
## $ Periode : chr [1:284] "Januari 2003" "Februari 2003" "Maret 2003" "April 2003" ...
## $ Data Inflasi: num [1:284] 8.68 7.6 7.17 7.62 7.15 6.98 6.27 6.51 6.33 6.48 ...
summary(inflasi)
## No Periode Data Inflasi
## Min. : 1.00 Length:284 Min. :-0.090
## 1st Qu.: 71.75 Class :character 1st Qu.: 3.042
## Median :142.50 Mode :character Median : 4.475
## Mean :142.50 Mean : 5.249
## 3rd Qu.:213.25 3rd Qu.: 6.678
## Max. :284.00 Max. :18.380
# 3. Pre-processing Data
# Mengecek apakah terdapat nilai kosong (Missing Values / NA)
jumlah_na <- sum(is.na(inflasi$`Data Inflasi`))
cat("Jumlah data NA (kosong):", jumlah_na, "\n")
## Jumlah data NA (kosong): 0
# 4. Model Time Series
# Menentukan start tahun 2003 bulan 1 (Januari)
ts_inflasi <- ts(
inflasi$`Data Inflasi`,
start = c(2003, 1),
frequency = 12
)
# 5. Exploratory Data Analysis (EDA) & Visualisasi
# Plot deret waktu dasar
ts.plot(ts_inflasi, main = "Plot Inflasi Bulanan",
ylab = "Inflasi (%)", xlab = "Tahun",
col = "blue", lwd = 2)
# Dekomposisi Data (melihat unsur trend, musiman dan acak)
dekomposisi <- decompose(ts_inflasi)
plot(dekomposisi)
# Plot ACF dan PACF
acf(ts_inflasi, main = "Plot ACF")
pacf(ts_inflasi, main = "Plot PACF")
# 6. Uji Stasioneritas (Augmented Dickey-Fuller Test)
# H0: Data tidak stasioner (mengandung akar unit)
# H1: Data stasioner
# Tolak H0 jika p-value < 0.05
adf.test(ts_inflasi)
## Warning in adf.test(ts_inflasi): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: ts_inflasi
## Dickey-Fuller = -4.4595, Lag order = 6, p-value = 0.01
## alternative hypothesis: stationary
# 7. Pemodelan SARIMA
# Mencari model terbaik berdasarkan nilai AIC terkecil
model_ts <- auto.arima(ts_inflasi, seasonal = TRUE)
summary(model_ts)
## Series: ts_inflasi
## ARIMA(2,1,2)(0,0,1)[12] with drift
##
## Coefficients:
## ar1 ar2 ma1 ma2 sma1 drift
## -0.5539 -0.6636 0.7001 0.6263 -0.8004 -0.0240
## s.e. 0.2614 0.1666 0.2678 0.2113 0.0445 0.0108
##
## sigma^2 = 0.5106: log likelihood = -309.76
## AIC=633.52 AICc=633.93 BIC=659.04
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.005939128 0.7056905 0.3732131 0.8590653 10.04494 0.1403115
## ACF1
## Training set 0.002853482
# 8. Uji Asumsi Residual (Model Diagnostik)
## a. Uji Autokorelasi Residual (Ljung-Box Test)
# H0: Tidak ada autokorelasi pada residual (Model sudah cukup baik)
# H1: Terdapat autokorelasi pada residual
# Tolak H0 jika p-value < 0.05
Box.test(model_ts$residuals, type = "Ljung-Box")
##
## Box-Ljung test
##
## data: model_ts$residuals
## X-squared = 0.0023369, df = 1, p-value = 0.9614
checkresiduals(model_ts)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(2,1,2)(0,0,1)[12] with drift
## Q* = 15.772, df = 19, p-value = 0.6724
##
## Model df: 5. Total lags used: 24
## b. Uji Normalitas Residual (Shapiro-Wilk Test)
# H0: Residual berdistribusi normal
# H1: Residual tidak berdistribusi normal
# Tolak H0 jika p-value < 0.05
shapiro.test(model_ts$residuals)
##
## Shapiro-Wilk normality test
##
## data: model_ts$residuals
## W = 0.64248, p-value < 2.2e-16
# 9. Forecasting menggunakan seluruh data
# Periode ramalan
h <- 4
# Ramalan SARIMA
ramalan_sarima <- forecast(model_ts, h)
print(ramalan_sarima)
## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## Sep 2026 3.180168 2.2644296 4.095907 1.7796669 4.580669
## Oct 2026 2.932595 1.5396769 4.325512 0.8023109 5.062878
## Nov 2026 2.910842 1.2296502 4.592035 0.3396807 5.482004
## Dec 2026 2.788644 0.8757752 4.701513 -0.1368365 5.714125
# Ramalan SMA
# Function Forecasting SMA
forecast_sma <- function(data, h, k) {
data_baru <- as.numeric(data)
ramalan <- numeric(h)
for (i in 1:h) {
ramalan[i] <- mean(tail(data_baru, k), na.rm = TRUE)
data_baru <- c(data_baru, ramalan[i])
}
return(ramalan)
}
# Forecasting 4 bulan ke depan dengan MA(3)
ramalan_sma <- forecast_sma(ts_inflasi, h, k = 3)
print(ramalan_sma)
## [1] 3.136667 3.068889 3.131852 3.112469
# Periode
nama_bulan <- c("Jan", "Feb", "Mar", "Apr", "May", "Jun",
"Jul", "Aug", "Sep", "Oct", "Nov", "Dec")
periode <- paste(nama_bulan[cycle(ramalan_sarima$mean)],
floor(time(ramalan_sarima$mean)))
# Tabel hasil peramalan
hasil_ramalan <- data.frame(
Periode = periode,
SARIMA = as.numeric(ramalan_sarima$mean),
SMA = as.numeric(ramalan_sma)
)
print(hasil_ramalan)
## Periode SARIMA SMA
## 1 Sep 2026 3.180168 3.136667
## 2 Oct 2026 2.932595 3.068889
## 3 Nov 2026 2.910842 3.131852
## 4 Dec 2026 2.788644 3.112469
# 10. Visualisasi Hasil Ramalan
# Visualisasi Metode SARIMA
plot(ramalan_sarima,
main = "Peramalan Inflasi 4 Bulan Kedepan",
ylab = "Inflasi (%)",
xlab = "Tahun")
# Visualisasi Peramalan 4 Bulan ke Depan
plot(
ramalan_sarima,
include = 12,
main = "Peramalan Inflasi 4 Bulan ke Depan",
ylab = "Inflasi (%)",
xlab = "",
fcol = "blue",
shadecols = c("#cce6ff", "#99ccff"),
xaxt = "n"
)
# Membuat label bulan dan tahun
waktu_plot <- c(tail(time(ts_inflasi), 12), time(ramalan_sarima$mean))
awal <- as.Date(sprintf("%d-%02d-01", floor(waktu_plot[1]), cycle(tail(ts_inflasi, 12))[1]))
label_waktu <- format(seq(awal, by = "month", length.out = length(waktu_plot)),"%b %Y")
axis(1, at = waktu_plot, labels = label_waktu, las = 2, cex.axis = 0.8)
grid()
# Visualisasi Metode SMA
sma <- rollmean(ts_inflasi, k = 3, align = "right", fill = NA)
# Menggabungkan SMA historis dan hasil ramalan
sma_gabungan <- ts(
c(as.numeric(sma), as.numeric(ramalan_sma)),
start = start(ts_inflasi),
frequency = frequency(ts_inflasi)
)
# Visualisasi data
ts.plot(
ts_inflasi, sma_gabungan,
col = c("black", "steelblue"),
lwd = c(1, 2),
main = "Peramalan Inflasi Menggunakan SMA(3)",
ylab = "Inflasi (%)",
xlab = "Tahun"
)
legend(
"topleft",
legend = c("Data Aktual", "SMA(3) dan Ramalan"),
col = c("black", "steelblue"),
lty = 1,
lwd = c(1, 2)
)
# Perbandingan Error SARIMA dan SMA(3)
# Membagi data latih dan data uji
h_uji <- 6
n <- length(ts_inflasi)
data_latih <- ts_inflasi[1:(n - h_uji)]
data_uji <- as.numeric(tail(ts_inflasi, h_uji))
# Ramalan SARIMA pada data uji
model_sarima_uji <- auto.arima(data_latih, seasonal = TRUE)
pred_sarima <- as.numeric(
forecast(model_sarima_uji, h = h_uji)$mean
)
# Ramalan SMA(3) pada data uji
pred_sma <- forecast_sma(data_latih, h = h_uji, k = 3)
# Fungsi menghitung error
hitung_error <- function(aktual, prediksi) {
c(
MAE = mean(abs(aktual - prediksi)),
RMSE = sqrt(mean((aktual - prediksi)^2)),
MAPE = mean(abs(
(aktual[aktual != 0] - prediksi[aktual != 0]) /
aktual[aktual != 0]
)) * 100
)
}
# Membandingkan error
hasil_error <- rbind(
SARIMA = hitung_error(data_uji, pred_sarima),
SMA3 = hitung_error(data_uji, pred_sma)
)
cat("Perbandingan Error SARIMA dan SMA(3):\n")
## Perbandingan Error SARIMA dan SMA(3):
print(round(hasil_error, 4))
## MAE RMSE MAPE
## SARIMA 1.8715 1.9010 63.3011
## SMA3 0.9415 1.0315 32.7726