future_forecast.R

suppressPackageStartupMessages({
  library(dplyr)
  library(readr)
})

table_dir <- "outputs/tables"
output_path <- file.path(table_dir, "future_forecast_summary.txt")

read_latest_actual_date <- function() {
  bus_path <- "data/processed/silver/ridership_bus_daily.csv"
  rail_path <- "data/processed/silver/ridership_rail_station_daily.csv"

  actual_dates <- c(
    read_csv(bus_path, show_col_types = FALSE) %>% summarise(latest_date = max(as.Date(date))) %>% pull(latest_date),
    read_csv(rail_path, show_col_types = FALSE) %>% summarise(latest_date = max(as.Date(date))) %>% pull(latest_date)
  )

  min(actual_dates, na.rm = TRUE)
}

read_forecast <- function(path, mode_name) {
  if (!file.exists(path)) {
    stop("Missing forecast file: ", path, call. = FALSE)
  }

  read_csv(path, show_col_types = FALSE) %>%
    mutate(
      date = as.Date(date),
      mode = mode_name
    ) %>%
    arrange(date) %>%
    slice_head(n = 30)
}

summarise_forecast <- function(data) {
  data %>%
    summarise(
      start_date = min(date),
      end_date = max(date),
      days = n(),
      total_ridership = sum(prediction, na.rm = TRUE),
      average_daily_ridership = mean(prediction, na.rm = TRUE),
      low_day = date[which.min(prediction)],
      low_day_ridership = min(prediction, na.rm = TRUE),
      high_day = date[which.max(prediction)],
      high_day_ridership = max(prediction, na.rm = TRUE),
      average_lower_80 = mean(lower_80, na.rm = TRUE),
      average_upper_80 = mean(upper_80, na.rm = TRUE),
      .groups = "drop"
    )
}

format_number <- function(x) {
  format(round(x), big.mark = ",", scientific = FALSE)
}

format_daily_lines <- function(data) {
  data %>%
    mutate(
      line = paste0(
        date,
        ": ",
        format_number(prediction),
        " [",
        format_number(lower_80),
        " to ",
        format_number(upper_80),
        "]"
      )
    ) %>%
    pull(line)
}

format_mode_block <- function(mode_name, summary_row, forecast_data) {
  c(
    paste0(mode_name, " forecast"),
    paste0("- Date range: ", summary_row$start_date, " to ", summary_row$end_date),
    paste0("- Forecast days: ", summary_row$days),
    paste0("- Total predicted ridership: ", format_number(summary_row$total_ridership)),
    paste0("- Average daily ridership: ", format_number(summary_row$average_daily_ridership)),
    paste0("- Lowest predicted day: ", summary_row$low_day, " (", format_number(summary_row$low_day_ridership), ")"),
    paste0("- Highest predicted day: ", summary_row$high_day, " (", format_number(summary_row$high_day_ridership), ")"),
    paste0(
      "- Average 80% range: ",
      format_number(summary_row$average_lower_80),
      " to ",
      format_number(summary_row$average_upper_80)
    ),
    "",
    format_daily_lines(forecast_data)
  )
}

last_true_data_date <- read_latest_actual_date()
bus_forecast <- read_forecast(file.path(table_dir, "bus_forecast_output.csv"), "Bus")
rail_forecast <- read_forecast(file.path(table_dir, "rail_system_forecast_output.csv"), "Rail")

bus_summary <- summarise_forecast(bus_forecast)
rail_summary <- summarise_forecast(rail_forecast)
combined_summary <- bind_rows(bus_forecast, rail_forecast) %>%
  summarise_forecast()

output_lines <- c(
  "WMATA 30-Day Future Forecast Summary",
  paste0("Generated: ", Sys.Date()),
  paste0("Last True Data Date: ", last_true_data_date),
  "",
  format_mode_block("Bus", bus_summary, bus_forecast),
  "",
  format_mode_block("Rail", rail_summary, rail_forecast),
  "",
  "Combined Bus + Rail forecast",
  paste0("- Date range: ", combined_summary$start_date, " to ", combined_summary$end_date),
  paste0("- Forecast days by mode: ", combined_summary$days / 2),
  paste0("- Total predicted ridership: ", format_number(combined_summary$total_ridership)),
  paste0("- Average daily ridership across both modes: ", format_number(combined_summary$total_ridership / (combined_summary$days / 2)))
)

dir.create(table_dir, recursive = TRUE, showWarnings = FALSE)
writeLines(output_lines, output_path)

message("Wrote ", output_path)