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

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

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

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

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

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

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


结果

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

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

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

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

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

更新:

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


新更新:

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

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

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

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

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

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


当前回答

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

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

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

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

其他回答

解决方案

数据文件:paw.txt。源代码:

from scipy import *
from operator import itemgetter

n = 5  # how many fingers are we looking for

d = loadtxt("paw.txt")
width, height = d.shape

# Create an array where every element is a sum of 2x2 squares.

fourSums = d[:-1,:-1] + d[1:,:-1] + d[1:,1:] + d[:-1,1:]

# Find positions of the fingers.

# Pair each sum with its position number (from 0 to width*height-1),

pairs = zip(arange(width*height), fourSums.flatten())

# Sort by descending sum value, filter overlapping squares

def drop_overlapping(pairs):
    no_overlaps = []
    def does_not_overlap(p1, p2):
        i1, i2 = p1[0], p2[0]
        r1, col1 = i1 / (width-1), i1 % (width-1)
        r2, col2 = i2 / (width-1), i2 % (width-1)
        return (max(abs(r1-r2),abs(col1-col2)) >= 2)
    for p in pairs:
        if all(map(lambda prev: does_not_overlap(p,prev), no_overlaps)):
            no_overlaps.append(p)
    return no_overlaps

pairs2 = drop_overlapping(sorted(pairs, key=itemgetter(1), reverse=True))

# Take the first n with the heighest values

positions = pairs2[:n]

# Print results

print d, "\n"

for i, val in positions:
    row = i / (width-1)
    column = i % (width-1)
    print "sum = %f @ %d,%d (%d)" % (val, row, column, i)
    print d[row:row+2,column:column+2], "\n"

输出不重叠方块。似乎选择的区域与您的示例中相同。

一些评论

棘手的部分是计算所有2x2平方的和。我猜你们都需要,所以可能会有重叠。我使用切片从原始2D数组中切割第一个/最后一个列和行,然后将它们全部重叠在一起并计算总和。

为了更好地理解它,对一个3x3数组进行成像:

>>> a = arange(9).reshape(3,3) ; a
array([[0, 1, 2],
       [3, 4, 5],
       [6, 7, 8]])

然后你可以把它切片:

>>> a[:-1,:-1]
array([[0, 1],
       [3, 4]])
>>> a[1:,:-1]
array([[3, 4],
       [6, 7]])
>>> a[:-1,1:]
array([[1, 2],
       [4, 5]])
>>> a[1:,1:]
array([[4, 5],
       [7, 8]])

现在想象一下,你把它们一个叠在另一个上面,然后在相同的位置上求和。这些和将与左上角位置相同的2x2平方的和完全相同:

>>> sums = a[:-1,:-1] + a[1:,:-1] + a[:-1,1:] + a[1:,1:]; sums
array([[ 8, 12],
       [20, 24]])

当你有2x2平方的和时,你可以用max来求最大值,或者排序,或者排序来求峰值。

为了记住峰值的位置,我将每个值(总和)与其在扁平数组中的序数位置(参见zip)耦合起来。然后在打印结果时再次计算行/列位置。

笔记

我允许2x2平方重叠。编辑后的版本过滤掉其中的一些,这样只有不重叠的正方形出现在结果中。

选择手指(一个想法)

另一个问题是如何从所有的峰中选择可能是手指的峰。我有个主意,可能有用也可能没用。我现在没有时间实现它,所以只是伪代码。

我注意到,如果前面的手指几乎保持在一个完美的圆圈上,后面的手指应该在这个圆圈里。此外,食指的间距也差不多相等。我们可以尝试使用这些启发式属性来检测手指。

伪代码:

select the top N finger candidates (not too many, 10 or 12)
consider all possible combinations of 5 out of N (use itertools.combinations)
for each combination of 5 fingers:
    for each finger out of 5:
        fit the best circle to the remaining 4
        => position of the center, radius
        check if the selected finger is inside of the circle
        check if the remaining four are evenly spread
        (for example, consider angles from the center of the circle)
        assign some cost (penalty) to this selection of 4 peaks + a rear finger
        (consider, probably weighted:
             circle fitting error,
             if the rear finger is inside,
             variance in the spreading of the front fingers,
             total intensity of 5 peaks)
choose a combination of 4 peaks + a rear peak with the lowest penalty

这是一种蛮力方法。如果N相对较小,那么我认为它是可行的。对于N=12,有C_12^5 = 792种组合,乘以5种选择后指的方法,所以每个爪子要评估3960种情况。

我不确定这是否回答了问题,但似乎你可以只寻找n个没有邻居的最高的山峰。

这是要点。注意,它是用Ruby编写的,但其思想应该很清楚。

require 'pp'

NUM_PEAKS = 5
NEIGHBOR_DISTANCE = 1

data = [[1,2,3,4,5],
        [2,6,4,4,6],
        [3,6,7,4,3],
       ]

def tuples(matrix)
  tuples = []
  matrix.each_with_index { |row, ri|
    row.each_with_index { |value, ci|
      tuples << [value, ri, ci]
    }
  }
  tuples
end

def neighbor?(t1, t2, distance = 1)
  [1,2].each { |axis|
    return false if (t1[axis] - t2[axis]).abs > distance
  }
  true
end

# convert the matrix into a sorted list of tuples (value, row, col), highest peaks first
sorted = tuples(data).sort_by { |tuple| tuple.first }.reverse

# the list of peaks that don't have neighbors
non_neighboring_peaks = []

sorted.each { |candidate|
  # always take the highest peak
  if non_neighboring_peaks.empty?
    non_neighboring_peaks << candidate
    puts "took the first peak: #{candidate}"
  else
    # check that this candidate doesn't have any accepted neighbors
    is_ok = true
    non_neighboring_peaks.each { |accepted|
      if neighbor?(candidate, accepted, NEIGHBOR_DISTANCE)
        is_ok = false
        break
      end
    }
    if is_ok
      non_neighboring_peaks << candidate
      puts "took #{candidate}"
    else
      puts "denied #{candidate}"
    end
  end
}

pp non_neighboring_peaks

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

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

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

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

如果你能够创建一些训练数据,那么使用神经网络可能是值得尝试的……但这需要手工标注许多样本。

物理学家对这个问题进行了深入的研究。在ROOT中有一个很好的实现。查看TSpectrum类(特别是TSpectrum2)和它们的文档。

引用:

M.Morhac et al.: Background elimination methods for multidimensional coincidence gamma-ray spectra. Nuclear Instruments and Methods in Physics Research A 401 (1997) 113-132. M.Morhac et al.: Efficient one- and two-dimensional Gold deconvolution and its application to gamma-ray spectra decomposition. Nuclear Instruments and Methods in Physics Research A 401 (1997) 385-408. M.Morhac et al.: Identification of peaks in multidimensional coincidence gamma-ray spectra. Nuclear Instruments and Methods in Research Physics A 443(2000), 108-125.

...对于那些没有订阅NIM的人:

Spectrum.doc SpectrumDec.ps.gz SpectrumSrc.ps.gz SpectrumBck.ps.gz