微信公众号搜"智元新知"关注
微信扫一扫可直接关注哦!

在 Scipy 曲线拟合中约束特定值

如何解决在 Scipy 曲线拟合中约束特定值

我有一个可能是非常基本的问题,但快速的谷歌搜索无法解决它。

所以我有一些实验数据需要用类似的方程来拟合

a * exp^{-x/t}

在需要更多组件的情况下,表达式为

a * exp^{-x/t1} + b * exp^{-x/t2} ... + n * exp^{-x/tn}

对于 n 个元素

现在我有以下代码

x = np.array([0.0001,0.0004,0.0006,0.0008,0.001,0.0015,0.002,0.004,0.006,0.008,0.01,0.05,0.1,0.2,0.5,0.6,0.8,1,1.5,2,4,6,8])
y1= np.array([5176350.00,5144208.69,4998297.04,4787100.79,4555731.93,4030741.17,3637802.79,2949911.45,2816472.26,2831962.09,2833262.53,2815205.34,2610685.14,3581566.94,1820610.74,2100882.80,1762737.50,1558251.40,997259.21,977892.00,518709.91,309594.88,186184.52])
y2 = np.array([441983.26,423371.31,399370.82,390603.58,378351.08,356511.93,349582.29,346425.39,351191.31,329363.40,325154.86,352906.21,333150.81,301613.81,94043.05,100885.77,86193.40,75548.26,27958.11,20262.68,27945.10])

    def fitcurve (x,a,b,t1,t2): 
        return a * np.exp(- x / t1) + b * np.exp(- x / t2) 
    
    popt,pcov = curve_fit(fitcurve,x,y)
    print('a = ',popt[0],'b = ',popt[1],'t1 = ',popt[2],'t2 = ',popt[3])
    
    
    plt.plot(x,y,'bo')
    plt.plot(x,fitcurve(x,*popt))

重要的是 a+b+...n = 等于 1。基本上是每个分量的百分比。理想情况下,我想绘制 1、2、3 和 4 个分量,看看哪些提供了更好的拟合

解决方法

恐怕您的数据无法用指数函数的简单总和拟合。你在图上画点是为了看看曲线的形状吗?

enter image description here

这看起来更像是一种逻辑函数(但不完全是逻辑函数),而不是指数的总和。

我可以提供一些建议来拟合指数总和(即使有关于系数总和的条件)。但这对您的数据没有用。当然,如果您有其他便于拟合指数总和的数据,我很乐意展示如何进行。

,

我不会进入模型拟合过程,但您可以做的是 argparse 可变数量的参数,然后尝试适应各种数量的指数。您可以利用 numpy 的广播功能来实现这一点。

编辑:您必须注意 argparse 中的元素数量。现在只有偶数有效。我让你来编辑那个部分(琐碎)。

目标

我们想要将 $$\sum_i^N a_i \exp(-b_i x)$$ 拟合为变量 $n$

输出:

enter image description here

实施:

from scipy import optimize,ndimage,interpolate
x = np.array([0.0001,0.0004,0.0006,0.0008,0.0010,0.0015,0.0020,0.0040,0.0060,0.0080,0.0100,0.0500,0.1000,0.2000,0.5000,0.6000,0.8000,1.0000,1.5000,2.0000,4.0000,6.0000,8.0000,10.0000])
y = np.array([416312.6500,387276.6400,364153.7600,350981.7000,336813.8800,314992.6100,310430.4600,318255.1700,318487.1700,291768.9700,276617.3000,305250.2100,272001.3500,260540.5600,173677.1900,155821.5500,151502.9700,83559.9000,256097.3600,20761.8400,1.0000])
# variable args fit
def fitcurve (x,*args): 
    args = np.array(args)
    half = len(args)//2
    y = args[:half] * np.exp(-x[:,None] * args[half:])
    return y.sum(-1)

# data seems to contain outlier?
# y = ndimage.median_filter(y,5)


popt,pcov = optimize.curve_fit(fitcurve,x,y,bounds = (0,np.inf),p0     = np.ones(6),# set variable size
                               maxfev = 1e3,)
fig,ax = plt.subplots()
ax.plot(x,'bo')
# ax.set_yscale('log')
ax.set_xscale('symlog')

ax.plot(x,fitcurve(x,*popt))
fig.show()

版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。