【问题标题】:Finding the maximum of a curve scipy寻找曲线scipy的最大值
【发布时间】:2015-04-30 22:16:39
【问题描述】:

我已将曲线拟合到一组数据点。我想知道如何找到我的曲线的最大值,然后我想注释那个点(我不想使用我的数据中的最大 y 值来做到这一点)。我无法准确地编写我的代码,但这是我的代码的基本布局。

import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

x = [1,2,3,4,5]
y = [1,4,16,4,1]

def f(x, p1, p2, p3):
    return p3*(p1/((x-p2)**2 + (p1/2)**2))   

p0 = (8, 16, 0.1) # guess perameters 
plt.plot(x,y,"ro")
popt, pcov = curve_fit(f, x, y, p0)
plt.plot(x, f(x, *popt)) 

还有没有办法找到峰宽?

我是否缺少可以执行此操作的简单内置函数?我可以区分函数并找到它为零的点吗?如果有怎么办?

【问题讨论】:

  • 如果y 是一个numpy 数组,你试过max(y)y.max() 吗?
  • 它是一个数组,但我被要求找到曲线的最大值而不是数组中的最大数据点。
  • 这与编程无关。那是纯数学。求一个函数的最大值,需要计算它的导数
  • 拟合函数后,您有一组参数 p1、p2、p3 来定义您的拟合。您的拟合最大值位于 (p2, 4*p3/p1)。你可以使用任何你想近似的答案,但这会给你准确的答案。
  • 我想我刚刚回答了这个问题:p1, p2, p3 = popt

标签: python matplotlib scipy max


【解决方案1】:

在您找到最佳参数以最大化您的功能后,您可以使用minimize_scalar(或scipy.optimize 中的其他方法之一)找到峰值。

请注意,在下面,我移动了x[2]=3.2,以便曲线的峰值不会落在数据点上,我们可以确定我们找到的是曲线的峰值,而不是数据。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit, minimize_scalar

x = [1,2,3.2,4,5]
y = [1,4,16,4,1]

def f(x, p1, p2, p3):
    return p3*(p1/((x-p2)**2 + (p1/2)**2))   

p0 = (8, 16, 0.1) # guess perameters 
plt.plot(x,y,"ro")
popt, pcov = curve_fit(f, x, y, p0)

# find the peak
fm = lambda x: -f(x, *popt)
r = minimize_scalar(fm, bounds=(1, 5))
print "maximum:", r["x"], f(r["x"], *popt)  #maximum: 2.99846874275 18.3928199902

x_curve = np.linspace(1, 5, 100)
plt.plot(x_curve, f(x_curve, *popt))
plt.plot(r['x'], f(r['x'], *popt), 'ko')
plt.show()

当然,与其优化函数,我们可以只计算一堆 x 值并接近:

x = np.linspace(1, 5, 10000)
y = f(x, *popt)
imax = np.argmax(y)
print imax, x[imax]     # 4996 2.99859985999

【讨论】:

  • 这正是我想要做的。谢谢!
【解决方案2】:

如果您不介意使用sympy,这很容易。假设您发布的代码已经运行:

import sympy

sym_x = sympy.symbols('x', real=True)
sym_f = f(sym_x, *popt)
sym_df = sym_f.diff()
solns = sympy.solve(sym_df)  # returns [3.0]

【讨论】:

  • 啊...确实有效,谢谢!您可能知道使用 scipy 的类似方式吗?
猜你喜欢
  • 2011-10-26
  • 1970-01-01
  • 2016-07-25
  • 2021-11-01
  • 2011-02-22
  • 2021-01-03
  • 2020-07-12
  • 1970-01-01
  • 1970-01-01
相关资源
最近更新 更多