for循环调用numpy 如何优化

如何用卷积的概念 或者其他方式,可以优化掉这两个for循环,加速算法



import numpy as np
data = np.zeros([500,500])
a = []
b = []
for i in range(50):
    d = data[i*10:i*10+10]
    for j in range(450):
        a .append(np.sum(d[:,j:j+50]))
    b.append(a)

举个例子,对于3*3卷积,将8领域值放入通道,然后矩阵乘法,再相加

通过一维卷积 可以优化调一个for循环,但是第一个for循环循环的就是维度不知道怎么优化好

import numpy as np
data = np.zeros([500500])
filter = np.array([1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1])
a = []
b = []
for i in range(500):
    d = data[i:i+1]
    a = np.convolve(d,filter,'valid')
    b.append(a)

省事的方法用numba

from numba import jit
data = np.zeros([500,500])
@jit(nopython=True)
def func(d):
    b = np.zeros((50,450))
    for i in range(50):
        for j in range(450):
            b[i,j] = np.sum(d[(i*10):(i*10+10),j:(j+50)])
    return b
b = func(data)

但是里面的求和其实可以进一步优化,因为有很多是重复求和

在回答问题之前,我想先和题主讨论一下:你期望得到的结果应该是一个50行450列的数组吧?而b却是50行22500列的数组,且每一行指向同一个对象。猜测题主的意图应该如下面代码所示,耗时约110ms。

>>> import numpy as np
>>> import time
>>> def A(data):
    t0 = time.time()
    a = list()
    for i in range(50):
        d = data[i*10:i*10+10]
        b = list()
        for j in range(450):
            b.append(np.sum(d[:,j:j+50]))
        a.append(b)
    t1 = time.time()
    print(t1-t0)

    
>>> data = np.arange(500*500).reshape(500,500)
>>> A(data)
0.11600136756896973

上面这段代码,只要改用数组切片,就可以大幅度提高效率:从110ms提升到12ms。

>>> def B(data):
    t0 = time.time()
    a = list()
    for i in range(0,50,10):
        b = list()
        for j in range(450):
            b.append(np.sum(data[i:i+10,j:j+50]))
        a.append(b)
    t1 = time.time()
    print(t1-t0)
    print(np.array(a).shape)

    
>>> B(data)
0.012001514434814453
(5, 450)

如果想要优化掉两层for循环,有些难度,不过最外层的循环是可以去掉的,速度可以进一步提升,耗时7ms。

>>> def C(data):
    t0 = time.time()
    d = np.sum(np.stack(np.vsplit(data, 50), axis=0), axis=0)
    result = np.zeros((50,450))
    for i in range(450):
        result[:, i] = np.sum(rows_sum[:,i:i+50], axis=1)
    t1 = time.time()
    print(t1-t0)
    print(result.shape)

    
>>> C(data)
0.007000923156738281
(50, 450)