最小化非线性方程组（指数积分）答案

【问题标题】：Minimize system of nonlinear equation (integral on exponent)最小化非线性方程组（指数积分）
【发布时间】：2019-08-26 10:34:20
【问题描述】：

一般：

我正在使用最大熵来查找正整数向量的分布，我可以估计均值和方差，并且我试图找到三个方程 a 和 b，

方程式：

integral(exp(a*x^2+bx+c) from (0 , infinity))-1
integral(xexp(ax^2+bx+c)from (0 , infinity))-均值
integral(x^2*exp(a*x^2+bx+c) from (0 , infinity))- mean^2 - var

([0,∞)之间的整数)

问题：

我正在尝试使用数值求解器，我使用了 sympy 的 fsolve 但我想我缺少一些知识。

我的代码：

import numpy as np
import sympy as sym
from scipy.optimize import *


def myFunction(x,*data):
    y = sym.symbols('y')
    m,v=data
    F = [0]*3
    x[0] = - abs(x[0])
    print(x)
    F[0] = (sym.integrate(sym.exp(x[0] * y ** 2 + x[1] * y + x[2]), (y, 0,sym.oo)) -1).evalf() 
    F[1] = (sym.integrate(y*sym.exp(x[0] * y ** 2 + x[1] * y + x[2]), (y, 0,sym.oo))-m).evalf()
    F[2] = (sym.integrate((y**2)*sym.exp(x[0] * y ** 2 + x[1] * y + x[2]), (y,0,sym.oo)) -v-m).evalf() 
    print(F)
    return F


data = (10,3.5) # mean and var for example
xGuess = [1, 1, 1]
z = fsolve(myFunction,xGuess,args = data)
print(z)

我的结果不太准确，有没有更好的解决方法？

积分(exp(a*x^2+bx+c))-1 = 5.67659292676884
积分(xexp(ax^2+bx+c))-均值 = -1.32123173796713
积分(x^2*exp(a*x^2+bx+c))-均值^2 - var = −2.20825624606312

谢谢

【问题讨论】：

这是Gaussian 分布吗？高斯公式（使用均值和方差）可以写成 exp(ax2 + bx +c)，例如 a=-1/(2*sigma2) ...跨度>
您好，感谢您的回答，我编辑了积分边界从 0 到无穷大的问题
看看scipy.optimize.minimize。如果您认为您的问题不是凸的，请查看global optimization 部分。

标签： python numpy entropy scipy-optimize

【解决方案1】：

我已经重写了用numpy 和lambdas（内联函数）替换sympy 的问题。另请注意，在您的问题陈述中，您用 $mean^2$ 减去第三个等式，但在您的代码中，您只减去 $mean$。

import numpy as np
from scipy.optimize import minimize
from scipy.integrate import quad


def myFunction(x,data):
    m,v=data
    F = np.zeros(3)  # use numpy array

    # use scipy.integrade.quad for integration of lambda functions
    # quad output is (result, error), so we just select the result value at the end
    F[0] = quad(lambda y: np.exp(x[0] * y ** 2 + x[1] * y + x[2]), 0, np.inf)[0] -1
    F[1] = quad(lambda y: y*np.exp(x[0] * y ** 2 + x[1] * y + x[2]), 0, np.inf)[0] -m
    F[2] = quad(lambda y: (y**2)*np.exp(x[0] * y ** 2 + x[1] * y + x[2]), 0, np.inf)[0] -v-m**2

    # minimize the squared error
    return np.sum(F**2)


data = (10,3.5) # mean and var for example
xGuess = [-1, 1, 1]
z = minimize(lambda x: myFunction(x, data), x0=xGuess,
             bounds=((None, 0), (None, None), (None, None)))  # use bounds for negative first coefficient
print(z)

# x: array([-0.99899311,  2.18819689,  1.85313181])

这看起来更合理吗？

【讨论】：