0

I am trying to make a curvfit using the lmfit library. My goal is to have the curve start at x=0, with f(0) = "a constant". I also want it's minima to be at x=1. However, when putting my constraints into the Model, f(1) > f(0).

F(x) = Aexp(Bx) + Cx2 + Dx + E
F’(x) = ABexp(B*x) + 2*Cx + D
F’’(x) = AB2 exp(Bx) +2C

F(0) = constant = A + E => E = constant - A
F’(1) = 0 = ABexp(B) + 2*C + D => D = - (2*C + BAexp(B))
F’’(1) > 0 => A*B2exp(B) +2*C > 0 => C > -A B2 *exp(B)/2

This gives us the following constrained equation:
F(x) = Cx(x-2) + A*(exp(Bx) - Bexp(B*x) -1) + constant

in addition to the following inequality which needs to be respected:
F’’(1) > 0 => A*B**2 *exp(B) + 2*C > 0 => C > -B2AeB/2

Image of the calculations in word for easier reading

def model_5(x, a, b, c):
    return  c*x*(x-2) + a*(np.exp(b*x) - x*b*np.exp(b) - 1) + 0.003591687375349475

model_5b = Model(model_5)

pars_5b = Parameters()
pars_5b.add('a', value = 0.007)
pars_5b.add('b', value = 0.05)
pars_5b.add('delta', value = 0.02, min=0, vary=True)
pars_5b.add('c', expr="delta - 0.5*a*(b**2)*exp(b) ")

model_5b.fit(y, x=X, params=pars_5b, method="trf")

# Output from above fit
Variables
name    value   initial value   min max vary    expression
a   9.9227e-05  0.007   -inf    inf True    
b   0.08377349  0.05    -inf    inf True    
delta   5.7477e-11  0.02    0.00000000  inf True    
c   -3.7855e-07 0.019990801377906712    -inf    inf False   delta - 0.5*a*(b**2)*exp(b)

# Create predictions with the fitted models values
pred_5b = model_5(X, 9.9227e-05, 0.08377349, -3.7855e-07)

# Check if f(1) equals the minima
pred_5b[1] == np.min(pred_5b) # returns False, should return true

When putting my constraints into the Model, f(1) > f(0), f(1) is not a minima. However, it should be when including the constraints. I therefore feel like I must be inputting the constraints wrongly. Can someone jelp guide me in the right direction?

The difference here, to other questions I have found on Stackoverflow is that the constraint includes an inequality with the parameters that are being optimized. That is C > -B2AeB/2

Mose Mehari Gebreselassie
  • 1
  • 1
  • 1
  • Possible almost duplicate of [How to fit a polynomial with some of the coefficients constrained?](https://stackoverflow.com/questions/48469889/how-to-fit-a-polynomial-with-some-of-the-coefficients-constrained) – pault Oct 28 '19 at 23:27
  • hanks for your reply pault. The difference here is that the constraint includes an inequality with the parameters that are being optimized. That is C > -B**2 *A*exp(B)/2 . In the question you found the constraint is expressed by numbers. – Mose Mehari Gebreselassie Oct 29 '19 at 00:12

1 Answers1

0

I think the constraint math is right (awesome!), and that the problem is at:

# Check if f(1) equals the minima
pred_5b[1] == np.min(pred_5b) # returns False, should return true

This is really only correct if x[1] = 1. That is, pred_5b[1] is the second element of pred_5b, not (necessarily) the value of pred_5b at x=1.

M Newville
  • 7,486
  • 2
  • 16
  • 29
  • Thanks for your reply M Newville. I can confirm that x[1] == 1. X is basically a vector of ages from 0-105. So x[z] = z, for all z's. Could I supply any additional information to make it easier for you to help out? I have added the math derivations and functions as an image in order to make it more readable.. – Mose Mehari Gebreselassie Oct 29 '19 at 08:54
  • Oh, well your constraints on the 1st and 2nd derivative really only ensure that x=1 will be a local minimum, not a global one. I think it will always work for x>=0, but for some x<0, the slope will change sign. – M Newville Oct 30 '19 at 01:11