<jupyter_start><jupyter_text>Определение параметров адаптивных моделей  Параметры адаптивных моделей ($\alpha, \beta, \gamma, \phi$) выбираются таким образом, чтобы ошибка при аппроксимации полученной моделью была как можно меньше, т.е. необходимо решить задачу оптимизации.  Оптимизация в SciPy Решать задачи оптимизации позволяет библиотека SciPy. Алгоритмы оптимизации лежат в модуле `scipy.optimize`. Найдем минимум следующей функции:$$f(x)=x^2+2x^4-3x.$$<jupyter_code>def fun(x):
    return x**2 + 2*x**4 - 3*x<jupyter_output><empty_output><jupyter_text>Воспользуемся функцией `minimize` из указанного модуля. В данной функции реализованы численные методы решения задач оптимизации. В качестве первого аргумента `minimize` принимает функцию, минимум которой необходимо найти, а в качестве второго аргумента &ndash; начальное значение переменных.<jupyter_code>from scipy.optimize import minimize
result = minimize(fun, 0)
result<jupyter_output><empty_output><jupyter_text>Значение целевой функции и оптимальное значение $x$:<jupyter_code>result.fun, result.x<jupyter_output><empty_output><jupyter_text>Функция нескольких переменных может быть передана в `minimize` двумя способами:<jupyter_code>def fun1(x1, x2):
    return x1**2 + 2*x1*x2 + x2**2
minimize(lambda x: fun1(x[0], x[1]), (2,2))
fun2 = lambda x: x[0]**2 + 2*x[0]*x[1] + x[1]**2
minimize(fun2, (2,2))<jupyter_output><empty_output><jupyter_text>Можно указать ограничения на значения переменных с помощью опций `constraints` и `bounds`.<jupyter_code>fun3 = lambda x: x[0] + x[1]
# x1 >= 0, x2 >= 0
bnds = ((0,None),(0,None))

# 3*x1 + 2*x2 >= 7
# x1 + 2*x2 >= 6
cons = ({'type': 'ineq', 'fun': lambda x: 3*x[0] + 2*x[1] - 7},
        {'type': 'ineq', 'fun': lambda x: x[0] + 2*x[1] - 6})
minimize(fun3, (2,2), bounds=bnds, constraints=cons)<jupyter_output><empty_output><jupyter_text>Задание 11. Определите функцию, которая принимает на вход временной ряд и параметры модели Хольта для линейного тренда (см. задание 2 с предыдущего занятия). Данная функция должна возвращать среднеквадратичную ошибку модели.2. Определите оптимальные параметры $\alpha$ и $\beta$ для ряда объемов пассажирских авиаперевозок в Австралии (файл ausair.csv). Постройте модель Хольта с найденными параметрами и получите прогноз на 10 лет вперед.<jupyter_code>import numpy as np
import pandas as pd

import matplotlib.pyplot as plt
%matplotlib inline
plt.rcParams["figure.figsize"] = (12,5)
def holt_method_MSE(data, alpha, beta):
    n = len(data)
    y, l, b = np.zeros(n), np.zeros(n), np.zeros(n)
    values = data.values.reshape(1, n)[0]
    y[0] = values[0]
    l[0] = values[0]
    b[0] = values[1] - values[0]
    for i in range(1, n):
        y[i] = l[i - 1] + b[i - 1]
        l[i] = alpha*values[i] + (1 - alpha)*(l[i - 1] + b[i - 1])
        b[i] = beta*(l[i] - l[i - 1]) + (1 - beta)*b[i - 1]
    
    mse = 1/n * sum((values[1:] - y[:-1])**2)
    return mse
passengers = pd.read_csv('data/ausair.csv', index_col='Year', parse_dates=True)
holt_method_MSE(passengers, 0.59, 0.43)
bnds = ((0,1),(0,1))
result = minimize(lambda x: holt_method_MSE(passengers, x[0], x[1]), (0,0), bounds=bnds)
alpha1, beta1 = result.x
result.fun, result.x
def holt_method(data, years, alpha, beta):
    n = len(data)
    y, l, b = np.zeros(n + years), np.zeros(n + years), np.zeros(n + years)
    values = data.values.reshape(1, n)[0]
    y[0] = values[0]
    l[0] = values[0]
    b[0] = values[1] - values[0]
    for i in range(1, n + years):
        y[i] = l[i - 1] + b[i - 1]
        if i >= n:
            base_values = y
        else:
            base_values = values
        l[i] = alpha*base_values[i] + (1 - alpha)*(l[i - 1] + b[i - 1])
        b[i] = beta*(l[i] - l[i - 1]) + (1 - beta)*b[i - 1]
    
    first_year = data.index[0]
    new_index = pd.date_range(first_year, periods=n+years, freq='YS')
    
    prediction = pd.DataFrame(y, columns=['Passengers'])
    prediction = prediction.set_index(new_index)
    prediction.index.name = 'Year'
    return prediction
plt.plot(passengers['Passengers'], label='Начальные данные')
plt.plot(holt_method(passengers, 10, alpha1, beta1)['Passengers'], label='Линейный тренд')
plt.legend()
plt.show()<jupyter_output><empty_output><jupyter_text>Задание 21. Для модели Хольта с аддитивным затухающим трендом (задание 3 с предыдущего занятия) проделайте те же шаги, что и в задании 1. В данном случае кроме параметров $\alpha$ и $\beta$ необходимо также определить оптимальный параметр $\phi$.2. Сравните на графике аппроксимации и прогнозы с помощью модели для линейного тренда и модели с аддитивным затухающим трендом.<jupyter_code>def halt_decay_method_MSE(data, alpha, beta, phi):
    n = len(data)
    y, l, b = np.zeros(n), np.zeros(n), np.zeros(n)
    values = data.values.reshape(1, n)[0]
    y[0] = values[0]
    l[0] = values[0]
    b[0] = values[1] - values[0]
    for i in range(1, n):
        y[i] = l[i - 1] + phi*b[i - 1]
        l[i] = alpha*values[i] + (1 - alpha)*(l[i - 1] + phi*b[i - 1])
        b[i] = beta*(l[i] - l[i - 1]) + (1 - beta)*phi*b[i - 1]
        
    mse = 1/n * sum((values[1:] - y[:-1])**2)
    return mse
halt_decay_method_MSE(passengers, 0.39, 0.10, 0.9)
bnds = ((0,1),(0,1),(0,1))
result = minimize(lambda x: halt_decay_method_MSE(passengers, x[0], x[1], x[2]), (0,0,0), bounds=bnds)
alpha2, beta2, phi2 = result.x
result.fun, result.x
def halt_decay_method(data, years, alpha, beta, phi):
    n = len(data)
    y, l, b = np.zeros(n + years), np.zeros(n + years), np.zeros(n + years)
    values = data.values.reshape(1, n)[0]
    y[0] = values[0]
    l[0] = values[0]
    b[0] = values[1] - values[0]
    for i in range(1, n + years):
        y[i] = l[i - 1] + phi*b[i - 1]
        if i >= n:
            base_values = y
        else:
            base_values = values
        l[i] = alpha*base_values[i] + (1 - alpha)*(l[i - 1] + phi*b[i - 1])
        b[i] = beta*(l[i] - l[i - 1]) + (1 - beta)*phi*b[i - 1]
    
    first_year = data.index[0]
    new_index = pd.date_range(first_year, periods=n+years, freq='YS')
    
    prediction = pd.DataFrame(y, columns=['Passengers'])
    prediction = prediction.set_index(new_index)
    prediction.index.name = 'Year'
    return prediction
plt.plot(halt_decay_method(passengers, 10, alpha2, beta2, phi2)['Passengers'], label='Аддитивный затухающий тренд')
plt.plot(holt_method(passengers, 10, alpha1, beta1)['Passengers'], label='Линейный тренд')
plt.legend()
plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.16009280742459397, 'fraction_numerical': 0.032946635730858466, 'mean_word_length': 3.861654135338346, 'pattern_counts': {'":': 0, '<': 29, '<?xml version=': 0, '>': 35, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '7', 'dir': 'jupyter-structured-clean-dedup', 'id': '215668', 'license': 'MIT', 'n_tokens_mistral': 2845, 'n_tokens_neox': 2931, 'n_words': 626, 'path': 'Forecasting-Methods/2. Adaptive methods/2. Adaptive forecasting methods. SciPy.ipynb', 'repo_name': 'cmlimm/uni-projects', 'size': '6465'}