3

I am trying to implement image convolution code using only NumPy, similar to what cv2.filter2D(...) does.

import numpy as np
import time

# kernal
H = np.array([[0,1,0],[1,-4,1],[0,1,0]])

# image
SEED = 23
img = np.random.RandomState(SEED).randint(10, size=(4, 4))

# shapes
Hi, Wi = img.shape
Hk, Wk = H.shape
hk = Hk//2
wk = Wk//2

# padding
new_img = np.pad(img, (hk, wk), 'constant', constant_values=0)
pHi, pWi = new_img.shape

print('img: ')
print(new_img)

print('kernal: ')
print(H)
print('\n')

# image convolution
##################################################
# Method 1
st = time.time()

out = np.zeros((Hi, Wi))
for i in range(hk, pHi-hk):
    for j in range(wk, pWi-wk):
        batch = new_img[i-hk:i+hk+1, j-wk:j+wk+1]
        out[i-hk][j-wk] = np.sum(batch*H)

print('out: ')
print(out)

print("[*] process time : %f" % (time.time()- st))
print('\n')
##################################################
# Method 2
st = time.time()
Hi, Wi = img.shape

out = np.zeros((Hi, Wi))

H_1d = H.ravel()
hl = len(H_1d)

for i in range(Wi):
    img_slice = new_img[:, i:Hk+i]
    for j in range(Hi):
        out[j][i] = np.sum(img_slice[j:j+Hk,:]*H)

    il = len(img_slice)
    h1 = 0
    h2 = hl
    
    while h2 <= il:
        index = h1//Hk
        # print(index)
        out[index][i] = np.sum(img_slice[h1:h2]*H_1d)
        h1 = h1 + Hk
        h2 = h2 + Hk

print('out: ')
print(out)

print("[*] process time : %f" % (time.time()- st))
print('\n')

##################################################

It works but I need a faster, more efficient implementation. In order to make algorithm faster I tried to find pattern in the convolution process looking at image here, but still couldn't figure out what should be done.

Could someone please help me to improve my current code to make it faster, using only NumPy?

Arty
  • 14,883
  • 6
  • 36
  • 69
abeparkk
  • 43
  • 1
  • 4
  • 4
    For a sliding one - https://stackoverflow.com/questions/30109068/implement-matlabs-im2col-sliding-in-python – Divakar Oct 29 '20 at 08:20
  • @abeparkk Implemented several solutions [in my answer](https://stackoverflow.com/a/64660822/941531), some a very-very fast, `2600x` times faster than original convolution in your question. Other solutions are compared in speed. – Arty Jan 10 '21 at 07:01

2 Answers2

5

I've implemented several very fast solutions. Depending on inputs sizes these solutions give speedup up to 2650x times compared to canonical simple convolution implementation!

Note: Not to read whole post/article and just to use boosted function straight away scroll to code below and copy just code of first conv2d_fast(...) function and use it in your code, it contains everything what's needed inside, also for short I copied just code of two best implementations to this code paste.

Maximal gained speedups achieved in my implementations for some of inputs: 1) conv2d_fast - 755x 2) conv2d_medium - 1830x 3) conv2d_fast_numba - 2650x.

Speedups of solutions for image size 460x512 and kernel size 3x3 (this kernel size is most commonly used in computer vision and neural networks for convolutions): 1) conv2d_fast - 52.5x 2) conv2d_medium - 126x 3) conv2d_fast_numba - 319x. On plots shown down below this reference image size and kernel size is signified by red point, other points are blue.

Speedup means how much faster is fast solution compared to base slow solution (canonical simplest convolution implementation). E.g. speedup of 20x means that if base slow solution spent 200ms inside function then fast solution spent 10ms.

Slow base reference function is called conv2d_slow(...) in my code. All timings/speedups and correctness of results are measured against this function. It implements canonical simplest algorithm of computing convolution, as shown in next formula/picture, taken from OpenCV filter2D(...) doc:

conv_formula

Algorithm conv2d_fast is pure NumPy-only, has no Python loops, it is fastest among NumPy-only versions.

conv2d_medium contains Python loops but is also very fast, it is slower than previous algorithm for some inputs (when kernel size is large and close to image size) but faster for others, it is 2x faster on average than conv2d_fast for common kernel sizes. The only reason not to choose this solution is because it is not pure-NumPy, because contains Python loops, also for large kernels it is significantly slower than conv2d_fast.

conv2d_fast_numba is Numba-based, Numba is JIT compiler that converts Python code to pure C++ code and compiles it to machine code, although it wasn't asked to implement non-NumPy solutions I still did this because Numba was created specifically to improve Numpy-based code hence it is closely related to NumPy. Numba solution is the fastest, it is 3-8x faster on average than previous fastest conv2d_medium.

All functions mentioned above need no extra preparation to use them (like importing some modules), just copy-paste single chosen function code and use it in your code.

My code needs installing next pip modules by running next command just once python -m pip install numpy numba matplotlib timerit. But main implementations of convolution conv2d_fast/conv2d_slow/conv2d_medium need just numpy module to be installed.

Try next code online!

# Needs: python -m pip install numpy numba matplotlib timerit

# ---------- Fastest NumPy-only solution ----------

def conv2d_fast(img, krn):
    import numpy as np
    
    is0, is1, ks0, ks1 = *img.shape, *krn.shape
    rs0, rs1 = is0 - ks0 + 1, is1 - ks1 + 1
    
    ix0 = np.arange(ks0)[:, None] + np.arange(rs0)[None, :]
    ix1 = np.arange(ks1)[:, None] + np.arange(rs1)[None, :]
    
    res = krn[:, None, :, None] * img[(ix0.ravel()[:, None], ix1.ravel()[None, :])].reshape(ks0, rs0, ks1, rs1)
    res = res.transpose(1, 3, 0, 2).reshape(rs0, rs1, -1).sum(axis = -1)
    
    return res
    
# ---------- Slowest NumPy+Python solution for reference (base of time measurement) ----------
    
def conv2d_slow(img, krn):
    import numpy as np

    is0, is1, ks0, ks1 = *img.shape, *krn.shape
    rs0, rs1 = is0 - ks0 + 1, is1 - ks1 + 1
    res = np.zeros((rs0, rs1), dtype = krn.dtype)
    
    for i in range(rs0):
        for j in range(rs1):
            res[i, j] = (krn * img[i : i + ks0, j : j + ks1]).sum()
            
    return res
    
# ---------- Medium NumPy+Python solution ----------
    
def conv2d_medium(img, krn):
    import numpy as np

    is0, is1, ks0, ks1 = *img.shape, *krn.shape
    rs0, rs1 = is0 - ks0 + 1, is1 - ks1 + 1
    res = np.zeros((rs0, rs1), dtype = krn.dtype)
    
    for k in range(ks0):
        for l in range(ks1):
            res[...] += krn[k, l] * img[k : k + rs0, l : l + rs1]
            
    return res
    
# ---------- Fastest of all, NumPy+Numba solution ----------
    
def conv2d_fast_numba(img, krn, *, state = {}):
    if 'f' not in state:
        import numpy as np
        
        def conv2d_fast_nm(img, krn):
            is0, is1, ks0, ks1 = *img.shape, *krn.shape
            rs0, rs1 = is0 - ks0 + 1, is1 - ks1 + 1
            res = np.zeros((rs0, rs1), dtype = krn.dtype)
            
            for k in range(ks0):
                for l in range(ks1):
                    for i in range(rs0):
                        for j in range(rs1):
                            res[i, j] += krn[k, l] * img[i + k, j + l]
            
            return res

        import numba
        state['f'] = numba.njit(cache = True, parallel = True, fastmath = True)(conv2d_fast_nm)
        
    return state['f'](img, krn)

# ---------- Testing correctness, measuring speed, drawing plots ----------

def test():    
    import math, matplotlib, matplotlib.pyplot as plt, numpy as np
    from timerit import Timerit
    Timerit._default_asciimode = True
    np.random.seed(0)
    
    for ifuncs, (fname, funcs) in enumerate([
        ('fast', (conv2d_slow, conv2d_fast)),
        ('medium', (conv2d_slow, conv2d_medium)),
        ('numba', (conv2d_slow, conv2d_fast_numba)),
        ('fast_bmed', (conv2d_medium, conv2d_fast)),
        ('numba_bmed', (conv2d_medium, conv2d_fast_numba)),
    ]):
        stats = []

        for krn_size in [1, 2, 3, 4, 6, 8, 12, 16, 24]:
            for img_size in [32, 64, 128, 256, 512]:
                ih, iw = max(1, math.floor(img_size * 0.9)), img_size
                kh, kw = max(1, math.floor(krn_size * 0.9)), krn_size
                if krn_size == 3:
                    kh = 3
                print(f'krn_size ({kh}, {kw}) img_size ({ih}, {iw})', flush = True)
                krn = np.random.uniform(-1., 1., (kh, kw))
                img = np.random.uniform(0., 255., (ih, iw))
                for ifn, f in enumerate(funcs):
                    print(f'{f.__name__}: ', end = '', flush = True)
                    work = ih * iw * kh * kw
                    tim = Timerit(num = min(20, math.ceil(2 ** 15 / work)) * 3, verbose = 1)
                    for i, t in enumerate(tim):
                        r = f(img, krn)
                    rt = tim.mean()
                    if ifn == 0:
                        bt, ba = rt, r
                    else:
                        assert np.allclose(ba, r)
                        print(f'speedup {round(bt / rt, 2)}x')
                        stats.append({
                            'img_size': img_size,
                            'krn_size': krn_size,
                            'speedup': bt / rt,
                        })

        stats = sorted(stats, key = lambda e: e['speedup'])

        x = np.arange(len(stats))
        y = np.array([e['speedup'] for e in stats])

        plt.rcParams['figure.figsize'] = (12.8, 7.2)

        for scale in ['linear', 'log']:
            plt.clf()
            plt.xlabel('iteration')
            plt.ylabel(f'speedup_{fname}_{scale}')
            plt.yscale(scale)
            plt.scatter(x, y, marker = '.', color = [('blue', 'red')[stats[ii]['krn_size'] == 3 and stats[ii]['img_size'] == 512] for ii in range(x.size)])
            for i in range(x.size):
                plt.annotate(
                    (f"k{str(stats[i]['krn_size']).zfill(2)}\ni{str(stats[i]['img_size']).zfill(3)}\n" +
                     f"x{round(stats[i]['speedup'], 2 if stats[i]['speedup'] < 10 else 1 if stats[i]['speedup'] < 100 else None)}"),
                    (x[i], y[i]), fontsize = 'xx-small',
                )
            plt.subplots_adjust(left = 0.07, right = 0.99, bottom = 0.08, top = 0.99)
            plt.xlim(left = -0.1)
            ymin, ymax = np.amin(y), np.amax(y)
            if scale == 'log':
                plt.ylim((ymin / 1.05, ymax * 1.8))
            elif scale == 'linear':
                plt.ylim((ymin - (ymax - ymin) * 0.02, ymax + (ymax - ymin) * 0.08))
                plt.yticks([ymin] + [e for e in plt.yticks()[0] if ymin + 10 ** -6 < e < ymax - 10 ** -6] + [ymax])
                #plt.gca().get_yaxis().set_major_formatter(matplotlib.ticker.FormatStrFormatter('%.1f'))
            plt.savefig(f'conv2d_numpy_{fname}_{scale}.png', dpi = 150)
            plt.show()

if __name__ == '__main__':
    test()

Next plots contain speedups of algorithms for different inputs sizes. x axis shows index of iteration, y axis shows speedup, points are sorted in order of ascending speedup. Each point is labeled kXX iYYY xZZZ, where XX is kernel width, YYY is image width, ZZZ is speedup of current algorithm for these kernel and image sizes. There are two versions of each plot, in linear scale (regular linear y axis) and log scale (y axis is scaled logarithmically). Each plot contains just one red point (the rest are blue), this red point shows most commonly used case, kernel 3x3 and image 460x512, this point should be considered as the most representable point of current graph.

Labels near points can look too small because StackOverflow shows plots at reduced resolution, open images in new browser tab and zoom in to see them in full 1920x1080 resolution.

  1. conv2d_fast linear:

fast_linear

  1. conv2d_fast log:

fast_log

  1. conv2d_medium linear:

medium_linear

  1. conv2d_medium log:

medium_log

  1. conv2d_fast_numba linear:

numba_linear

  1. conv2d_fast_numba log:

numba_log

  1. conv2d_fast compared to conv2d_medium, linear:

fast_bmedium_linear

  1. conv2d_fast compared to conv2d_medium, log:

fast_bmedium_log

  1. conv2d_fast_numba compared to conv2d_medium, linear:

numba_bmedium_linear

  1. conv2d_fast_numba compared to conv2d_medium, log:

numba_bmedium_log

Console output of code above:

krn_size (1, 1) img_size (28, 32)
conv2d_slow: Timed best=15.037 ms, mean=15.409 +- 0.2 ms
conv2d_fast: Timed best=100.701 us, mean=102.412 +- 1.4 us
speedup 150.46x
krn_size (1, 1) img_size (57, 64)
conv2d_slow: Timed best=61.679 ms, mean=65.878 +- 6.0 ms
conv2d_fast: Timed best=151.541 us, mean=192.332 +- 28.5 us
speedup 342.52x
krn_size (1, 1) img_size (115, 128)
conv2d_slow: Timed best=249.408 ms, mean=252.172 +- 2.4 ms
conv2d_fast: Timed best=328.012 us, mean=333.878 +- 6.6 us
speedup 755.28x
krn_size (1, 1) img_size (230, 256)
conv2d_slow: Timed best=1.013 s, mean=1.013 +- 0.0 s
conv2d_fast: Timed best=1.706 ms, mean=1.706 +- 0.0 ms
speedup 593.72x
krn_size (1, 1) img_size (460, 512)
conv2d_slow: Timed best=4.101 s, mean=4.101 +- 0.0 s
conv2d_fast: Timed best=10.182 ms, mean=10.182 +- 0.0 ms
speedup 402.72x
krn_size (1, 2) img_size (28, 32)
conv2d_slow: Timed best=14.418 ms, mean=15.877 +- 1.8 ms
conv2d_fast: Timed best=132.476 us, mean=135.178 +- 2.1 us
speedup 117.45x
krn_size (1, 2) img_size (57, 64)
conv2d_slow: Timed best=60.057 ms, mean=61.136 +- 1.3 ms
conv2d_fast: Timed best=231.222 us, mean=238.359 +- 11.6 us
speedup 256.49x
krn_size (1, 2) img_size (115, 128)
conv2d_slow: Timed best=264.794 ms, mean=268.452 +- 3.7 ms
conv2d_fast: Timed best=674.112 us, mean=729.595 +- 55.5 us
speedup 367.95x
krn_size (1, 2) img_size (230, 256)
conv2d_slow: Timed best=999.722 ms, mean=999.722 +- 0.0 ms
conv2d_fast: Timed best=4.910 ms, mean=4.910 +- 0.0 ms
speedup 203.61x
krn_size (1, 2) img_size (460, 512)
conv2d_slow: Timed best=4.360 s, mean=4.360 +- 0.0 s
conv2d_fast: Timed best=18.148 ms, mean=18.148 +- 0.0 ms
speedup 240.24x
krn_size (3, 3) img_size (28, 32)
conv2d_slow: Timed best=16.771 ms, mean=17.395 +- 0.7 ms
conv2d_fast: Timed best=256.153 us, mean=268.080 +- 13.5 us
speedup 64.89x
krn_size (3, 3) img_size (57, 64)
conv2d_slow: Timed best=72.964 ms, mean=72.964 +- 0.0 ms
conv2d_fast: Timed best=863.783 us, mean=863.783 +- 0.0 us
speedup 84.47x
krn_size (3, 3) img_size (115, 128)
conv2d_slow: Timed best=302.496 ms, mean=302.496 +- 0.0 ms
conv2d_fast: Timed best=6.490 ms, mean=6.490 +- 0.0 ms
speedup 46.61x
krn_size (3, 3) img_size (230, 256)
conv2d_slow: Timed best=1.281 s, mean=1.281 +- 0.0 s
conv2d_fast: Timed best=24.801 ms, mean=24.801 +- 0.0 ms
speedup 51.66x
krn_size (3, 3) img_size (460, 512)
conv2d_slow: Timed best=5.134 s, mean=5.134 +- 0.0 s
conv2d_fast: Timed best=97.726 ms, mean=97.726 +- 0.0 ms
speedup 52.53x
krn_size (3, 4) img_size (28, 32)
conv2d_slow: Timed best=16.275 ms, mean=16.351 +- 0.1 ms
conv2d_fast: Timed best=286.950 us, mean=295.749 +- 11.9 us
speedup 55.29x
krn_size (3, 4) img_size (57, 64)
conv2d_slow: Timed best=72.170 ms, mean=72.170 +- 0.0 ms
conv2d_fast: Timed best=1.824 ms, mean=1.824 +- 0.0 ms
speedup 39.57x
krn_size (3, 4) img_size (115, 128)
conv2d_slow: Timed best=305.169 ms, mean=305.169 +- 0.0 ms
conv2d_fast: Timed best=8.462 ms, mean=8.462 +- 0.0 ms
speedup 36.06x
krn_size (3, 4) img_size (230, 256)
conv2d_slow: Timed best=1.245 s, mean=1.245 +- 0.0 s
conv2d_fast: Timed best=31.527 ms, mean=31.527 +- 0.0 ms
speedup 39.5x
krn_size (3, 4) img_size (460, 512)
conv2d_slow: Timed best=5.262 s, mean=5.262 +- 0.0 s
conv2d_fast: Timed best=128.232 ms, mean=128.232 +- 0.0 ms
speedup 41.03x
krn_size (5, 6) img_size (28, 32)
conv2d_slow: Timed best=14.060 ms, mean=14.507 +- 0.4 ms
conv2d_fast: Timed best=469.288 us, mean=478.087 +- 8.8 us
speedup 30.34x
krn_size (5, 6) img_size (57, 64)
conv2d_slow: Timed best=67.638 ms, mean=67.638 +- 0.0 ms
conv2d_fast: Timed best=3.542 ms, mean=3.542 +- 0.0 ms
speedup 19.1x
krn_size (5, 6) img_size (115, 128)
conv2d_slow: Timed best=299.806 ms, mean=299.806 +- 0.0 ms
conv2d_fast: Timed best=18.730 ms, mean=18.730 +- 0.0 ms
speedup 16.01x
krn_size (5, 6) img_size (230, 256)
conv2d_slow: Timed best=1.294 s, mean=1.294 +- 0.0 s
conv2d_fast: Timed best=77.809 ms, mean=77.809 +- 0.0 ms
speedup 16.64x
krn_size (5, 6) img_size (460, 512)
conv2d_slow: Timed best=5.336 s, mean=5.336 +- 0.0 s
conv2d_fast: Timed best=317.518 ms, mean=317.518 +- 0.0 ms
speedup 16.81x
krn_size (7, 8) img_size (28, 32)
conv2d_slow: Timed best=12.076 ms, mean=12.076 +- 0.0 ms
conv2d_fast: Timed best=861.827 us, mean=861.827 +- 0.0 us
speedup 14.01x
krn_size (7, 8) img_size (57, 64)
conv2d_slow: Timed best=64.761 ms, mean=64.761 +- 0.0 ms
conv2d_fast: Timed best=8.889 ms, mean=8.889 +- 0.0 ms
speedup 7.29x
krn_size (7, 8) img_size (115, 128)
conv2d_slow: Timed best=293.768 ms, mean=293.768 +- 0.0 ms
conv2d_fast: Timed best=32.908 ms, mean=32.908 +- 0.0 ms
speedup 8.93x
krn_size (7, 8) img_size (230, 256)
conv2d_slow: Timed best=1.245 s, mean=1.245 +- 0.0 s
conv2d_fast: Timed best=139.752 ms, mean=139.752 +- 0.0 ms
speedup 8.91x
krn_size (7, 8) img_size (460, 512)
conv2d_slow: Timed best=5.535 s, mean=5.535 +- 0.0 s
conv2d_fast: Timed best=599.906 ms, mean=599.906 +- 0.0 ms
speedup 9.23x
krn_size (10, 12) img_size (28, 32)
conv2d_slow: Timed best=8.776 ms, mean=8.776 +- 0.0 ms
conv2d_fast: Timed best=1.801 ms, mean=1.801 +- 0.0 ms
speedup 4.87x
krn_size (10, 12) img_size (57, 64)
conv2d_slow: Timed best=55.155 ms, mean=55.155 +- 0.0 ms
conv2d_fast: Timed best=13.861 ms, mean=13.861 +- 0.0 ms
speedup 3.98x
krn_size (10, 12) img_size (115, 128)
conv2d_slow: Timed best=275.688 ms, mean=275.688 +- 0.0 ms
conv2d_fast: Timed best=65.345 ms, mean=65.345 +- 0.0 ms
speedup 4.22x
krn_size (10, 12) img_size (230, 256)
conv2d_slow: Timed best=1.215 s, mean=1.215 +- 0.0 s
conv2d_fast: Timed best=319.263 ms, mean=319.263 +- 0.0 ms
speedup 3.8x
krn_size (10, 12) img_size (460, 512)
conv2d_slow: Timed best=5.413 s, mean=5.413 +- 0.0 s
conv2d_fast: Timed best=1.360 s, mean=1.360 +- 0.0 s
speedup 3.98x
krn_size (14, 16) img_size (28, 32)
conv2d_slow: Timed best=6.660 ms, mean=6.660 +- 0.0 ms
conv2d_fast: Timed best=2.498 ms, mean=2.498 +- 0.0 ms
speedup 2.67x
krn_size (14, 16) img_size (57, 64)
conv2d_slow: Timed best=49.958 ms, mean=49.958 +- 0.0 ms
conv2d_fast: Timed best=24.663 ms, mean=24.663 +- 0.0 ms
speedup 2.03x
krn_size (14, 16) img_size (115, 128)
conv2d_slow: Timed best=273.521 ms, mean=273.521 +- 0.0 ms
conv2d_fast: Timed best=120.138 ms, mean=120.138 +- 0.0 ms
speedup 2.28x
krn_size (14, 16) img_size (230, 256)
conv2d_slow: Timed best=1.329 s, mean=1.329 +- 0.0 s
conv2d_fast: Timed best=537.025 ms, mean=537.025 +- 0.0 ms
speedup 2.47x
krn_size (14, 16) img_size (460, 512)
conv2d_slow: Timed best=5.615 s, mean=5.615 +- 0.0 s
conv2d_fast: Timed best=2.623 s, mean=2.623 +- 0.0 s
speedup 2.14x
krn_size (21, 24) img_size (28, 32)
conv2d_slow: Timed best=2.203 ms, mean=2.203 +- 0.0 ms
conv2d_fast: Timed best=955.684 us, mean=955.684 +- 0.0 us
speedup 2.3x
krn_size (21, 24) img_size (57, 64)
conv2d_slow: Timed best=41.180 ms, mean=41.180 +- 0.0 ms
conv2d_fast: Timed best=40.526 ms, mean=40.526 +- 0.0 ms
speedup 1.02x
krn_size (21, 24) img_size (115, 128)
conv2d_slow: Timed best=288.575 ms, mean=288.575 +- 0.0 ms
conv2d_fast: Timed best=277.724 ms, mean=277.724 +- 0.0 ms
speedup 1.04x
krn_size (21, 24) img_size (230, 256)
conv2d_slow: Timed best=1.334 s, mean=1.334 +- 0.0 s
conv2d_fast: Timed best=1.288 s, mean=1.288 +- 0.0 s
speedup 1.04x
krn_size (21, 24) img_size (460, 512)
conv2d_slow: Timed best=6.351 s, mean=6.351 +- 0.0 s
conv2d_fast: Timed best=5.996 s, mean=5.996 +- 0.0 s
speedup 1.06x
krn_size (1, 1) img_size (28, 32)
conv2d_slow: Timed best=15.644 ms, mean=16.388 +- 0.9 ms
conv2d_medium: Timed best=34.708 us, mean=36.174 +- 1.8 us
speedup 453.03x
krn_size (1, 1) img_size (57, 64)
conv2d_slow: Timed best=62.969 ms, mean=77.952 +- 15.0 ms
conv2d_medium: Timed best=54.261 us, mean=62.300 +- 4.8 us
speedup 1251.24x
krn_size (1, 1) img_size (115, 128)
conv2d_slow: Timed best=285.944 ms, mean=487.718 +- 148.2 ms
conv2d_medium: Timed best=246.865 us, mean=266.581 +- 14.8 us
speedup 1829.53x
krn_size (1, 1) img_size (230, 256)
conv2d_slow: Timed best=1.474 s, mean=1.474 +- 0.0 s
conv2d_medium: Timed best=1.400 ms, mean=1.400 +- 0.0 ms
speedup 1052.98x
krn_size (1, 1) img_size (460, 512)
conv2d_slow: Timed best=4.434 s, mean=4.434 +- 0.0 s
conv2d_medium: Timed best=5.222 ms, mean=5.222 +- 0.0 ms
speedup 849.22x
krn_size (1, 2) img_size (28, 32)
conv2d_slow: Timed best=14.759 ms, mean=15.430 +- 1.2 ms
conv2d_medium: Timed best=69.904 us, mean=77.803 +- 6.8 us
speedup 198.33x
krn_size (1, 2) img_size (57, 64)
conv2d_slow: Timed best=60.825 ms, mean=62.482 +- 1.8 ms
conv2d_medium: Timed best=138.342 us, mean=139.711 +- 1.3 us
speedup 447.22x
krn_size (1, 2) img_size (115, 128)
conv2d_slow: Timed best=261.207 ms, mean=263.777 +- 2.6 ms
conv2d_medium: Timed best=343.167 us, mean=352.699 +- 9.5 us
speedup 747.88x
krn_size (1, 2) img_size (230, 256)
conv2d_slow: Timed best=1.031 s, mean=1.031 +- 0.0 s
conv2d_medium: Timed best=1.887 ms, mean=1.887 +- 0.0 ms
speedup 546.21x
krn_size (1, 2) img_size (460, 512)
conv2d_slow: Timed best=4.126 s, mean=4.126 +- 0.0 s
conv2d_medium: Timed best=10.341 ms, mean=10.341 +- 0.0 ms
speedup 398.96x
krn_size (3, 3) img_size (28, 32)
conv2d_slow: Timed best=16.776 ms, mean=16.975 +- 0.2 ms
conv2d_medium: Timed best=283.528 us, mean=287.829 +- 6.3 us
speedup 58.98x
krn_size (3, 3) img_size (57, 64)
conv2d_slow: Timed best=73.292 ms, mean=73.292 +- 0.0 ms
conv2d_medium: Timed best=636.472 us, mean=636.472 +- 0.0 us
speedup 115.15x
krn_size (3, 3) img_size (115, 128)
conv2d_slow: Timed best=308.051 ms, mean=308.051 +- 0.0 ms
conv2d_medium: Timed best=1.198 ms, mean=1.198 +- 0.0 ms
speedup 257.21x
krn_size (3, 3) img_size (230, 256)
conv2d_slow: Timed best=1.277 s, mean=1.277 +- 0.0 s
conv2d_medium: Timed best=7.335 ms, mean=7.335 +- 0.0 ms
speedup 174.1x
krn_size (3, 3) img_size (460, 512)
conv2d_slow: Timed best=5.306 s, mean=5.306 +- 0.0 s
conv2d_medium: Timed best=42.131 ms, mean=42.131 +- 0.0 ms
speedup 125.94x
krn_size (3, 4) img_size (28, 32)
conv2d_slow: Timed best=16.628 ms, mean=18.073 +- 2.1 ms
conv2d_medium: Timed best=378.363 us, mean=388.629 +- 9.5 us
speedup 46.51x
krn_size (3, 4) img_size (57, 64)
conv2d_slow: Timed best=71.848 ms, mean=71.848 +- 0.0 ms
conv2d_medium: Timed best=570.966 us, mean=570.966 +- 0.0 us
speedup 125.84x
krn_size (3, 4) img_size (115, 128)
conv2d_slow: Timed best=307.969 ms, mean=307.969 +- 0.0 ms
conv2d_medium: Timed best=2.292 ms, mean=2.292 +- 0.0 ms
speedup 134.36x
krn_size (3, 4) img_size (230, 256)
conv2d_slow: Timed best=1.419 s, mean=1.419 +- 0.0 s
conv2d_medium: Timed best=9.847 ms, mean=9.847 +- 0.0 ms
speedup 144.14x
krn_size (3, 4) img_size (460, 512)
conv2d_slow: Timed best=5.808 s, mean=5.808 +- 0.0 s
conv2d_medium: Timed best=73.120 ms, mean=73.120 +- 0.0 ms
speedup 79.44x
krn_size (5, 6) img_size (28, 32)
conv2d_slow: Timed best=14.524 ms, mean=14.705 +- 0.2 ms
conv2d_medium: Timed best=918.044 us, mean=930.998 +- 13.0 us
speedup 15.79x
krn_size (5, 6) img_size (57, 64)
conv2d_slow: Timed best=71.519 ms, mean=71.519 +- 0.0 ms
conv2d_medium: Timed best=1.415 ms, mean=1.415 +- 0.0 ms
speedup 50.54x
krn_size (5, 6) img_size (115, 128)
conv2d_slow: Timed best=299.480 ms, mean=299.480 +- 0.0 ms
conv2d_medium: Timed best=3.706 ms, mean=3.706 +- 0.0 ms
speedup 80.8x
krn_size (5, 6) img_size (230, 256)
conv2d_slow: Timed best=1.647 s, mean=1.647 +- 0.0 s
conv2d_medium: Timed best=38.446 ms, mean=38.446 +- 0.0 ms
speedup 42.83x
krn_size (5, 6) img_size (460, 512)
conv2d_slow: Timed best=5.496 s, mean=5.496 +- 0.0 s
conv2d_medium: Timed best=136.823 ms, mean=136.823 +- 0.0 ms
speedup 40.17x
krn_size (7, 8) img_size (28, 32)
conv2d_slow: Timed best=16.223 ms, mean=16.223 +- 0.0 ms
conv2d_medium: Timed best=2.033 ms, mean=2.033 +- 0.0 ms
speedup 7.98x
krn_size (7, 8) img_size (57, 64)
conv2d_slow: Timed best=82.617 ms, mean=82.617 +- 0.0 ms
conv2d_medium: Timed best=3.206 ms, mean=3.206 +- 0.0 ms
speedup 25.77x
krn_size (7, 8) img_size (115, 128)
conv2d_slow: Timed best=304.944 ms, mean=304.944 +- 0.0 ms
conv2d_medium: Timed best=6.716 ms, mean=6.716 +- 0.0 ms
speedup 45.4x
krn_size (7, 8) img_size (230, 256)
conv2d_slow: Timed best=1.295 s, mean=1.295 +- 0.0 s
conv2d_medium: Timed best=40.545 ms, mean=40.545 +- 0.0 ms
speedup 31.94x
krn_size (7, 8) img_size (460, 512)
conv2d_slow: Timed best=5.400 s, mean=5.400 +- 0.0 s
conv2d_medium: Timed best=260.807 ms, mean=260.807 +- 0.0 ms
speedup 20.71x
krn_size (10, 12) img_size (28, 32)
conv2d_slow: Timed best=9.135 ms, mean=9.135 +- 0.0 ms
conv2d_medium: Timed best=3.154 ms, mean=3.154 +- 0.0 ms
speedup 2.9x
krn_size (10, 12) img_size (57, 64)
conv2d_slow: Timed best=57.130 ms, mean=57.130 +- 0.0 ms
conv2d_medium: Timed best=4.768 ms, mean=4.768 +- 0.0 ms
speedup 11.98x
krn_size (10, 12) img_size (115, 128)
conv2d_slow: Timed best=277.271 ms, mean=277.271 +- 0.0 ms
conv2d_medium: Timed best=12.841 ms, mean=12.841 +- 0.0 ms
speedup 21.59x
krn_size (10, 12) img_size (230, 256)
conv2d_slow: Timed best=1.221 s, mean=1.221 +- 0.0 s
conv2d_medium: Timed best=85.651 ms, mean=85.651 +- 0.0 ms
speedup 14.25x
krn_size (10, 12) img_size (460, 512)
conv2d_slow: Timed best=5.475 s, mean=5.475 +- 0.0 s
conv2d_medium: Timed best=552.321 ms, mean=552.321 +- 0.0 ms
speedup 9.91x
krn_size (14, 16) img_size (28, 32)
conv2d_slow: Timed best=7.463 ms, mean=7.463 +- 0.0 ms
conv2d_medium: Timed best=7.304 ms, mean=7.304 +- 0.0 ms
speedup 1.02x
krn_size (14, 16) img_size (57, 64)
conv2d_slow: Timed best=63.592 ms, mean=63.592 +- 0.0 ms
conv2d_medium: Timed best=11.538 ms, mean=11.538 +- 0.0 ms
speedup 5.51x
krn_size (14, 16) img_size (115, 128)
conv2d_slow: Timed best=268.824 ms, mean=268.824 +- 0.0 ms
conv2d_medium: Timed best=34.761 ms, mean=34.761 +- 0.0 ms
speedup 7.73x
krn_size (14, 16) img_size (230, 256)
conv2d_slow: Timed best=1.398 s, mean=1.398 +- 0.0 s
conv2d_medium: Timed best=153.959 ms, mean=153.959 +- 0.0 ms
speedup 9.08x
krn_size (14, 16) img_size (460, 512)
conv2d_slow: Timed best=6.059 s, mean=6.059 +- 0.0 s
conv2d_medium: Timed best=1.063 s, mean=1.063 +- 0.0 s
speedup 5.7x
krn_size (21, 24) img_size (28, 32)
conv2d_slow: Timed best=1.928 ms, mean=1.928 +- 0.0 ms
conv2d_medium: Timed best=10.410 ms, mean=10.410 +- 0.0 ms
speedup 0.19x
krn_size (21, 24) img_size (57, 64)
conv2d_slow: Timed best=43.078 ms, mean=43.078 +- 0.0 ms
conv2d_medium: Timed best=21.525 ms, mean=21.525 +- 0.0 ms
speedup 2.0x
krn_size (21, 24) img_size (115, 128)
conv2d_slow: Timed best=265.293 ms, mean=265.293 +- 0.0 ms
conv2d_medium: Timed best=47.004 ms, mean=47.004 +- 0.0 ms
speedup 5.64x
krn_size (21, 24) img_size (230, 256)
conv2d_slow: Timed best=1.384 s, mean=1.384 +- 0.0 s
conv2d_medium: Timed best=335.177 ms, mean=335.177 +- 0.0 ms
speedup 4.13x
krn_size (21, 24) img_size (460, 512)
conv2d_slow: Timed best=6.308 s, mean=6.308 +- 0.0 s
conv2d_medium: Timed best=2.224 s, mean=2.224 +- 0.0 s
speedup 2.84x
krn_size (1, 1) img_size (28, 32)
conv2d_slow: Timed best=16.040 ms, mean=16.287 +- 0.3 ms
conv2d_fast_numba: Timed best=31.285 us, mean=33.168 +- 1.6 us
speedup 491.05x
krn_size (1, 1) img_size (57, 64)
conv2d_slow: Timed best=63.740 ms, mean=65.235 +- 1.0 ms
conv2d_fast_numba: Timed best=41.062 us, mean=42.801 +- 4.1 us
speedup 1524.16x
krn_size (1, 1) img_size (115, 128)
conv2d_slow: Timed best=259.237 ms, mean=261.639 +- 2.9 ms
conv2d_fast_numba: Timed best=88.969 us, mean=98.746 +- 12.5 us
speedup 2649.63x
krn_size (1, 1) img_size (230, 256)
conv2d_slow: Timed best=1.046 s, mean=1.046 +- 0.0 s
conv2d_fast_numba: Timed best=406.715 us, mean=406.715 +- 0.0 us
speedup 2571.19x
krn_size (1, 1) img_size (460, 512)
conv2d_slow: Timed best=4.335 s, mean=4.335 +- 0.0 s
conv2d_fast_numba: Timed best=2.957 ms, mean=2.957 +- 0.0 ms
speedup 1465.91x
krn_size (1, 2) img_size (28, 32)
conv2d_slow: Timed best=14.989 ms, mean=15.260 +- 0.2 ms
conv2d_fast_numba: Timed best=34.707 us, mean=35.840 +- 0.8 us
speedup 425.79x
krn_size (1, 2) img_size (57, 64)
conv2d_slow: Timed best=60.887 ms, mean=61.745 +- 0.7 ms
conv2d_fast_numba: Timed best=46.440 us, mean=49.959 +- 5.4 us
speedup 1235.9x
krn_size (1, 2) img_size (115, 128)
conv2d_slow: Timed best=258.905 ms, mean=259.000 +- 0.1 ms
conv2d_fast_numba: Timed best=123.188 us, mean=139.320 +- 16.1 us
speedup 1859.04x
krn_size (1, 2) img_size (230, 256)
conv2d_slow: Timed best=1.040 s, mean=1.040 +- 0.0 s
conv2d_fast_numba: Timed best=547.991 us, mean=547.991 +- 0.0 us
speedup 1897.39x
krn_size (1, 2) img_size (460, 512)
conv2d_slow: Timed best=4.180 s, mean=4.180 +- 0.0 s
conv2d_fast_numba: Timed best=4.705 ms, mean=4.705 +- 0.0 ms
speedup 888.42x
krn_size (3, 3) img_size (28, 32)
conv2d_slow: Timed best=17.002 ms, mean=17.359 +- 0.3 ms
conv2d_fast_numba: Timed best=44.973 us, mean=45.853 +- 1.3 us
speedup 378.58x
krn_size (3, 3) img_size (57, 64)
conv2d_slow: Timed best=73.428 ms, mean=73.428 +- 0.0 ms
conv2d_fast_numba: Timed best=101.190 us, mean=101.190 +- 0.0 us
speedup 725.65x
krn_size (3, 3) img_size (115, 128)
conv2d_slow: Timed best=317.337 ms, mean=317.337 +- 0.0 ms
conv2d_fast_numba: Timed best=326.058 us, mean=326.058 +- 0.0 us
speedup 973.25x
krn_size (3, 3) img_size (230, 256)
conv2d_slow: Timed best=1.311 s, mean=1.311 +- 0.0 s
conv2d_fast_numba: Timed best=1.428 ms, mean=1.428 +- 0.0 ms
speedup 917.72x
krn_size (3, 3) img_size (460, 512)
conv2d_slow: Timed best=5.461 s, mean=5.461 +- 0.0 s
conv2d_fast_numba: Timed best=17.109 ms, mean=17.109 +- 0.0 ms
speedup 319.22x
......................
Arty
  • 14,883
  • 6
  • 36
  • 69
1

Here is a faster method using strides (note that view_as_windows uses numpy strides under the hood. If you have to strictly use numpy, simply use strides from numpy package. You can refer to @Divakar's answer for equivalent implementation of view_as_windows in numpy):

from skimage.util import view_as_windows
out = view_as_windows(new_img, (Hk,Wk))
out = np.einsum('ijkl,kl->ij',out,H)

output:

out: 
[[  0  -5 -10 -19]
 [-10  -7  -2 -18]
 [  5  -7  11   7]
 [-12  -9  11 -18]]

timing for your input:

#Method 1
[*] process time : 0.000615

#Method 2
[*] process time : 0.000475

#proposed method in this answer
[*] process time : 0.000348

timing for larger arrays using benchit:

(Method 2 has errors for larger arrays)

Method 3 is orders of magnitude faster than mehtod 1.

enter image description here

Ehsan
  • 12,072
  • 2
  • 20
  • 33