Clever Geek Handbook

Minimize objective function using limfit.minimize in Python

I have a problem with the package minification procedure lmfit.minimize

. In fact, I could not create the correct target function for my problem.

Defining the problem

  • My function:, yn = a_11*x1**2 + a_12*x2**2 + ... + a_m*xn**2

    where xn

    - unknowns, a_m

    are coefficients. n = 1..N, m = 1..M
  • In my case, N=5

    for x1,..,x5

    and M=3

    for y1, y2, y3

    .

I need to find optimum: x1, x2,...,x5

so he can satisfyy

My question is:

  • Error: ValueError: operands could not be broadcast together with shapes (3,) (3,5)

    .
  • Have I created the object function correctly for my problem in Python?

My code:

import numpy as np
from lmfit import Parameters, minimize

def func(x,a):
    return np.dot(a, x**2)

def residual(pars, a, y):
    vals = pars.valuesdict()
    x = vals['x']
    model = func(x,a)
    return y - model

def main():
    # simple one: a(M,N) = a(3,5)
    a = np.array([ [ 0, 0, 1, 1, 1 ],
                   [ 1, 0, 1, 0, 1 ],
                   [ 0, 1, 0, 1, 0 ] ])
    # true values of x
    x_true = np.array([10, 13, 5, 8, 40])
    # data without noise
    y = func(x_true,a)

    #************************************
    # Apriori x0
    x0 = np.array([2, 3, 1, 4, 20])
    fit_params = Parameters()
    fit_params.add('x', value=x0)

    out = minimize(residual, fit_params, args=(a, y))
    print out
if __name__ == '__main__':
    main()
+3


source to share


1 answer


Using directly scipy.optimize.minimize()

, the code below solves this problem. Note that with more points yn

you will get the same result as x_true

, otherwise there will be multiple solutions. You can minimize the impact of poor optimization by adding borders (see the parameter bounds

used below).

import numpy as np
from scipy.optimize import minimize

def residual(x, a, y):
    s = ((y - a.dot(x**2))**2).sum()
    return s

def main():
    M = 3
    N = 5
    a = np.random.random((M, N))
    x_true = np.array([10, 13, 5, 8, 40])
    y = a.dot(x_true**2)

    x0 = np.array([2, 3, 1, 4, 20])
    bounds = [[0, None] for x in x0]
    out = minimize(residual, x0=x0, args=(a, y), method='L-BFGS-B', bounds=bounds)
    print(out.x)


If M>=N

you can also use scipy.optimize.leastsq

for this task:

import numpy as np
from scipy.optimize import leastsq

def residual(x, a, y):
    return y - a.dot(x**2)

def main():
    M = 5
    N = 5
    a = np.random.random((M, N))
    x_true = np.array([10, 13, 5, 8, 40])
    y = a.dot(x_true**2)

    x0 = np.array([2, 3, 1, 4, 20])
    out = leastsq(residual, x0=x0, args=(a, y))
    print(out[0])
+3


source






More articles:


All Articles