在Python中计算直方图峰值

Question

在Python中，如何计算直方图的峰值？

我尝试了以下代码：

import numpy as np
from scipy.signal import argrelextrema

data = [0, 1, 2, 3, 4, 0, 1, 2, 3, 4, 0, 1, 2, 3, 4, 1, 2, 3, 4,

        5, 6, 7, 8, 9, 5, 6, 7, 8, 9, 5, 6, 7, 8, 9,

        12,

        15, 16, 17, 18, 19, 15, 16, 17, 18, 

        19, 20, 21, 22, 23, 24,]

h = np.histogram(data, bins=[0, 5, 10, 15, 20, 25])
hData = h[0]
peaks = argrelextrema(hData, np.greater)

但结果是：

(array([3]),)

我希望它能够找到第0列和第3列的峰值。

请注意，这些峰值跨越了超过1个列。我不希望它将跨越多个列的峰值视为额外的峰值。

我也可以接受其他获取峰值的方法。

注意：

>>> h[0]
array([19, 15,  1, 10,  5])
>>> 

- brian

---

argrelextrema 返回局部最大值和最小值 - 数组中只有索引三是局部最大值。你如何定义峰值？一旦定义了要求并找出逻辑，您可以编写自己的解决方案。 - wwii

在 使用Numpy在1D numpy数组中查找局部最大值/最小值 的问题中，被接受的答案应该可以帮助你入门。 - wwii

@wwii 我不确定如何用言语描述，但我希望如果这是一个连续函数，导数为零，斜率趋近于正值。 - brian

关于端点怎么样？ - wwii

我希望端点是包含在内的。如果这是一个连续函数，最终我想找到所有的最大值。 - brian

请查看我之前评论中的链接。 - wwii

Answer 1

在计算拓扑学中，持久同调的形式提供了一个定义“峰”的方法，似乎可以满足您的需求。在一维的情况下，峰由以下图中的蓝色条表示： 该算法的描述可参见Stack Overflow answer中的一个峰值检测问题。

有趣的是，这种方法不仅识别了峰值，而且自然地量化了“重要性”。

在此博客文章中提供了上述答案的简单高效实现（与排序数字一样快速）和源材料。

- S. Huber

---

哦，也许你可以为那些不熟悉计算拓扑的人添加一个使用注释/示例？比如说，我有一些列表/数组，我从中绘制我的直方图-我该如何具体使用你的示例实现。我注意到get_persistent_homology()当输入被排序或未排序时返回不同的结果。另外，这个函数究竟返回什么？是一个由“Peak”元素组成的列表，其中每个“Peak”都有左/右/出生/死亡字段...我怎样从中获取真正的“Peak”呢？ - maxschlepzig

Answer 2

尝试使用findpeaks库。

pip install findpeaks

# Your input data:
data = [0, 1, 2, 3, 4, 0, 1, 2, 3, 4, 0, 1, 2, 3, 4, 1, 2, 3, 4, 5, 6, 7, 8, 9, 5, 6, 7, 8, 9, 5, 6, 7, 8, 9, 12, 15, 16, 17, 18, 19, 15, 16, 17, 18,  19, 20, 21, 22, 23, 24,]

# import library
from findpeaks import findpeaks

# Find some peaks using the smoothing parameter.
fp = findpeaks(lookahead=1, interpolate=10)
# fit
results = fp.fit(data)
# Make plot
fp.plot()

# Results with respect to original input data.
results['df']

# Results based on interpolated smoothed data.
results['df_interp']

- erdogant

Answer 3

我写了一个简单的函数。

def find_peaks(a):
  x = np.array(a)
  max = np.max(x)
  length = len(a)
  ret = []
  for i in range(length):
      ispeak = True
      if i-1 > 0:
          ispeak &= (x[i] > 1.8 * x[i-1])
      if i+1 < length:
          ispeak &= (x[i] > 1.8 * x[i+1])
    
      ispeak &= (x[i] > 0.05 * max)
      if ispeak:
          ret.append(i)
  return ret

我将峰值定义为大于邻居的180%且大于最大值的5%。当然，您可以根据需要调整这些数值，以找到最适合您问题的设置。

- Andrea Mauro