fix(models): improve the SES implementation

README.md

@@ -170,10 +170,10 @@ Uses a weighted average of all past observations where the weights decrease expo
 
 |Model | Point Forecast | Probabilistic Forecast | Insample fitted values | Probabilistic fitted values |Exogenous features|
 |:------|:-------------:|:----------------------:|:---------------------:|:----------------------------:|:----------------:|
-|[SimpleExponentialSmoothing](https://nixtlaverse.nixtla.io/statsforecast/src/core/models.html#simpleexponentialsmoothing)|✅|||||
-|[SimpleExponentialSmoothingOptimized](https://nixtlaverse.nixtla.io/statsforecast/src/core/models.html#simpleexponentialsmoothingoptimized)|✅|||||
-|[SeasonalExponentialSmoothing](https://nixtlaverse.nixtla.io/statsforecast/src/core/models.html#seasonalexponentialsmoothing)|✅|||||
-|[SeasonalExponentialSmoothingOptimized](https://nixtlaverse.nixtla.io/statsforecast/src/core/models.html#seasonalexponentialsmoothingoptimized)|✅|||||
+|[SimpleExponentialSmoothing](https://nixtlaverse.nixtla.io/statsforecast/src/core/models.html#simpleexponentialsmoothing)|✅||✅|||
+|[SimpleExponentialSmoothingOptimized](https://nixtlaverse.nixtla.io/statsforecast/src/core/models.html#simpleexponentialsmoothingoptimized)|✅||✅|||
+|[SeasonalExponentialSmoothing](https://nixtlaverse.nixtla.io/statsforecast/src/core/models.html#seasonalexponentialsmoothing)|✅||✅|||
+|[SeasonalExponentialSmoothingOptimized](https://nixtlaverse.nixtla.io/statsforecast/src/core/models.html#seasonalexponentialsmoothingoptimized)|✅||✅|||
 |[Holt](https://nixtlaverse.nixtla.io/statsforecast/src/core/models.html#holt)|✅|✅|✅|✅||
 |[HoltWinters](https://nixtlaverse.nixtla.io/statsforecast/src/core/models.html#holtwinters)|✅|✅|✅|✅||
 

nbs/src/core/models.ipynb

@@ -58,11 +58,11 @@
     "#| export\n",
     "import warnings\n",
     "from math import trunc\n",
-    "from typing import Any, Dict, List, Optional, Sequence, Tuple, Union\n",
+    "from typing import Any, Dict, List, Optional, Tuple, Union\n",
     "\n",
     "import numpy as np\n",
     "from numba import njit\n",
-    "from scipy.optimize import minimize\n",
+    "from scipy.optimize import minimize_scalar\n",
     "from scipy.special import inv_boxcox\n",
     "\n",
     "from statsforecast.arima import (\n",
@@ -3655,19 +3655,20 @@
     "    This function returns the one step ahead prediction\n",
     "    as well as the mean squared error of the fit.\n",
     "    \"\"\"\n",
-    "    smoothed = x[0]\n",
-    "    n = x.size\n",
-    "    mse = 0.\n",
-    "    fitted = np.full(n, np.nan, dtype=x.dtype)\n",
+    "    n = len(x)\n",
+    "    complement = 1 - alpha\n",
+    "    fitted = np.empty_like(x)\n",
+    "    fitted[0] = x[0]\n",
+    "    mse = 0.0\n",
+    "    j = 0\n",
     "\n",
     "    for i in range(1, n):\n",
-    "        smoothed = (alpha * x[i - 1] + (1 - alpha) * smoothed).item()\n",
-    "        error = x[i] - smoothed\n",
-    "        mse += error * error\n",
-    "        fitted[i] = smoothed\n",
+    "        fitted[i] = alpha * x[j] + complement * fitted[j]\n",
+    "        mse += (x[i] - fitted[i]) ** 2\n",
+    "        j += 1\n",
     "\n",
     "    mse /= n\n",
-    "    forecast = alpha * x[-1] + (1 - alpha) * smoothed\n",
+    "    forecast = alpha * x[j] + complement * fitted[j]\n",
     "    return forecast, mse, fitted\n",
     "\n",
     "\n",
@@ -3701,16 +3702,14 @@
     "\n",
     "def _optimized_ses_forecast(\n",
     "        x: np.ndarray,\n",
-    "        bounds: Sequence[Tuple[float, float]] = [(0.1, 0.3)]\n",
+    "        bounds: Tuple[float, float] = (0.1, 0.3)\n",
     "    ) -> Tuple[float, np.ndarray]:\n",
     "    r\"\"\"Searches for the optimal alpha and computes SES one step forecast.\"\"\"\n",
-    "    alpha = minimize(\n",
+    "    alpha = minimize_scalar(\n",
     "        fun=_ses_mse,\n",
-    "        x0=(0,),\n",
-    "        args=(x,),\n",
     "        bounds=bounds,\n",
-    "        method='L-BFGS-B'\n",
-    "    ).x[0]\n",
+    "        args=(x,),\n",
+    "    ).x\n",
     "    forecast, fitted = _ses_forecast(x, alpha)\n",
     "    return forecast, fitted\n",
     "\n",
@@ -3930,7 +3929,7 @@
     "#in the `ses` function\n",
     "np.testing.assert_allclose(\n",
     "    ses.predict_in_sample()['fitted'][[0, 1, -1]], \n",
-    "    np.array([np.nan, 118 - 6., 432 + 31.447525])\n",
+    "    np.array([112, 118 - 6., 432 + 31.447525])\n",
     ")"
    ]
   },
@@ -4057,7 +4056,7 @@
     "        h: int, # forecasting horizon\n",
     "        fitted: bool, # fitted values\n",
     "    ):\n",
-    "    fcst_, fitted_vals = _optimized_ses_forecast(y, [(0.01, 0.99)])\n",
+    "    fcst_, fitted_vals = _optimized_ses_forecast(y, (0.01, 0.99))\n",
     "    mean = _repeat_val(val=fcst_, h=h)\n",
     "    fcst = {'mean': mean}\n",
     "    if fitted:\n",
@@ -4746,7 +4745,7 @@
     "    fitted_vals = np.full_like(y, np.nan)\n",
     "    for i in range(season_length):\n",
     "        init_idx = (i + n % season_length)\n",
-    "        season_vals[i], fitted_vals[init_idx::season_length] = _optimized_ses_forecast(y[init_idx::season_length], [(0.01, 0.99)])\n",
+    "        season_vals[i], fitted_vals[init_idx::season_length] = _optimized_ses_forecast(y[init_idx::season_length], (0.01, 0.99))\n",
     "    out = _repeat_val_seas(season_vals=season_vals, h=h)\n",
     "    fcst = {'mean': out}\n",
     "    if fitted:\n",
@@ -13027,9 +13026,9 @@
  ],
  "metadata": {
   "kernelspec": {
-   "display_name": "python3",
+   "display_name": "statsforecast",
    "language": "python",
-   "name": "python3"
+   "name": "statsforecast"
   }
  },
  "nbformat": 4,

python/statsforecast/models.py

@@ -12,11 +12,11 @@
 # %% ../../nbs/src/core/models.ipynb 4
 import warnings
 from math import trunc
-from typing import Any, Dict, List, Optional, Sequence, Tuple, Union
+from typing import Any, Dict, List, Optional, Tuple, Union
 
 import numpy as np
 from numba import njit
-from scipy.optimize import minimize
+from scipy.optimize import minimize_scalar
 from scipy.special import inv_boxcox
 
 from statsforecast.arima import (
@@ -2136,19 +2136,20 @@ def _ses_fcst_mse(x: np.ndarray, alpha: float) -> Tuple[float, float, np.ndarray
     This function returns the one step ahead prediction
     as well as the mean squared error of the fit.
     """
-    smoothed = x[0]
-    n = x.size
+    n = len(x)
+    complement = 1 - alpha
+    fitted = np.empty_like(x)
+    fitted[0] = x[0]
     mse = 0.0
-    fitted = np.full(n, np.nan, dtype=x.dtype)
+    j = 0
 
     for i in range(1, n):
-        smoothed = (alpha * x[i - 1] + (1 - alpha) * smoothed).item()
-        error = x[i] - smoothed
-        mse += error * error
-        fitted[i] = smoothed
+        fitted[i] = alpha * x[j] + complement * fitted[j]
+        mse += (x[i] - fitted[i]) ** 2
+        j += 1
 
     mse /= n
-    forecast = alpha * x[-1] + (1 - alpha) * smoothed
+    forecast = alpha * x[j] + complement * fitted[j]
     return forecast, mse, fitted
 
 
@@ -2181,12 +2182,14 @@ def _probability(x: np.ndarray) -> np.ndarray:
 
 
 def _optimized_ses_forecast(
-    x: np.ndarray, bounds: Sequence[Tuple[float, float]] = [(0.1, 0.3)]
+    x: np.ndarray, bounds: Tuple[float, float] = (0.1, 0.3)
 ) -> Tuple[float, np.ndarray]:
     r"""Searches for the optimal alpha and computes SES one step forecast."""
-    alpha = minimize(
-        fun=_ses_mse, x0=(0,), args=(x,), bounds=bounds, method="L-BFGS-B"
-    ).x[0]
+    alpha = minimize_scalar(
+        fun=_ses_mse,
+        bounds=bounds,
+        args=(x,),
+    ).x
     forecast, fitted = _ses_forecast(x, alpha)
     return forecast, fitted
 
@@ -2380,7 +2383,7 @@ def _ses_optimized(
     h: int,  # forecasting horizon
     fitted: bool,  # fitted values
 ):
-    fcst_, fitted_vals = _optimized_ses_forecast(y, [(0.01, 0.99)])
+    fcst_, fitted_vals = _optimized_ses_forecast(y, (0.01, 0.99))
     mean = _repeat_val(val=fcst_, h=h)
     fcst = {"mean": mean}
     if fitted:
@@ -2759,7 +2762,7 @@ def _seasonal_ses_optimized(
     for i in range(season_length):
         init_idx = i + n % season_length
         season_vals[i], fitted_vals[init_idx::season_length] = _optimized_ses_forecast(
-            y[init_idx::season_length], [(0.01, 0.99)]
+            y[init_idx::season_length], (0.01, 0.99)
         )
     out = _repeat_val_seas(season_vals=season_vals, h=h)
     fcst = {"mean": out}