Serielle vektorielle Subtraktion in Python oder "subtraktive Faltung"?

Question

Zunächst möchte ich mich für den unklaren Titel der Frage entschuldigen: Der Grund ist, dass ich den mathematischen Prozess bei der Arbeit nicht identifizieren konnte. Hier

ist die Situation, in Kürze:

1. Ich habe zwei Vektoren, f1 und f2, unterschiedlich lang.
2. Ich möchte den kleinsten quadratischen Abstand zwischen f1 und f2, elementweise berechnen.

Hier ist, wie ich vorgehen,

import numpy as np 

def distLstSq(f1,f2): 
    "Return the least square distance between f1 and f2 per unit element" 
    return np.sum(np.square(np.subtract(f1,f2)))/len(f1) 

f1 = np.arange(100) 
f2 = np.random.random_integers(1,100,5) 

nBuf = len(f2) 
dist = np.empty(len(f1)-nBuf) 
for i in range(len(f1)-nBuf): 
    temp = f1[i:i+nBuf] 
    dist[i] = distLstSq(temp,f2) 

jedoch aufgrund des großen Vektor f1 (aus einer Datei von 4MB generiert), ich frage mich, ob es nicht mehr elegant pythonic Lösung war, Verwenden Sie weniger CPU-Ressourcen in kürzerer Zeit. Möglicherweise eine Art "subtraktive Faltung", wobei f2 elementweise über f1 gleitet, mit der Subtraktionsoperation bei jedem Schritt.

Vielen Dank für Ihre Eingabe!

Bertrand

Answer 1

Zuerst möchte ich, dass die Anzahl der Begriffe hinweisen, in dist nicht len(f1)-nBuf ist aber len(f1)-nBuf+1. Auf diese Weise überlappt der kürzere Vektor vollständig mit dem längeren.

Ignorieren der Division durch len(f1), die durch eine Konstante nur Skalierung ist, können Sie das folgende für jede Scheibe von f2 Berechnung:

Ich denke also, die Anzahl der Operationen mit einigen reduzieren Vorberechnung. Und wir können auch np.convolve verwenden, um die Kreuzterme zu finden.

f1_squared = f1**2 
f2_squared_sum = np.sum(f2**2) 
nBuf = len(f2) 
cross_terms = -2*np.convolve(f1, f2[::-1], "valid") 
# reverse f2 to get what we want. 
# "valid" returns where vectors completely overlap 
squared_distance = [f2_squared_sum + np.sum(f1_squared[i:i+nBuf]) + cross_terms[i] 
        for i in xrange(len(cross_terms))] 
mean_squared_distance = np.array(squared_distance)/nBuf 

Ihre Version:

nBuf = len(f2) 
dist = np.empty(len(f1)-nBuf+1) 
for i in xrange(len(f1)-nBuf+1): 
    temp = f1[i:i+nBuf] 
    dist[i] = distLstSq(temp,f2) 

Basierend auf timeit.timeit, meine Version ist 30-50% schneller.

Die Leistung kann weiter verbessert werden, indem np.sum(f1_squared[i:i+nBuf]) optimiert wird, da dies wiederholte Operationen erfordert. Es sollte einen Weg geben, um die Summe zu teilen.

Auch ich denke, scipy.signal.fftconvolve schneller als np.convolve sein kann, aber das hängt von der Länge des kürzeren Vektor (see here)