二维阵列的峰值检测

我在帮助一家兽医诊所测量狗爪下的压力。我使用Python进行数据分析，现在我正试图将爪子划分为(解剖学上的)子区域。

我为每个爪子制作了一个2D数组，其中包括爪子随时间加载的每个传感器的最大值。这里有一个单爪的例子，我使用Excel绘制我想要“检测”的区域。传感器周围是2 * 2的方框带有局部最大值，它们加起来的和最大。

所以我尝试了一些实验，并决定简单地寻找每一列和每一行的最大值(由于爪子的形状，不能只看一个方向)。这似乎能很好地“检测”到不同脚趾的位置，但也能标记出相邻的传感器。

那么告诉Python哪些最大值是我想要的最好方法是什么呢?

注意:2x2的方块不能重叠，因为它们必须是分开的脚趾!

此外，我选择了2x2作为方便，任何更高级的解决方案都是受欢迎的，但我只是一个人类运动科学家，所以我既不是真正的程序员也不是数学家，所以请保持“简单”。

下面是一个可以用np.loadtxt加载的版本

结果

所以我尝试了@jextee的解决方案(见下面的结果)。正如你所看到的，它对前爪很有效，但对后腿就不那么有效了。

更具体地说，它无法识别第四个脚趾的小峰。这显然是固有的事实，即循环从上向下看最低值，而不考虑它在哪里。

有没有人知道如何调整@jextee的算法，让它也能找到第四个脚趾?

因为我还没有处理其他的试验，所以我无法提供其他的样品。但我之前给出的数据是每个爪子的平均值。这个文件是一个数组，其中包含9只爪子的最大数据，按照它们与盘子接触的顺序排列。

这张照片显示了它们在空间上是如何分布在板块上的。

更新:

我已经为任何感兴趣的人建立了一个博客，我已经设置了一个OneDrive与所有的原始测量。因此，对于任何要求更多数据的人:给你更多力量!

新更新:

所以在我得到关于爪子检测和爪子分类的问题的帮助后，我终于能够检查每个爪子的脚趾检测!事实证明，除了像我自己例子中那样大的爪子，它在其他地方都不太管用。当然，事后看来，是我自己的错，我这么武断地选择了2x2。

这里有一个很好的例子，它错在哪里:指甲被识别为脚趾，而“脚跟”太宽了，它被识别了两次!

爪子太大了，所以采用没有重叠的2x2尺寸，会导致一些脚趾被检测两次。反过来，在小型犬中，它经常找不到第5个脚趾，我怀疑这是由于2x2面积太大造成的。

在对我所有的测量都尝试了当前的解决方案后，我得出了一个惊人的结论，几乎所有的小狗都没有发现第5个脚趾，而在超过50%的大型狗的冲击中，它会发现更多!

显然我需要改变它。我自己的想法是改变社区的大小，让小狗住的地方更小，大狗住的地方更大。但是generate_binary_structure不允许我改变数组的大小。

因此，我希望其他人有更好的建议来定位脚趾，也许是脚趾面积与爪子大小的比例?

当前回答

谢谢你的原始数据。我在火车上，这是我到的最远的地方(我的站就要到了)。我用regexps按摩了你的txt文件，并将其放入一个html页面，使用一些javascript进行可视化。我在这里分享它是因为一些人，比如我自己，可能会发现它比python更容易被破解。

我认为一个很好的方法是尺度和旋转不变，我的下一步将是研究高斯的混合。(每个爪垫是高斯分布的中心)。

    <html>
<head>
    <script type="text/javascript" src="http://vis.stanford.edu/protovis/protovis-r3.2.js"></script> 
    <script type="text/javascript">
    var heatmap = [[[0,0,0,0,0,0,0,4,4,0,0,0,0],
[0,0,0,0,0,7,14,22,18,7,0,0,0],
[0,0,0,0,11,40,65,43,18,7,0,0,0],
[0,0,0,0,14,61,72,32,7,4,11,14,4],
[0,7,14,11,7,22,25,11,4,14,65,72,14],
[4,29,79,54,14,7,4,11,18,29,79,83,18],
[0,18,54,32,18,43,36,29,61,76,25,18,4],
[0,4,7,7,25,90,79,36,79,90,22,0,0],
[0,0,0,0,11,47,40,14,29,36,7,0,0],
[0,0,0,0,4,7,7,4,4,4,0,0,0]
],[
[0,0,0,4,4,0,0,0,0,0,0,0,0],
[0,0,11,18,18,7,0,0,0,0,0,0,0],
[0,4,29,47,29,7,0,4,4,0,0,0,0],
[0,0,11,29,29,7,7,22,25,7,0,0,0],
[0,0,0,4,4,4,14,61,83,22,0,0,0],
[4,7,4,4,4,4,14,32,25,7,0,0,0],
[4,11,7,14,25,25,47,79,32,4,0,0,0],
[0,4,4,22,58,40,29,86,36,4,0,0,0],
[0,0,0,7,18,14,7,18,7,0,0,0,0],
[0,0,0,0,4,4,0,0,0,0,0,0,0],
],[
[0,0,0,4,11,11,7,4,0,0,0,0,0],
[0,0,0,4,22,36,32,22,11,4,0,0,0],
[4,11,7,4,11,29,54,50,22,4,0,0,0],
[11,58,43,11,4,11,25,22,11,11,18,7,0],
[11,50,43,18,11,4,4,7,18,61,86,29,4],
[0,11,18,54,58,25,32,50,32,47,54,14,0],
[0,0,14,72,76,40,86,101,32,11,7,4,0],
[0,0,4,22,22,18,47,65,18,0,0,0,0],
[0,0,0,0,4,4,7,11,4,0,0,0,0],
],[
[0,0,0,0,4,4,4,0,0,0,0,0,0],
[0,0,0,4,14,14,18,7,0,0,0,0,0],
[0,0,0,4,14,40,54,22,4,0,0,0,0],
[0,7,11,4,11,32,36,11,0,0,0,0,0],
[4,29,36,11,4,7,7,4,4,0,0,0,0],
[4,25,32,18,7,4,4,4,14,7,0,0,0],
[0,7,36,58,29,14,22,14,18,11,0,0,0],
[0,11,50,68,32,40,61,18,4,4,0,0,0],
[0,4,11,18,18,43,32,7,0,0,0,0,0],
[0,0,0,0,4,7,4,0,0,0,0,0,0],
],[
[0,0,0,0,0,0,4,7,4,0,0,0,0],
[0,0,0,0,4,18,25,32,25,7,0,0,0],
[0,0,0,4,18,65,68,29,11,0,0,0,0],
[0,4,4,4,18,65,54,18,4,7,14,11,0],
[4,22,36,14,4,14,11,7,7,29,79,47,7],
[7,54,76,36,18,14,11,36,40,32,72,36,4],
[4,11,18,18,61,79,36,54,97,40,14,7,0],
[0,0,0,11,58,101,40,47,108,50,7,0,0],
[0,0,0,4,11,25,7,11,22,11,0,0,0],
[0,0,0,0,0,4,0,0,0,0,0,0,0],
],[
[0,0,4,7,4,0,0,0,0,0,0,0,0],
[0,0,11,22,14,4,0,4,0,0,0,0,0],
[0,0,7,18,14,4,4,14,18,4,0,0,0],
[0,4,0,4,4,0,4,32,54,18,0,0,0],
[4,11,7,4,7,7,18,29,22,4,0,0,0],
[7,18,7,22,40,25,50,76,25,4,0,0,0],
[0,4,4,22,61,32,25,54,18,0,0,0,0],
[0,0,0,4,11,7,4,11,4,0,0,0,0],
],[
[0,0,0,0,7,14,11,4,0,0,0,0,0],
[0,0,0,4,18,43,50,32,14,4,0,0,0],
[0,4,11,4,7,29,61,65,43,11,0,0,0],
[4,18,54,25,7,11,32,40,25,7,11,4,0],
[4,36,86,40,11,7,7,7,7,25,58,25,4],
[0,7,18,25,65,40,18,25,22,22,47,18,0],
[0,0,4,32,79,47,43,86,54,11,7,4,0],
[0,0,0,14,32,14,25,61,40,7,0,0,0],
[0,0,0,0,4,4,4,11,7,0,0,0,0],
],[
[0,0,0,0,4,7,11,4,0,0,0,0,0],
[0,4,4,0,4,11,18,11,0,0,0,0,0],
[4,11,11,4,0,4,4,4,0,0,0,0,0],
[4,18,14,7,4,0,0,4,7,7,0,0,0],
[0,7,18,29,14,11,11,7,18,18,4,0,0],
[0,11,43,50,29,43,40,11,4,4,0,0,0],
[0,4,18,25,22,54,40,7,0,0,0,0,0],
[0,0,4,4,4,11,7,0,0,0,0,0,0],
],[
[0,0,0,0,0,7,7,7,7,0,0,0,0],
[0,0,0,0,7,32,32,18,4,0,0,0,0],
[0,0,0,0,11,54,40,14,4,4,22,11,0],
[0,7,14,11,4,14,11,4,4,25,94,50,7],
[4,25,65,43,11,7,4,7,22,25,54,36,7],
[0,7,25,22,29,58,32,25,72,61,14,7,0],
[0,0,4,4,40,115,68,29,83,72,11,0,0],
[0,0,0,0,11,29,18,7,18,14,4,0,0],
[0,0,0,0,0,4,0,0,0,0,0,0,0],
]
];
</script>
</head>
<body>
    <script type="text/javascript+protovis">    
    for (var a=0; a < heatmap.length; a++) {
    var w = heatmap[a][0].length,
    h = heatmap[a].length;
var vis = new pv.Panel()
    .width(w * 6)
    .height(h * 6)
    .strokeStyle("#aaa")
    .lineWidth(4)
    .antialias(true);
vis.add(pv.Image)
    .imageWidth(w)
    .imageHeight(h)
    .image(pv.Scale.linear()
        .domain(0, 99, 100)
        .range("#000", "#fff", '#ff0a0a')
        .by(function(i, j) heatmap[a][j][i]));
vis.render();
}
</script>
  </body>
</html>

2010-09-11 01:07:08

其他回答

有趣的问题。我想尝试的解决方案如下。

应用一个低通滤波器，比如用二维高斯掩码卷积。这将给你一堆(可能，但不一定是浮点数)值。使用每个爪垫(或脚趾)的已知近似半径执行2D非最大抑制。

这将为你提供最大的职位，而不是有多个候选人。澄清一下，第1步中蒙版的半径也应该与第2步中使用的半径相似。这个半径是可以选择的，或者兽医可以事先明确测量它(它会随着年龄/品种等而变化)。

建议的一些解决方案(均值漂移、神经网络等)可能在某种程度上是可行的，但它们过于复杂，可能并不理想。

2010-09-10 22:39:09

大概的轮廓……

你可能想要使用连接组件算法来隔离每个爪子区域。wiki对此有一个很好的描述(有一些代码):http://en.wikipedia.org/wiki/Connected_Component_Labeling

你必须决定是使用4连通性还是8连通性。就我个人而言，对于大多数问题，我更喜欢6连通性。无论如何，一旦你将每个“爪印”分离为一个连通区域，就很容易遍历该区域并找到最大值。一旦你找到了最大值，你就可以迭代地扩大区域，直到达到预定的阈值，以确定它是给定的“脚趾”。

这里有一个微妙的问题是，一旦你开始使用计算机视觉技术来识别某个东西是右/左/前/后爪，你开始观察单个脚趾，你就必须开始考虑旋转、倾斜和平移。这是通过分析所谓的“矩”来完成的。在视觉应用中有几个不同的时刻需要考虑:

中心矩:平移不变量归一化矩:缩放和平移不变量胡矩:平移、尺度和旋转不变量

更多关于moments的信息可以在wiki上搜索“image moments”。

2010-09-10 19:05:25

我用局部最大滤波器检测了峰值。下面是你的第一个4个爪子数据集的结果:

我还在9个爪子的第二个数据集上运行了它，效果也很好。

你可以这样做:

import numpy as np
from scipy.ndimage.filters import maximum_filter
from scipy.ndimage.morphology import generate_binary_structure, binary_erosion
import matplotlib.pyplot as pp

#for some reason I had to reshape. Numpy ignored the shape header.
paws_data = np.loadtxt("paws.txt").reshape(4,11,14)

#getting a list of images
paws = [p.squeeze() for p in np.vsplit(paws_data,4)]


def detect_peaks(image):
    """
    Takes an image and detect the peaks usingthe local maximum filter.
    Returns a boolean mask of the peaks (i.e. 1 when
    the pixel's value is the neighborhood maximum, 0 otherwise)
    """

    # define an 8-connected neighborhood
    neighborhood = generate_binary_structure(2,2)

    #apply the local maximum filter; all pixel of maximal value 
    #in their neighborhood are set to 1
    local_max = maximum_filter(image, footprint=neighborhood)==image
    #local_max is a mask that contains the peaks we are 
    #looking for, but also the background.
    #In order to isolate the peaks we must remove the background from the mask.

    #we create the mask of the background
    background = (image==0)

    #a little technicality: we must erode the background in order to 
    #successfully subtract it form local_max, otherwise a line will 
    #appear along the background border (artifact of the local maximum filter)
    eroded_background = binary_erosion(background, structure=neighborhood, border_value=1)

    #we obtain the final mask, containing only peaks, 
    #by removing the background from the local_max mask (xor operation)
    detected_peaks = local_max ^ eroded_background

    return detected_peaks


#applying the detection and plotting results
for i, paw in enumerate(paws):
    detected_peaks = detect_peaks(paw)
    pp.subplot(4,2,(2*i+1))
    pp.imshow(paw)
    pp.subplot(4,2,(2*i+2) )
    pp.imshow(detected_peaks)

pp.show()

之后你所需要做的就是在蒙版上使用scipy. nmage .measurements.label来标记所有不同的对象。然后你就可以单独和他们玩了。

注意，该方法工作得很好，因为背景没有噪声。如果是的话，你会在背景中检测到一堆其他不想要的峰。另一个重要因素是社区的大小。如果峰值大小发生变化，您将需要调整它(应该保持大致成比例)。

2010-09-11 03:38:07

这是一个图像配准问题。总的策略是:

有一个已知的例子，或者有一些先验的数据。将数据与示例相匹配，或将示例与数据相匹配。如果您的数据在一开始就大致对齐，这是有帮助的。

这里有一个粗略而现成的方法，“可能起作用的最愚蠢的方法”:

从五个脚趾坐标开始，大致在你期望的位置。使用每一个，迭代地爬到山顶。即给定当前位置，移动到最大邻近像素，如果它的值大于当前像素。当你的脚趾坐标停止移动时停止。

为了解决方向问题，你可以为基本方向设置8个左右的初始设置(北，东北等)。单独运行每一个，并丢弃任何两个或多个脚趾最终位于同一像素的结果。我会再思考一下这个问题，但这种事情在图像处理中还在研究中——没有正确的答案!

稍微复杂一点的想法:(加权)k -均值聚类。没那么糟。

从五个脚趾坐标开始，但现在这些是“集群中心”。

然后迭代直到收敛:

将每个像素分配到最近的集群(只需为每个集群制作一个列表)。计算每个簇的质心。对于每个聚类，这是:Sum(坐标*强度值)/Sum(坐标) 移动每个星团到新的质心。

这种方法几乎肯定会得到更好的结果，你可以得到每个簇的质量，这可能有助于识别脚趾。

(同样，您已经预先指定了集群的数量。对于聚类，您必须以一种或另一种方式指定密度:要么选择集群的数量，在这种情况下合适，要么选择集群半径，看看最终有多少。后者的一个例子是mean-shift。)

很抱歉缺少实现细节或其他细节。我可以把它写出来，但我有最后期限了。如果到下周还没有效果，请告诉我，我会试试的。

2010-09-10 14:54:44

使用持久同源分析你的数据集，我得到了以下结果(点击放大):

这是这个SO答案中描述的峰值检测方法的2d版本。上图仅显示了按持久性排序的0维持久性同源类。

我确实使用scipy.misc.imresize()将原始数据集提升了2倍。但是，请注意，我确实将四个爪子视为一个数据集;把它分成四份会使问题更简单。

方法。这背后的思想非常简单:考虑函数的函数图，该函数为每个像素分配其级别。它是这样的:

现在考虑高度为255的水位不断下降到较低的水位。在局部极大值时，岛屿弹出(出生)。两个岛屿在鞍点合并;我们认为低岛与高岛合并(死亡)。所谓的持久性图(0维同调类，我们的岛屿)描述了所有岛屿的死亡-出生值:

岛屿的持续存在就是出生水平和死亡水平之间的差异;点到灰色主对角线的垂直距离。这个数字通过减少持久性来标记岛屿。

第一张图片显示了这些岛屿的诞生地点。该方法不仅给出了局部极大值，而且通过上述持久性量化了它们的“意义”。这样就会过滤掉所有持久性过低的岛屿。但是，在您的示例中，每个岛(即每个局部最大值)都是您要寻找的峰值。

Python代码可以在这里找到。

2017-11-08 21:41:43

二维阵列的峰值检测

推荐文章

最新文章

标签