logo标签列表

当迭代次数高达40,000时，如何加速图像分析中的“for循环”

pythonperformancefor-loopnumpyequality
11

这段代码的前提条件很长，所以我会尽力概括。WB/RG/BYColor是基础图像，FIDO是该基础图像的覆盖层，应用于它。S_wb/rg/by是最终的输出图像。WB/RG/BYColor与FIDO的大小相同。

对于FIDO中的每个独特元素，我们想计算该区域在基础图片中的平均颜色。下面的代码可以实现此目的，但是由于numFIDOs非常大（高达40,000），因此这需要很长时间。

这些平均值是针对三个独立的RGB通道计算的。

sX=200
sY=200
S_wb = np.zeros((sX, sY))
S_rg = np.zeros((sX, sY))
S_by = np.zeros((sX, sY))
uniqueFIDOs, unique_counts = np.unique(FIDO, return_counts=True) 
numFIDOs = uniqueFIDOs.shape  
for i in np.arange(0,numFIDOs[0]):
    Lookup = FIDO==uniqueFIDOs[i]
    # Get average of color signals for this FIDO
    S_wb[Lookup] = np.sum(WBColor[Lookup])/unique_counts[i]
    S_rg[Lookup] = np.sum(RGColor[Lookup])/unique_counts[i]
    S_by[Lookup] = np.sum(BYColor[Lookup])/unique_counts[i]

这需要大约7.89秒才能运行，不算太长，但会包含在另一个循环中，所以会逐渐累积！

我已经尝试了向量化（如下所示），但我无法做到。

FIDOsize = unique_counts[0:numFIDOs[0]:1]
Lookup = FIDO ==uniqueFIDOs[0:numFIDOs[0]:1]
S_wb[Lookup] = np.sum(WBColor[Lookup])/FIDOsize
S_rg[Lookup] = np.sum(RGColor[Lookup])/FIDOsize
S_by[Lookup] = np.sum(BYColor[Lookup])/FIDOsize

数组大小匹配错误

您可能希望使您的代码自包含并添加运行时测量。 - cel
能否附加变量？为了使其自包含，我必须将许多行代码复制到这里。 - William Baker Morrison
一个明显的事情是，你可以将7个缓慢的查找变成一个：indexarr = FIDO==uniqueFIDOs[i]。然后只需执行S[indexarr] = np.sum(Color[indexarr])/FIDOsize。你也可以通过使用uniqueFIDOs, unique_counts = np.unique(FIDO, return_counts=True)来简化FIDOsize计算，然后只需FIDOsize = unique_counts[i]。 - wflynny
当涉及到（时间关键的）图像分析时，类似pyCUDA这样的工具可能值得一试 - 至少对于某些应用情况来说是这样。 - Marco13
5个回答
5

根据我的计时，这个方法比你原来的方法快了大约10倍。我使用以下数组进行测试：

import numpy as np

sX=200
sY=200

FIDO = np.random.randint(0, sX*sY, (sX, sY))
WBColor = np.random.randint(0, sX*sY, (sX, sY))
RGColor = np.random.randint(0, sX*sY, (sX, sY))
BYColor = np.random.randint(0, sX*sY, (sX, sY))

这是我计时的部分：
import collections

colors = {'wb': WBColor, 'rg': RGColor, 'by': BYColor}
planes = colors.keys()
S = {plane: np.zeros((sX, sY)) for plane in planes}

for plane in planes:
    counts = collections.defaultdict(int)
    sums = collections.defaultdict(int)
    for (i, j), f in np.ndenumerate(FIDO):
        counts[f] += 1
        sums[f] += colors[plane][i, j]
    for (i, j), f in np.ndenumerate(FIDO):
        S[plane][i, j] = sums[f]/counts[f]

可能是因为虽然Python中的循环较慢，但是这种方法遍历数据次数较少。

请注意，如果FIDO中的唯一值较少，原始版本会变得更快。对于大多数情况，这需要大致相同的时间。

你可以使用 collections.Counter 替代 collections.defaultdict(int)。 - Bakuriu
defaultdict更快。此外，对于这个应用程序，大部分时间实际上是花在遍历数据上的，因此我希望能够更新总和而不仅仅是计数。 - chthonicdaemon
在Python3中，Counter比defaultdict(int)快两倍。 - Bakuriu
我的 MacBook Air 上的 Python 3.5：Counter：530 毫秒，defaultdict(int)：420 毫秒。 - chthonicdaemon
4

正如@lejlot之前建议的那样，该代码很难向量化。除非您事先知道哪些像素属于每个FIDO，否则不能并行运行。我不知道您是否将FIDO称为超像素，但我通常处理这种问题，并且迄今为止找到的最佳解决方案如下：

  • Flatten the data:

    data = data.reshape(-1, 3)
labels = FIDO.copy()

    Here data is your (Width, Height, 3) image, rather than the separate 3 vectors that you have. It gets flattened to (Width * Height, 3).

  • Relabel FIDO to 0..N-1 range, where N=num unique FIDO:

    from skimage.segmentation import relabel_sequential

labels = relabel_sequential(labels)[0]
labels -= labels.min()

    The above, from scikit-image, transforms your FIDO array to the [0, N-1] range, which is much easier to work with later.

  • Last, code in cython a simple function to calculate the mean for each of the FIDO;s (as they are ordered from 0 to N, you can do it in a 1D array with length N):

    def fmeans(double[:, ::1] data, long[::1] labels, long nsp):
    cdef long n,  N = labels.shape[0]
    cdef int K = data.shape[1]
    cdef double[:, ::1] F = np.zeros((nsp, K), np.float64)
    cdef int[::1] sizes = np.zeros(nsp, np.int32)
    cdef long l, b
    cdef double t

    for n in range(N):
        l = labels[n]
        sizes[l] += 1

        for z in range(K):
            t = data[n, z]
            F[l, z] += t

    for n in range(nsp):
        for z in range(K):
            F[n, z] /= sizes[n]

return np.asarray(F)
您可以稍后调用该函数（一旦使用cython编译），操作非常简单：

mean_colors = fmeans(data, labels.flatten(), labels.max()+1) # labels.max()+1 == N

然后可以恢复平均颜色的图像如下所示：
mean_img = mean_colors[labels]

如果您不想使用cython进行编码，scikit-image还通过使用图形结构和networkx提供了绑定，但速度较慢： http://scikit-image.org/docs/dev/auto_examples/plot_rag_mean_color.html 上面的示例包含了获取每个超像素平均颜色作为labels1（您的FIDO）的函数调用。
注意：使用cython方法更快，因为它只迭代一次图像而不是迭代唯一FIDO N的数量并扫描图像（大小M = Width x Height）。 因此，计算成本的顺序为O(M+N)，而不是原始方法的O(M*N)。
示例测试：
import numpy as np
from skimage.segmentation import relabel_sequential

sX=200
sY=200

FIDO = np.random.randint(0, sX*sY, (sX, sY))
data = np.random.rand(sX, sY, 3) # Your image

压平和重新标记：

data = data.reshape(-1, 3)
labels = relabel_sequential(FIDO)[0]
labels -= labels.min()

获取均值：

>>> %timeit color_means = fmeans(data, labels.flatten(), labels.max()+1)
1000 loops, best of 3: 520 µs per loop

对于一个200x200像素的图像，需要0.5毫秒（半毫秒）：

print labels.max()+1 # --> 25787 unique FIDO
print color_means.shape # --> (25287, 3), the mean color of each FIDO

您可以使用智能索引来恢复平均颜色的图像：

mean_image = color_means[labels]
print mean_image.shape # --> (200, 200, 3)

我怀疑你无法通过原始的Python方法达到那种速度（或者至少我没有找到如何做到）。

很遗憾，Skimage 0.10中的“relabel_sequential”在我的系统（Linux7.0）上似乎不可用。 - William Baker Morrison
relabel_sequential 应该在 scikit-image 0.10.x 中可用。要安装它，您可以在 Linux 终端中键入 pip install scikit-image，或者如果您使用 anaconda，则可以键入 conda install scikit-image。但是，如果您仍然需要一个纯 Python 方法，@B.M. 的 pandas 选项似乎是一个合适的选择。 - Imanol Luengo
@WBM 输入 easy_install --upgrade pip 来升级 pip （你可能需要在之前加上sudo）之后，再尝试上面的命令来重新安装 scikit-image。 - Imanol Luengo
@imalunego 要求已经满足（使用--upgrade进行升级）：scikit-image在/usr/lib/pymodules/python2.7中 /usr/local/lib/python2.7/dist-packages/pip-7.1.2-py2.7.egg/pip/_vendor/requests/packages/urllib3/util/ssl_.py:90: InsecurePlatformWarning: 真正的SSLContext对象不可用。这会阻止urllib3适当地配置SSL，并可能导致某些SSL连接失败。有关更多信息，请参见https://urllib3.readthedocs.org/en/latest/security.html#insecureplatformwarning。 InsecurePlatformWarning - William Baker Morrison
尝试按照命令建议进行升级，即pip install --upgrade scikit-image。 - Imanol Luengo
显示剩余3条评论
4

你的代码不够优化，因为你在FIDO每个区域都要扫描整张图片。更好的方法是将每个区域的像素分组并先计算平均值。 pandas提供了这种计算的好工具（只用一个通道）。然后你可以在这些区域上应用平均值：

import numpy as np
import pandas as pd     
sX=200
sY=200
Nreg=sX*sY
WBColor=np.random.randint(0,256,(sX,sY))
FIDO=np.random.randint(0,Nreg,(sX,sY))


def oldloop():
    S_wb = np.zeros((sX, sY))
    uniqueFIDOs, unique_counts = np.unique(FIDO, return_counts=True) 
    numFIDOs = uniqueFIDOs.shape 
    for i in np.arange(0,numFIDOs[0]):
        Lookup = FIDO==uniqueFIDOs[i]
        S_wb[Lookup] = np.sum(WBColor[Lookup])/unique_counts[i]
    return S_wb

def newloop():
    index=pd.Index(FIDO.flatten(),name='region')
    means= pd.DataFrame(WBColor.flatten(),index).groupby(level='region').mean()
    lookup=np.zeros(Nreg)
    lookup[means.index]=means.values
    return lookup[FIDO]

在这种情况下，速度大约快了200倍：
In [32]: np.allclose(oldloop(),newloop())
Out[32]: True

In [33]: %timeit -n1 oldloop()
1 loops, best of 3: 3.92 s per loop

In [34]: %timeit -n100 newloop()
100 loops, best of 3: 20.5 ms per loop

编辑

另一种很酷的现代方法是使用numba。您编写基本的Python代码，几乎以C的速度运行：

from numba import jit

@jit
def numbaloops():
    counts=np.zeros(Nreg)
    sums=np.zeros(Nreg)
    S = np.empty((sX, sY))
    for x in range(sX):
        for y in range(sY):
            region=FIDO[x,y]
            value=WBColor[x,y]
            counts[region]+=1
            sums[region]+=value
    for x in range(sX):
        for y in range(sY):
            region=FIDO[x,y]
            S[x,y]=sums[region]/counts[region]
    return S

您现在的速度大约快了4000倍：

In [45]: np.allclose(oldloop(),numbaloops())
Out[45]: True

In [46]: %timeit -n1000 numbaloops()
1000 loops, best of 3: 1.06 ms per loop
你可以将Numba解决方案中的最后两个循环合并。 - chthonicdaemon
@chthonicdaemon：谢谢。我会考虑的，现在只有两个循环，可以再提高2倍的效率。 - B. M.
@B. M. 我在使用 Pandas 方法时遇到了这个错误： >>> means=pd.DataFrame(WBColor.flatten(),index).groupby(level='region').mean() .....code cut out..... File "/usr/lib/pymodules/python2.7/pandas/core/groupby.py", line 1055, in _get_grouper raise ValueError('level > 0 only valid with MultiIndex') ValueError: level > 0 only valid with MultiIndex - William Baker Morrison
@B. M.现在正在使用pandas 0.17.0进行工作，但我仍然收到这个错误。 - William Baker Morrison
你尝试过 lookup[pdmeans.index.values]=pdmeans.values 吗？ - B. M.
显示剩余3条评论
3

这已经在Scipy中实现了，所以你可以这样做：

from scipy.ndimage.measurements import mean as labeled_mean

labels = np.arange(FIDO.max()+1, dtype=int)
S_wb = labeled_mean(WBColor, FIDO, labels)[FIDO]
S_rg = labeled_mean(RGColor, FIDO, labels)[FIDO]
S_by = labeled_mean(BYColor, FIDO, labels)[FIDO]

假设FIDO包含相对较小的整数，那么可以通过np.unique(FIDO, return_inverse=True)进行转换。如果不是这种情况，那么请您自行处理。
当图片大小为200x200且FIDO包含从零到40,000的随机整数时，这段简单的代码比原来的代码快了1000倍。
有没有想法可以避免这个错误？ValueError: operands could not be broadcast together with shapes (38890) (40000) - William Baker Morrison
@WBM - 我很难重现这个错误。你能否提供一些细节，比如这个错误发生在哪一行，涉及到的数组形状是什么？为了测试，我使用了以下代码：FIDO = np.random.randint(40000, size=(200, 200)) 和 WBColor = np.random.rand(200, 200)。 - user2379410
FIDO是200x200，WBColor是200x200（-1到1）。这是我的当前错误：S_wb = labeled_mean(WBColor, FIDO, labels)[FIDO] IndexError: 数组用作索引必须是整数（或布尔）类型。 - William Baker Morrison
@WBM - 如果FIDO已经包含整数，您可以执行FIDO = FIDO.astype(int)，然后它应该可以工作。如果FIDO中的元素不是整数，您可以先执行notused，FIDO_int = np.unique(FIDO, return_inverse=True)，然后执行FIDO_int.shape = (200, 200)，并使用FIDO_int代替FIDO。 - user2379410
3
简而言之：Python 中的循环速度较慢。你可以采取以下措施之一：
  • 向量化（你尝试过，但声称“无法使用”），那你指的是什么意思？如果有可能进行向量化，则始终适用
  • 切换到 Cython，并将迭代器值声明为 int

上述两种方法都基于将瓶颈循环转换为 C 循环。

这个函数看起来很难向量化，所以你可能只能选择第二个选项。 - lejlot
网页内容由stack overflow 提供, 点击上面的
可以查看英文原文，
原文链接