Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Welcome To Ask or Share your Answers For Others

Categories

0 votes
898 views
in Technique[技术] by (71.8m points)

numpy - Fitting with constraints on derivative Python

While trying to create an optimization algorithm, I had to put constraints on the curve fitting of my set.

Here is my problem, I have an array :

Z = [10.3, 10, 10.2, ...]
L = [0, 20, 40, ...]

I need to find a function that fits Z with condition on slope which is the derivative of the function I'm looking for.

Suppose f is my function, f should fit Z and have a condition on f its derivative, it shouldnt exceed a special value.

Are there any libraries in python that can help me achieve this task ?

See Question&Answers more detail:os

与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome To Ask or Share your Answers For Others

1 Reply

0 votes
by (71.8m points)

The COBYLA minimzer can handle such problems. In the following example a polynomial of degree 3 is fitted with the constraint that the derivative is positive everywhere.

from matplotlib import pylab as plt

import numpy as np
from scipy.optimize import minimize

def func(x, pars):
    a,b,c,d=pars
    return a*x**3+b*x**2+c*x+d

x = np.linspace(-4,9,60)
y = func(x, (.3,-1.8,1,2))
y += np.random.normal(size=60, scale=4.0)

def resid(pars):
    return ((y-func(x,pars))**2).sum()

def constr(pars):
    return np.gradient(func(x,pars))

con1 = {'type': 'ineq', 'fun': constr}
res = minimize(resid, [.3,-1,1,1], method='cobyla', options={'maxiter':50000}, constraints=con1)
print res

f=plt.figure(figsize=(10,4))
ax1 = f.add_subplot(121)
ax2 = f.add_subplot(122)

ax1.plot(x,y,'ro',label='data')
ax1.plot(x,func(x,res.x),label='fit')
ax1.legend(loc=0) 
ax2.plot(x,constr(res.x),label='slope')
ax2.legend(loc=0)
plt.show()

sample data and fit


与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
OGeek|极客中国-欢迎来到极客的世界,一个免费开放的程序员编程交流平台!开放,进步,分享!让技术改变生活,让极客改变未来! Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Click Here to Ask a Question

...