R/long_term_future.R

In oRaklE: Multi-Horizon Electricity Demand Forecasting in High Resolution

#' Long-term trend predictions for future years
#'
#' This function extends the long-term trend predictions generated by \code{\link{long_term_lm}} until a specified future year.
#' The unknown macro-economic covariates are either obtained from the WEO Outlook with \code{\link{long_term_future_data}} or can be supplied manually.
#' The function also produces and saves visualizations of the actual and the predicted demand over the training, test, and future periods.
#' @param longterm_future_macro_data Dataframe. Generated by \code{\link{long_term_future_data}}
#' @param data_directory The path to the directory where the data will be saved and where the function will look for
#' the long-term models from \code{\link{long_term_lm}}. The default is set to a temporary directory.
#' @param model_list A list with the models from \code{\link{long_term_lm}}. Only needs to be specified if the models
#' are not in the data directory.
#' @param verbose A boolean value indicating if you want the generated plots to be shown (set to TRUE if yes).
#' @return A list with the extended initial dataframe with the future predictions for each of the 3 best long term models and one plot with the respective results
#' per model.
#' \describe{
#'   \item{longterm_future_predictions}{A dataframe with the input data and additional columns for test_set_steps and for best three models longterm_model_predictions1, longterm_model_predictions2 and longterm_model_predictions3.}
#'   \item{logterm_future_plots}{A list with the respective plots for each model.}
#' }
#' @export
#'
#' @seealso See also function \code{\link{mid_term_future}} and \code{\link{short_term_future}} for the other prediction models and \code{\link{long_term_future_data}} for the covariate download.
#'
#' @examples
#' example_longterm_future_predictions <- long_term_future(example_longterm_future_macro_data)
long_term_future <- function(longterm_future_macro_data, data_directory = tempdir(), model_list = NULL, verbose = FALSE) {
  if ("example" %in% colnames(longterm_future_macro_data)) {
    if (unique(longterm_future_macro_data$example) == TRUE) {
      variables <- c("GNI", "industrial_value_added", "rural_population")
      f <- stats::as.formula(paste("avg_hourly_demand", paste(variables, collapse = " + "),
        sep = " ~ "
      ))

      model <- stats::lm(f, data = longterm_future_macro_data[1:14])
      LT <- stats::predict(model, longterm_future_macro_data)
      expected_LT <- c(
        54429.04, 55113.21, 55527.12, 55770.34, 56881.23, 56116.27,
        55496.54, 54588.04, 54757.23, 54511.24, 54499.45, 54548.50,
        54067.50, 52758.44, 50637.02, 53057.05, 49703.69, 48458.31,
        48048.71, 46018.16, 44717.06, 43486.50, 41999.38
      )

      LT <- round(LT, 2)
      expected_LT <- round(expected_LT, 2)
      if (all.equal(unname(LT), expected_LT)) {
        return(oRaklE::example_longterm_predictions)
      } else {
        stop("The example in long_term_future() failed. Please contact the package maintainer at schwenzer@europa-uni.de")
      }
    }
  }

  if (grepl("Rtmp", data_directory)) {
    message(paste(
      "\nThis function will try to save the plots and find the long-term models in a folder called", unique(longterm_future_macro_data$country),
      "\nin the current data directory:", data_directory
    ))
    message("\nIf the long-term models are not found, a list with the models has to be passed in the *model_list* argument.")

    message("\nPlease choose an option:")
    message("\n1: Keep it as a tempdir")
    message(paste("2: Save data and look for the models in the current working directory (", getwd(), ")", sep = ""))
    message("3: Set the directory manually\n")

    choice <- readline(prompt = "Enter the option number (1, 2, or 3): ")


    if (choice == "1") {
      message("\nData will be saved in a temporary directory and cleaned up when R is shut down.\n")
    } else if (choice == "2") {
      data_directory <- getwd()
      message(paste0("\nData will be saved in the current working directory in ", data_directory, "/", unique(longterm_future_macro_data$country), "/data"))
      message("\nYou can specify the *data_directory* parameter in the following functions as ", data_directory)
    } else if (choice == "3") {
      new_dir <- readline(prompt = "Enter the full path of the directory where you want to save the data: ")
      data_directory <- new_dir
      if (!dir.exists(data_directory)) {
        stop("The specified data_directory does not exist: ", data_directory, "\nPlease run the function again.")
      }
      message("\nData will be saved in the specified directory: ", data_directory, "/", unique(longterm_future_macro_data$country), "/data")
    } else {
      message("Invalid input. Keeping the temporary directory.\nData will be cleaned up when R is shut down.\n")
    }
  } else {
    if (!dir.exists(data_directory)) {
      stop("The specified data_directory does not exist: ", data_directory, "\nPlease run the function again.")
    }
    message("\nData will be saved in the specified working directory in ", data_directory, "/", unique(longterm_future_macro_data$country), "/data")
  }

  new_row_start <- min(which(is.na(longterm_future_macro_data$avg_hourly_demand)))

  if (inherits(model_list, "list")) {
    message("Taking the models specified in model_list.")
    i <- 1
    for (m in model_list) {
      prediction_column <- which(colnames(longterm_future_macro_data) == paste0("longterm_model_predictions", i))
      longterm_future_macro_data[new_row_start:nrow(longterm_future_macro_data), prediction_column] <- stats::predict(m, newdata = longterm_future_macro_data)[new_row_start:nrow(longterm_future_macro_data)]
      i <- i + 1
    }
  } else {
    model_path <- paste0(data_directory, "/", unique(longterm_future_macro_data$country), "/models/longterm/best_lm_model1.Rdata")

    if (file.exists(model_path)) {
      for (i in 1:3) {
        model_path <- paste0(data_directory, "/", unique(longterm_future_macro_data$country), "/models/longterm/best_lm_model", i, ".Rdata")

        loaded_model <- load(model_path)
        best_lm_model <- get(loaded_model)
        prediction_column <- which(colnames(longterm_future_macro_data) == paste0("longterm_model_predictions", i))
        longterm_future_macro_data[new_row_start:nrow(longterm_future_macro_data), prediction_column] <- stats::predict(best_lm_model, newdata = longterm_future_macro_data)[new_row_start:nrow(longterm_future_macro_data)]
      }
    } else {
      stop("\nPlease either specify the base path where the country data is saved (e.g. the current working directory or supply a list with models for the *model_list* variable.")
    }
  }

  country <- unique(longterm_future_macro_data$country)
  if (!file.exists(paste0(data_directory, "/", country))) {
    dir.create(paste0(data_directory, "/", country))
  }
  if (!file.exists(paste0(data_directory, "/", country, "/data"))) {
    dir.create(paste0(data_directory, "/", country, "/data"))
  }
  if (!file.exists(paste0(data_directory, "/", country, "/plots"))) {
    dir.create(paste0(data_directory, "/", country, "/plots"))
  }

  utils::write.csv(longterm_future_macro_data, paste0(data_directory, "/", unique(longterm_future_macro_data$country), "/data/longterm_future_predictions.csv"), row.names = FALSE)

  intercept <- longterm_future_macro_data$year[(new_row_start - 1)] - unique(longterm_future_macro_data$test_set_steps)
  training_text_index <- min(longterm_future_macro_data$year, na.rm = T) + ((intercept - min(longterm_future_macro_data$year, na.rm = T)) / 2)
  test_set_end <- longterm_future_macro_data$year[(new_row_start - 1)]
  max_value <- max(c(
    max(longterm_future_macro_data$avg_hourly_demand, na.rm = T), max(longterm_future_macro_data$longterm_model_predictions1, na.rm = T),
    max(longterm_future_macro_data$longterm_model_predictions2, na.rm = T),
    max(longterm_future_macro_data$longterm_model_predictions3, na.rm = T)
  ))
  future_set <- length(longterm_future_macro_data$year[new_row_start:nrow(longterm_future_macro_data)])

  suppressWarnings(
    lt_plot <- ggplot(longterm_future_macro_data) +
      geom_line(aes(longterm_future_macro_data$year, longterm_future_macro_data$avg_hourly_demand, color = "actual"), lwd = 1) +
      geom_line(aes(longterm_future_macro_data$year, longterm_future_macro_data$longterm_model_predictions1, color = "Model1")) +
      geom_line(aes(longterm_future_macro_data$year, longterm_future_macro_data$longterm_model_predictions2, color = "Model2")) +
      geom_line(aes(longterm_future_macro_data$year, longterm_future_macro_data$longterm_model_predictions3, color = "Model3")) +
      geom_vline(xintercept = intercept, linetype = 2) +
      geom_vline(xintercept = test_set_end, linetype = 3) +
      ggthemes::theme_foundation(base_size = 14, base_family = "sans") +
      xlab("\nYear") +
      ylab("Avg Hourly Demand p. Year\n [MW]\n") +
      ggtitle(paste("Long Term Model Results -", unique(longterm_future_macro_data$country), "\n")) +
      theme(
        plot.title = element_text(
          face = "bold",
          size = rel(1.2), hjust = 0.5
        ),
        plot.subtitle = element_text(size = rel(1), hjust = 0.5),
        text = element_text(),
        panel.background = element_rect(colour = NA),
        plot.background = element_rect(colour = NA),
        panel.border = element_rect(colour = NA),
        axis.title = element_text(face = "bold", size = rel(1)),
        axis.title.y = element_text(angle = 90, vjust = 2),
        axis.title.x = element_text(vjust = -0.2),
        axis.text = element_text(),
        axis.line.x = element_line(colour = "black"),
        axis.line.y = element_line(colour = "black"),
        axis.ticks = element_line(),
        panel.grid.major = element_line(colour = "#f0f0f0"),
        panel.grid.minor = element_blank(),
        legend.key = element_rect(colour = NA),
        legend.position = "bottom",
        legend.direction = "horizontal",
        legend.key.size = unit(0.2, "cm"),
        plot.margin = unit(c(10, 5, 5, 5), "mm"),
        strip.background = element_rect(colour = "#f0f0f0", fill = "#f0f0f0"),
        strip.text = element_text(face = "bold")
      ) +
      theme(legend.title = element_blank()) +
      guides(color = guide_legend(override.aes = list(linewidth = 2))) +
      annotate("text", x = training_text_index, y = (max_value + max_value * 0.02), label = "Training", size = 4, hjust = 0.5, vjust = 0) +
      annotate("text", x = (intercept + unique(longterm_future_macro_data$test_set_steps) / 2), y = (max_value + max_value * 0.02), label = "Test", size = 4, hjust = 0.5, vjust = 0) +
      annotate("text", x = (longterm_future_macro_data$year[new_row_start] + future_set / 2), y = (max_value + max_value * 0.02), label = "Unknown", size = 4, hjust = 0.5, vjust = 0)
  )

  if (verbose == FALSE) {
    message("\nVerbose is set to FALSE. Set to TRUE if you want to see the generated plots automatically. The plots are saved in the output under *plots* and in the plots folder in ", data_directory)
  } else {
    suppressWarnings(
      print(lt_plot)
    )
  }

  suppressWarnings(
    lt_plot2 <- ggplot(longterm_future_macro_data) +
      geom_line(aes(longterm_future_macro_data$year, longterm_future_macro_data$avg_hourly_demand, color = "actual"), lwd = 1) +
      geom_line(aes(longterm_future_macro_data$year, longterm_future_macro_data$longterm_model_predictions1, color = "Model1")) +
      geom_line(aes(longterm_future_macro_data$year, longterm_future_macro_data$longterm_model_predictions2, color = "Model2")) +
      geom_line(aes(longterm_future_macro_data$year, longterm_future_macro_data$longterm_model_predictions3, color = "Model3")) +
      xlab("\nYear") +
      ylab("Avg Hourly Demand p. Year\n [MW]\n") +
      geom_vline(xintercept = intercept, linetype = 2) +
      geom_vline(xintercept = test_set_end, linetype = 3) +
      ggthemes::theme_foundation(base_size = 14, base_family = "sans") +
      ggtitle(paste("Long Term Model Results -", unique(longterm_future_macro_data$country), "\n")) +
      theme(
        plot.title = element_text(
          face = "bold",
          size = rel(1.2), hjust = 0.5
        ),
        plot.subtitle = element_text(size = rel(1), hjust = 0.5),
        text = element_text(),
        panel.background = element_rect(colour = NA),
        plot.background = element_rect(colour = NA),
        panel.border = element_rect(colour = NA),
        axis.title = element_text(face = "bold", size = rel(1)),
        axis.title.y = element_text(angle = 90, vjust = 2),
        axis.title.x = element_text(vjust = -0.2),
        axis.text = element_text(),
        axis.line.x = element_line(colour = "black"),
        axis.line.y = element_line(colour = "black"),
        axis.ticks = element_line(),
        panel.grid.major = element_line(colour = "#f0f0f0"),
        panel.grid.minor = element_blank(),
        legend.key = element_rect(colour = NA),
        legend.position = "bottom",
        legend.direction = "horizontal",
        legend.key.size = unit(0.2, "cm"),
        plot.margin = unit(c(10, 5, 5, 5), "mm"),
        strip.background = element_rect(colour = "#f0f0f0", fill = "#f0f0f0"),
        strip.text = element_text(face = "bold")
      ) +
      theme(legend.title = element_blank()) +
      theme(axis.title = element_text(size = 23)) +
      theme(legend.text = element_text(size = 23)) +
      theme(axis.text = element_text(size = 20)) +
      theme(plot.title = element_text(size = 26)) +
      guides(color = guide_legend(override.aes = list(linewidth = 2))) +
      theme(legend.title = element_blank()) +
      guides(color = guide_legend(override.aes = list(linewidth = 2))) +
      annotate("text", x = training_text_index, y = (max_value + max_value * 0.02), label = "Training", size = 4, hjust = 0.5, vjust = 0) +
      annotate("text", x = (intercept + unique(longterm_future_macro_data$test_set_steps) / 2), y = (max_value + max_value * 0.02), label = "Test", size = 4, hjust = 0.5, vjust = 0) +
      annotate("text", x = (longterm_future_macro_data$year[new_row_start] + future_set / 2), y = (max_value + max_value * 0.02), label = "Unknown", size = 4, hjust = 0.5, vjust = 0)
  )
  suppressWarnings(
    ggsave(filename = paste0(data_directory, "/", unique(longterm_future_macro_data$country), "/plots/Long_term_results_future.png"), plot = lt_plot2, width = 12, height = 8)
  )

  longterm_predictions_future <- longterm_future_macro_data
  all_plots <- list(
    longterm_future_plot = lt_plot
  )
  return(list("longterm_future_predictions" = longterm_predictions_future, "logterm_future_plots" = all_plots))
}