59

我们使用数据采集卡从设备中获取读数，该设备将信号增加到峰值，然后回落到原始值附近。为了找到峰值，我们目前在数组中搜索最高读数，并使用索引来确定我们计算中使用的峰值时间。

如果最大值是我们正在寻找的峰值，则此方法效果很好，但如果设备无法正常工作，我们可以看到第二个峰值，该峰值可能高于初始峰值。我们在 90 秒的时间内从 16 个设备中每秒读取 10 个读数。

我最初的想法是循环检查读数以查看前一个点和下一个点是否小于当前点以找到一个峰值并构建一个峰值数组。也许我们应该查看当前位置两侧的多个点的平均值，以考虑系统中的噪声。这是最好的方法还是有更好的技术？

我们确实使用 LabVIEW，我查看了LAVA 论坛，有很多有趣的例子。这是我们测试软件的一部分，我们正在努力避免使用过多的非标准 VI 库，因此我希望获得有关所涉及的流程/算法而不是特定代码的反馈。

4

9 回答 9

86

有很多很多经典的峰值检测方法，其中任何一种都可能有效。您必须特别了解限制数据质量的因素。以下是基本说明：

  1. 在数据中的任意两点之间，(x(0), y(0))和(x(n), y(n))相加并称之为（“旅行”）并将（“上升”）设置y(i + 1) - y(i)为适当的小。 表示一个峰值。如果不太可能由于噪声导致大行程或噪声围绕基本曲线形状对称分布，则此方法可以正常工作。对于您的应用，接受分数高于给定阈值的最早峰值，或分析每次上升的行程曲线值以获得更多有趣的属性。0 <= i < nTRy(n) - y(0) + kkT/R > 1

  2. 使用匹配过滤器对与标准峰形的相似度进行评分（本质上，使用针对某个形状的归一化点积来获得相似度的余弦度量）

  3. 根据标准峰形去卷积并检查高值（尽管我经常发现 2 对简单仪器输出的噪声不太敏感）。

  4. 平滑数据并检查等距点的三元组 where, ifx0 < x1 < x2, y1 > 0.5 * (y0 + y2)或检查欧几里得距离，如下所示： D((x0, y0), (x1, y1)) + D((x1, y1), (x2, y2)) > D((x0, y0),(x2, y2))，这依赖于三角形不等式。使用简单的比率将再次为您提供评分机制。

  5. 为您的数据拟合一个非常简单的 2 高斯混合模型（例如，Numerical Recipes 有一个很好的现成代码块）。取较早的高峰。这将正确处理重叠峰。

  6. 在数据中找到与简单高斯、柯西、泊松或你有什么关系的曲线的最佳匹配。在广泛的范围内评估此曲线，并在记录其峰值位置后从数据副本中减去它。重复。取模型参数（可能是标准偏差，但某些应用程序可能关心峰度或其他特征）满足某些标准的最早峰值。当从数据中减去峰值时，请注意留下的伪影。最佳匹配可能由上面#2 中建议的匹配评分类型确定。

我已经完成了您之前正在做的事情：在 DNA 序列数据中找到峰值，在从测量曲线估计的导数中找到峰值，以及在直方图中找到峰值。

我鼓励您仔细注意适当的基线。维纳滤波或其他滤波或简单的直方图分析通常是在存在噪声的情况下进行基线的简单方法。

最后，如果您的数据通常很嘈杂，并且您从卡中获取数据作为未引用的单端输出（甚至是引用的，只是没有差分），并且如果您对每个数据点的大量观察结果进行平均，请尝试对这些数据进行排序观察并丢弃第一个和最后一个四分位数并平均剩余的四分位数。有许多这样的异常值消除策略非常有用。

于 2008-09-04T18:07:31.570 回答
10

您可以尝试信号平均，即对于每个点，将周围 3 个或更多点的值平均。如果噪声信号很大，那么即使这样也可能无济于事。

我意识到这与语言无关，但猜测您使用的是 LabView，LabView 附带了许多预打包的信号处理 VI，您可以使用它们来进行平滑和降噪。NI 论坛是在这类事情上获得更专业帮助的好地方。

于 2008-08-06T11:12:48.890 回答
6

已经对这个问题进行了一些详细的研究。

ROOT（核/粒子物理分析工具）TSpectrum*类中有一组非常最新的实现。该代码适用于一维到三维数据。

ROOT 源代码可用，因此您可以根据需要获取此实现。

从TSpectrum类文档：

此类中使用的算法已在以下参考文献中发布：

[1] M.Morhac 等人：多维符合伽马射线光谱的背景消除方法。物理研究中的核仪器和方法 A 401 (1997) 113-132。

[2] M.Morhac 等人：高效的一维和二维金反卷积及其在伽马射线光谱分解中的应用。物理研究中的核仪器和方法 A 401 (1997) 385-408。

[3] M.Morhac 等人：识别多维符合伽马射线光谱中的峰值。研究物理学中的核仪器和方法 A 443（2000），108-125。

对于那些没有 NIM 在线订阅的人，这些论文是从课程文档中链接的。

所做的简短版本是直方图展平以消除噪声，然后在展平直方图中通过蛮力检测局部最大值。

于 2008-08-23T14:28:23.130 回答
6

我想为这个线程贡献一个我自己开发的算法：

它基于分散原理：如果一个新的数据点与某个移动均值相差给定 x 个标准差，则算法发出信号（也称为z-score）。该算法非常稳健，因为它构造了单独的移动均值和偏差，因此信号不会破坏阈值。因此，无论先前信号的数量如何，都以大致相同的精度识别未来信号。该算法需要 3 个输入lag = the lag of the moving window：threshold = the z-score at which the algorithm signals和influence = the influence (between 0 and 1) of new signals on the mean and standard deviation。例如，a lagof 5 将使用最后 5 个观测值来平滑数据。threshold如果数据点与移动平均值相差 3.5 个标准差，则A值为 3.5。而influence0.5 的 an 给出了一半的信号正常数据点的影响。同样，influence为 0 的 an 完全忽略重新计算新阈值的信号：因此 0 的影响是最稳健的选择。

它的工作原理如下：

伪代码

# Let y be a vector of timeseries data of at least length lag+2
# Let mean() be a function that calculates the mean
# Let std() be a function that calculates the standard deviaton
# Let absolute() be the absolute value function

# Settings (the ones below are examples: choose what is best for your data)
set lag to 5;          # lag 5 for the smoothing functions
set threshold to 3.5;  # 3.5 standard deviations for signal
set influence to 0.5;  # between 0 and 1, where 1 is normal influence, 0.5 is half

# Initialise variables
set signals to vector 0,...,0 of length of y;   # Initialise signal results
set filteredY to y(1,...,lag)                   # Initialise filtered series
set avgFilter to null;                          # Initialise average filter
set stdFilter to null;                          # Initialise std. filter
set avgFilter(lag) to mean(y(1,...,lag));       # Initialise first value
set stdFilter(lag) to std(y(1,...,lag));        # Initialise first value

for i=lag+1,...,t do
  if absolute(y(i) - avgFilter(i-1)) > threshold*stdFilter(i-1) then
    if y(i) > avgFilter(i-1)
      set signals(i) to +1;                     # Positive signal
    else
      set signals(i) to -1;                     # Negative signal
    end
    # Adjust the filters
    set filteredY(i) to influence*y(i) + (1-influence)*filteredY(i-1);
    set avgFilter(i) to mean(filteredY(i-lag,i),lag);
    set stdFilter(i) to std(filteredY(i-lag,i),lag);
  else
    set signals(i) to 0;                        # No signal
    # Adjust the filters
    set filteredY(i) to y(i);
    set avgFilter(i) to mean(filteredY(i-lag,i),lag);
    set stdFilter(i) to std(filteredY(i-lag,i),lag);
  end
end

演示

鲁棒阈值算法的演示

> 有关更多信息，请参阅原始答案

于 2015-11-03T21:46:15.990 回答
4

这种方法基本上来自大卫马尔的书《愿景》

用峰的预期宽度高斯模糊您的信号。这样可以消除噪声尖峰，并且您的相位数据完好无损。

然后边缘检测（LOG会做）

那么你的边缘就是特征的边缘（比如山峰）。在边缘之间寻找峰，按大小对峰进行排序，就完成了。

我对此使用了变体，它们工作得很好。

于 2008-09-02T01:00:13.913 回答
2

我认为您希望将您的信号与预期的示例信号进行交叉关联。但是，我学习信号处理已经很长时间了，即便如此，我也没有太在意。

于 2008-08-06T11:38:05.393 回答
0

您可以将一些标准偏差应用于您的逻辑并注意超过 x% 的峰值。

于 2008-08-06T11:17:29.240 回答
0

我不太了解仪器，所以这可能完全不切实际，但话又说回来，这可能是一个有用的不同方向。如果你知道读数是如何失败的，并且在这些失败的情况下峰值之间有一定的间隔，为什么不在每个间隔进行梯度下降。如果下降带您回到之前搜索过的区域，您可以放弃它。根据采样表面的形状，这也可以帮助您比搜索更快地找到峰。

于 2008-08-06T11:38:37.607 回答
0

期望的峰值和不需要的第二个峰值之间是否存在质的差异？如果两个峰值都是“尖锐的”——即持续时间很短——当在频域中查看信号时（通过进行 FFT），您将在大多数频段获得能量。但是，如果“好”峰值可靠地在“坏”峰值中不存在的频率处存在能量，或者反之亦然，那么您可能能够以这种方式自动区分它们。

于 2008-09-02T12:33:16.163 回答