Решение нелинейной системы сильно зависит от начального предположения
У меня есть система нелинейных уравнений
- xF - xA * exp(-(wAA * xA + wBB * xB)) = 0
- xF - xB * exp(-(wBB * xB + wAB * xA)) = 0
- хА + хВ + хF -1 = 0
где переменные — это xA, xB и xF, а wAA, wBB и wAB — заданные параметры. Переменные неотрицательны и находятся в интервале [0,1] Вот минимальный пример кода
import math
from scipy.optimize import least_squares
coupleInteractions = [10000., 10., 0.]
def equations(p):
xA, xB, xF = p
wAA, wBB, wAB = coupleInteractions
eq1 = xF - xA * math.exp(-(wAA * xA + wBB * xB))
eq2 = xF - xB * math.exp(-(wBB * xB + wAB * xA))
eq3 = xA + xB + xF -1
return (eq1, eq2, eq3)
guess = [0., 1., 0.]
bounds = ([0., 0., 0.], [1., 1., 1.])
root = least_squares(equations, x0=guess, bounds=bounds, xtol=1.e-15, gtol = 1.e-15, ftol = 1.e-15, loss = "cauchy")
print(root)
print(root.x)
Я хочу протестировать код в некоторых предельных случаях. Например, если все параметры равны нулю, wAA = wAB = wBB = 0, у меня есть простая линейная система с xA = xB = xF = 1,/3. В этом случае код отлично работает вне зависимости от догадок.
Я хочу проверить, как в опубликованном коде, случай wAA >> wBB и wAA >> wAB. С математической точки зрения я получаю xA = 1, xB = 0 и xF = 0. Но с функцией least_square я получаю неправильный ответ. Я намеренно использую "неправильную" догадку для проверки метода, потому что в промежуточных случаях не хочу, чтобы результаты сильно зависели от первоначального выбора. Есть ли ошибка в моем коде?
Я выполнил новый тест, чтобы увидеть, может ли глобальная оптимизация работать. Я определил функцию:
def scalarEquations(p):
xA, xB, xF = p
wAA, wBB, wAB = coupleInteractions
eq1 = xF - xA * math.exp(-(wAA * xA + wBB * xB))
eq2 = xF - xB * math.exp(-(wBB * xB + wAB * xA))
eq3 = xA + xB + xF -1
return (eq1**2 + eq2**2 + eq3 **2)
Он возвращает сумму квадратов для 3 уравнений. Затем я выполнил подход грубой силы.
from scipy.optimize import brute
rranges = (slice(0,1, 0.05), slice(0,1, 0.005), slice(0,1, 0.005))
t = brute(scalarEquations, rranges)
print(t)
Я получаю такое же «неправильное» решение. Я знаю, что глобальный оптимизатор должен найти абсолютный минимум. Я не понимаю.
Ответы 2
Ваш код выглядит нормально. Похоже, что существует локальный минимум около [0,1,0] для wAA big и wBB больше e, и оптимизация сходится к нему.
вы, вероятно, не можете. Если есть локальный минимум, вы можете ожидать, что любой оптимизатор сойдется к нему, когда он будет запущен достаточно близко. (И «близко» может быть довольно далеко)
Хорошо спасибо. Мне интересно, могу ли я использовать глобальный оптимизатор.
Наверное, я нашел решение. Я применил к обеим сторонам логарифмической функции моей системы уравнений. Кроме того, я определяю новую функцию потерь с этой идеей. Я меняю нижнюю границу на 1.e-15, чтобы избежать ошибки домена. Наконец, я выполняю грубую оптимизацию. Вот мой новый кусок кода
def scalarEquations(p):
xA, xB, xF = p
wAA, wBB, wAB = coupleInteractions
eq1 = math.log(xF) - math.log(xA * math.exp(-(wAA * xA + wAB * xB)))
eq2 = math.log(xF) - math.log(xB * math.exp(-(wBB * xB + wAB * xA)))
eq3 = math.log(xA + xB + xF) -math.log(1.)
return (eq1**2 + eq2**2 + eq3 **2)
rranges = ((1.e-15, 1.), (1.e-15, 1.), (1.e-15, 1.) )
t = brute(scalarEquations, rranges, Ns = 100)
print(t)
Как этого избежать?