如何解决LMFIT:使用多项式模型时约束输出
我正在使用 LMFIT 将分段多项式拟合到正弦波的第一象限。 我希望能够对多项式输出添加约束 - 而不是对其参数。 例如,我想确保输出 >= 0 和
我知道使用 np.polyfit 可能会更好,但最终我想添加更多非线性约束,并且 LMFIT 框架更加灵活。
import numpy as np
from lmfit.models import LinearModel
#split sine wave in 4 segments with 1024 points
nseg = 4
frac = 2**10
npoints = nseg*frac
xfrac = np.linspace(0,1,num=frac,endpoint=False)
x = np.linspace(0,num=npoints,endpoint=False)
y = np.sin(x*np.pi/2)
yseg = np.reshape(y,(nseg,frac))
mod = LinearModel()
coeff = []
bestfit = []
for i in range(nseg):
pars = mod.guess(yseg[i],x=xfrac)
out = mod.fit(yseg[i],pars,x=xfrac)
coeff.append([out.best_values['slope'],out.best_values['intercept']])
bestfit.append(out.best_fit)
bestfit = np.reshape(bestfit,(1,npoints))[0]
解决方法
事实证明这是通过在参数本身上添加约束来完成的,这些约束变成了对模型输出的正确约束。 使用自定义模型进行线性插值,可以按如下方式完成:
def func(x,c0,c1):
return c0 + c1*x
pmodel = Model(func)
params = Parameters()
params.add('c0')
params.add('clip',value=0,max=1.0,vary=True)
params.add('c1',expr='clip-c0')
,
一种选择可能是使用样条曲线。 一个快速而肮脏的方法,只是为了展示这个想法,可能看起来像这样:
import matplotlib.pyplot as plt
import numpy as np
## quich and dirty spline function
def l_spline(x,abc ):
if isinstance( x,( list,tuple,np.ndarray ) ):
out = [ l_spline( elem,abc ) for elem in x]
else:
a,b,c = abc
if x < a:
f = lambda t: 0
elif x < b:
f = lambda t: ( t - a ) / ( b - a )
elif x < c:
f = lambda t: -( t - c ) / (c - b )
else:
f = lambda t: 0
out = f(x)
return out
### test data
xl = np.linspace( 0,4,150 )
sl = np.fromiter( ( np.sin( elem ) for elem in xl ),np.float )
### test splines with manual double knots on first and last
yl = dict()
yl[0] = l_spline( xl,( 0,.4 ) )
for i in range(1,10 ):
yl[i] = l_spline( xl,( (i - 1 ) * 0.4,i * 0.4,(i + 1 ) * 0.4 ) )
yl[10] = l_spline( xl,( 3.6,4 ) )
## This is the most simple linear least square for the coefficients
AT = list()
for i in range( 11 ):
AT.append( yl[i] )
AT = np.array( AT )
A = np.transpose( AT )
U = np.dot( AT,A )
UI = np.linalg.inv( U )
K = np.dot( UI,AT )
v = np.dot( K,sl )
## adding up the weigthed sum
out = np.zeros( len( sl ) )
for a,l in zip( v,AT ):
out += a * l
### plotting
fig = plt.figure()
ax = fig.add_subplot( 1,1,1 )
ax.plot( xl,sl,ls=':' )
for i in range( 11 ):
ax.plot( xl,yl[i] )
ax.plot( xl,out,color='k')
plt.show()
代替简单的线性优化,可以使用更复杂的函数来确保没有参数大于 1。这会自动确保函数不会超过 1。可以通过设置相应的 b 样条来建立不动点为固定值,即不拟合其参数。
版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。