df <- read.csv2('wisman3.csv', stringsAsFactors = FALSE, sep = '\t')
## Warning in read.table(file = file, header = header, sep = sep, quote = quote, :
## incomplete final line found by readTableHeader on 'wisman3.csv'
head(df)
library(tidyr)
## Warning: package 'tidyr' was built under R version 4.5.2
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
df_long <- df %>%
  pivot_longer(cols = -X, names_to = "Month", values_to = "Value")

time_values <- df_long$Value
ts_data <- ts(time_values, start = c(2024, 1), frequency = 12)
plot(ts_data, type = "l", main = "Time Series Data", ylab = "Value", xlab = "Time")

# Decompose
decomp <- decompose(ts_data, type="multiplicative")

plot(decomp)

pure_noise_manual <- ts_data / (decomp$trend * decomp$seasonal)
plot(pure_noise_manual, type = "l", main = "Pure Noise (Manual)", ylab = "Value", xlab = "Time")

seasonal_add <- decomp$seasonal
trend_add <- decomp$trend
residual_add <- decomp$random

deseasonalized_add <- ts_data - seasonal_add #- trend_add

plot(deseasonalized_add, type = "l", main = "Deseasonalized Data (Additive)", ylab = "Value", xlab = "Time")

h <- 12  # forecast horizon

last_value <- tail(deseasonalized_add, 1)

naive_forecast <- rep(last_value, h)

# Reapply seasonal component
season_future <- tail(seasonal_add, 12)

final_forecast_add <- naive_forecast + season_future
plot(final_forecast_add, type = "l", main = "Final Forecast (Additive)", ylab = "Value", xlab = "Time")

moving_average <- function(x, k) {
  stats::filter(x, rep(1/k, k), sides = 1)
}
k <- 3

ma_values <- moving_average(deseasonalized_add, k)

# Use last MA value for forecast
last_ma <- tail(na.omit(ma_values), 1)

ma_forecast <- rep(last_ma, h)

# Add seasonality back
final_ma_forecast_add <- ma_forecast + season_future
ts.plot(ts_data, col = "black", xlim = c(2024, 2028))

# Add forecast
lines(ts(final_forecast_add, start = c(2026,1), frequency = 12),
      col = "blue")

lines(ts(final_ma_forecast_add, start = c(2026,1), frequency = 12),
      col = "red")

legend("topleft",
       legend = c("Historis", "Naive", "Moving Avg"),
       col = c("black", "blue", "red"),
       lty = 1)

seasonal_naive <- rep(tail(ts_data, 6), length.out = h) + seasonal_add[0:6]

forecast <- numeric(h)
temp <- deseasonalized_add

for(i in 1:h){
  next_val <- mean(tail(temp, k))
  forecast[i] <- next_val
  temp <- c(temp, next_val)
}

final_forecast_add2 <- forecast + season_future
ts.plot(ts_data, col = "black", xlim = c(2024, 2027))
lines(ts(final_forecast_add2, start = c(2026,1), frequency = 12),
      col = "red")
lines(ts(seasonal_naive, start = c(2026,1), frequency = 12),
      col = "blue")
legend("topleft",
       legend = c("Historis", "Naive", "Moving Avg"),
       col = c("black", "blue", "red"),
       lty = 1)