Глава 4. ОПТИМИЗАЦИЯ

4.1. Градиентный спуск

\[ \boldsymbol{\theta}^{(t+1)} = \boldsymbol{\theta}^{(t)} - \eta \nabla J(\boldsymbol{\theta}^{(t)}) \]

Объяснение: Градиентный спуск — это алгоритм оптимизации, который итеративно обновляет параметры в направлении отрицательного градиента, чтобы минимизировать функцию стоимости \(J(\boldsymbol{\theta})\). \(\eta\) - скорость обучения.

Пример: Для \(J(\theta) = \theta^2\) и \(\eta = 0.1\). \(\nabla J(\theta) = 2\theta\). Обновление: \(\theta^{(t+1)} = \theta^{(t)} - 0.1 \cdot (2\theta^{(t)}) = \theta^{(t)} - 0.2\theta^{(t)} = 0.8\theta^{(t)}\).

Реализация:


import numpy as np

# gradient - функция, вычисляющая градиент J по theta
# theta - начальные параметры (numpy array)
# eta - скорость обучения
# steps - количество итераций
def gradient_descent(gradient, theta, eta, steps):
  """Выполняет градиентный спуск."""
  theta = np.copy(theta) # Копируем, чтобы не изменять исходный массив
  for _ in range(steps):
    grad = gradient(theta)
    theta -= eta * grad
  return theta

# Пример использования:
def cost_gradient(theta): # Градиент для J(theta) = theta^2 (одномерный случай)
    return 2 * theta

initial_theta = np.array([10.0])
learning_rate = 0.1
num_steps = 50

final_theta = gradient_descent(cost_gradient, initial_theta, learning_rate, num_steps)
print(f"Theta after gradient descent: {final_theta}") # Должен быть близок к 0

4.2. Стохастический градиентный спуск (SGD)

\[ \boldsymbol{\theta}^{(t+1)} = \boldsymbol{\theta}^{(t)} - \eta \nabla J_i(\boldsymbol{\theta}^{(t)}) \]

Объяснение: SGD вычисляет градиенты по отдельным точкам данных (или небольшим подвыборкам - мини-батчам), обновляя параметры чаще. Широко используется в МО из-за своей эффективности на больших наборах данных.

Пример: Для \(J_i(\theta) = (\theta - y_i)^2\) (например, ошибка для одного примера), обновление основано на градиенте \(\nabla J_i = 2(\theta - y_i)\) для одной точки данных \(y_i\) на каждой итерации.

Реализация:


import numpy as np
import random

# gradient - функция, вычисляющая градиент J_i по theta для одного примера data_i
# theta - начальные параметры (numpy array)
# eta - скорость обучения
# data - набор данных (например, список пар (x_i, y_i))
# steps - количество итераций (эпох или обновлений)
def stochastic_gradient_descent(gradient, theta, eta, data, steps):
  """Выполняет стохастический градиентный спуск."""
  theta = np.copy(theta)
  num_samples = len(data)
  for _ in range(steps):
    # Выбираем случайный пример
    random_index = random.randrange(num_samples)
    data_i = data[random_index]
    # Вычисляем градиент по этому примеру
    grad = gradient(theta, data_i)
    theta -= eta * grad
  return theta

# Пример использования (упрощенный, только y):
# J_i(theta) = (theta - y_i)^2 => grad = 2 * (theta - y_i)
def simple_sgd_gradient(theta, y_i):
    return 2 * (theta - y_i)

# Предположим, у нас есть данные y_values
y_values = np.array([1.0, 1.5, 3.0, 2.5])
initial_theta = np.array([0.0]) # Находим среднее значение через SGD
learning_rate = 0.01
num_steps = 1000 # Больше шагов, т.к. обновления шумные

final_theta_sgd = stochastic_gradient_descent(simple_sgd_gradient, initial_theta, learning_rate, y_values, num_steps)
print(f"Theta after SGD: {final_theta_sgd}") # Должен быть близок к среднему y_values (2.0)
print(f"Actual mean: {np.mean(y_values)}")

4.3. Градиентный спуск с моментом

\[ \mathbf{v}^{(t+1)} = \beta \mathbf{v}^{(t)} + (1-\beta)\nabla J(\boldsymbol{\theta}^{(t)}) \quad (\text{или } \mathbf{v}^{(t+1)} = \beta \mathbf{v}^{(t)} + \nabla J(\boldsymbol{\theta}^{(t)})) \] \[ \boldsymbol{\theta}^{(t+1)} = \boldsymbol{\theta}^{(t)} - \eta \mathbf{v}^{(t+1)} \] (Примечание: Формула в OCR \(v^{(t+1)} = \beta v^{(t)} - \eta \nabla J(\theta^{(t)}), \theta^{(t+1)} = \theta^{(t)} + v^{(t+1)}\) немного отличается, но описывает тот же принцип накопления скорости. \(\beta\) - коэффициент момента.)

Объяснение: Метод моментов добавляет экспоненциально взвешенное скользящее среднее прошлых градиентов к текущему обновлению, улучшая скорость сходимости и стабильность.

Пример: При \(\beta = 0.9\), \(\eta = 0.1\) обновление скорости сглаживает колебания градиентного спуска.

Реализация (согласно OCR):


import numpy as np

# gradient - функция градиента
# theta - начальные параметры
# eta - скорость обучения
# beta - коэффициент момента
# steps - количество шагов
def momentum_gradient_descent(gradient, theta, eta, beta, steps):
  """Градиентный спуск с моментом."""
  theta = np.copy(theta)
  v = np.zeros_like(theta) # Инициализируем скорость нулями
  for _ in range(steps):
    grad = gradient(theta)
    v = beta * v - eta * grad # Обновляем скорость (согласно OCR)
    theta += v                # Обновляем параметры (согласно OCR)
    # Альтернативная, более частая формулировка:
    # v = beta * v + (1 - beta) * grad # Скользящее среднее градиента
    # theta -= eta * v # Обновление параметрами со скоростью
  return theta

# Пример использования:
def cost_gradient(theta): return 2 * theta
initial_theta = np.array([10.0])
learning_rate = 0.1
momentum_beta = 0.9
num_steps = 50

final_theta_mom = momentum_gradient_descent(cost_gradient, initial_theta, learning_rate, momentum_beta, num_steps)
print(f"Theta after momentum GD: {final_theta_mom}")

4.4. Ускоренный градиент Нестерова (NAG)

\[ \mathbf{v}^{(t+1)} = \beta \mathbf{v}^{(t)} - \eta \nabla J(\boldsymbol{\theta}^{(t)} + \beta \mathbf{v}^{(t)}) \] \[ \boldsymbol{\theta}^{(t+1)} = \boldsymbol{\theta}^{(t)} + \mathbf{v}^{(t+1)} \]

Объяснение: NAG улучшает метод моментов, вычисляя градиент в "предсказанной" точке (lookahead position), что приводит к более точным обновлениям.

Пример: При \(\beta = 0.9\) NAG предвосхищает будущее направление, уменьшая перескоки в колебательных сценариях.

Реализация:


import numpy as np

def nesterov_gradient_descent(gradient, theta, eta, beta, steps):
  """Ускоренный градиент Нестерова."""
  theta = np.copy(theta)
  v = np.zeros_like(theta)
  for _ in range(steps):
    # Вычисляем градиент в предсказанной точке
    lookahead_theta = theta + beta * v
    grad = gradient(lookahead_theta)
    # Обновляем скорость
    v = beta * v - eta * grad
    # Обновляем параметры
    theta += v
  return theta

# Пример использования:
def cost_gradient(theta): return 2 * theta
initial_theta = np.array([10.0])
learning_rate = 0.1
nesterov_beta = 0.9
num_steps = 50

final_theta_nag = nesterov_gradient_descent(cost_gradient, initial_theta, learning_rate, nesterov_beta, num_steps)
print(f"Theta after NAG: {final_theta_nag}")

4.5. RMSProp

\[ s^{(t+1)} = \beta s^{(t)} + (1 - \beta)[\nabla J(\boldsymbol{\theta}^{(t)})]^2 \] \[ \boldsymbol{\theta}^{(t+1)} = \boldsymbol{\theta}^{(t)} - \frac{\eta}{\sqrt{s^{(t+1)}} + \epsilon} \nabla J(\boldsymbol{\theta}^{(t)}) \] (Квадрат градиента и корень из s понимаются поэлементно. \(\epsilon\) - малое число для стабильности.)

Объяснение: RMSProp масштабирует скорость обучения для каждого параметра на основе скользящего среднего квадратов градиентов, улучшая сходимость для невыпуклых задач.

Пример: При \(\beta = 0.9\) RMSProp адаптирует размер шага для каждого параметра, стабилизируя обновления.

Реализация:


import numpy as np

def rmsprop(gradient, theta, eta, beta, epsilon, steps):
  """RMSProp оптимизация."""
  theta = np.copy(theta)
  s = np.zeros_like(theta) # Скользящее среднее квадратов градиентов
  for _ in range(steps):
    grad = gradient(theta)
    s = beta * s + (1 - beta) * grad**2
    theta -= eta / (np.sqrt(s) + epsilon) * grad
  return theta

# Пример использования:
def cost_gradient(theta): return 2 * theta
initial_theta = np.array([10.0])
learning_rate = 0.1 # Может потребоваться другая скорость обучения для RMSProp
beta_rms = 0.9
epsilon_rms = 1e-8
num_steps = 100 # Может потребоваться больше шагов

final_theta_rms = rmsprop(cost_gradient, initial_theta, learning_rate, beta_rms, epsilon_rms, num_steps)
print(f"Theta after RMSProp: {final_theta_rms}")

4.6. Оптимизация Adam

\[ m^{(t+1)} = \beta_1 m^{(t)} + (1 - \beta_1)\nabla J(\boldsymbol{\theta}^{(t)}) \] \[ s^{(t+1)} = \beta_2 s^{(t)} + (1 - \beta_2)[\nabla J(\boldsymbol{\theta}^{(t)})]^2 \] \[ \hat{m} = \frac{m^{(t+1)}}{1 - \beta_1^t}, \quad \hat{s} = \frac{s^{(t+1)}}{1 - \beta_2^t} \] \[ \boldsymbol{\theta}^{(t+1)} = \boldsymbol{\theta}^{(t)} - \frac{\eta}{\sqrt{\hat{s}} + \epsilon} \hat{m} \] (Операции с градиентом и \(s\) поэлементные. \(t\) - номер итерации, начиная с 1.)

Объяснение: Adam объединяет идеи моментов и RMSProp, адаптируя размеры шагов и сглаживая обновления. Это один из самых популярных алгоритмов оптимизации в МО.

Пример: При \(\beta_1 = 0.9, \beta_2 = 0.999\) Adam балансирует момент и масштабирование для каждого параметра.

Реализация:


import numpy as np

def adam(gradient, theta, eta, beta1, beta2, epsilon, steps):
  """Adam оптимизация."""
  theta = np.copy(theta)
  m = np.zeros_like(theta) # Первый момент (момент)
  s = np.zeros_like(theta) # Второй момент (RMSProp)
  t = 0
  for _ in range(steps):
    t += 1 # Счетчик итераций (начинается с 1 для коррекции смещения)
    grad = gradient(theta)
    m = beta1 * m + (1 - beta1) * grad
    s = beta2 * s + (1 - beta2) * grad**2
    # Коррекция смещения
    m_hat = m / (1 - beta1**t)
    s_hat = s / (1 - beta2**t)
    # Обновление параметров
    theta -= eta / (np.sqrt(s_hat) + epsilon) * m_hat
  return theta

# Пример использования:
def cost_gradient(theta): return 2 * theta
initial_theta = np.array([10.0])
learning_rate = 0.1 # Adam часто менее чувствителен к нач. скорости обучения
beta1_adam = 0.9
beta2_adam = 0.999
epsilon_adam = 1e-8
num_steps = 100

final_theta_adam = adam(cost_gradient, initial_theta, learning_rate, beta1_adam, beta2_adam, epsilon_adam, num_steps)
print(f"Theta after Adam: {final_theta_adam}")

4.7. Регуляризованная целевая функция оптимизации

\[ J_{\text{reg}}(\boldsymbol{\theta}) = J(\boldsymbol{\theta}) + \lambda R(\boldsymbol{\theta}) \]

Объяснение: Регуляризация штрафует сложность модели для предотвращения переобучения. Распространенные регуляризаторы включают нормы L1 (лассо) и L2 (гребневая регрессия).

Пример: Для \(R(\boldsymbol{\theta}) = \|\boldsymbol{\theta}\|_2^2 = \sum \theta_i^2\) (L2 регуляризация), \(J_{\text{reg}}(\boldsymbol{\theta}) = J(\boldsymbol{\theta}) + \lambda \|\boldsymbol{\theta}\|_2^2\).

Реализация:


# loss - функция основной потери J(theta)
# theta - параметры
# reg - функция регуляризатора R(theta)
# lam - коэффициент регуляризации lambda
def regularized_objective(loss, theta, reg, lam):
  """Вычисляет регуляризованную целевую функцию."""
  return loss(theta) + lam * reg(theta)

# Пример использования:
def mse_loss(theta): # Пример основной функции потерь
    # Предположим, theta - это один параметр, а целевое значение y=0
    y_true = 0
    y_pred = theta # Простейшая модель
    return (y_pred - y_true)**2

def l2_regularizer(theta): # L2 регуляризатор
    return np.sum(theta**2)

theta_val = np.array([2.0])
lambda_reg = 0.1

total_loss = regularized_objective(mse_loss, theta_val, l2_regularizer, lambda_reg)
# J = (2-0)^2 = 4. R = 2^2 = 4. J_reg = 4 + 0.1 * 4 = 4.4
print(f"Regularized loss: {total_loss}")

4.8. Затухание скорости обучения (Learning Rate Decay)

\[ \eta_t = \frac{\eta_0}{1 + \gamma t} \quad \text{(пример: обратное временное затухание)} \] (Другие варианты: экспоненциальное \(\eta_t = \eta_0 e^{-\gamma t}\), ступенчатое и т.д.)

Объяснение: Затухание скорости обучения постепенно уменьшает скорость обучения, чтобы улучшить стабильность сходимости по мере продвижения обучения.

Пример: Для \(\eta_0 = 0.1\), \(\gamma = 0.01\), на шаге \(t = 10\), \(\eta_{10} = 0.1 / (1 + 0.01 \cdot 10) = 0.1 / 1.1 \approx 0.0909\).

Реализация (обратное временное затухание):


def learning_rate_decay(eta0, gamma, t):
  """Вычисляет затухающую скорость обучения (обратное время)."""
  return eta0 / (1 + gamma * t)

# Пример
initial_lr = 0.1
decay_rate = 0.01
step = 10
current_lr = learning_rate_decay(initial_lr, decay_rate, step)
print(f"Learning rate at step {step}: {current_lr}")

4.9. Ограничение градиента (Gradient Clipping)

\[ \mathbf{g} = \text{clip}(\mathbf{g}, -\tau, \tau) \quad \text{(по значению)} \] \[ \mathbf{g} = \frac{\tau}{\|\mathbf{g}\|} \mathbf{g} \quad \text{if } \|\mathbf{g}\| > \tau \quad \text{(по норме)} \]

Объяснение: Ограничение градиента лимитирует величину градиента для предотвращения "взрыва градиентов" в глубоких нейронных сетях.

Пример: При \(\tau = 1.0\) ограничить градиенты диапазоном [-1, 1] (по значению) или ограничить норму вектора градиента значением 1.0 (по норме).

Реализация (ограничение по значению):


import numpy as np

# grad - вектор градиента
# tau - порог ограничения
def gradient_clipping_value(grad, tau):
  """Ограничивает градиент по значению."""
  return np.clip(grad, -tau, tau)

# Реализация (ограничение по норме):
def gradient_clipping_norm(grad, tau):
    """Ограничивает градиент по норме."""
    norm = np.linalg.norm(grad)
    if norm > tau:
        return grad * (tau / norm)
    else:
        return grad

# Пример
gradient_vector = np.array([0.5, -1.5, 2.0])
clip_threshold = 1.0

clipped_by_value = gradient_clipping_value(gradient_vector, clip_threshold)
clipped_by_norm = gradient_clipping_norm(gradient_vector, clip_threshold)

print(f"Original gradient: {gradient_vector}")
print(f"Clipped by value (-{clip_threshold}, {clip_threshold}): {clipped_by_value}") # [ 0.5, -1. ,  1. ]
print(f"Clipped by norm ({clip_threshold}): {clipped_by_norm}") # Вектор с нормой 1.0

4.10. Мини-батчевый градиентный спуск (Minibatch Gradient Descent)

\[ \boldsymbol{\theta}^{(t+1)} = \boldsymbol{\theta}^{(t)} - \eta \nabla J_{B_t}(\boldsymbol{\theta}^{(t)}) \] (где \(J_{B_t}\) - функция стоимости, вычисленная по мини-батчу \(B_t\))

Объяснение: Мини-батчевый градиентный спуск вычисляет обновления, используя небольшие случайные подмножества данных (мини-батчи), балансируя шум SGD и стабильность батчевого градиентного спуска.

Пример: Использовать размер мини-батча \(B = 32\) для вычисления обновлений по меньшим подмножествам данных.

Реализация:


import numpy as np
import random

# gradient - функция, вычисляющая градиент J по theta для батча данных
# theta - начальные параметры
# eta - скорость обучения
# data - полный набор данных
# batch_size - размер мини-батча
# steps - количество шагов/обновлений
def minibatch_gradient_descent(gradient, theta, eta, data, batch_size, steps):
  """Мини-батчевый градиентный спуск."""
  theta = np.copy(theta)
  num_samples = len(data)
  indices = np.arange(num_samples)

  for _ in range(steps):
    # Перемешиваем данные и берем батч (или выбираем случайно без перемешивания)
    np.random.shuffle(indices)
    batch_indices = indices[:batch_size]
    # В реальных библиотеках часто используется data[batch_indices]
    # Здесь предполагаем, что data - список, и gradient работает с подсписком
    batch_data = [data[i] for i in batch_indices]

    # Вычисляем градиент по батчу
    grad = gradient(theta, batch_data)
    theta -= eta * grad
  return theta

# Пример использования (очень упрощенный):
# J_B(theta) = 1/|B| sum_{i in B} (theta - y_i)^2 => grad = 2/|B| sum_{i in B} (theta - y_i)
def simple_minibatch_gradient(theta, batch_y):
    return 2 * np.mean(theta - batch_y)

y_values = np.array([1.0, 1.5, 3.0, 2.5, 1.8, 2.2])
initial_theta = np.array([0.0])
learning_rate = 0.05
batch_sz = 2
num_steps = 500

final_theta_mini = minibatch_gradient_descent(simple_minibatch_gradient, initial_theta, learning_rate, y_values, batch_sz, num_steps)
print(f"Theta after Minibatch GD: {final_theta_mini}") # Должен быть близок к среднему y_values (~2.0)
print(f"Actual mean: {np.mean(y_values)}")

4.11. Координатный спуск

\[ \theta_j^{(t+1)} = \theta_j^{(t)} - \eta \frac{\partial J(\boldsymbol{\theta})}{\partial \theta_j} \quad (\text{обновление для одной координаты } j) \] (Часто \(\eta\) адаптируется или используется точная минимизация по координате)

Объяснение: Координатный спуск оптимизирует один параметр за раз, циклически проходя по всем параметрам до сходимости. Эффективен для задач высокой размерности.

Пример: Минимизировать \(J(\theta_1, \theta_2) = (\theta_1 - 1)^2 + (\theta_2 - 2)^2\) путем поочередного обновления \(\theta_1\) (при фиксированной \(\theta_2\)) и \(\theta_2\) (при фиксированной \(\theta_1\)).

Реализация:


import numpy as np

# gradient - функция, возвращающая частную производную dJ/d(theta_j)
# theta - начальные параметры
# eta - скорость обучения (может не использоваться, если решается точно)
# steps - количество полных циклов по координатам
def coordinate_descent(gradient_partial, theta, eta, steps):
    """Координатный спуск."""
    theta = np.copy(theta)
    num_params = len(theta)
    for _ in range(steps):
        for j in range(num_params):
            # Вычисляем частную производную по j-й координате
            grad_j = gradient_partial(theta, j)
            # Обновляем j-ю координату
            theta[j] -= eta * grad_j
            # Примечание: в некоторых вариантах eta=0 и решается одномерная
            # задача минимизации J по theta_j при фиксированных остальных theta_k
    return theta

# Пример использования: J(t1, t2) = (t1-1)^2 + (t2-2)^2
# dJ/dt1 = 2*(t1-1), dJ/dt2 = 2*(t2-2)
def partial_grad_example(theta, j):
    if j == 0: # dJ/d(theta_1)
        return 2 * (theta[0] - 1)
    elif j == 1: # dJ/d(theta_2)
        return 2 * (theta[1] - 2)
    else:
        return 0

initial_theta = np.array([5.0, -5.0])
learning_rate = 0.1
num_cycles = 100

final_theta_coord = coordinate_descent(partial_grad_example, initial_theta, learning_rate, num_cycles)
print(f"Theta after Coordinate Descent: {final_theta_coord}") # Должен быть близок к [1.0, 2.0]

4.12. Регуляризация Elastic Net

\[ J_{\text{reg}}(\boldsymbol{\theta}) = J(\boldsymbol{\theta}) + \lambda_1 \|\boldsymbol{\theta}\|_1 + \lambda_2 \|\boldsymbol{\theta}\|_2^2 \] (где \(\|\boldsymbol{\theta}\|_1 = \sum |\theta_i|\) - L1 норма, \(\|\boldsymbol{\theta}\|_2^2 = \sum \theta_i^2\) - квадрат L2 нормы)

Объяснение: Elastic Net сочетает L1 и L2 регуляризацию для обработки разреженности и мультиколлинеарности. Часто используется в задачах регрессии.

Пример: Для \(\lambda_1 = 0.1\), \(\lambda_2 = 0.2\) и \(J(\boldsymbol{\theta}) = \|\mathbf{y} - X\boldsymbol{\theta}\|_2^2\), вычислить регуляризованную целевую функцию.

Реализация (вычисление значения функции):


import numpy as np

# loss - функция основной потери J(theta)
# theta - параметры
# lam1 - коэффициент L1 регуляризации
# lam2 - коэффициент L2 регуляризации
def elastic_net_objective(loss, theta, lam1, lam2):
  """Вычисляет целевую функцию Elastic Net."""
  l1_term = np.sum(np.abs(theta))
  l2_term = np.sum(theta**2) # Иногда используется 0.5 * lam2 * sum(theta**2)
  return loss(theta) + lam1 * l1_term + lam2 * l2_term

# Пример использования (с MSE из предыдущего примера):
def mse_loss_vector(theta):
    # Пример: предсказание y=0 вектором theta
    y_true = np.zeros_like(theta)
    y_pred = theta # Модель y=theta
    return np.mean((y_pred - y_true)**2) # Используем MSE

theta_val = np.array([2.0, -1.0])
lambda1 = 0.1
lambda2 = 0.2

total_loss_en = elastic_net_objective(mse_loss_vector, theta_val, lambda1, lambda2)
# loss = mean((2-0)^2 + (-1-0)^2) = mean(4+1)=2.5
# l1 = |2|+|-1| = 3
# l2 = 2^2 + (-1)^2 = 5
# J_reg = 2.5 + 0.1 * 3 + 0.2 * 5 = 2.5 + 0.3 + 1.0 = 3.8
print(f"Elastic Net objective value: {total_loss_en}")

4.13. Оптимизация Adagrad

\[ G^{(t)}_{jj} = \sum_{i=1}^{t} [\nabla J(\boldsymbol{\theta}^{(i)})_j]^2 \] \[ \theta_j^{(t+1)} = \theta_j^{(t)} - \frac{\eta}{\sqrt{G^{(t)}_{jj}} + \epsilon} \nabla J(\boldsymbol{\theta}^{(t)})_j \] (Обновление для j-го параметра. \(G\) - диагональная матрица сумм квадратов градиентов.)

Объяснение: Adagrad адаптирует скорость обучения для каждого параметра на основе истории градиентов, улучшая производительность на разреженных данных.

Пример: При \(\eta = 0.1\) адаптивно масштабировать обновления для различных признаков.

Реализация:


import numpy as np

def adagrad(gradient, theta, eta, epsilon, steps):
  """Adagrad оптимизация."""
  theta = np.copy(theta)
  G = np.zeros_like(theta) # Сумма квадратов градиентов по каждому параметру
  for _ in range(steps):
    grad = gradient(theta)
    G += grad**2
    theta -= eta / (np.sqrt(G) + epsilon) * grad
  return theta

# Пример использования:
def cost_gradient(theta): return 2 * theta
initial_theta = np.array([10.0])
learning_rate = 0.5 # Adagrad может требовать большей начальной eta
epsilon_adagrad = 1e-8
num_steps = 100

final_theta_adagrad = adagrad(cost_gradient, initial_theta, learning_rate, epsilon_adagrad, num_steps)
print(f"Theta after Adagrad: {final_theta_adagrad}")

4.14. Оптимизация AdamW

\[ \boldsymbol{\theta}^{(t+1)}_{\text{Adam}} = \boldsymbol{\theta}^{(t)} - \text{AdamStep}(\nabla J(\boldsymbol{\theta}^{(t)})) \] \[ \boldsymbol{\theta}^{(t+1)} = \boldsymbol{\theta}^{(t+1)}_{\text{Adam}} - \eta \lambda \boldsymbol{\theta}^{(t)} \] (AdamW разделяет обновление Adam и затухание весов \(\lambda\).)

Объяснение: AdamW модифицирует Adam, отделяя затухание весов (L2 регуляризацию) от обновлений градиента, улучшая регуляризацию и обобщение в моделях МО.

Пример: При \(\lambda = 0.01\) регуляризовать веса наряду с адаптивными скоростями обучения.

Реализация:


import numpy as np

# gradient, theta, eta, beta1, beta2, epsilon, steps - как в Adam
# lam - коэффициент затухания весов (weight decay)
def adamw(gradient, theta, eta, beta1, beta2, lam, epsilon, steps):
  """AdamW оптимизация."""
  theta = np.copy(theta)
  m = np.zeros_like(theta)
  s = np.zeros_like(theta)
  t = 0
  for _ in range(steps):
    t += 1
    grad = gradient(theta)
    m = beta1 * m + (1 - beta1) * grad
    s = beta2 * s + (1 - beta2) * grad**2
    m_hat = m / (1 - beta1**t)
    s_hat = s / (1 - beta2**t)
    # Обновление Adam + Затухание весов
    adam_update = eta / (np.sqrt(s_hat) + epsilon) * m_hat
    weight_decay = eta * lam * theta # Применяем затухание к текущим весам
    theta -= (adam_update + weight_decay)
    # Примечание: иногда затухание применяется иначе: theta = theta * (1 - eta * lam) - adam_update
  return theta

# Пример использования:
def cost_gradient(theta): return 2 * theta
initial_theta = np.array([10.0])
learning_rate = 0.1
beta1_adamw = 0.9
beta2_adamw = 0.999
lambda_adamw = 0.01 # Коэффициент затухания
epsilon_adamw = 1e-8
num_steps = 100

final_theta_adamw = adamw(cost_gradient, initial_theta, learning_rate, beta1_adamw, beta2_adamw, lambda_adamw, epsilon_adamw, num_steps)
print(f"Theta after AdamW: {final_theta_adamw}")

4.15. Метод моментов "Тяжелый шарик" (Momentum "Heavy Ball")

\[ \boldsymbol{\theta}^{(t+1)} = \boldsymbol{\theta}^{(t)} + \beta(\boldsymbol{\theta}^{(t)} - \boldsymbol{\theta}^{(t-1)}) - \eta \nabla J(\boldsymbol{\theta}^{(t)}) \]

Объяснение: Этот вариант метода моментов включает инерционный член для улучшения скорости сходимости для сильно выпуклых задач.

Пример: При \(\beta = 0.9\) "тяжелый шарик" ускоряет градиентный спуск.

Реализация:


import numpy as np

def heavy_ball(gradient, theta, eta, beta, steps):
  """Метод моментов 'Тяжелый шарик'."""
  theta = np.copy(theta)
  prev_theta = theta.copy() # Theta на шаге t-1
  # Или инициализировать иначе, например, theta_minus_1 = theta - eta * gradient(theta)
  theta_history = [prev_theta, theta] # Храним два последних значения

  for _ in range(steps):
    grad = gradient(theta_history[-1]) # Градиент в текущей точке theta(t)
    # Вычисляем следующее theta
    current_theta = theta_history[-1]
    prev_theta = theta_history[-2]
    next_theta = current_theta + beta * (current_theta - prev_theta) - eta * grad

    # Обновляем историю
    theta_history.append(next_theta)
    # Оставляем только два последних значения
    if len(theta_history) > 2:
        theta_history.pop(0)

  return theta_history[-1] # Возвращаем последнее вычисленное значение

# Пример использования:
def cost_gradient(theta): return 2 * theta
initial_theta = np.array([10.0])
learning_rate = 0.1
heavy_ball_beta = 0.9 # Коэффициент инерции
num_steps = 50

final_theta_hb = heavy_ball(cost_gradient, initial_theta, learning_rate, heavy_ball_beta, num_steps)
print(f"Theta after Heavy Ball: {final_theta_hb}")

4.16. Проекция / Проекционный градиентный спуск

\[ \boldsymbol{\theta}^{(t+1)} = \text{Proj}_C (\boldsymbol{\theta}^{(t)} - \eta \nabla J(\boldsymbol{\theta}^{(t)})) \] (где \(\text{Proj}_C(\mathbf{z})\) - проекция точки \(\mathbf{z}\) на допустимое множество C)

Объяснение: Проекционный градиентный спуск гарантирует, что обновления остаются в пределах допустимого множества C, часто используется для оптимизации с ограничениями.

Пример: Для \(C = \{\boldsymbol{\theta} : \|\boldsymbol{\theta}\|_2 \le 1\}\), проецировать \(\boldsymbol{\theta}\) на единичный шар после каждого шага градиентного спуска.

Реализация:


import numpy as np

# gradient - функция градиента
# theta - начальные параметры
# eta - скорость обучения
# projection - функция, выполняющая проекцию на множество C
# steps - количество шагов
def projected_gradient_descent(gradient, theta, eta, projection, steps):
  """Проекционный градиентный спуск."""
  theta = np.copy(theta)
  for _ in range(steps):
    # Шаг градиентного спуска
    grad = gradient(theta)
    theta_temp = theta - eta * grad
    # Шаг проекции
    theta = projection(theta_temp)
  return theta

# Пример проекции на L2-шар радиуса R=1
def project_onto_l2_ball(v, radius=1.0):
    norm = np.linalg.norm(v)
    if norm > radius:
        return v * (radius / norm)
    else:
        return v

# Пример использования:
def cost_gradient(theta): return 2 * theta # Пример градиента
initial_theta = np.array([10.0, -5.0])
learning_rate = 0.1
num_steps = 50
radius_limit = 1.0 # Ограничение на норму

final_theta_proj = projected_gradient_descent(cost_gradient, initial_theta, learning_rate,
                                              lambda v: project_onto_l2_ball(v, radius_limit),
                                              num_steps)
print(f"Theta after Projected GD (L2 norm <= {radius_limit}): {final_theta_proj}")
print(f"Norm of final theta: {np.linalg.norm(final_theta_proj)}") # Должна быть <= radius_limit

4.17. Метод Ньютона

\[ \boldsymbol{\theta}^{(t+1)} = \boldsymbol{\theta}^{(t)} - [H(\boldsymbol{\theta}^{(t)})]^{-1} \nabla J(\boldsymbol{\theta}^{(t)}) \] (где \(H(\boldsymbol{\theta}^{(t)})\) - матрица Гессе функции J в точке \(\boldsymbol{\theta}^{(t)}\))

Объяснение: Метод Ньютона использует информацию второго порядка через Гессиан для улучшения сходимости, особенно для квадратичных функций стоимости.

Пример: Для \(J(\theta) = \theta^2\). \(\nabla J = 2\theta\), \(H = 2\). Обновление: \(\theta^{(t+1)} = \theta^{(t)} - [2]^{-1} (2\theta^{(t)}) = \theta^{(t)} - \frac{1}{2} (2\theta^{(t)}) = \theta^{(t)} - \theta^{(t)} = 0\). Сходимость за один шаг для квадратичной функции.

Реализация:


import numpy as np

# gradient - функция градиента
# hessian - функция, вычисляющая матрицу Гессе
# theta - начальные параметры
# steps - количество шагов
# eta - скорость обучения (часто eta=1 для чистого метода Ньютона)
def newtons_method(gradient, hessian, theta, steps, eta=1.0):
  """Метод Ньютона."""
  theta = np.copy(theta)
  for _ in range(steps):
    grad = gradient(theta)
    hess = hessian(theta)
    try:
      # Решаем систему H * delta = -grad вместо вычисления обратной матрицы
      delta = np.linalg.solve(hess, -grad)
      # Или вычисляем обратную матрицу (менее стабильно и эффективно)
      # hess_inv = np.linalg.inv(hess)
      # delta = -np.dot(hess_inv, grad)
    except np.linalg.LinAlgError:
      print("Предупреждение: Гессиан вырожден. Переход к шагу градиентного спуска.")
      # Можно откатиться к градиентному спуску или использовать регуляризацию
      delta = -eta * grad # Простой откат

    # Обновляем параметры
    theta += eta * delta # eta используется для демпфированного метода Ньютона
  return theta

# Пример использования: J(t1, t2) = (t1-1)^2 + (t2-2)^2
def grad_example(theta): return np.array([2*(theta[0]-1), 2*(theta[1]-2)])
def hess_example(theta): return np.array([[2.0, 0.0], [0.0, 2.0]])

initial_theta = np.array([5.0, -5.0])
num_steps = 1 # Достаточно одного шага для этой функции

final_theta_newton = newtons_method(grad_example, hess_example, initial_theta, num_steps)
print(f"Theta after Newton's method: {final_theta_newton}") # Должен быть [1.0, 2.0]

4.18. Проксимальный градиентный метод

\[ \boldsymbol{\theta}^{(t+1)} = \text{prox}_{\eta R}(\boldsymbol{\theta}^{(t)} - \eta \nabla J(\boldsymbol{\theta}^{(t)})) \] (где \(J\) - гладкая часть, \(R\) - возможно, негладкая часть (регуляризатор), \(\text{prox}_{\eta R}(\mathbf{z}) = \arg\min_{\mathbf{x}} \left( R(\mathbf{x}) + \frac{1}{2\eta}\|\mathbf{x} - \mathbf{z}\|_2^2 \right)\) - проксимальный оператор)

Объяснение: Проксимальный градиентный метод обобщает градиентный спуск для работы с негладкими членами регуляризации, такими как L1 норма.

Пример: Для \(R(\boldsymbol{\theta}) = \lambda \|\boldsymbol{\theta}\|_1\), \(\text{prox}_{\eta R}\) вычисляет операцию мягкого порога (soft thresholding) для каждого параметра.

Реализация:


import numpy as np

# gradient - градиент гладкой части J
# theta - начальные параметры
# eta - скорость обучения
# prox - функция проксимального оператора prox_{eta*R}
# steps - количество шагов
def proximal_gradient(gradient, theta, eta, prox, steps):
  """Проксимальный градиентный метод."""
  theta = np.copy(theta)
  for _ in range(steps):
    # Шаг градиентного спуска по гладкой части
    grad = gradient(theta)
    theta_temp = theta - eta * grad
    # Шаг проксимального оператора для негладкой части
    theta = prox(theta_temp, eta) # Передаем eta, если prox зависит от нее
  return theta

# Пример prox для L1: Soft Thresholding
def soft_thresholding(x, threshold):
    return np.sign(x) * np.maximum(0, np.abs(x) - threshold)

# prox_{eta * lambda * ||.||_1} (z)
def prox_l1(z, eta, lam):
    return soft_thresholding(z, eta * lam)

# Пример использования: J(theta) = 0.5 * ||theta - a||^2, R(theta) = lam * ||theta||_1
def grad_smooth(theta, a): return theta - a
target_a = np.array([1.5, -0.5, 0.0, 2.0, -1.2])
lambda_l1 = 0.3
initial_theta = np.zeros_like(target_a)
learning_rate = 0.1
num_steps = 100

# Создаем функцию prox с фиксированным lambda
prox_operator = lambda z, eta: prox_l1(z, eta, lambda_l1)
# Создаем функцию градиента с фиксированным a
smooth_grad = lambda theta: grad_smooth(theta, target_a)

final_theta_prox = proximal_gradient(smooth_grad, initial_theta, learning_rate, prox_operator, num_steps)
print(f"Theta after Proximal Gradient (L1): {final_theta_prox}")
# Результат должен быть разреженным и близким к target_a с поправкой на L1

4.19. Проксимальный градиент с L1 (ISTA)

\[ \boldsymbol{\theta}^{(t+1)} = \text{soft}(\boldsymbol{\theta}^{(t)} - \eta \nabla J(\boldsymbol{\theta}^{(t)}), \eta \lambda) \] (где \(\text{soft}(\mathbf{z}, \tau)_i = \text{sign}(z_i) \max(0, |z_i| - \tau)\) - оператор мягкого порога)

Объяснение: Итеративный алгоритм мягкого порога (Iterative Shrinkage-Thresholding Algorithm - ISTA) применяет операцию мягкого порога для обновления параметров в задачах разреженной оптимизации (например, Лассо).

Пример: Для \(J(\boldsymbol{\theta}) = \frac{1}{2}\|\mathbf{y} - X\boldsymbol{\theta}\|_2^2\) и \(R(\boldsymbol{\theta}) = \lambda\|\boldsymbol{\theta}\|_1\), применить мягкий порог к каждому \(\theta_i\).

Реализация:


import numpy as np

# gradient - градиент гладкой части J
# theta - начальные параметры
# eta - скорость обучения
# lam - коэффициент L1 регуляризации lambda
# steps - количество шагов
def ista(gradient, theta, eta, lam, steps):
  """ISTA - Итеративный алгоритм мягкого порога."""
  theta = np.copy(theta)

  def soft_thresholding(x, threshold):
      return np.sign(x) * np.maximum(0, np.abs(x) - threshold)

  for _ in range(steps):
    grad = gradient(theta)
    theta_temp = theta - eta * grad
    theta = soft_thresholding(theta_temp, eta * lam) # Применяем мягкий порог
  return theta

# Пример использования: Лассо регрессия
# J(theta) = 0.5 * ||y - X*theta||^2 => grad J = X^T (X*theta - y)
np.random.seed(0)
n_samples, n_features = 50, 100
X = np.random.randn(n_samples, n_features)
true_coef = np.zeros(n_features)
true_coef[0:10] = np.random.randn(10) # Только первые 10 признаков важны
y = X @ true_coef + np.random.randn(n_samples) * 0.1

def lasso_gradient(theta, X_data, y_data):
    return X_data.T @ (X_data @ theta - y_data)

initial_theta = np.zeros(n_features)
# Скорость обучения для ISTA часто выбирается как 1/L, где L - константа Липшица градиента
# Для лассо L = ||X^T X||_2 (спектральная норма)
L = np.linalg.norm(X.T @ X, ord=2)
learning_rate = 1.0 / L
lambda_lasso = 0.1 # Коэффициент регуляризации
num_steps = 200

# Создаем функцию градиента с фиксированными X, y
grad_func = lambda theta: lasso_gradient(theta, X, y)

final_theta_ista = ista(grad_func, initial_theta, learning_rate, lambda_lasso, num_steps)
print(f"Number of non-zero coeffs in ISTA: {np.sum(np.abs(final_theta_ista) > 1e-6)}")
print(f"First 15 coeffs: {final_theta_ista[:15]}")
# Должны быть ненулевые значения среди первых 10 и нули/малые значения дальше

4.20. Метод штрафов

\[ J_{\text{penalty}}(\boldsymbol{\theta}, \mu) = J(\boldsymbol{\theta}) + \mu \sum_i [h_i(\boldsymbol{\theta})]^2 \quad (+ \mu \sum_j [\max(0, g_j(\boldsymbol{\theta}))]^2) \] (для ограничений \(h_i(\boldsymbol{\theta})=0\) и \(g_j(\boldsymbol{\theta})\le 0\). \(\mu > 0\) - параметр штрафа.) (Формула в OCR \(J(\theta) + \frac{1}{\mu} h(\theta)^2\) немного отличается нотацией/масштабированием.)

Объяснение: Метод штрафов решает задачи оптимизации с ограничениями путем добавления штрафа за нарушение ограничений к целевой функции.

Пример: Для ограничения \(h(\theta) = \|\theta\|_2^2 - 1 = 0\), штрафовать отклонения от единичной сферы. \(J_{\text{penalty}} = J(\theta) + \mu (\|\theta\|_2^2 - 1)^2\).

Реализация (вычисление значения):


# loss - основная функция потерь J(theta)
# theta - параметры
# penalty_constraints - список функций-ограничений [h1(theta), h2(theta), ...]
# mu - коэффициент штрафа
def penalty_method_objective(loss, theta, penalty_constraints, mu):
    """Вычисляет целевую функцию метода штрафов."""
    penalty_term = 0
    for h in penalty_constraints:
        penalty_term += h(theta)**2
    # Добавить сюда обработку неравенств g(theta)<=0, если нужно
    # penalty_term += np.sum(np.maximum(0, g(theta))**2)
    return loss(theta) + mu * penalty_term

# Пример использования: минимизировать J=theta^2 при ограничении h=theta-1=0
def loss_simple(theta): return theta[0]**2
def constraint_h(theta): return theta[0] - 1.0

theta_val = np.array([1.5])
mu_penalty = 10.0 # Коэффициент штрафа

objective_val = penalty_method_objective(loss_simple, theta_val, [constraint_h], mu_penalty)
# J = 1.5^2 = 2.25. h = 1.5 - 1 = 0.5. J_pen = 2.25 + 10 * (0.5)^2 = 2.25 + 10 * 0.25 = 2.25 + 2.5 = 4.75
print(f"Penalty objective value: {objective_val}")

4.21. Метод дополненного Лагранжиана

\[ L_A(\boldsymbol{\theta}, \boldsymbol{\lambda}, \mu) = J(\boldsymbol{\theta}) + \sum_i \lambda_i h_i(\boldsymbol{\theta}) + \frac{\mu}{2} \sum_i [h_i(\boldsymbol{\theta})]^2 \] (для ограничений равенства \(h_i(\boldsymbol{\theta})=0\). \(\lambda_i\) - множители Лагранжа.)

Объяснение: Метод дополненного Лагранжиана сочетает подходы Лагранжиана и штрафов для решения задач оптимизации с ограничениями. Он чередует обновление параметров и множителей Лагранжа.

Пример: Для \(J(\boldsymbol{\theta}) = \|\boldsymbol{\theta}\|_2^2\) и \(h(\boldsymbol{\theta}) = \|\boldsymbol{\theta}\|_1 - 1 = 0\), вычислять обновления для \(\boldsymbol{\theta}\), \(\lambda\) и \(\mu\).

Реализация (один шаг обновления):


import numpy as np
from scipy.optimize import minimize # Используем для минимизации по theta

# loss - функция J(theta)
# constraints_h - список функций ограничений hi(theta)=0
# theta - текущие параметры
# lam - текущие множители Лагранжа (lambda)
# mu - текущий параметр штрафа
def augmented_lagrangian_step(loss, constraints_h, theta, lam, mu):
    """Один шаг метода дополненного Лагранжиана (обновление theta)."""

    # Определяем функцию дополненного Лагранжиана для минимизации по theta
    def objective(t):
        cost = loss(t)
        penalty = 0
        lagrange = 0
        for i in range(len(constraints_h)):
            h_i = constraints_h[i](t)
            lagrange += lam[i] * h_i
            penalty += (mu / 2.0) * (h_i**2)
        return cost + lagrange + penalty

    # Определяем градиент (если можем, иначе scipy будет его оценивать)
    # def objective_grad(t):
    #     grad_loss = ... # Градиент J(t)
    #     grad_penalty = 0
    #     grad_lagrange = 0
    #     for i in range(len(constraints_h)):
    #         h_i = constraints_h[i](t)
    #         grad_h_i = ... # Градиент hi(t)
    #         grad_lagrange += lam[i] * grad_h_i
    #         grad_penalty += mu * h_i * grad_h_i
    #     return grad_loss + grad_lagrange + grad_penalty

    # Минимизируем по theta
    result = minimize(objective, theta, method='L-BFGS-B') # Пример метода
    new_theta = result.x

    # Обновление множителей Лагранжа (lambda)
    new_lam = np.copy(lam)
    for i in range(len(constraints_h)):
        new_lam[i] += mu * constraints_h[i](new_theta)

    # Обновление mu (опционально, часто увеличивают)
    new_mu = mu # или mu * growth_factor

    return new_theta, new_lam, new_mu

# Пример использования (упрощенный): J=theta^2, h=theta-1=0
def loss_simple(theta): return theta[0]**2
def constraint_h(theta): return theta[0] - 1.0

theta_init = np.array([1.5])
lam_init = np.array([0.0]) # Нач. множитель Лагранжа
mu_init = 1.0

# Выполняем один шаг
theta_1, lam_1, mu_1 = augmented_lagrangian_step(loss_simple, [constraint_h], theta_init, lam_init, mu_init)

print(f"Theta after 1 step: {theta_1}")
print(f"Lambda after 1 step: {lam_1}")
# Для сходимости нужно много шагов и правильное обновление mu

4.22. Метод двойственного подъема (Dual Ascent)

\[ \boldsymbol{\theta}^{(t+1)} = \arg \min_{\boldsymbol{\theta}} L(\boldsymbol{\theta}, \boldsymbol{\lambda}^{(t)}) \] \[ \boldsymbol{\lambda}^{(t+1)} = \boldsymbol{\lambda}^{(t)} + \eta \, h(\boldsymbol{\theta}^{(t+1)}) \quad (\text{для } h(\boldsymbol{\theta})=0) \] (где \(L(\boldsymbol{\theta}, \boldsymbol{\lambda}) = J(\boldsymbol{\theta}) + \boldsymbol{\lambda}^T h(\boldsymbol{\theta})\) - Лагранжиан)

Объяснение: Метод двойственного подъема оптимизирует двойственную задачу для оптимизации с ограничениями путем итеративного обновления множителей Лагранжа.

Пример: Для ограничения \(h(\boldsymbol{\theta}) = \|\boldsymbol{\theta}\|_1 - 1 = 0\), обновлять \(\lambda\) на основе нарушения ограничения.

Реализация (один шаг):


import numpy as np
from scipy.optimize import minimize

# loss - функция J(theta)
# constraint_h - функция ограничения h(theta)=0 (векторная или скалярная)
# theta - текущие параметры
# lam - текущий множитель Лагранжа
# eta - скорость обучения для lambda
def dual_ascent_step(loss, constraint_h, theta, lam, eta_dual):
    """Один шаг метода двойственного подъема."""

    # Определяем Лагранжиан
    def lagrangian(t, current_lambda):
        # Убедимся, что constraint_h возвращает массив numpy
        h_val = np.asarray(constraint_h(t))
        # Используем скалярное произведение для lambda^T h(theta)
        return loss(t) + np.dot(current_lambda, h_val)

    # Минимизируем Лагранжиан по theta при фиксированном lambda
    result = minimize(lambda t: lagrangian(t, lam), theta, method='L-BFGS-B')
    new_theta = result.x

    # Обновляем lambda (двойственный подъем)
    h_at_new_theta = np.asarray(constraint_h(new_theta))
    new_lam = lam + eta_dual * h_at_new_theta

    return new_theta, new_lam

# Пример: J=0.5*theta^2, h=theta-1=0
def loss_qp(theta): return 0.5 * theta[0]**2
def constraint_qp(theta): return np.array([theta[0] - 1.0])

theta_init = np.array([1.5])
lam_init = np.array([0.0]) # Нач. множитель Лагранжа
dual_learning_rate = 0.1

# Выполняем один шаг
theta_1, lam_1 = dual_ascent_step(loss_qp, constraint_qp, theta_init, lam_init, dual_learning_rate)

print(f"Theta after 1 dual ascent step: {theta_1}") # Должен стремиться к 1
print(f"Lambda after 1 dual ascent step: {lam_1}") # Должен стремиться к -1 (т.к. grad J + lambda grad h = 0 => theta + lambda*1 = 0 => 1 + lambda = 0)

4.23. Метод доверительной области (Trust Region Method)

\[ \boldsymbol{\theta}^{(t+1)} = \boldsymbol{\theta}^{(t)} + \Delta^{(t)} \] \[ \Delta^{(t)} = \arg \min_{\|\Delta\| \le \Delta_{\max}} \{ m(\Delta) = J(\boldsymbol{\theta}^{(t)}) + \nabla J(\boldsymbol{\theta}^{(t)})^T \Delta + \frac{1}{2} \Delta^T H(\boldsymbol{\theta}^{(t)}) \Delta \} \] (Минимизация квадратичной модели \(m(\Delta)\) внутри шара - доверительной области)

Объяснение: Метод доверительной области ограничивает размер шага областью, где квадратичная аппроксимация функции стоимости считается действительной, обеспечивая стабильность.

Пример: Для \(J(\boldsymbol{\theta}) = \|\mathbf{y} - X\boldsymbol{\theta}\|_2^2\), вычислять шаги \(\Delta\), ограниченные \(\|\Delta\| \le \Delta_{\max}\).

Реализация (упрощенный пример шага):

Решение подзадачи доверительной области является нетривиальным. Часто используются итерационные методы, такие как метод сопряженных градиентов (Steihaug) или метод 'dogleg'. Ниже приведен очень упрощенный вариант, который просто масштабирует шаг Ньютона, если он выходит за пределы области (это не совсем корректный метод доверительной области, а скорее его эмуляция).


import numpy as np

# loss - функция потерь
# gradient - функция градиента
# hessian - функция гессиана
# theta - текущие параметры
# delta_max - максимальный радиус доверительной области
# steps - количество шагов
def trust_region_simplified(gradient, hessian, theta, delta_max, steps):
    """Упрощенная эмуляция метода доверительной области."""
    theta = np.copy(theta)
    for _ in range(steps):
        grad = gradient(theta)
        hess = hessian(theta)

        try:
            # Попытка найти шаг Ньютона (решение H*delta = -grad)
            delta = np.linalg.solve(hess, -grad)
        except np.linalg.LinAlgError:
            # Если Гессиан плохой, используем шаг градиента
            print("Warning: Hessian issue, using gradient step.")
            delta = -grad # Нужна адаптивная eta здесь
            # Подбираем eta так, чтобы шаг был в пределах delta_max
            grad_norm = np.linalg.norm(grad)
            if grad_norm > 0:
                delta = - (delta_max / grad_norm) * grad

        # Проверяем норму шага
        delta_norm = np.linalg.norm(delta)
        if delta_norm > delta_max:
            # Масштабируем шаг, чтобы он поместился в доверительную область
            # (Это не совсем верно для TR, но просто для примера)
            delta = delta * (delta_max / delta_norm)

        # Обновляем параметры
        # В настоящем TR методе нужно еще проверить,
        # уменьшилась ли функция на самом деле (сравнить с предсказанием модели m(delta))
        # и адаптировать delta_max.
        theta += delta
    return theta

# Пример использования: J(t1, t2) = (t1-1)^2 + (t2-2)^2
def grad_example(theta): return np.array([2*(theta[0]-1), 2*(theta[1]-2)])
def hess_example(theta): return np.array([[2.0, 0.0], [0.0, 2.0]])

initial_theta = np.array([5.0, -5.0])
max_radius = 0.5 # Маленькая доверительная область
num_steps = 50

final_theta_tr = trust_region_simplified(grad_example, hess_example, initial_theta, max_radius, num_steps)
print(f"Theta after simplified Trust Region: {final_theta_tr}") # Будет сходиться медленнее из-за малого радиуса

4.24. Метод барьеров (Внутренней точки)

\[ \min_{\boldsymbol{\theta}} J_{\text{barrier}}(\boldsymbol{\theta}) = J(\boldsymbol{\theta}) - \mu \sum_i \ln(-g_i(\boldsymbol{\theta})) \] (для ограничений-неравенств \(g_i(\boldsymbol{\theta}) \le 0\). \(\mu > 0\) - параметр барьера, уменьшающийся к 0.)

Объяснение: Метод барьеров решает задачи оптимизации с ограничениями-неравенствами, добавляя логарифмический барьер, штрафующий приближение к границе допустимой области. Обновления остаются внутри допустимой области.

Пример: Для ограничения \(g(\theta) = \|\theta\|_1 - 1 \le 0\), использовать \(-\mu \ln(1 - \|\theta\|_1)\) как барьерный член.

Реализация (один шаг):


import numpy as np
from scipy.optimize import minimize

# loss - основная функция J(theta)
# constraints_g - список функций ограничений gi(theta) <= 0
# theta - текущие параметры (должны быть строго внутри допустимой области)
# mu - текущий параметр барьера
def barrier_method_step(loss, constraints_g, theta, mu):
    """Один шаг метода барьеров (минимизация барьерной функции)."""

    def barrier_objective(t):
        cost = loss(t)
        barrier_term = 0
        for g in constraints_g:
            g_val = g(t)
            if g_val >= 0: # Точка вне или на границе - барьер бесконечен
                return np.inf
            barrier_term -= mu * np.log(-g_val) # Логарифмический барьер
        return cost + barrier_term

    # Минимизируем барьерную функцию
    # Начальная точка theta должна быть допустимой (все g_i < 0)
    result = minimize(barrier_objective, theta, method='L-BFGS-B')
    new_theta = result.x

    # Уменьшаем mu для следующего шага
    new_mu = mu * 0.9 # Пример фактора уменьшения

    return new_theta, new_mu

# Пример: J=theta^2, g=theta-1 <= 0
def loss_simple(theta): return theta[0]**2
def constraint_g(theta): return theta[0] - 1.0

theta_init = np.array([0.5]) # Начальная точка внутри области theta < 1
mu_init = 1.0

# Выполняем один шаг
theta_1, mu_1 = barrier_method_step(loss_simple, [constraint_g], theta_init, mu_init)

print(f"Theta after 1 barrier step: {theta_1}") # Должен быть ближе к 0, но < 1
print(f"New mu: {mu_1}")
# Для сходимости к истинному решению (theta=0) нужно много шагов с уменьшающимся mu

4.25. Имитация отжига (Simulated Annealing)

\[ P(\text{принять } \Delta E) = \begin{cases} 1, & \text{if } \Delta E < 0 \\ e^{-\Delta E / T}, & \text{if } \Delta E \ge 0 \end{cases} \] (где \(\Delta E\) - изменение энергии/стоимости, T - температура)

Объяснение: Имитация отжига — это вероятностный алгоритм оптимизации, вдохновленный отжигом в металлургии. Он исследует пространство решений, принимая с некоторой вероятностью худшие решения, чтобы избежать локальных минимумов.

Пример: Минимизировать \(J(\theta) = \theta^2\) с начальной температурой \(T = 1\), постепенно ее понижая.

Реализация:


import numpy as np
import math

# loss - функция потерь (энергии)
# theta - начальное состояние (параметры)
# T - начальная температура
# cooling_rate - коэффициент охлаждения (e.g., 0.95-0.99)
# steps - количество шагов
# neighbor_func - функция генерации соседнего состояния (опционально)
def simulated_annealing(loss, theta, T, cooling_rate, steps, neighbor_func=None):
    """Имитация отжига."""
    current_theta = np.copy(theta)
    current_loss = loss(current_theta)
    best_theta = current_theta
    best_loss = current_loss
    current_T = T

    for _ in range(steps):
        # Генерируем соседнее состояние
        if neighbor_func:
            new_theta = neighbor_func(current_theta)
        else:
            # Простой пример: добавляем случайный шум
            noise = np.random.uniform(-0.5, 0.5, size=current_theta.shape) # Масштаб шума можно адаптировать
            new_theta = current_theta + noise

        new_loss = loss(new_theta)
        delta_E = new_loss - current_loss

        # Решаем, принять ли новое состояние
        if delta_E < 0 or np.random.rand() < math.exp(-delta_E / current_T):
            current_theta = new_theta
            current_loss = new_loss

            # Обновляем лучшее найденное решение
            if current_loss < best_loss:
                best_theta = current_theta
                best_loss = current_loss

        # Охлаждаем температуру
        current_T *= cooling_rate

    return best_theta, best_loss


# Пример использования: J(theta) = theta^2
def cost_func(theta): return theta[0]**2

initial_theta = np.array([10.0])
initial_T = 10.0
cool_rate = 0.99
num_steps = 1000

best_theta_sa, min_loss_sa = simulated_annealing(cost_func, initial_theta, initial_T, cool_rate, num_steps)
print(f"Best theta found by SA: {best_theta_sa}") # Должен быть близок к 0
print(f"Minimum loss found by SA: {min_loss_sa}")
Предыдущая глава    Следующая глава

Другие статьи по этой теме: