Skip to content

Instantly share code, notes, and snippets.

@andrequeiroz

andrequeiroz/holtwinters.py

Last active February 3, 2025 04:35
Show Gist options
  • Save andrequeiroz/5888967 to your computer and use it in GitHub Desktop.
Save andrequeiroz/5888967 to your computer and use it in GitHub Desktop.
Download ZIP

Revisions

  1. andrequeiroz revised this gist Aug 5, 2015. 1 changed file with 22 additions and 0 deletions.
    22 changes: 22 additions & 0 deletions holtwinters.py
    View file
    Original file line number Diff line number Diff line change
    @@ -1,3 +1,25 @@
    #The MIT License (MIT)
    #
    #Copyright (c) 2015 Andre Queiroz
    #
    #Permission is hereby granted, free of charge, to any person obtaining a copy
    #of this software and associated documentation files (the "Software"), to deal
    #in the Software without restriction, including without limitation the rights
    #to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
    #copies of the Software, and to permit persons to whom the Software is
    #furnished to do so, subject to the following conditions:
    #
    #The above copyright notice and this permission notice shall be included in
    #all copies or substantial portions of the Software.
    #
    #THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
    #IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
    #FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
    #AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
    #LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
    #OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
    #THE SOFTWARE.
    #
    # Holt-Winters algorithms to forecasting
    # Coded in Python 2 by: Andre Queiroz
    # Description: This module contains three exponential smoothing algorithms. They are Holt's linear trend method and Holt-Winters seasonal methods (additive and multiplicative).
  2. andrequeiroz revised this gist Jun 12, 2015. 1 changed file with 2 additions and 1 deletion.
    3 changes: 2 additions & 1 deletion holtwinters.py
    View file
    Original file line number Diff line number Diff line change
    @@ -5,6 +5,7 @@
    # Hyndman, R. J.; Athanasopoulos, G. (2013) Forecasting: principles and practice. http://otexts.com/fpp/. Accessed on 07/03/2013.
    # Byrd, R. H.; Lu, P.; Nocedal, J. A Limited Memory Algorithm for Bound Constrained Optimization, (1995), SIAM Journal on Scientific and Statistical Computing, 16, 5, pp. 1190-1208.

    from __future__ import division
    from sys import exit
    from math import sqrt
    from numpy import array
    @@ -58,7 +59,7 @@ def RMSE(params, *args):
    a.append(alpha * (Y[i] / s[i]) + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    s.append(gamma * (Y[i] / (a[i] + b[i])) + (1 - gamma) * s[i])
    y.append(a[i + 1] + b[i + 1] + s[i + 1])
    y.append((a[i + 1] + b[i + 1]) * s[i + 1])

    else:

  3. andrequeiroz revised this gist Jul 5, 2013. 1 changed file with 48 additions and 49 deletions.
    97 changes: 48 additions & 49 deletions holtwinters.py
    View file
    Original file line number Diff line number Diff line change
    @@ -5,68 +5,64 @@
    # Hyndman, R. J.; Athanasopoulos, G. (2013) Forecasting: principles and practice. http://otexts.com/fpp/. Accessed on 07/03/2013.
    # Byrd, R. H.; Lu, P.; Nocedal, J. A Limited Memory Algorithm for Bound Constrained Optimization, (1995), SIAM Journal on Scientific and Statistical Computing, 16, 5, pp. 1190-1208.

    from sys import exit
    from math import sqrt
    from numpy import array
    from scipy.optimize import fmin_l_bfgs_b

    def rmse_linear(params, *args):
    def RMSE(params, *args):

    alpha, beta = params
    Y = args[0]
    a = [Y[0]]
    b = [Y[1] - Y[0]]
    y = [a[0] + b[0]]
    type = args[1]
    rmse = 0

    for i in range(len(Y)):
    if type == 'linear':

    a.append(alpha * Y[i] + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    y.append(a[i + 1] + b[i + 1])

    rmse = sqrt(sum([(m - n) ** 2 for m, n in zip(Y, y[:-1])]) / len(Y))
    alpha, beta = params
    a = [Y[0]]
    b = [Y[1] - Y[0]]
    y = [a[0] + b[0]]

    return rmse

    def rmse_additive(params, *args):
    for i in range(len(Y)):

    alpha, beta, gamma = params
    Y = args[0]
    m = args[1]
    a = [sum(Y[0:m]) / float(m)]
    b = [(sum(Y[m:2 * m]) - sum(Y[0:m])) / m ** 2]
    s = [Y[i] - a[0] for i in range(m)]
    y = [a[0] + b[0] + s[0]]
    rmse = 0
    a.append(alpha * Y[i] + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    y.append(a[i + 1] + b[i + 1])

    for i in range(len(Y)):
    else:

    a.append(alpha * (Y[i] - s[i]) + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    s.append(gamma * (Y[i] - a[i] - b[i]) + (1 - gamma) * s[i])
    y.append(a[i + 1] + b[i + 1] + s[i + 1])

    rmse = sqrt(sum([(m - n) ** 2 for m, n in zip(Y, y[:-1])]) / len(Y))
    alpha, beta, gamma = params
    m = args[2]
    a = [sum(Y[0:m]) / float(m)]
    b = [(sum(Y[m:2 * m]) - sum(Y[0:m])) / m ** 2]

    return rmse

    def rmse_multiplicative(params, *args):
    if type == 'additive':

    alpha, beta, gamma = params
    Y = args[0]
    m = args[1]
    a = [sum(Y[0:m]) / float(m)]
    b = [(sum(Y[m:2 * m]) - sum(Y[0:m])) / m ** 2]
    s = [Y[i] / a[0] for i in range(m)]
    y = [(a[0] + b[0]) * s[0]]
    rmse = 0
    s = [Y[i] - a[0] for i in range(m)]
    y = [a[0] + b[0] + s[0]]

    for i in range(len(Y)):
    for i in range(len(Y)):

    a.append(alpha * (Y[i] / s[i]) + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    s.append(gamma * (Y[i] / (a[i] + b[i])) + (1 - gamma) * s[i])
    y.append(a[i + 1] + b[i + 1] + s[i + 1])
    a.append(alpha * (Y[i] - s[i]) + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    s.append(gamma * (Y[i] - a[i] - b[i]) + (1 - gamma) * s[i])
    y.append(a[i + 1] + b[i + 1] + s[i + 1])

    elif type == 'multiplicative':

    s = [Y[i] / a[0] for i in range(m)]
    y = [(a[0] + b[0]) * s[0]]

    for i in range(len(Y)):

    a.append(alpha * (Y[i] / s[i]) + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    s.append(gamma * (Y[i] / (a[i] + b[i])) + (1 - gamma) * s[i])
    y.append(a[i + 1] + b[i + 1] + s[i + 1])

    else:

    exit('Type must be either linear, additive or multiplicative')

    rmse = sqrt(sum([(m - n) ** 2 for m, n in zip(Y, y[:-1])]) / len(Y))

    @@ -80,8 +76,9 @@ def linear(x, fc, alpha = None, beta = None):

    initial_values = array([0.3, 0.1])
    boundaries = [(0, 1), (0, 1)]
    type = 'linear'

    parameters = fmin_l_bfgs_b(rmse_linear, x0 = initial_values, args = (Y,), bounds = boundaries, approx_grad = True)
    parameters = fmin_l_bfgs_b(RMSE, x0 = initial_values, args = (Y, type), bounds = boundaries, approx_grad = True)
    alpha, beta = parameters[0]

    a = [Y[0]]
    @@ -110,8 +107,9 @@ def additive(x, m, fc, alpha = None, beta = None, gamma = None):

    initial_values = array([0.3, 0.1, 0.1])
    boundaries = [(0, 1), (0, 1), (0, 1)]
    type = 'additive'

    parameters = fmin_l_bfgs_b(rmse_additive, x0 = initial_values, args = (Y, m), bounds = boundaries, approx_grad = True)
    parameters = fmin_l_bfgs_b(RMSE, x0 = initial_values, args = (Y, type, m), bounds = boundaries, approx_grad = True)
    alpha, beta, gamma = parameters[0]

    a = [sum(Y[0:m]) / float(m)]
    @@ -140,10 +138,11 @@ def multiplicative(x, m, fc, alpha = None, beta = None, gamma = None):

    if (alpha == None or beta == None or gamma == None):

    initial_values = array([0.3, 0.1, 0.1])
    initial_values = array([0.0, 1.0, 0.0])
    boundaries = [(0, 1), (0, 1), (0, 1)]
    type = 'multiplicative'

    parameters = fmin_l_bfgs_b(rmse_multiplicative, x0 = initial_values, args = (Y, m), bounds = boundaries, approx_grad = True)
    parameters = fmin_l_bfgs_b(RMSE, x0 = initial_values, args = (Y, type, m), bounds = boundaries, approx_grad = True)
    alpha, beta, gamma = parameters[0]

    a = [sum(Y[0:m]) / float(m)]
  4. andrequeiroz revised this gist Jul 4, 2013. 1 changed file with 139 additions and 30 deletions.
    169 changes: 139 additions & 30 deletions holtwinters.py
    View file
    Original file line number Diff line number Diff line change
    @@ -1,58 +1,167 @@
    # Holt-Winters algorithms to forecasting
    # Coded in Python 2 by: Andre Queiroz
    # Description: This module contains three exponential smoothing algorithms. They are Holt's linear trend method and Holt-Winters seasonal methods (additive and multiplicative).
    # References:
    # Brockwell, P. J.; Davis, R. A. Introduction to Time Series and Forecasting. New York: Springer-Verlang, 2002. 434 p.
    # Hyndman, R. J.; Athanasopoulos, G. (2013) Forecasting: principles and practice. http://otexts.com/fpp/. Accessed on 07/03/2013.
    # Byrd, R. H.; Lu, P.; Nocedal, J. A Limited Memory Algorithm for Bound Constrained Optimization, (1995), SIAM Journal on Scientific and Statistical Computing, 16, 5, pp. 1190-1208.

    from math import sqrt
    from numpy import array
    from scipy.optimize import fmin_l_bfgs_b

    def mse_non_seasonal(params, *args):
    def rmse_linear(params, *args):

    alpha, beta = params
    Y = args[0]
    a = [Y[1]]
    a = [Y[0]]
    b = [Y[1] - Y[0]]
    y = [a[0] + b[0]]
    e = []
    mse = 0
    rmse = 0

    for i in range(2, len(Y)):
    for i in range(len(Y)):

    a.append(alpha * Y[i] + (1 - alpha) * (a[i - 2] + b[i - 2]))
    b.append(beta * (a[i - 1] - a[i - 2]) + (1 - beta) * b[i - 2])
    y.append(a[i - 1] + b[i - 1])
    e.append((Y[i] - y[i - 2]) ** 2)
    a.append(alpha * Y[i] + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    y.append(a[i + 1] + b[i + 1])

    mse = sqrt(sum([(m - n) ** 2 for m, n in zip(Y[2:], y)]) / len(Y[2:]))
    rmse = sqrt(sum([(m - n) ** 2 for m, n in zip(Y, y[:-1])]) / len(Y))

    return mse
    return rmse

    def rmse_additive(params, *args):

    def non_seasonal(x, fc):
    alpha, beta, gamma = params
    Y = args[0]
    m = args[1]
    a = [sum(Y[0:m]) / float(m)]
    b = [(sum(Y[m:2 * m]) - sum(Y[0:m])) / m ** 2]
    s = [Y[i] - a[0] for i in range(m)]
    y = [a[0] + b[0] + s[0]]
    rmse = 0

    for i in range(len(Y)):

    a.append(alpha * (Y[i] - s[i]) + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    s.append(gamma * (Y[i] - a[i] - b[i]) + (1 - gamma) * s[i])
    y.append(a[i + 1] + b[i + 1] + s[i + 1])

    rmse = sqrt(sum([(m - n) ** 2 for m, n in zip(Y, y[:-1])]) / len(Y))

    return rmse

    def rmse_multiplicative(params, *args):

    alpha, beta, gamma = params
    Y = args[0]
    m = args[1]
    a = [sum(Y[0:m]) / float(m)]
    b = [(sum(Y[m:2 * m]) - sum(Y[0:m])) / m ** 2]
    s = [Y[i] / a[0] for i in range(m)]
    y = [(a[0] + b[0]) * s[0]]
    rmse = 0

    for i in range(len(Y)):

    a.append(alpha * (Y[i] / s[i]) + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    s.append(gamma * (Y[i] / (a[i] + b[i])) + (1 - gamma) * s[i])
    y.append(a[i + 1] + b[i + 1] + s[i + 1])

    rmse = sqrt(sum([(m - n) ** 2 for m, n in zip(Y, y[:-1])]) / len(Y))

    valores_iniciais = array([0.0, 0.0])
    limites = [(0, 1), (0, 1)]
    return rmse

    parametros = fmin_l_bfgs_b(mse_non_seasonal, x0 = valores_iniciais, args = (x,), bounds = limites, approx_grad = True)
    def linear(x, fc, alpha = None, beta = None):

    alpha, beta = parametros[0]
    a = [x[1]]
    b = [x[1] - x[0]]
    Y = x[:]

    if (alpha == None or beta == None):

    initial_values = array([0.3, 0.1])
    boundaries = [(0, 1), (0, 1)]

    parameters = fmin_l_bfgs_b(rmse_linear, x0 = initial_values, args = (Y,), bounds = boundaries, approx_grad = True)
    alpha, beta = parameters[0]

    a = [Y[0]]
    b = [Y[1] - Y[0]]
    y = [a[0] + b[0]]
    e = []
    mse = 0
    rmse = 0

    for i in range(len(Y) + fc):

    if i == len(Y):
    Y.append(a[-1] + b[-1])

    a.append(alpha * Y[i] + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    y.append(a[i + 1] + b[i + 1])

    rmse = sqrt(sum([(m - n) ** 2 for m, n in zip(Y[:-fc], y[:-fc - 1])]) / len(Y[:-fc]))

    return Y[-fc:], alpha, beta, rmse

    def additive(x, m, fc, alpha = None, beta = None, gamma = None):

    Y = x[:]

    if (alpha == None or beta == None or gamma == None):

    initial_values = array([0.3, 0.1, 0.1])
    boundaries = [(0, 1), (0, 1), (0, 1)]

    parameters = fmin_l_bfgs_b(rmse_additive, x0 = initial_values, args = (Y, m), bounds = boundaries, approx_grad = True)
    alpha, beta, gamma = parameters[0]

    a = [sum(Y[0:m]) / float(m)]
    b = [(sum(Y[m:2 * m]) - sum(Y[0:m])) / m ** 2]
    s = [Y[i] - a[0] for i in range(m)]
    y = [a[0] + b[0] + s[0]]
    rmse = 0

    for i in range(len(Y) + fc):

    if i == len(Y):
    Y.append(a[-1] + b[-1] + s[-m])

    a.append(alpha * (Y[i] - s[i]) + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    s.append(gamma * (Y[i] - a[i] - b[i]) + (1 - gamma) * s[i])
    y.append(a[i + 1] + b[i + 1] + s[i + 1])

    rmse = sqrt(sum([(m - n) ** 2 for m, n in zip(Y[:-fc], y[:-fc - 1])]) / len(Y[:-fc]))

    return Y[-fc:], alpha, beta, gamma, rmse

    def multiplicative(x, m, fc, alpha = None, beta = None, gamma = None):

    Y = x[:]

    if (alpha == None or beta == None or gamma == None):

    initial_values = array([0.3, 0.1, 0.1])
    boundaries = [(0, 1), (0, 1), (0, 1)]

    parameters = fmin_l_bfgs_b(rmse_multiplicative, x0 = initial_values, args = (Y, m), bounds = boundaries, approx_grad = True)
    alpha, beta, gamma = parameters[0]

    a = [sum(Y[0:m]) / float(m)]
    b = [(sum(Y[m:2 * m]) - sum(Y[0:m])) / m ** 2]
    s = [Y[i] / a[0] for i in range(m)]
    y = [(a[0] + b[0]) * s[0]]
    rmse = 0

    for i in range(2, len(x) + fc):
    for i in range(len(Y) + fc):

    if i == len(x):
    x.append(a[-1] + b[-1])
    if i == len(Y):
    Y.append((a[-1] + b[-1]) * s[-m])

    a.append(alpha * x[i] + (1 - alpha) * (a[i - 2] + b[i - 2]))
    b.append(beta * (a[i - 1] - a[i - 2]) + (1 - beta) * b[i - 2])
    y.append(a[i - 1] + b[i - 1])
    e.append((x[i] - y[i - 2]) ** 2)
    a.append(alpha * (Y[i] / s[i]) + (1 - alpha) * (a[i] + b[i]))
    b.append(beta * (a[i + 1] - a[i]) + (1 - beta) * b[i])
    s.append(gamma * (Y[i] / (a[i] + b[i])) + (1 - gamma) * s[i])
    y.append((a[i + 1] + b[i + 1]) * s[i + 1])

    mse = sqrt(sum([(m - n) ** 2 for m, n in zip(x[2:-fc], y[:-fc - 1])]) / len(x[2:-fc]))
    rmse = sqrt(sum([(m - n) ** 2 for m, n in zip(Y[:-fc], y[:-fc - 1])]) / len(Y[:-fc]))

    return x[-fc:], alpha, beta, mse
    return Y[-fc:], alpha, beta, gamma, rmse
  5. andrequeiroz revised this gist Jul 1, 2013. 1 changed file with 30 additions and 25 deletions.
    55 changes: 30 additions & 25 deletions holtwinters.py
    View file
    Original file line number Diff line number Diff line change
    @@ -1,37 +1,42 @@
    # Holt-Winters algorithms to forecasting
    # Coded in Python 2 by: Andre Queiroz
    # Reference: Brockwell, P. J.; Davis, R. A. Introduction to Time Series and Forecasting. New York: Springer-Verlang, 2002. 434 p.
    # References:
    # Brockwell, P. J.; Davis, R. A. Introduction to Time Series and Forecasting. New York: Springer-Verlang, 2002. 434 p.
    # Byrd, R. H.; Lu, P.; Nocedal, J. A Limited Memory Algorithm for Bound Constrained Optimization, (1995), SIAM Journal on Scientific and Statistical Computing, 16, 5, pp. 1190-1208.

    from sys import maxint
    from math import sqrt
    from numpy import array
    from scipy.optimize import fmin_l_bfgs_b

    def non_seasonal(x, fc):
    def mse_non_seasonal(params, *args):

    alpha_opt = 0
    beta_opt = 0
    erro = maxint
    alpha, beta = params
    Y = args[0]
    a = [Y[1]]
    b = [Y[1] - Y[0]]
    y = [a[0] + b[0]]
    e = []
    mse = 0

    for alpha in range(101):
    for i in range(2, len(Y)):

    for beta in range(101):
    a.append(alpha * Y[i] + (1 - alpha) * (a[i - 2] + b[i - 2]))
    b.append(beta * (a[i - 1] - a[i - 2]) + (1 - beta) * b[i - 2])
    y.append(a[i - 1] + b[i - 1])
    e.append((Y[i] - y[i - 2]) ** 2)

    mse = sqrt(sum([(m - n) ** 2 for m, n in zip(Y[2:], y)]) / len(Y[2:]))

    a = [x[1]]
    b = [x[1] - x[0]]
    y = [a[0] + b[0]]
    e = []
    return mse

    for i in range(2, len(x)):
    def non_seasonal(x, fc):

    a.append(alpha / float(100) * x[i] + (1 - alpha / float(100)) * (a[i - 2] + b[i - 2]))
    b.append(beta / float(100) * (a[i - 1] - a[i - 2]) + (1 - beta / float(100)) * b[i - 2])
    y.append(a[i - 1] + b[i - 1])
    e.append((x[i] - y[i - 2]) ** 2)
    valores_iniciais = array([0.0, 0.0])
    limites = [(0, 1), (0, 1)]

    if sum(e) < erro:
    erro = sum(e)
    alpha_opt = alpha / float(100)
    beta_opt = beta / float(100)
    parametros = fmin_l_bfgs_b(mse_non_seasonal, x0 = valores_iniciais, args = (x,), bounds = limites, approx_grad = True)

    alpha, beta = parametros[0]
    a = [x[1]]
    b = [x[1] - x[0]]
    y = [a[0] + b[0]]
    @@ -43,11 +48,11 @@ def non_seasonal(x, fc):
    if i == len(x):
    x.append(a[-1] + b[-1])

    a.append(alpha_opt * x[i] + (1 - alpha_opt) * (a[i - 2] + b[i - 2]))
    b.append(beta_opt * (a[i - 1] - a[i - 2]) + (1 - beta_opt) * b[i - 2])
    a.append(alpha * x[i] + (1 - alpha) * (a[i - 2] + b[i - 2]))
    b.append(beta * (a[i - 1] - a[i - 2]) + (1 - beta) * b[i - 2])
    y.append(a[i - 1] + b[i - 1])
    e.append((x[i] - y[i - 2]) ** 2)

    mse = sqrt(sum([(m - n) ** 2 for m, n in zip(x[2:-fc], y[:-fc - 1])]) / len(x[2:-fc]))

    return x[-fc:], alpha_opt, beta_opt, mse
    return x[-fc:], alpha, beta, mse
  6. andrequeiroz created this gist Jun 28, 2013.
    53 changes: 53 additions & 0 deletions holtwinters.py
    View file
    Original file line number Diff line number Diff line change
    @@ -0,0 +1,53 @@
    # Holt-Winters algorithms to forecasting
    # Coded in Python 2 by: Andre Queiroz
    # Reference: Brockwell, P. J.; Davis, R. A. Introduction to Time Series and Forecasting. New York: Springer-Verlang, 2002. 434 p.

    from sys import maxint
    from math import sqrt

    def non_seasonal(x, fc):

    alpha_opt = 0
    beta_opt = 0
    erro = maxint

    for alpha in range(101):

    for beta in range(101):

    a = [x[1]]
    b = [x[1] - x[0]]
    y = [a[0] + b[0]]
    e = []

    for i in range(2, len(x)):

    a.append(alpha / float(100) * x[i] + (1 - alpha / float(100)) * (a[i - 2] + b[i - 2]))
    b.append(beta / float(100) * (a[i - 1] - a[i - 2]) + (1 - beta / float(100)) * b[i - 2])
    y.append(a[i - 1] + b[i - 1])
    e.append((x[i] - y[i - 2]) ** 2)

    if sum(e) < erro:
    erro = sum(e)
    alpha_opt = alpha / float(100)
    beta_opt = beta / float(100)

    a = [x[1]]
    b = [x[1] - x[0]]
    y = [a[0] + b[0]]
    e = []
    mse = 0

    for i in range(2, len(x) + fc):

    if i == len(x):
    x.append(a[-1] + b[-1])

    a.append(alpha_opt * x[i] + (1 - alpha_opt) * (a[i - 2] + b[i - 2]))
    b.append(beta_opt * (a[i - 1] - a[i - 2]) + (1 - beta_opt) * b[i - 2])
    y.append(a[i - 1] + b[i - 1])
    e.append((x[i] - y[i - 2]) ** 2)

    mse = sqrt(sum([(m - n) ** 2 for m, n in zip(x[2:-fc], y[:-fc - 1])]) / len(x[2:-fc]))

    return x[-fc:], alpha_opt, beta_opt, mse