Metodo ottimizzato per il calcolo della distanza del coseno in Python

Question

Ho scritto un metodo per calcolare la distanza del coseno tra due array:

def cosine_distance(a, b):
    if len(a) != len(b):
        return False
    numerator = 0
    denoma = 0
    denomb = 0
    for i in range(len(a)):
        numerator += a[i]*b[i]
        denoma += abs(a[i])**2
        denomb += abs(b[i])**2
    result = 1 - numerator / (sqrt(denoma)*sqrt(denomb))
    return result

Eseguirlo può essere molto lento su un array di grandi dimensioni. Esiste una versione ottimizzata di questo metodo che sarebbe più veloce?

Aggiornamento: ad oggi ho provato tutti i suggerimenti, incluso Scipy. Ecco la versione da battere, che include suggerimenti di Mike e Steve:

def cosine_distance(a, b):
    if len(a) != len(b):
        raise ValueError, "a and b must be same length" #Steve
    numerator = 0
    denoma = 0
    denomb = 0
    for i in range(len(a)):       #Mike's optimizations:
        ai = a[i]             #only calculate once
        bi = b[i]
        numerator += ai*bi    #faster than exponent (barely)
        denoma += ai*ai       #strip abs() since it's squaring
        denomb += bi*bi
    result = 1 - numerator / (sqrt(denoma)*sqrt(denomb))
    return result

Answer 1

Se puoi usare SciPy, puoi usare cosine da spatial.distance :

http://docs.scipy.org/doc/scipy/reference/spatial.distance.html

Se non puoi usare SciPy, potresti provare a ottenere un piccolo aumento di velocità riscrivendo il tuo Python (EDIT: ma non ha funzionato come pensavo, vedi sotto).

from itertools import izip
from math import sqrt

def cosine_distance(a, b):
    if len(a) != len(b):
        raise ValueError, "a and b must be same length"
    numerator = sum(tup[0] * tup[1] for tup in izip(a,b))
    denoma = sum(avalue ** 2 for avalue in a)
    denomb = sum(bvalue ** 2 for bvalue in b)
    result = 1 - numerator / (sqrt(denoma)*sqrt(denomb))
    return result

È meglio sollevare un'eccezione quando le lunghezze di aeb non corrispondono.

Usando le espressioni del generatore all'interno delle chiamate a sum () puoi calcolare i tuoi valori con la maggior parte del lavoro svolto dal codice C all'interno di Python. Questo dovrebbe essere più veloce dell'uso di un ciclo per .

Non ho cronometrato, quindi non riesco a indovinare quanto potrebbe essere più veloce. Ma il codice SciPy è quasi certamente scritto in C o C ++ e dovrebbe essere il più veloce possibile.

Se stai facendo bioinformatica in Python, dovresti comunque usare SciPy.

EDIT: Darius Bacon ha cronometrato il mio codice e lo ha trovato più lento. Quindi ho cronometrato il mio codice e ... sì, è più lento. La lezione per tutti: quando stai cercando di accelerare le cose, non indovinare, misurare.

Sono sconcertato dal motivo per cui il mio tentativo di mettere più lavoro sugli interni C di Python è più lento. L'ho provato per elenchi di lunghezza 1000 ed era ancora più lento.

Non posso passare altro tempo a cercare di hackerare Python in modo intelligente. Se hai bisogno di più velocità, ti suggerisco di provare SciPy.

EDIT: ho appena testato a mano, senza tempo. Trovo che in breve aeb, il vecchio codice è più veloce; per lunghi aeb, il nuovo codice è più veloce; in entrambi i casi la differenza non è grande. (Ora mi chiedo se posso fidarmi del timeit sul mio computer Windows; voglio provare di nuovo questo test su Linux.) Non cambierei il codice di lavoro per cercare di farlo più velocemente. E ancora una volta ti esorto a provare SciPy. : -)

Answer 2

(inizialmente pensavo) non lo accelererai molto senza irrompere in C (come intorpidimento o scipy) o cambiare ciò che calcoli. Ma ecco come lo proverei, comunque:

from itertools import imap
from math import sqrt
from operator import mul

def cosine_distance(a, b):
    assert len(a) == len(b)
    return 1 - (sum(imap(mul, a, b))
                / sqrt(sum(imap(mul, a, a))
                       * sum(imap(mul, b, b))))

È circa due volte più veloce in Python 2.6 con array di 500k elementi. (Dopo aver cambiato la mappa in imap, seguendo Jarret Hardie.)

Ecco una versione ottimizzata del codice rivisto del poster originale:

from itertools import izip

def cosine_distance(a, b):
    assert len(a) == len(b)
    ab_sum, a_sum, b_sum = 0, 0, 0
    for ai, bi in izip(a, b):
        ab_sum += ai * bi
        a_sum += ai * ai
        b_sum += bi * bi
    return 1 - ab_sum / sqrt(a_sum * b_sum)

È brutto, ma esce più velocemente. . .

Modifica: e prova Psyco ! Accelera la versione finale di un altro fattore di 4. Come potrei dimenticare?

Answer 3

Non è necessario prendere abs () di a [i] e b [i] se lo quadrate.

Memorizza a [i] e b [i] in variabili temporanee, per evitare di eseguire l'indicizzazione più di una volta. Forse il compilatore può ottimizzarlo, ma forse no.

Controlla nell'operatore ** 2 . Lo sta semplificando in un moltiplicare o sta usando una funzione di potenza generale (log - moltiplica per 2 - antilog)?

Non fare sqrt due volte (anche se il costo è piccolo). Esegui sqrt (denoma * denomb) .

Answer 4

Simile alla risposta di Darius Bacon, ho giocato con l'operatore e gli itertools per produrre una risposta più veloce. Quanto segue sembra essere 1/3 più veloce su un array di 500 elementi secondo timeit:

from math import sqrt
from itertools import imap
from operator import mul

def op_cosine(a, b):
    dot_prod = sum(imap(mul, a, b))
    a_veclen = sqrt(sum(i ** 2 for i in a))
    b_veclen = sqrt(sum(i ** 2 for i in b))

    return 1 - dot_prod / (a_veclen * b_veclen)

Answer 5

Questo è più veloce per array di circa 1000+ elementi.

from numpy import array
def cosine_distance(a, b):
    a=array(a)
    b=array(b)
    numerator=(a*b).sum()
    denoma=(a*a).sum()
    denomb=(b*b).sum()
    result = 1 - numerator / sqrt(denoma*denomb)
    return result

Answer 6

L'uso del codice C all'interno di SciPy vince molto per gli array di input lunghi. L'uso di Python semplice e diretto vince per array di input brevi; Il codice basato su izip () di Darius Bacon è stato testato al meglio. Pertanto, la soluzione definitiva è decidere quale utilizzare in fase di esecuzione, in base alla lunghezza degli array di input:

from scipy.spatial.distance import cosine as scipy_cos_dist

from itertools import izip
from math import sqrt

def cosine_distance(a, b):
    len_a = len(a)
    assert len_a == len(b)
    if len_a > 200:  # 200 is a magic value found by benchmark
        return scipy_cos_dist(a, b)
    # function below is basically just Darius Bacon's code
    ab_sum = a_sum = b_sum = 0
    for ai, bi in izip(a, b):
        ab_sum += ai * bi
        a_sum += ai * ai
        b_sum += bi * bi
    return 1 - ab_sum / sqrt(a_sum * b_sum)

Ho realizzato un cablaggio di prova che ha testato le funzioni con input di lunghezza diversa e ho scoperto che intorno alla lunghezza 200 la funzione SciPy ha iniziato a vincere. Più grandi sono gli array di input, più grande vince. Per matrici di lunghezza molto breve, diciamo lunghezza 3, vince il codice più semplice. Questa funzione aggiunge una piccola quantità di spese generali per decidere in che modo farlo, quindi lo fa nel modo migliore.

Se sei interessato, ecco il cablaggio di prova:

from darius2 import cosine_distance as fn_darius2
fn_darius2.__name__ = "fn_darius2"

from ult import cosine_distance as fn_ult
fn_ult.__name__ = "fn_ult"

from scipy.spatial.distance import cosine as fn_scipy
fn_scipy.__name__ = "fn_scipy"

import random
import time

lst_fn = [fn_darius2, fn_scipy, fn_ult]

def run_test(fn, lst0, lst1, test_len):
    start = time.time()
    for _ in xrange(test_len):
        fn(lst0, lst1)
    end = time.time()
    return end - start

for data_len in range(50, 500, 10):
    a = [random.random() for _ in xrange(data_len)]
    b = [random.random() for _ in xrange(data_len)]
    print "len(a) ==", len(a)
    test_len = 10**3
    for fn in lst_fn:
        n = fn.__name__
        r = fn(a, b)
        t = run_test(fn, a, b, test_len)
        print "%s:\t%f seconds, result %f" % (n, t, r)

Answer 7

def cd(a,b):
    if(len(a)!=len(b)):
        raise ValueError, "a and b must be the same length"
    rn = range(len(a))
    adb = sum([a[k]*b[k] for k in rn])
    nma = sqrt(sum([a[k]*a[k] for k in rn]))
    nmb = sqrt(sum([b[k]*b[k] for k in rn]))

    result = 1 - adb / (nma*nmb)
    return result

Answer 8

La tua soluzione aggiornata ha ancora due radici quadrate. Puoi ridurlo a uno sostituendo la riga sqrt con:

  

risultato = 1 - numeratore /   (Sqrt (denoma * denomb))

Una moltiplicazione è in genere un po 'più veloce di una sqrt. Potrebbe non sembrare molto come viene chiamato solo una volta nella funzione, ma sembra che tu stia calcolando molte distanze del coseno, quindi il miglioramento si sommerà.

Il tuo codice sembra che dovrebbe essere maturo per le ottimizzazioni vettoriali. Quindi, se il supporto cross-platofrm non è un problema e si desidera accelerarlo ulteriormente, è possibile codificare il codice della distanza del coseno in C e assicurarsi che il compilatore stia vettorizzando in modo aggressivo il codice risultante (anche Pentium II è in grado di vettorializzare in virgola mobile )