---
title: "Revisi MPDW"
author: "Dalilah Husna - G1401221003"
date: "2025-04-13"
output: html_document
---



## LOAD LIBRARY

```{r}
library(ggplot2)
library(tsibble)
library(tseries)
library(MASS)
library(forecast)
library(TSA)
library(TTR)
library(aTSA)
library(graphics)
```

## READ DATA

```{r}
library(readxl)
data <- read_xlsx("C:\\Users\\USER\\OneDrive - apps.ipb.ac.id\\SEM 5 Husna\\Metode Peramalan Deret Waktu\\Tugas Kelompok\\mpdwlengkap.xlsx")
data.ts <- ts(data$jumlah)
data
```

## PLOT

```{r}
plot(data.ts, xlab = "Tanggal", ylab = "Jumlah Kunjungan")
```

## DATA LATIH

```{r}
datatrain<-data$jumlah[1:400]
train.ts<-ts(datatrain)
plot.ts(train.ts, lty=1, xlab="Waktu", ylab="Jumlah Kunjungan", main="Plot Data Train")
```

## DATA TEST

```{r}
datatest<-data$jumlah[401:500]
test.ts<-ts(datatest)
plot.ts(test.ts, lty=1, xlab="Waktu", ylab="Jumlah Kunjungan", main="Plot Data Test")
```

## UJI STASIONERITAS DATA

### PLOT ACF

```{r}
acf(train.ts)
```

### PLOT BOXCOX
```{r}
index <- seq(1:400)
bc = boxcox(train.ts~index, lambda = seq(-2,2,by=0.01))
lambda <- bc$x[which.max(bc$y)]
lambda
bc$x[bc$y > max(bc$y) - 1/2 * qchisq(.95,1)]
```

## TRANSFORMASI BOXCOX

```{r}
transformed <- BoxCox(data.ts, 0.14)
transformed.ts <- ts(transformed)
```

```{r}
trt <- transformed.ts[1:400]
trt.ts <- ts(trt)
tet <- transformed.ts[401:500]
tet.ts <- ts(tet)
```

```{r}
acf(trt.ts)
```

## DIFFERENCING

```{r}
trtdiff <- diff(trt.ts, differences = 1)
plot.ts(trtdiff)
```

```{r}
acf(trtdiff)
```

```{r}
pacf(trtdiff)
```

```{r}
eacf(trtdiff)
```
MA(1) ARIMA(0,1,1) ARIMA(0,1,2) ARIMA(1,1,1) ARIMA(1,1,2)

## PENDUGAAN PARAMETER MODEL TENTATIF

### MA(1)

```{r}
model1 <- Arima(trtdiff, order = c(0,0,1), method = "ML")
summary(model1)
lmtest::coeftest(model1)
```
### ARIMA(0,1,1)

```{r}
model2 <- Arima(trtdiff, order = c(0,0,1), method = "ML")
summary(model2)
lmtest::coeftest(model2)
```
### ARIMA(0,1,2)

```{r}
model3 <- Arima(trtdiff, order = c(0,0,2), method = "ML")
summary(model3)
lmtest::coeftest(model3)
```
### ARIMA(1,1,2)

```{r}
model4 <- Arima(trtdiff, order = c(1,0,2), method = "ML")
summary(model4)
lmtest::coeftest(model4)
```
### ARIMA(1,1,1)

```{r}
model5 <- Arima(trtdiff, order = c(1,0,1), method = "ML")
summary(model5)
lmtest::coeftest(model5)
```

## UJI DIAGNOSTIK

### EKSPLORASI SISAAN

```{r}
sisaanmodel <- model4$residuals
par(mfrow=c(2,2))
qqnorm(sisaanmodel)
qqline(sisaanmodel, col = "blue", lwd = 2)
plot(c(1:length(sisaanmodel)),sisaanmodel)
acf(sisaanmodel)
pacf(sisaanmodel)
```
## UJI FORMAL

```{r}
#1) Sisaan Menyebar Normal
ks.test(sisaanmodel,"pnorm") 
```
```{r}
#2) Sisaan saling bebas/tidak ada autokorelasi 
Box.test(sisaanmodel, type = "Ljung")  #tak tolak H0 > sisaan saling bebas
```
```{r}
#3) Sisaan homogen 
Box.test((sisaanmodel)^2, type = "Ljung")  #tak tolak H0 > sisaan homogen
```
```{r}
#4) Nilai tengah sisaan sama dengan nol 
t.test(sisaanmodel, mu = 0, conf.level = 0.95)  #tak tolak h0 > nilai tengah sisaan sama dengan 0
```

## OVERFITTING
ARIMA(1,1,2) jadi ARIMA(1,1,3) ARIMA(2,1,2)

```{r}
model1a=Arima(trtdiff, order=c(1,0,3),method="ML")
summary(model1a) 
lmtest::coeftest(model1a) 
```
```{r}
model1b=Arima(trtdiff, order=c(2,0,2),method="ML")
summary(model1b) 
lmtest::coeftest(model1b) 
```
ARIMA(2,1,2)

## SISAAN MODEL TERBAIK

```{r}
sisaanmodell <- model1b$residuals
par(mfrow=c(2,2))
qqnorm(sisaanmodell)
qqline(sisaanmodell, col = "black", lwd = 2)
plot(c(1:length(sisaanmodell)), sisaanmodell,
     type = "p",              # Plot tipe point
     main = "Residual vs Order", 
     xlab = "Order", 
     ylab = "Residual", 
     col = "black")            # Warna residual
acf(sisaanmodell)
pacf(sisaanmodell)
```

```{r}
#1) Sisaan Menyebar Normal
ks.test(sisaanmodell,"pnorm") 
```
```{r}
#2) Sisaan saling bebas/tidak ada autokorelasi 
Box.test(sisaanmodell, type = "Ljung")  #tak tolak H0 > sisaan saling bebas
```
```{r}
#3) Sisaan homogen 
Box.test((sisaanmodell)^2, type = "Ljung")  #tak tolak H0 > sisaan homogen
```
```{r}
#4) Nilai tengah sisaan sama dengan nol 
t.test(sisaanmodell, mu = 0, conf.level = 0.95)  #tak tolak h0 > nilai tengah sisaan sama dengan 0
```
## PERAMALAN
```{r}
ramalan <- forecast::forecast(model1b, h = 100) 
ramalan
dataramalan <- ramalan$mean
plot(ramalan)
```
```{r}
tra <- InvBoxCox(ramalan$mean, 0.14)
```

```{r}
traa <- diffinv(tra, differences = 1)
```

```{r}
perbandingan.da<-matrix(data=c(head(test.ts, n=100), tra[-1]),
                     nrow = 100, ncol = 2)
colnames(perbandingan.da)<-c("Aktual","Hasil Forecast")
perbandingan.da
```
```{r}
accuracy(ts(traa[-1]), head(test.ts, n=100))
```

## UJI EFEK ARCH GARCH
```{r}
library(FinTS)
```

```{r}
cobaefek <- residuals(model1b)
arch_test <- ArchTest(cobaefek, lags = 15)
print(arch_test)
```
GARCH

## TAMBAH EFEK ARCH GARCH
```{r}
library(rugarch)
```

```{r}
spec1 <- ugarchspec(variance.model = list(model = "sGARCH", garchOrder = c(1,1)),
                   mean.model = list(armaOrder = c(2,2), include.mean = TRUE),
                   distribution.model = "norm")  

garch_fit1 <- ugarchfit(spec = spec1, data = trtdiff)
garch_fit1 
```

ARCH(1) : tidk signifikan
ARCH(2) : 2.5658, signifikan semua
GARCH(1,1) : 2.4897, signifikan semua
GARCH(1,2) : tidk signifikan
GARCH(2,1) : tidk signifikan
GARCH(2,2) : tidk signifikan
GARCH(3,1) -(3,5) : tidk signifikan
GARCH(3,2) : tidk signifikan
GARCH(3,3) : 
```{r}
# Ambil residual standar dari model GARCH
garch_residuals1 <- residuals(garch_fit1, standardize = TRUE)

# Melakukan uji ARCH pada residual standar
ArchTest(garch_residuals1)
```
```{r}
forc <- ugarchforecast(garch_fit1, data = trtdiff, n.ahead = 100)
ramalan <- forc@forecast$seriesFor[,1]
pt_1 <- train.ts #nilai akhir data latih
hasil.forc <- diffinv(ramalan, differences = 1) + pt_1
hasil.forcc <- InvBoxCox(hasil.forc, 0.14)
perbandingan <- data.frame("Aktual"= train.ts,
                           "Ramalan" = hasil.forcc)
perbandingan
```




```{r}
forc <- ugarchforecast(garch_fit1, data = trtdiff, n.ahead = 100, n.roll = 0)
forc
```
```{r}
plot(forc, which= 1)
```

```{r}
forc <- ugarchforecast(garch_fit1, data = trtdiff, n.ahead = 100)
ramalan <- forc@forecast$seriesFor[,1]
pt_1 <- train.ts #nilai akhir data latih
hasil.forc <- diffinv(ramalan, differences = 1) + pt_1
hasil.forcc <- InvBoxCox(hasil.forc, 0.14)
perbandingan <- data.frame("Aktual"= train.ts,
                           "Ramalan" = hasil.forcc)
perbandingan
```


```{r}
forecast_and_mape <- function(horizon, model, data_test, last_value, lambda) {
  forc <- ugarchforecast(model, data = trtdiff, n.ahead = horizon)
  ramalan <- forc@forecast$seriesFor[,1]

  # Transformasi balik differencing (dengan nilai awal dari trt.ts)
  hasil_diffinv <- diffinv(ramalan, differences = 1, xi = tail(last_value, 1))[-1]

  # Transformasi balik Box-Cox
  hasil_invboxcox <- InvBoxCox(hasil_diffinv, lambda)

  # Ambil aktual dari test data sebanyak 'horizon'
  actual <- head(data_test, horizon)

  # Hitung MAPE
  mape_val <- mape(actual, hasil_invboxcox[1:horizon])

  # Tampilkan hasil
  list(
    horizon = horizon,
    forecast = hasil_invboxcox[1:horizon],
    actual = actual,
    mape = mape_val
  )
}


```


```{r}
# Nilai terakhir dari trt.ts (bukan train.ts)
last_value <- trt.ts
lambda <- 0.14

hasil_7 <- forecast_and_mape(7, garch_fit1, test.ts, last_value, lambda)
hasil_30 <- forecast_and_mape(30, garch_fit1, test.ts, last_value, lambda)
hasil_60 <- forecast_and_mape(60, garch_fit1, test.ts, last_value, lambda)
hasil_100 <- forecast_and_mape(100, garch_fit1, test.ts, last_value, lambda)

# Tampilkan MAPEs
cat("MAPE 7 periode:", hasil_7$mape, "\n")
cat("MAPE 30 periode:", hasil_30$mape, "\n")
cat("MAPE 60 periode:", hasil_60$mape, "\n")
cat("MAPE 100 periode:", hasil_100$mape, "\n")

```
```{r}
hasil_7
```

```{r}
pred1 = c(14.96, 15.47, 15.68, 15.83, 15.97, 16.1, 16.23)
pred2 = c(15.79,16.84,17.56,18.04,18.39,18.67,18.92)
aktual7 = c(12,  14,  18,  17,  14,  19,  18)
errorgarch = abs(aktual7 - pred1)
errorgarch
errorlstm = abs(aktual7- pred2)
errorlstm
```
```{r}
pred2 = c(15.79,16.84,17.56,18.04,18.39,18.67,18.92)

```

```{r}
t.test(errorlstm, errorgarch, paired = TRUE)
```