diff --git a/src/TimeSeriesInterface.jl b/src/TimeSeriesInterface.jl
index 9f9374f..f0a17d2 100644
--- a/src/TimeSeriesInterface.jl
+++ b/src/TimeSeriesInterface.jl
@@ -7,8 +7,8 @@ using RecipesBase
 using Statistics
 
 include("timeseries.jl")
-include("forecast.jl")
 include("models.jl")
+include("forecast.jl")
 include("FileFormats/FileFormats.jl")
 
 end
\ No newline at end of file
diff --git a/src/forecast.jl b/src/forecast.jl
index b45c91d..a4f7d7b 100644
--- a/src/forecast.jl
+++ b/src/forecast.jl
@@ -3,6 +3,7 @@ export PointForecast, ScenariosForecast, QuantilesForecast
 export PointForecastMetrics, ScenariosForecastMetrics
 
 export forecast_metrics
+export cross_validation
 
 abstract type Forecast end
 abstract type ProbabilisticForecast <: Forecast end
@@ -12,7 +13,7 @@ abstract type ProbabilisticForecast <: Forecast end
 
 Define the results of point forecasts.    
 """
-mutable struct PointForecast{T<:Real} <: Forecast
+mutable struct PointForecast{T <: Real} <: Forecast
     name::String
     timestamps::Vector{DateTime}
     forecast::Vector{T}
@@ -21,7 +22,7 @@ mutable struct PointForecast{T<:Real} <: Forecast
         name::String,
         timestamps::Vector{DateTime},
         forecast::Vector{T},
-    ) where {T<:Real}
+    ) where {T <: Real}
 
         if length(timestamps) != length(forecast)
             throw(DimensionMismatch("timestamps and forecast do not have the same length."))
@@ -44,13 +45,13 @@ function PointForecast(
     name::String,
     timestamps::Vector{Date},
     forecast::Vector{T},
-) where {T<:Real}
+) where {T <: Real}
 
     return PointForecast(name, DateTime.(timestamps), forecast)
 end
 
 ## Evaluation Metrics for Point Forecast
-struct PointForecastMetrics{T<:Real}
+struct PointForecastMetrics{T <: Real}
     errors::Vector{T}
     absolute_percentage_errors::Vector{T}
 end
@@ -82,12 +83,44 @@ error(real::Vector{T}, forecast::Vector{T}) where {T} = real .- forecast
 absolute_percentage_error(real::Vector{T}, forecast::Vector{T}) where {T} =
     abs.(error(real, forecast) ./ real)
 
+"""
+    cross_validation(fit_input::FitInput{T}, fit_function::Function, predict_function::Function, metric_function::Function,
+                            min_history::Int, horizon::Int)
+
+Function that receives all avaiable data in fit_input, a fit_function and a prediction_funtion,
+a validation metric (eg. mape), a min_history as a starting point and a prediction horizon.
+Returns the result of the cross validation separeted by lead time.
+"""
+function cross_validation(fit_input::FitInput{T},
+                            fit_function::Function, 
+                            predict_function::Function,
+                            metric_function::Function,
+                            min_history::Int, 
+                            horizon::Int) where T
+    n = length(fit_input.dependent[1].timestamps)
+    metric_matrix = Matrix{Float64}(undef, length(min_history:(n - horizon)), horizon)
+    for i = min_history:(n - horizon)
+        training_dependent = map(x -> TimeSeries(x.name, x.timestamps[1:i], x.vals[1:i]), fit_input.dependent)
+        training_exogenous = map(x -> TimeSeries(x.name, x.timestamps[1:i], x.vals[1:i]), fit_input.exogenous)
+        training_fit_input = FitInput(fit_input.parameters, training_dependent, training_exogenous)
+        validation_fit_result = fit_function(training_fit_input)
+        validation_timestamps_forecast = fit_input.dependent[1].timestamps[(i + 1):(i + horizon)]
+        validation_exogenous_forecast = map(x -> TimeSeries(x.name, x.timestamps[(i + 1):(i + horizon)], x.vals[(i + 1):(i + horizon)]), fit_input.exogenous)
+        validation_simulate_input = SimulateInput(training_fit_input, validation_timestamps_forecast,
+                                                    validation_exogenous_forecast, validation_fit_result)
+        point_forecast = predict_function(validation_simulate_input)
+        observed = map(x -> TimeSeries(x.name, x.timestamps[(i + 1):(i + horizon)], x.vals[(i + 1):(i + horizon)]), fit_input.dependent)
+        metric_matrix[(i - min_history + 1), :] .= metric_function(observed[1].vals, point_forecast.forecast)
+    end
+    return mean(metric_matrix, dims=1)
+end
+
 """
     ScenariosForecast
 
 Define the probabilistic forecast results calculated from scenarios.   
 """
-mutable struct ScenariosForecast{T<:Real} <: ProbabilisticForecast
+mutable struct ScenariosForecast{T <: Real} <: ProbabilisticForecast
     name::String
     timestamps::Vector{DateTime}
     scenarios::Matrix{T}
@@ -100,7 +133,7 @@ mutable struct ScenariosForecast{T<:Real} <: ProbabilisticForecast
         scenarios::Matrix{T},
         quantiles_probabilities::Vector{T},
         quantiles::Matrix{T},
-    ) where {T<:Real}
+    ) where {T <: Real}
 
         if length(timestamps) != size(scenarios, 1)
             throw(DimensionMismatch("timestamps and scenarios do not have the same length."))
@@ -133,7 +166,7 @@ function ScenariosForecast(
     scenarios::Matrix{T},
     quantiles_probabilities::Vector{T},
     quantiles::Matrix{T},
-) where {T<:Real}
+) where {T <: Real}
 
     return ScenariosForecast(
         name,
@@ -179,7 +212,7 @@ function forecast_metrics(
 end
 
 function get_quantiles(quantile_probs::Vector{T}, scenarios::Matrix{T}) where {T}
-    quantiles = mapslices(x -> quantile(x, quantile_probs), scenarios; dims = 2)
+    quantiles = mapslices(x -> quantile(x, quantile_probs), scenarios; dims=2)
     return quantiles
 end
 
@@ -298,7 +331,7 @@ function mean_crps(forecast::Vector{ScenariosForecastMetrics})
     for (i, forec) in enumerate(forecast)
         m_crps[:, i] = forec.crps
     end
-    return vec(mean(m_crps, dims = 2))
+    return vec(mean(m_crps, dims=2))
 end
 
 
@@ -307,7 +340,7 @@ end
 
 Define the probabilistic forecast results calculated from distributions.      
 """
-mutable struct QuantilesForecast{T<:Real} <: ProbabilisticForecast
+mutable struct QuantilesForecast{T <: Real} <: ProbabilisticForecast
     name::String
     timestamps::Vector{DateTime}
     quantiles_probabilities::Vector{T}
@@ -318,7 +351,7 @@ mutable struct QuantilesForecast{T<:Real} <: ProbabilisticForecast
         timestamps::Vector{DateTime},
         quantiles_probabilities::Vector{T},
         quantiles::Matrix{T},
-    ) where {T<:Real}
+    ) where {T <: Real}
 
         if length(timestamps) != size(quantiles, 1)
             throw(DimensionMismatch("timestamps and quantiles do not have the same length."))
@@ -346,7 +379,7 @@ function QuantilesForecast(
     timestamps::Vector{Date},
     quantiles_probabilities::Vector{T},
     quantiles::Matrix{T},
-) where {T<:Real}
+) where {T <: Real}
 
     return QuantilesForecast(
         name,