Matlab vs Python 2D convolution performance

Question

I like to prototype algorithms in Matlab, but I have the requirement of putting them on a server that also runs quite a bit of Python code. Hence I quickly converted the code to Python and compared the two. The Matlab implementation runs ~1000 times faster (from timing function calls - no profiling). Anyone know off hand why the performance of Python is so slow?

Matlab

% init random data
w = 800;
h = 1200;
hmap = zeros(w,h);

npts = 250;
for i=1:npts
    hmap(randi(w),randi(h)) = hmap(randi(w),randi(h))+1;
end

% Params
disksize = 251;
nBreaks = 25;
saturation = .9;
floorthresh =.05;

fh = fspecial('gaussian', disksize, disksize/7);
hmap = conv2(hmap, fh, 'same');

% Scaling, paritioning etc
hmap = hmap/(max(max(hmap)));
hmap(hmap<floorthresh) = 0;
hmap = round(nBreaks * hmap)/nBreaks;
hmap = hmap * (1/saturation);

% Show the image
imshow(hmap, [0,1])
colormap('jet')

Python

import numpy as np
from scipy.signal import convolve2d as conv2

# Test data parameters
w = 800
h = 1200
npts = 250

# generate data
xvals = np.random.randint(w, size=npts)
yvals = np.random.randint(h, size=npts)

# Heatmap parameters
gaussianSize = 250
nbreaks = 25

# Preliminary function definitions
def populateMat(w, h, xvals, yvals):
    container = np.zeros((w,h))

    for idx in range(0,xvals.size):
        x = xvals[idx]
        y = yvals[idx]
        container[x,y] += 1

    return container

def makeGaussian(size, fwhm):

    x = np.arange(0, size, 1, float)
    y = x[:,np.newaxis]

    x0 = y0 = size // 2

    return np.exp(-4*np.log(2) * ((x-x0)**2 + (y-y0)**2) / fwhm**2)


# Create the data matrix
dmat = populateMat(w,h,xvals,yvals)
h = makeGaussian(gaussianSize, fwhm=gaussianSize/2)

# Convolve
dmat2 = conv2(dmat, h, mode='same')

# Scaling etc
dmat2 = dmat2 / dmat2.max()
dmat2 = np.round(nbreaks*dmat2)/nbreaks
# Show
imshow(dmat2)

Use numpy for the convolution. Have a look at stackoverflow.com/questions/16121269/… — kezzos
– kezzos, Commented Nov 12, 2014 at 9:48
are you saying that the implementation in scipy.signal.convolve2d is the culprit? Will investigate. — Joe
– Joe, Commented Nov 12, 2014 at 9:50
I read the post and it seems to say that the only reason to use numpy is to avoid the scipy dependency. I don't have that restriction. Tried ndimage.convolve but this throws a MemoryError. Will keep trying the differnt options suggested in that thread. — Joe
– Joe, Commented Nov 12, 2014 at 9:59
My guess is that matlab uses an optimized algorithm for the Gaussian filter, exploiting separability of the the filter and possibly using a recursive approximation. Try using separability under Python. — user1196549
– user1196549, Commented Nov 12, 2014 at 10:11
This is trivial: define two 1D Gaussian, horizontal and vertical and apply them in sequence. — user1196549
– user1196549, Commented Nov 12, 2014 at 10:13

Joe · Accepted Answer · 2014-11-12 10:29:12Z

Ok, problem solved for me thanks to suggestion from @Yves Daust's comments;

The filter scipy.ndimage.filters.gaussian_filter utilises the separability of the kernel and reduces the running time to within a single order of magnitude of the matlab implementation.

import numpy as np
from scipy.ndimage.filters import gaussian_filter as gaussian


# Test data parameters
w = 800
h = 1200
npts = 250

# generate data
xvals = np.random.randint(w, size=npts)
yvals = np.random.randint(h, size=npts)

# Heatmap parameters
gaussianSize = 250
nbreaks = 25

# Preliminary function definitions
def populateMat(w, h, xvals, yvals):
    container = np.zeros((w,h))

    for idx in range(0,xvals.size):
        x = xvals[idx]
        y = yvals[idx]
        container[x,y] += 1

    return container


# Create the data matrix
dmat = populateMat(w,h,xvals,yvals)

# Convolve
dmat2 = gaussian(dmat, gaussianSize/7)

# Scaling etc
dmat2 = dmat2 / dmat2.max()
dmat2 = np.round(nbreaks*dmat2)/nbreaks

# Show
imshow(dmat2)

Collectives™ on Stack Overflow

Matlab vs Python 2D convolution performance

Matlab

Python

1 Answer 1

Comments

Your Answer

Linked

Hot Network Questions

Collectives™ on Stack Overflow

Matlab

Python

1 Answer 1

Comments

Your Answer

Sign up or log in

Post as a guest

Linked

Related