format files

Author: andre_ramos

Project.toml

@@ -8,7 +8,6 @@ Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
 GLMNet = "8d5ece8b-de18-5317-b113-243142960cc6"
 HTTP = "cd3eb016-35fb-5094-929b-558a96fad6f3"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
-PyCall = "438e738f-606a-5dbb-bf0a-cddfbfd45ab0"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"

README.md

@@ -240,6 +240,14 @@ To reproduce the paper results, run the following experiments:
 
 Evaluates SSL and SS (StateSpaceModels) benchmark models on M4 competition dataset across all granularities (Monthly, Quarterly, Daily, Hourly, Weekly, Yearly).
 
+**Before running:** Add PyCall to your Julia environment:
+```julia
+using Pkg
+Pkg.add("PyCall")
+```
+
+The script also requires Python packages (`statsmodels`, `numpy`) for the SS benchmark evaluation.
+
 ```shell
 julia paper_tests/m4_test/m4_test.jl
 ```
@@ -254,8 +262,6 @@ Results are saved in:
 - `paper_tests/m4_test/results_SSL/` - SSL model results by granularity
 - `paper_tests/m4_test/metrics_results/` - Summary metrics and benchmark results
 
-Note: Requires PyCall and Python packages (`statsmodels`, `numpy`) for the SS benchmark evaluation.
-
 ### Simulation Parameter Study
 
 Compares SSL vs Kalman filter on simulated data.

paper_tests/m4_test/evaluate_model.jl

@@ -25,7 +25,7 @@ function evaluate_SSL(
     else
         ξ_threshold = 0
     end
-    
+
     model = StateSpaceLearning.StructuralModel(
         normalized_y;
         level="stochastic",
@@ -38,7 +38,7 @@ function evaluate_SSL(
         ξ_threshold=ξ_threshold,
         ζ_threshold=param["ζ_threshold"],
         ω_threshold=param["ω_threshold"],
-        ϕ_threshold=param["ϕ_threshold"]
+        ϕ_threshold=param["ϕ_threshold"],
     )
 
     if selection == "split"
@@ -78,26 +78,25 @@ function evaluate_SSL(
     return initialization_df, results_df
 end
 
-
 function evaluate_SS(input::Dict, m::Int, H::Int, frequency::Int)
     """
     Evaluate statsmodels UnobservedComponents model using PyCall
     Similar to evaluate_ss in m4_test.py
     Requires PyCall to be imported in the calling scope
     """
-    
+
     y_train = input["train"]
     y_test = input["test"]
 
     y_train_py = PyObject(y_train)
     y_test_py = PyObject(y_test)
     H_py = PyObject(H)
     frequency_py = PyObject(frequency)
-    
+
     py"""
     import statsmodels.api as sm
     import numpy as np
-    
+
     def evaluate_ss_py(y_train, y_test, H, frequency):
         # Prepare model components
         if frequency > 1:
@@ -146,21 +145,17 @@ function evaluate_SS(input::Dict, m::Int, H::Int, frequency::Int)
         
         return forecast_values, simulation
     """
-    
-    prediction, scenarios = py"evaluate_ss_py"(
-        y_train_py, y_test_py, 
-        H_py, frequency_py
-    )
-    
+
+    prediction, scenarios = py"evaluate_ss_py"(y_train_py, y_test_py, H_py, frequency_py)
+
     # Convert back to Julia arrays
     prediction = Vector{Float64}(prediction)
     scenarios = Matrix{Float64}(scenarios)
-    
+
     # Calculate metrics
-    mase = MASE(y_train, y_test, prediction; m = m)
+    mase = MASE(y_train, y_test, prediction; m=m)
     smape = sMAPE(y_test, prediction)
     crps = CRPS(scenarios, y_test)
-    
+
     return DataFrame([[mase], [smape], [crps]], [:MASE, :sMAPE, :CRPS])
-    
-end
+end

paper_tests/m4_test/m4_test.jl

@@ -34,21 +34,112 @@ include("metrics.jl")
 include("evaluate_model.jl")
 include("prepare_data.jl")
 
-complete_dict_vec = Dict("M" => build_train_test_dict(df_train, df_test;name="M"),
-                         "D" => build_train_test_dict(df_train_daily, df_test_daily;name="D"),
-                         "Q" => build_train_test_dict(df_train_quarterly, df_test_quarterly;name="Q"),
-                         "H" => build_train_test_dict(df_train_hourly, df_test_hourly;name="H"),
-                         "W" => build_train_test_dict(df_train_weekly, df_test_weekly;name="W"),
-                         "Y" => build_train_test_dict(df_train_yearly, df_test_yearly;name="Y"),
+complete_dict_vec = Dict(
+    "M" => build_train_test_dict(df_train, df_test; name="M"),
+    "D" => build_train_test_dict(df_train_daily, df_test_daily; name="D"),
+    "Q" => build_train_test_dict(df_train_quarterly, df_test_quarterly; name="Q"),
+    "H" => build_train_test_dict(df_train_hourly, df_test_hourly; name="H"),
+    "W" => build_train_test_dict(df_train_weekly, df_test_weekly; name="W"),
+    "Y" => build_train_test_dict(df_train_yearly, df_test_yearly; name="Y"),
 )
 
 parameters_dict = Dict(
-    "M" => Dict("freq_seasonal" => 12, "cycle_period" => 0, "m" => 12, "ξ_threshold" => 1, "ζ_threshold" => 12, "ω_threshold" => 12, "ϕ_threshold" => 0, "seasonal" => "stochastic", "cycle" => "none", "sample_size" => 60, "name" => "MONTH", "H" => 18, "NAIVE_sMAPE" => 14.427, "NAIVE_MASE" => 1.063),
-    "Q" => Dict("freq_seasonal" => 4, "cycle_period" => 0, "m" => 4, "ξ_threshold" => 1, "ζ_threshold" => 4, "ω_threshold" => 4, "ϕ_threshold" => 0, "seasonal" => "stochastic", "cycle" => "none", "sample_size" => "all", "name" => "QUARTERLY", "H" => 8, "NAIVE_sMAPE" => 11.012, "NAIVE_MASE" => 1.371),
-    "D" => Dict("freq_seasonal" => 1, "cycle_period" => 0, "m" => 1, "ξ_threshold" => 1, "ζ_threshold" => 7, "ω_threshold" => 1, "ϕ_threshold" => 0, "seasonal" => "none", "cycle" => "none", "sample_size" => 90, "name" => "DAILY", "H" => 14, "NAIVE_sMAPE" => 3.045, "NAIVE_MASE" => 3.278),
-    "W" => Dict("freq_seasonal" => 1, "cycle_period" => 0, "m" => 1, "ξ_threshold" => 1, "ζ_threshold" => 4, "ω_threshold" => 1, "ϕ_threshold" => 0, "seasonal" => "none", "cycle" => "none", "sample_size" => 104, "name" => "WEEKLY", "H" => 13, "NAIVE_sMAPE" => 9.161, "NAIVE_MASE" => 2.777),
-    "Y" => Dict("freq_seasonal" => 1, "cycle_period" => 0, "m" => 1, "ξ_threshold" => 1, "ζ_threshold" => 2, "ω_threshold" => 1, "ϕ_threshold" => 0, "seasonal" => "none", "cycle" => "none", "sample_size" => "all", "name" => "YEARLY", "H" => 6, "NAIVE_sMAPE" => 16.342, "NAIVE_MASE" => 3.974),
-    "H" => Dict("freq_seasonal" => 168, "cycle_period" => [24], "m" => 24, "ξ_threshold" => 1, "ζ_threshold" => 168, "ω_threshold" => 168, "ϕ_threshold" => 12, "seasonal" => "stochastic", "cycle" => "stochastic", "sample_size" => 720, "name" => "HOURLY", "H" => 48, "NAIVE_sMAPE" => 18.383, "NAIVE_MASE" => 2.395),
+    "M" => Dict(
+        "freq_seasonal" => 12,
+        "cycle_period" => 0,
+        "m" => 12,
+        "ξ_threshold" => 1,
+        "ζ_threshold" => 12,
+        "ω_threshold" => 12,
+        "ϕ_threshold" => 0,
+        "seasonal" => "stochastic",
+        "cycle" => "none",
+        "sample_size" => 60,
+        "name" => "MONTH",
+        "H" => 18,
+        "NAIVE_sMAPE" => 14.427,
+        "NAIVE_MASE" => 1.063,
+    ),
+    "Q" => Dict(
+        "freq_seasonal" => 4,
+        "cycle_period" => 0,
+        "m" => 4,
+        "ξ_threshold" => 1,
+        "ζ_threshold" => 4,
+        "ω_threshold" => 4,
+        "ϕ_threshold" => 0,
+        "seasonal" => "stochastic",
+        "cycle" => "none",
+        "sample_size" => "all",
+        "name" => "QUARTERLY",
+        "H" => 8,
+        "NAIVE_sMAPE" => 11.012,
+        "NAIVE_MASE" => 1.371,
+    ),
+    "D" => Dict(
+        "freq_seasonal" => 1,
+        "cycle_period" => 0,
+        "m" => 1,
+        "ξ_threshold" => 1,
+        "ζ_threshold" => 7,
+        "ω_threshold" => 1,
+        "ϕ_threshold" => 0,
+        "seasonal" => "none",
+        "cycle" => "none",
+        "sample_size" => 90,
+        "name" => "DAILY",
+        "H" => 14,
+        "NAIVE_sMAPE" => 3.045,
+        "NAIVE_MASE" => 3.278,
+    ),
+    "W" => Dict(
+        "freq_seasonal" => 1,
+        "cycle_period" => 0,
+        "m" => 1,
+        "ξ_threshold" => 1,
+        "ζ_threshold" => 4,
+        "ω_threshold" => 1,
+        "ϕ_threshold" => 0,
+        "seasonal" => "none",
+        "cycle" => "none",
+        "sample_size" => 104,
+        "name" => "WEEKLY",
+        "H" => 13,
+        "NAIVE_sMAPE" => 9.161,
+        "NAIVE_MASE" => 2.777,
+    ),
+    "Y" => Dict(
+        "freq_seasonal" => 1,
+        "cycle_period" => 0,
+        "m" => 1,
+        "ξ_threshold" => 1,
+        "ζ_threshold" => 2,
+        "ω_threshold" => 1,
+        "ϕ_threshold" => 0,
+        "seasonal" => "none",
+        "cycle" => "none",
+        "sample_size" => "all",
+        "name" => "YEARLY",
+        "H" => 6,
+        "NAIVE_sMAPE" => 16.342,
+        "NAIVE_MASE" => 3.974,
+    ),
+    "H" => Dict(
+        "freq_seasonal" => 168,
+        "cycle_period" => [24],
+        "m" => 24,
+        "ξ_threshold" => 1,
+        "ζ_threshold" => 168,
+        "ω_threshold" => 168,
+        "ϕ_threshold" => 12,
+        "seasonal" => "stochastic",
+        "cycle" => "stochastic",
+        "sample_size" => 720,
+        "name" => "HOURLY",
+        "H" => 48,
+        "NAIVE_sMAPE" => 18.383,
+        "NAIVE_MASE" => 2.395,
+    ),
 )
 
 # Function to append results to CSV file
@@ -70,7 +161,7 @@ function run_config(
     selection::String,
     information_criteria::String,
     α::AbstractFloat,
-    param::Dict
+    param::Dict,
 )
     results_df = DataFrame()
     initialization_df = DataFrame()
@@ -93,7 +184,7 @@ function run_config(
             α,
             selection,
             information_criteria,
-            param
+            param,
         )
 
         if i % clear_df_number == 0 || i == length(dict_vec)
@@ -177,7 +268,10 @@ function run_benchmark_model(dict_vec::Vector, param::Dict)
         "CRPS" => [crps],
         "Median CRPS" => [median_crps],
     )
-    CSV.write("paper_tests/m4_test/metrics_results/BENCHMARK_$(param["name"])_SUMMARY.csv", summary_df)
+    CSV.write(
+        "paper_tests/m4_test/metrics_results/BENCHMARK_$(param["name"])_SUMMARY.csv",
+        summary_df,
+    )
     return results_df
 end
 
@@ -189,21 +283,30 @@ function main()
             for selection in ["split", "fixed_alpha"]
                 if selection == "fixed_alpha"
                     information_criteria = "aic"
-                    alpha_set = [0.1, 0.3, 0.5, 0.7, 0.9, 1.0] 
-                else 
+                    alpha_set = [0.1, 0.3, 0.5, 0.7, 0.9, 1.0]
+                else
                     information_criteria = "aic"
                     alpha_set = [-1.0]
                 end
                 for α in alpha_set
                     @info "Running configuration: Param=$(gran), Outlier=$(outlier), IC=$(information_criteria), α=$(α)"
                     results_table = run_config(
-                        complete_dict_vec[gran], results_table, outlier, selection, information_criteria, α, parameters_dict[gran]
+                        complete_dict_vec[gran],
+                        results_table,
+                        outlier,
+                        selection,
+                        information_criteria,
+                        α,
+                        parameters_dict[gran],
                     )
                 end
             end
         end
         filename = parameters_dict[gran]["name"]
-        CSV.write("paper_tests/m4_test/metrics_results/SSL_$(filename)_METRICS_RESULTS.csv", results_table)
+        CSV.write(
+            "paper_tests/m4_test/metrics_results/SSL_$(filename)_METRICS_RESULTS.csv",
+            results_table,
+        )
         @info "Running benchmark model for $(parameters_dict[gran]["name"])"
         run_benchmark_model(complete_dict_vec[gran], parameters_dict[gran])
     end

paper_tests/simulation_param/evaluate_models.jl

@@ -19,30 +19,21 @@ function align_components!(
 end
 
 function get_SSL_results(
-    y_train::Vector{Fl},
-    s::Int,
-    μ_true, ν_true, γ_true,
-    inf_criteria::String
+    y_train::Vector{Fl}, s::Int, μ_true, ν_true, γ_true, inf_criteria::String
 ) where {Fl<:AbstractFloat}
-
     min_y = minimum(y_train)
     max_y = maximum(y_train)
     normy = (y_train .- min_y) ./ (max_y - min_y)
 
     model = StateSpaceLearning.StructuralModel(
-        normy;
-        freq_seasonal=s,
-        outlier=false,
-        ξ_threshold=0,
-        ζ_threshold=1,
-        ω_threshold=1
+        normy; freq_seasonal=s, outlier=false, ξ_threshold=0, ζ_threshold=1, ω_threshold=1
     )
     StateSpaceLearning.fit!(
         model;
         information_criteria=inf_criteria,
         ϵ=0.05,
         penalize_initial_states=true,
-        α = 0.1
+        α=0.1,
     )
 
     μ_hat = model.output.decomposition["trend"] .* (max_y - min_y)
@@ -52,4 +43,4 @@ function get_SSL_results(
     align_components!(μ_hat, ν_hat, γ_hat, μ_true, ν_true, γ_true)
 
     return μ_hat, ν_hat, γ_hat
-end
+end