Untitled

#ARIMA Model

library(tseries)
library(dplyr)
library(forecast)

# Assuming 'dda' is your data frame and 'p' is the price column
# Convert price to a time series object
price_ts <- ts(dda$p, start = c(2017, 2), frequency = 4)

# First difference and seasonal difference
dy <- dda$p - lag(dda$p, n=1)
qy <- dda$p - lag(dda$p, n=4)
dy <- na.omit(dy)  # Remove NAs generated by lagging
qy <- na.omit(qy)  # Remove initial NAs from seasonal differencing

# Visualize the differenced data
windows()
par(mfrow=c(2,2))
plot(dy, main="(1-B)Y", col="blue", type="l", lwd=2, xlab="Index", ylab="Change in Price")
grid(nx = NULL, ny = NULL, col = "gray", lty = "dotted")
plot(dy, main="(1-B)Y", col="blue", type="l", lwd=2, xlab="Index (Quarterly)", ylab="Change in Price")
grid(nx = NULL, ny = NULL, col = "gray", lty = "dotted")
acf(dy, col="red", lag.max=10)
acf(qy, col="red", lag.max=10)

# Double differencing for combined seasonal and first differences
dqy <- dy - lag(dy, n=4)
dqy <- na.omit(dqy)

# Visualize the doubly differenced data
windows()
layout(matrix(c(1,1,2,3), nrow=2, ncol=2, byrow=T))
plot(dqy, type="l", col="blue", lwd=2, main="(1-B)(1-B^4)Y")
acf(dqy, col="red", lag.max=10)
pacf(dqy, col="red", lag.max=10)

# Split the data into training and testing sets
train_end <- floor(0.85 * length(price_ts))
train_data <- window(price_ts, end = c(2017 + (train_end - 1) %/% 4, (train_end - 1) %% 4 + 1))
test_data <- window(price_ts, start = c(2017 + train_end %/% 4, train_end %% 4 + 1))

# Fit ARIMA model on training data
m0 <- arima(train_data, order = c(1,0,0), seasonal = list(order = c(0,0,0), period = 4))
summary(m0)


Box.test(m0$residuals, lag=20, type="Ljung-Box")
windows()
acf(m0$residuals,lag.max=10,col="red")

acf(residuals(m0), main="ACF of Residuals")
Box.test(residuals(m0), type="Ljung-Box")


# Forecast future values using the model and evaluate on the test data
future_forecast <- forecast(m0, h = length(test_data))

# Determine the maximum value for y-axis limit
windows()
plot(future_forecast, main = "ARIMA Forecast", ylim = c(0, 200), cex.axis = 2, cex.main = 2, cex.lab = 2)

lines(test_data, col = "red", type = "o")  # Overlay test data for comparison


#Second Arima model
# Fit ARIMA model on training data
m1 <- arima(train_data, order = c(2,0,0), seasonal = list(order = c(0,0,0), period = 4))
summary(m1)

Box.test(m1$residuals, lag=20, type="Ljung-Box")
windows()
acf(m1$residuals,lag.max=10,col="red")

acf(residuals(m1), main="ACF of Residuals")
Box.test(residuals(m1), type="Ljung-Box")


# Forecast future values using the model and evaluate on the test data
future_forecast2 <- forecast(m1, h = length(test_data))

# Determine the maximum value for y-axis limit
windows()
plot(future_forecast2, main = "ARIMA Forecast", ylim = c(0, 200), cex.axis = 2, cex.main = 2, cex.lab = 2)

lines(test_data, col = "red", type = "o")  # Overlay test data for comparison