R/01_seasonal_fit.R

In weatherMRJD: Weather Analysis and Markov Regime Switching Jump Diffusion Models

#' Fit Fourier Seasonal Baseline Trend
#'
#' Fits a truncated Fourier series to daily atmospheric temperature observations.
#'
#' @param temp Numeric vector of observed daily average temperatures.
#' @param day_of_year Numeric vector of day-of-year indices (1 to 365/366).
#' @return A list containing the fitted OLS model, fitted values, and isolated residuals.
#' @examples
#' doy <- rep(1:365, times = 2)
#' temp <- 15 + 10 * sin(2 * pi * doy / 365) + stats::rnorm(730, sd = 2)
#' fit <- fit_seasonal_trend(temp, doy)
#' print(fit$rmse)
#' @export
fit_seasonal_trend <- function(temp, day_of_year) {
  sin_term <- sin(2 * pi * day_of_year / 365.0)
  cos_term <- cos(2 * pi * day_of_year / 365.0)
  time_idx <- seq_along(temp)
  
  ols_fit <- stats::lm(temp ~ time_idx + sin_term + cos_term)
  fitted_vals <- stats::fitted(ols_fit)
  residuals_vals <- stats::residuals(ols_fit)
  
  list(
    model = ols_fit,
    fitted = fitted_vals,
    residuals = residuals_vals,
    rmse = sqrt(mean(residuals_vals^2)),
    mae = mean(abs(residuals_vals))
  )
}