Функция минимизации scipy.optimize не дает правильного ответа

Я пытаюсь решить проблему минимизации, используя функцию минимизации Scipy. Целевая функция - это просто отношение двух многомерных нормальных распределений с разными средними значениями и дисперсией. Я надеюсь найти максимум функции g_func, что эквивалентно нахождению минимума функции g_optimization. Также я добавил ограничение x [0] = 0. Здесь x - вектор с 8 элементами. Целевая функция g_optimization следующая:

import numpy as np
from scipy.optimize import minimize
# Set up mean and variance for two MVN distributions
n_trait = 8
sigma = np.full((n_trait, n_trait),0.0005)
np.fill_diagonal(sigma,0.005)
omega = np.full((n_trait, n_trait),0.0000236)
np.fill_diagonal(omega,0.0486)
sigma_pos = np.linalg.inv(np.linalg.inv(sigma)+np.linalg.inv(omega))
mu_pos = np.array([-0.01288244,0.08732091,0.01049617,0.0860966,0.10055626,0.07952922,0.04363669,-0.0061975])
mu_pri = 0
sigma_pri = omega

#objective function
def g_func(beta,mu_sim_pos):
    g1 = ((np.linalg.det(sigma_pri))**(1/2))/((np.linalg.det(sigma_pos))**(1/2))
    g2 = (-1/2)*np.linalg.multi_dot([np.transpose(beta-mu_sim_pos),np.linalg.inv(sigma_pos),beta-mu_sim_pos])
    g3 = (1/2)*np.linalg.multi_dot([np.transpose(beta-mu_pri),np.linalg.inv(sigma_pri),beta-mu_pri])
    g = g1*np.exp(g2+g3)
    return g
def g_optimization(beta,mu_sim_pos):
    return -1*g_func(beta,mu_sim_pos)

#optimization
start_point = np.full(8,0)
cons = ({'type': 'eq',
         'fun' : lambda x: np.array([x[0]])})

anws = minimize (g_optimization, [start_point], args=(mu_pos),
                 constraints=cons, options={'maxiter': 50}, tol=0.001)
anws

Оптимизация останавливается после двух итераций, и минимальное значение, которое дает функция, равно 0 в точке np.array ([0,10.32837891, -1.62396508,10.13790152,12.38752653,9.11615259,3.53201544, -4.22115517]). Это не может быть правдой, потому что даже если мы подключим начальную точку np.zeros (8) к функции g_optimization, получим результат -657.0041125829354, что меньше 0. Таким образом, предлагаемое решение определенно не является минимальным.

g_optimization(np.zeros(8),mu_pos) #gives solution of -657.0041125829354

Я не уверен, где я ошибся.


person GarfieldH    schedule 21.08.2020    source источник


Ответы (1)


Я бы попробовал другой решатель. Например L-BFGS-B работает хорошо. Вы можете просмотреть все варианты здесь.

anws = minimize (g_optimization, [start_point], args=(mu_pos), method='L-BFGS-B',
                 constraints=cons, options={'maxiter': 50}, tol=0.001)

print(anws)

# success: True
# message: b'CONVERGENCE: REL_REDUCTION_OF_F_<=_FACTR*EPSMCH'
# fun: -21688.00879938617
# x: array([-0.0101048, 0.09937778, 0.01543875, 0.0980401, 0.11383878, 0.09086455, 0.05164822, -0.00280081])

РЕДАКТИРОВАТЬ: L-BFGS-B не может обрабатывать общие ограничения h(x)=0, только ограничивающие рамки для переменных:

Границы переменных для методов L-BFGS-B, TNC, SLSQP, Powell и trust-constr. Есть два способа указать границы: Класс Instance of Bounds. Последовательность пар (min, max) для каждого элемента в x. Нет используется для указания отсутствия границы.

В вашем случае вы должны определить 8 пар нижних и верхних пределов. Для x [0] вы должны сделать жесткую привязку, поскольку метод не может обрабатывать x_low == x_high.

bounds = [(None, None)] * 8
bounds[0] = (0, 0.00001)
anws = minimize (g_optimization, [start_point], args=(mu_pos), method='L-BFGS-B', bounds=bounds,
                 options={'maxiter': 50}, tol=0.001)
# fun: -21467.48153792194
# x: array([0., 0.10039832, 0.01641271, 0.0990599, 0.11486735, 0.09188037, 0.05264228, -0.00183697])

Другой вариант - исключить значение x [0] из задачи оптимизации:

def g_optimization(beta,mu_sim_pos):
    beta2 = np.empty(8)
    beta2[0] = 0
    beta2[1:] = beta
    return -1*g_func(beta2, mu_sim_pos)

start_point = np.zeros(7)  # exclude x[0]

anws = minimize(g_optimization, [start_point], args=(mu_pos), method='L-BFGS-B',
                 options={'maxiter': 50}, tol=0.001)
# fun: -21467.47686079844
# x: array([0.10041797, 0.01648995, 0.09908046, 0.11487707, 0.09190585, 0.05269467, -0.00174722])
#          ^ missing x[0]
person scleronomic    schedule 21.08.2020
comment
Спасибо за вашу помощь. Однако я хочу добавить к этой задаче ограничение x [0] = 0. Решение, предоставляемое L-BFGS-B, по-видимому, не удовлетворяет этому ограничению. Вы что-нибудь об этом знаете? Спасибо еще раз. - person GarfieldH; 21.08.2020
comment
Кажется, что я могу использовать связанный аргумент для 'L-BFGS-B' вместо ограничения, но я впервые добавляю этот аргумент, и я не совсем уверен, как добавить только связанный с x [0]. Не могли бы вы мне подсказать? Спасибо еще раз. - person GarfieldH; 21.08.2020
comment
Еще раз спасибо. Это очень помогает! - person GarfieldH; 21.08.2020