first tests for Neighbors2 after debug; TODO: some missing forecasts
[talweg.git] / pkg / R / F_Neighbors.R
index 600c5c8..5b2c899 100644 (file)
@@ -33,21 +33,25 @@ NeighborsForecaster = R6::R6Class("NeighborsForecaster",
                        # Indices of similar days for cross-validation; TODO: 45 = magic number
                        sdays = getSimilarDaysIndices(today, limit=45, same_season=FALSE)
 
+                       cv_days = intersect(fdays,sdays)
+                       # Limit to 20 most recent matching days (TODO: 20 == magic number)
+                       cv_days = sort(cv_days,decreasing=TRUE)[1:min(20,length(cv_days))]
+
                        # Function to optimize h : h |--> sum of prediction errors on last 45 "similar" days
                        errorOnLastNdays = function(h, kernel, simtype)
                        {
                                error = 0
                                nb_jours = 0
-                               for (i in intersect(fdays,sdays))
+                               for (i in seq_along(cv_days))
                                {
                                        # mix_strategy is never used here (simtype != "mix"), therefore left blank
                                        prediction = private$.predictShapeAux(data,
-                                               fdays, i, horizon, h, kernel, simtype, FALSE)
+                                               fdays, cv_days[i], horizon, h, kernel, simtype, FALSE)
                                        if (!is.na(prediction[1]))
                                        {
                                                nb_jours = nb_jours + 1
                                                error = error +
-                                                       mean((data$getCenteredSerie(i+1)[1:horizon] - prediction)^2)
+                                                       mean((data$getCenteredSerie(cv_days[i]+1)[1:horizon] - prediction)^2)
                                        }
                                }
                                return (error / nb_jours)
@@ -96,14 +100,11 @@ NeighborsForecaster = R6::R6Class("NeighborsForecaster",
                                h_endo = ifelse(simtype=="mix", h[1], h)
 
                                # Distances from last observed day to days in the past
-                               distances2 = rep(NA, length(fdays))
-                               for (i in seq_along(fdays))
-                               {
-                                       delta = data$getCenteredSerie(today) - data$getCenteredSerie(fdays[i])
-                                       # Require at least half of non-NA common values to compute the distance
-                                       if ( !any( is.na(delta) ) )
-                                               distances2[i] = mean(delta^2)
-                               Centered}
+                               serieToday = data$getSerie(today)
+                               distances2 = sapply(fdays, function(i) {
+                                       delta = serieToday - data$getSerie(i)
+                                       mean(delta^2)
+                               })
 
                                sd_dist = sd(distances2)
                                if (sd_dist < .Machine$double.eps)
@@ -136,12 +137,10 @@ NeighborsForecaster = R6::R6Class("NeighborsForecaster",
                                sigma_inv = solve(sigma) #TODO: use pseudo-inverse if needed?
 
                                # Distances from last observed day to days in the past
-                               distances2 = rep(NA, nrow(M)-1)
-                               for (i in 2:nrow(M))
-                               {
-                                       delta = M[1,] - M[i,]
-                                       distances2[i-1] = delta %*% sigma_inv %*% delta
-                               }
+                               distances2 = sapply(seq_along(fdays), function(i) {
+                                       delta = M[1,] - M[i+1,]
+                                       delta %*% sigma_inv %*% delta
+                               })
 
                                sd_dist = sd(distances2)
                                if (sd_dist < .Machine$double.eps)