Scipy curve fit (optimization) - vectorizing a conditional to identify threshold using a custom function

user_4383027 0
Viewed 203

I'm trying to use scipy curve_fit to capture the value of a0 parameter. As of now, it is not changing (always comes out as 1):

X = [[1,2,3],[4,5,6]]
def func(X, a0, c):
     x1 = X[0]; x2 = X[1]
     a = x1*x2
     result = np.where(a(a<a0), -c*(a + np.sqrt(x2)), -c*x1)
     return result

Popt, Cov = scipy.curve_fit(func, X, y) 
a0, c = Popt
Predicted = func(X, a0, c) # a0 and c are constants

I get the values for c, which is a scalar, without any problem. I can't explain why a0 (also a scalar) is always 1, and I am not sure how to fix it. I did see elsewhere on SO that np.where can be used the way I have used it here, but apparently not for curve_fit function. Maybe I need to use a different method of optimization, and I'd like some pointers to do this using scipy methods.

Edit: I tried the construct suggested by Brad, but that's not it.

scikit-learn curve-fitting scipy optimization python
1 Answers

Updated!

This should work. note that the a variable is a vector in this example of length 3 because it is computed by the element wise multiplication of the first and second elements of X which is a 2x3 matrix. Therefore a0 can either be a scalar or a vector of length 3 and c can also be a scalar or a vector of length 3.

import numpy as np


X = np.array([[1, 2, 3], [4, 5, 6]])

a0 = np.array([8,25,400])
# a0 = 2


# Code works whether C is scalar or a matrix since it can be broadcast to matrix a below.
# c = 3 # Uncomment this for scalar
c = np.array([8, 12, 2000])  # Element wise


def func(X, a0, c):
    x = X[0]
    y = X[1]
    a = x * y
    print(a.shape)
    result = np.where(a < a0, c * (a + np.sqrt(y)), c * x)
    return result


func(X, a0, c)

This is a minimum amount of code that works. Notice I removed the y>0 and defined a to be the same size as c. Now you get the correct insertions because the first parameter of np.where is now the same size as the second and third parameters. Before (x<a) & (y>0) always evaluated to True or False and that is a scalar in this context. If a was a N dimensional array you would have received a ValueError because the operands could not be broadcast together

import numpy as np


c = np.array([[22,34],[33,480]])


def func(X, a):
     x = X[0]; y = X[1]
     return np.where(c[(x<a)], -c*(a + np.sqrt(y)), -c*x)


X = [25, 600]
a = np.array([[2,14],[33,22]])

func(X,a)

This also works if c is a constant and a was the array you wanted manipulated

import numpy as np


c = 2


def func(X, a):
     x = X[0]; y = X[1]
     return np.where(a[(x<a)], -c*(a + np.sqrt(y)), -c*x)


X = [25, 600]
a = np.array([[2,14],[33,22]])

func(X,a)
user_87551020
Related