如何统计图像中的虫子数量？

Question

我有一张图片，想用Python来计算上面有多少个虫子（就是那些连续的颜色或灰色块）。我该怎么做呢？

到目前为止，我看过ImageChops、SciPy和PIL，但不太确定该用哪个……

我觉得可以用ndimage.gaussian_filter()，然后再用scipy.ndimage.measurements.label()，但我还不太清楚怎么用后者来计算我在高斯模糊后的图像中的蓝点……它看起来像这样

---

好的，

根据上面的图片，我现在有了这段代码：

#! /usr/bin/python

import numpy as np
import scipy
import pylab
import pymorph
import mahotas
from PIL import Image
import PIL.ImageOps
from scipy import ndimage

image = Image.open('bugs.jpg')   
inverted_image = PIL.ImageOps.invert(image)    
inverted_image.save('in_bugs.jpg')
dna = mahotas.imread('in_bugs.jpg')
#pylab.imshow(dna)
pylab.gray()
#pylab.show()
T = mahotas.thresholding.otsu(dna)
pylab.imshow(dna > T)
#pylab.show()
dnaf = ndimage.gaussian_filter(dna, 8)
T = mahotas.thresholding.otsu(dnaf)
pylab.imshow(dnaf > T)
#pylab.show()
labeled,nr_objects = ndimage.label(dnaf > T)
print nr_objects
pylab.imshow(labeled)
pylab.jet()
pylab.show()

问题是，这段代码给我的结果是5，虽然还不错，但我希望能更准确一点，我想看到2。那我该怎么做呢？在应用高斯滤波之前模糊一下图片会有帮助吗？

谢谢你的帮助！

Ron

Answer 1

这是一个计算机视觉的任务，你可以用一个很流行的工具叫做OpenCV来解决。

你可能想要用一些形态学操作来处理图像，比如开运算，这样可以去掉一些杂乱的噪声。然后，你可以数一下那些面积超过某个阈值的连通区域有多少个。

---

这里有一些可以参考的资源：

• 物体识别（维基百科）
• 计算机视觉（维基百科）
• OpenCV食谱

Answer 2

你的高斯滤波已经差不多可以用了，但你考虑的半径比实际需要的要大得多。比如，我们可以拿半径为15的卷积核来举个例子。下面是我们得到的效果：

这里有两个明显的谷（其实在图上看起来像山峰），而经过滤波后的图像直方图显示，大部分数据现在都接近可能的最大值。

如果只看直方图的一部分，我们可以更清楚地看到我们感兴趣的数据：那些较暗的区域。

所以，简单设置一个阈值在 0.5，我们得到的结果（正好对应虫子的位置）是：

根据你实现的方式（或者你使用的库的实现方式），这个阈值可能会有所不同。但通过观察直方图，你可以找到一个合适的阈值。如果你不想通过直方图来猜测这个阈值，那么你需要在高斯滤波之外对图像进行预处理。这样做之后，你的图像会变得简单得多，像Otsu方法这样的算法就能自动找到你想要的阈值。进行形态学闭合操作，然后再进行形态学开操作，最后用Otsu进行二值化，得到的结果是：

这些形状更接近最初的样子，因为我们没有依赖线性低通滤波器，这种滤波器最多只能模糊轮廓。

编辑：

由于问题现在包含了一些代码，我觉得有必要解释一下为什么使用Otsu的方法是错误的，因为代码的实现方式不对。Otsu的阈值方法实际上是期望数据是双峰的，但正如上面的直方图所示，这里并不是这种情况。Otsu会给出一个接近右侧巨大峰值的阈值，而0.5这个好的点则远离那里。为了复制这个答案中展示的第一个结果，这里有一些基本代码：

import sys
import numpy
from PIL import Image
from scipy.ndimage import gaussian_filter, label

img = Image.open(sys.argv[1]).convert('L')
im = numpy.array(img)

im_g = gaussian_filter(im, 3)
im_norm = (im_g - im_g.min()) / (float(im_g.max()) - im_g.min())
im_norm[im_norm < 0.5] = 0
im_norm[im_norm >= 0.5] = 1

result = 255 - (im_norm * 255).astype(numpy.uint8)
print u"Objects: %d" % label(result)[1]

Image.fromarray(result).save(sys.argv[2])

注意，这段代码使用了 sigma = 3（而最初使用的是7.5）作为高斯核，scipy内部会构建一个半径大4倍的窗口。对于这张特定的图像，sigma的范围从2到10都能得到同样的效果——检测到2个物体。