1
$\begingroup$

I have a curve fitting problem of the form: $$ \textbf{y} = f(\textbf{x}, a,b,c,d) + \varepsilon $$ $$ f(x, a,b,c,d) = \frac{b}{e^{x\cdot a}+c}+d $$ with the constraint \begin{equation} \begin{aligned} \arg\min_{a,b,c,d} \quad & \sum_{i=1}^{n}{(f(\text{x})_i-\text{y}_i)^2}\\ \textrm{s.t.} \quad & -c-\frac{b}{d}-1 = 0\\ \quad & a = 1 \\ \end{aligned} \end{equation}

where

  • $\textbf{x}, \textbf{y} \in \mathbb{R}^n$. $\textbf{x}, \textbf{y}$ are the data given by an experiment
  • $a \in \mathbb{R}$
  • $b \in \mathbb{R}$
  • $c \in \mathbb{R}$
  • $d \in \mathbb{R}$ Note: $a = 1$ because this parameter is being measured by a measurement in an experience and is, therefore fixed.

The simulated data and the following code are as follows:

import cvxpy as cp
import numpy as np
import matplotlib.pyplot as plt
n = 20
np.random.seed(1)
def sigm01(x,a,b,c,d):
  return b/(np.exp(x*a)+c)+d

start = -10
end = -start
x = np.arange(start,end,0.01)
a = 1
b = -2
c = 1
d = 1
y = sigm01(x,a,b,c,d) + np.random.randn(len(x))/10
#plt.plot(x,y)
#plt.show()

I would like to use an optimization package (such as scipy, cvxpy or Convex.jl).

How could I setup the problem?

$\endgroup$
9
  • $\begingroup$ So a,b,c are fixed and you're minimizing with respect to x? $\endgroup$
    – Richard
    Dec 25 '20 at 16:49
  • $\begingroup$ thanks for your comment @Richard. No $x$ is given. ($\textbf{x}$ and $\textbf{y}$ are the data). And I want to find the optimal $a,b,c,d$ to minimize $\sum_{i=1}^{n}{(\text{x}_i-\text{y}_i)^2}$. It is a curve fitting problem $\endgroup$
    – ecjb
    Dec 25 '20 at 16:56
  • $\begingroup$ Could you edit your question to clarify? d should be part of the argmin and the phrase "a=1 in the generated data" seems confusing. $\endgroup$
    – Richard
    Dec 25 '20 at 17:16
  • 1
    $\begingroup$ Is there any reason why you have that weird cost function $||\mathbf{x}-\mathbf{y}||^2$ instead of least square? I mean are you sure you are not looking for $||\mathbf{y}-f(\mathbf{x})||^2$? Cause I have a hard time understanding why you want to minimize the distance between measured $\mathbf{y}$ and your features $\mathbf{x}$. $\endgroup$ Dec 25 '20 at 19:18
  • 1
    $\begingroup$ Those constraints are not basically constraints in a classical definition. One of them is simply a given parameter for $a$ which you can replace $a$ with 1 everywhere. The other one just replaces $c$ with $b$ and $d$ everywhere and then finally you just fit $\mathbf{y}$ to $f(\mathbf{x},b,d)$. Simply the scipy’s example given in their documentation works perfectly fine here. $\endgroup$ Dec 25 '20 at 23:04
2
$\begingroup$

With your constraints of $a = 1$ and $c = -\frac{b}{d} - 1$, your $f$ looks like this:

$$\hat{\mathbf{y}} = f(\mathbf{x},b,d) = \frac{bd}{d\exp{(\mathbf{x})}-b - d} + d$$

You're trying to solve a nonlinear least square problem basically:

$$b,d = \arg \min \sum_{i=1}^{N} (f(x_{i})-y_{i})^{2}$$

It's a job for scipy.optimize.curve_fit:

import numpy as np
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt

def func(x,b,d):
    return ((b*d)/(d*np.exp(x)-b-d))+d

xdata = np.linspace(-10, 10, 10000)
y = func(xdata, -2, 1)
np.random.seed(1729)
y_noise = 0.2 * np.random.normal(size=xdata.size)
ydata = y + y_noise

popt, pcov = curve_fit(func, xdata, ydata, bounds=([-3, 0], [-1, 2.]))

plt.plot(xdata, ydata, 'b-', label='data')
plt.plot(xdata, func(xdata, *popt), 'r-',
         label='fit: b=%5.3f, d=%5.3f' % tuple(popt))
plt.xlabel('x')
plt.ylabel('y')
plt.legend(loc='best')
plt.savefig('curve_fit.png')
plt.show()

The final fit looks like this:

enter image description here

$\endgroup$
2
  • 1
    $\begingroup$ Happy Holidays! I hope things are well with you. $\endgroup$
    – Richard
    Dec 26 '20 at 17:56
  • 1
    $\begingroup$ @Richard Happy holidays and thanks for your kind message! $\endgroup$ Dec 26 '20 at 18:24

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.