Python/SciPy的寻峰算法

Question

我可以通过寻找一阶导数的零交叉点或其他东西来自己编写一些东西，但它似乎是一个足够常见的函数，可以包含在标准库中。有人知道吗？

我的特殊应用是一个二维数组，但通常用于查找fft中的峰值等

具体来说，在这类问题中，存在多个强峰值，然后有许多较小的“峰值”，这些都是由噪声引起的，应该忽略不计。这些只是示例；不是我的实际数据：

一维峰值：

二维峰值：

峰值查找算法将找到这些峰值的位置（不只是它们的值），理想情况下将找到真正的样本间峰值，而不仅仅是具有最大值的索引，可能使用quadratic interpolation或其他方法。

通常你只关心几个强峰，所以选择它们要么是因为它们高于某个阈值，要么是因为它们是按振幅排列的有序列表的第一个n峰。

就像我说的，我自己知道怎么写这样的东西。我只是问是否有一个预先存在的功能或包是已知的工作良好。

更新：

我translated a MATLAB script并且它在一维情况下工作正常，但是可能更好。

更新更新：

一维情况下的sixtenbecreated a better version。

Answer 1

我正在研究一个类似的问题，我发现一些最好的参考文献来自化学（从质量规范数据中的峰值发现）。要对峰值查找算法进行彻底的检查，请阅读this。这是我所遇到的最清晰的峰值查找技术评论之一。（小波是在噪声数据中寻找此类峰值的最佳方法。）。

你的山峰看起来很清晰，没有隐藏在噪音中。既然如此，我建议使用光滑的savtizky-golay导数来找到峰值（如果你只是区分上面的数据，你会有一堆误报）。这是一种非常有效的技术，而且很容易实现（您确实需要一个矩阵类w/basic操作）。如果你找到第一个S-G导数的过零点，我想你会很高兴的。

Answer 2

我不认为你要找的是西皮提供的。在这种情况下，我会自己编写代码。

scipy.interpolate中的样条插值和平滑非常好，可能有助于拟合峰值，然后找到其最大值的位置。

Answer 3

函数^{}，顾名思义，对此很有用。但了解其参数width、threshold、distance以及最重要的参数prominence是获得良好峰提取的重要条件。

根据我的测试和文档，突出度的概念是“有用的概念”，它可以保留好的峰值，并丢弃有噪声的峰值。

什么是(topographic) prominence？它是从山顶到任何更高地形所需的最低下降高度，如图所示：

想法是：

The higher the prominence, the more "important" the peak is.

测试：

我特意用了一个（有噪声的）频率变化的正弦波，因为它显示出许多困难。我们可以看到，width参数在这里不是很有用，因为如果设置的最小值width太高，那么它将无法跟踪高频部分非常接近的峰值。如果将width设置得太低，则信号左侧会出现许多不需要的峰值。与distance相同的问题。threshold只与直接邻居进行比较，这在这里是不有用的。prominence是最好的解决方案。请注意，您可以组合这些参数中的许多！

代码：

import numpy as np
import matplotlib.pyplot as plt 
from scipy.signal import find_peaks

x = np.sin(2*np.pi*(2**np.linspace(2,10,1000))*np.arange(1000)/48000) + np.random.normal(0, 1, 1000) * 0.15
peaks, _ = find_peaks(x, distance=20)
peaks2, _ = find_peaks(x, prominence=1)      # BEST!
peaks3, _ = find_peaks(x, width=20)
peaks4, _ = find_peaks(x, threshold=0.4)     # Required vertical distance to its direct neighbouring samples, pretty useless
plt.subplot(2, 2, 1)
plt.plot(peaks, x[peaks], "xr"); plt.plot(x); plt.legend(['distance'])
plt.subplot(2, 2, 2)
plt.plot(peaks2, x[peaks2], "ob"); plt.plot(x); plt.legend(['prominence'])
plt.subplot(2, 2, 3)
plt.plot(peaks3, x[peaks3], "vg"); plt.plot(x); plt.legend(['width'])
plt.subplot(2, 2, 4)
plt.plot(peaks4, x[peaks4], "xk"); plt.plot(x); plt.legend(['threshold'])
plt.show()