Ist die Faltung in Numpy langsamer als in Matlab?

Question

Faltung in Matlab scheint doppelt so schnell zu sein wie Faltung in Numpy.

Python-Code (dauert 19 Sekunden auf meinem Rechner):

import numpy as np 
from scipy import ndimage 
import time 

img = np.ones((512,512,512)) 
kernel = np.ones((5,5,5))/125 

start_time = time.time() 
ndimage.convolve(img,kernel,mode='constant') 
print "Numpy execution took ", (time.time() - start_time), "seconds" 

Matlab-Code (dauert 8,7 Sekunden auf meinem Rechner):

img = ones(512,512,512); 
kernel = ones(5,5,5)/125; 
tic 
convn(img, kernel, 'same'); 
toc 

Die zwei identische Ergebnisse.

Gibt es eine Möglichkeit, Numpy zu verbessern, um Matlabs Leistung hier zu erreichen oder zu übertreffen?

Interessanterweise ist dieser Faktor oder ~ 2 Unterschied in den Laufzeiten konsistent bei vielen Eingabegrößen.

Answer 1

Keine Antwort; nur einen Kommentar:

Bevor Sie die Leistung vergleichen zu können, müssen Sie die beiden Funktionen das gleiche Ergebnis stellen Sie sicher zurückkommen:

Wenn Matlab convn das gleiche Ergebnis wie Oktaves convn zurückgibt, dann convn anders als ndimage.convolve:

octave> convn(ones(3,3), ones(2,2)) 
ans = 

    1 2 2 1 
    2 4 4 2 
    2 4 4 2 
    1 2 2 1 

In [99]: ndimage.convolve(np.ones((3,3)), np.ones((2,2))) 
Out[99]: 
array([[ 4., 4., 4.], 
     [ 4., 4., 4.], 
     [ 4., 4., 4.]]) 

ndimage.convolve hat andere Modi, 'reflect','constant','nearest','mirror', 'wrap', aber keine von diesem Spiel convn ‚s default ("full") Verhalten.

---

Für 2D-Arrays scipy.signal.convolve2d ist schneller als scipy.signal.convolve.

Für 3D-Arrays, scheint scipy.signal.convolve das gleiche Verhalten wie convn(A,B) zu haben: standardmäßig

octave> x = convn(ones(3,3,3), ones(2,2,2)) 
x = 

ans(:,:,1) = 

    1 2 2 1 
    2 4 4 2 
    2 4 4 2 
    1 2 2 1 

ans(:,:,2) = 

    2 4 4 2 
    4 8 8 4 
    4 8 8 4 
    2 4 4 2 

ans(:,:,3) = 

    2 4 4 2 
    4 8 8 4 
    4 8 8 4 
    2 4 4 2 

ans(:,:,4) = 

    1 2 2 1 
    2 4 4 2 
    2 4 4 2 
    1 2 2 1 

In [109]: signal.convolve(np.ones((3,3,3), dtype='uint8'), np.ones((2,2,2), dtype='uint8')) 
Out[109]: 
array([[[1, 2, 2, 1], 
     [2, 4, 4, 2], 
     [2, 4, 4, 2], 
     [1, 2, 2, 1]], 

     [[2, 4, 4, 2], 
     [4, 8, 8, 4], 
     [4, 8, 8, 4], 
     [2, 4, 4, 2]], 

     [[2, 4, 4, 2], 
     [4, 8, 8, 4], 
     [4, 8, 8, 4], 
     [2, 4, 4, 2]], 

     [[1, 2, 2, 1], 
     [2, 4, 4, 2], 
     [2, 4, 4, 2], 
     [1, 2, 2, 1]]], dtype=uint8) 

---

Beachten Sie, dass np.ones((n,m,p)) schafft eine Schwimmer Array. Matlab ones(n,m,p) scheint ein Array von Ints zu erstellen. Um einen guten Vergleich zu machen, sollten Sie versuchen, den dtype der numpy Arrays dem Typ der Matlab-Matrizen anzupassen.

Answer 2

Welche genaue Operation machst du? Es gibt eine Reihe von Optimierungen, die ndimage bietet, wenn Sie keine allgemeine N-d-Faltung benötigen.

Zum Beispiel Ihre aktuelle Operation:

img = np.ones((512,512,512)) 
kernel = np.ones((5,5,5))/125 
result = ndimage.convolve(img, kernel) 

entspricht:

img = np.ones((512,512,512)) 
result = ndimage.uniform_filter(img, 5) 

Allerdings gibt es einen großen Unterschied in der Ausführungsgeschwindigkeit:

In [22]: %timeit ndimage.convolve(img, kernel) 
1 loops, best of 3: 25.9 s per loop 

In [23]: %timeit ndimage.uniform_filter(img, 5) 
1 loops, best of 3: 8.69 s per loop 

Der Unterschied verursacht wird durch uniform_filter Preforming eine 1-d-Faltung entlang jeder Achse, anstelle einer generischen 3 D Faltung.

In Fällen, in denen der Kernel symmetrisch ist, können Sie diese Vereinfachungen vornehmen und eine signifikante Beschleunigung erzielen.

Ich bin mir nicht sicher über convn, aber oft Matlab-Funktionen machen diese Art von Optimierungen hinter den Kulissen, wenn Ihre Eingabedaten bestimmte Kriterien erfüllen. Scipy verwendet häufiger einen Algorithmus pro Funktion und erwartet, dass der Benutzer weiß, welcher in welchem ​​Fall auszuwählen ist.

---

Sie erwähnten einen Filter "Laplace of Gaussian". Ich bekomme hier immer eine Berührung mit der Terminologie.

Ich denke, was Sie in Bezug auf ndimage Funktionen wollen, ist entweder scipy.ndimage.gaussian_laplace oder scipy.ndimage.gaussian_filter mit order=2 (die Filter durch die zweite Ableitung des Gaußschen).

In jedem Fall vereinfachen beide die Operation bis zu einer 1-d Faltung über jede Achse, was eine signifikante (2-3x) Beschleunigung geben sollte.