0

I want to use scipy.optimize.minimize on a function with an array type variable input. This is what I hope to do.

I have the following signal,

import numpy as np
time = np.linspace(0, 1, 501)
data = np.cos(2 * np.pi * 4 * time) + np.cos(2 * np.pi * 9 * time) + np.cos(2 * np.pi * 20 * time)
noise = np.sqrt(1 / 25) * np.random.randn(501)
signal = data + noise 

and I am hoping to find a curve fit for this signal. Since I created the data myself, I know that a sum of cosine functions will work for this. So the function that I hope to optimize is the following:

def cos_sum(x, P):

    assert isinstance(P, np.ndarray)
    assert P.shape[0] == P.shape[1]

    sums = []
    for param in P:
        a, b, c = param
        sums.append(a * np.cos(b * (x - c)))
    sums = np.array(sums)

    return np.sum(sums, axis=0) 

In order to use minimize to find the correct parameters, I create this residual function.

def resid(params, x):
    assert isinstance(params, np.ndarray)
    fit = cos_sum(x, params)
    residual = np.sqrt(np.mean(np.abs(fit - signal)) ** 2)
    return residual

I now need to create a guess. Since I already know how the signal was created, I made the following guess:

guess_A = np.random.normal(1, .2, size=3)
guess_B = 2 * np.pi * np.array([4, 9, 20], dtype=float)
guess_C = np.random.normal(0, .2, size=3)
guess = np.array([guess_A, guess_B, guess_C]).T

However, I am unable to run the following:

from scipy.optimize import minimize
optimization = minimize(resid, guess, args=(time))

and I receive the following error:

Traceback (most recent call last):
  File "/Users/nickeisenberg/GitRepos/Python_Misc/Misc/minimize_curvefit_vector_variables.py", line 70, in <module>
    optimization = minimize(resid, guess, args=(time))
  File "/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/scipy/optimize/_minimize.py", line 676, in minimize
    res = _minimize_bfgs(fun, x0, args, jac, callback, **options)
  File "/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/scipy/optimize/_optimize.py", line 1296, in _minimize_bfgs
    sf = _prepare_scalar_function(fun, x0, jac, args=args, epsilon=eps,
  File "/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/scipy/optimize/_optimize.py", line 263, in _prepare_scalar_function
    sf = ScalarFunction(fun, x0, args, grad, hess,
  File "/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/scipy/optimize/_differentiable_functions.py", line 158, in __init__
    self._update_fun()
  File "/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/scipy/optimize/_differentiable_functions.py", line 251, in _update_fun
    self._update_fun_impl()
  File "/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/scipy/optimize/_differentiable_functions.py", line 155, in update_fun
    self.f = fun_wrapped(self.x)
  File "/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/scipy/optimize/_differentiable_functions.py", line 137, in fun_wrapped
    fx = fun(np.copy(x), *args)
  File "/Users/nickeisenberg/GitRepos/Python_Misc/Misc/minimize_curvefit_vector_variables.py", line 53, in resid
    fit = cos_sum(x, params)
  File "/Users/nickeisenberg/GitRepos/Python_Misc/Misc/minimize_curvefit_vector_variables.py", line 30, in cos_sum
    assert P.shape[0] == P.shape[1]
IndexError: tuple index out of range

Is this possible to do?

1
  • What do you mean "I am unable"? Do you get some error? If so add it to the question. Commented Oct 17, 2022 at 16:05

1 Answer 1

1

The problem is that minimize treats parameters as 1D array. Adding print just before failing assert shows, that it reshaped guess array to 1D array.

Changing cos_sum so that it accepts 1D array of params fixes this problem. Here is the code.

def cos_sum(x, P):
    assert isinstance(P, np.ndarray)

    sums = []
    for i in range(0, P.shape[0], 3):
        a, b, c = P[i], P[i+1], P[i+2]
        sums.append(a * np.cos(b * (x - c)))
    sums = np.array(sums)

    return np.sum(sums, axis=0) 

Result:

optimization = minimize(resid, guess, args=(time))
print(optimization.x)
[  0.96805816  25.1679919    0.25020317   1.00261543  56.44511497
   0.77872223   1.00966167 125.71622787   0.55004217]
plt.figure(figsize=(8, 4))
plt.plot(data)
plt.plot([cos_sum(t, optimization.x) for t in time], 'C1--')
plt.show()

enter image description here

Sign up to request clarification or add additional context in comments.

Comments

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.