如何在 python中安装一个阶跃函数

我有一个关于使用像 curve_fit 这样的 scipy 例程拟合阶跃函数的问题。

我有一个关于使用像 curve_fit 这样的 scipy 例程拟合阶跃函数的问题。

import numpy as np
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt
xobs=np.linspace(0,10,100)
yl=np.random.rand(50); yr=np.random.rand(50)+100
yobs=np.concatenate((yl,yr),axis=0)
def model(x,rf,T1,T2):
#1:   x=np.vectorize(x)
    if x<rf:
        ret= T1
    else:
        ret= T2
    return ret
#2: model=np.vectorize(model)
popt, pcov = curve_fit(model, xobs, yobs, [40.,0.,100.])

上面写着

ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()

如果我添加 # 1 或 # 2 它运行,但并不真正适合数据:

OptimizeWarning: Covariance of the parameters could not be estimated     category=OptimizeWarning)
[ 40.          50.51182064  50.51182064] [[ inf  inf  inf]
[ inf  inf  inf]
[ inf  inf  inf]]

有人知道如何解决这个问题吗?THX

2

下面是我所做的,我保留了xobsyobs

import numpy as np
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt
xobs=np.linspace(0,10,100)
yl=np.random.rand(50); yr=np.random.rand(50)+100
yobs=np.concatenate((yl,yr),axis=0)

现在,必须生成 Heaviside 函数。为了给您这个函数的概述,请考虑 Heaviside 函数的半最大值约定:

Heaviside

在 Python 中,这相当于:def f(x): return 0.5 * (np.sign(x) + 1)

一个示例图将是:

xval = sorted(np.concatenate([np.linspace(-5,5,100),[0]])) # includes x = 0
yval = f(xval)
plt.plot(xval,yval,'ko-')
plt.ylim(-0.1,1.1)
plt.xlabel('x',size=18)
plt.ylabel('H(x)',size=20)

Heaviside2

现在,绘制xobsyobs给出:

plt.plot(xobs,yobs,'ko-')
plt.ylim(-10,110)
plt.xlabel('xobs',size=18)
plt.ylabel('yobs',size=20)

xobs_yobs

请注意,比较两个数字,第二个图移动了 5 个单位,最大值从 1.0 增加到 100。我推断第二个图的函数可以表示如下:

Hfit

或在 Python 中:(0.5 * (np.sign(x-5) + 1) * 100 = 50 * (np.sign(x-5) + 1)

结合这些图可以得到(其中Fit表示上述拟合函数)

Heaviside_combined

现在,假设您不知道这个正确的拟合函数是如何产生的,则创建一个广义拟合函数:def f(x,a,b,c): return a * (np.sign(x-b) + c),其中理论上为a = 50b = 5c = 1

进行估算:

popt,pcov=curve_fit(f,xobs,yobs,bounds=([49,4.75,0],[50,5,2])).

现在,bounds = ([lower bound of each parameter (a,b,c)],[upper bound of each parameter])。从技术上讲,这意味着 49 & lt;a& lt;50 、 4.75 & lt;b& lt;5 和 0 & lt;c& lt;2。

Here are MY results for popt and pcov: Results

pcov表示 popt 的估计协方差,对角线提供了参数估计[Source]的方差。

结果表明,参数估计值pcov接近理论值。

基本上,广义 Heaviside 函数可以表示为:a * (np.sign(x-b) + c)

以下是将生成参数估计和相应协方差的代码:

import numpy as np
from scipy.optimize import curve_fit
xobs = np.linspace(0,10,100)
yl = np.random.rand(50); yr=np.random.rand(50)+100
yobs = np.concatenate((yl,yr),axis=0)
def f(x,a,b,c): return a * (np.sign(x-b) + c) # Heaviside fitting function
popt, pcov = curve_fit(f,xobs,yobs,bounds=([49,4.75,0],[50,5,2]))
print 'popt = %s' % popt
print 'pcov = \n %s' % pcov

最后,请注意poptpcov的估计值会有所不同。

1

这个问题很老,但如果它对其他人有用:Heaviside 函数在步骤中不可微分,这会导致最小化问题。

在我的情况下,拟合重边函数总是失败。

x = np.linspace(0,10,101)
y = np.heaviside((x-5), 0.) 
def sigmoid(x, x0,b):
    return scipy.special.expit((x-x0)*b)
args, cov = optim.curve_fit(sigmoid, x, y)
plt.scatter(x,y)
plt.plot(x, sigmoid(x, *args))
print(args)
> 
[  5.05006427 532.21427701]

fit step data with a sigmoid

本站系公益性非盈利分享网址,本文来自用户投稿,不代表码文网立场,如若转载,请注明出处

(751)
抓取 gdelt数据时出现属性错误
上一篇
如何使用SeleniumPython提取 NewEgg上的图形卡发布的所有标题文本
下一篇

相关推荐

发表评论

登录 后才能评论

评论列表(67条)