Método optimizado para calcular la distancia del coseno en Python

Question

Escribí un método para calcular la distancia del coseno entre dos matrices:

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

Ejecutarlo puede ser muy lento en una gran matriz. ¿Existe una versión optimizada de este método que se ejecute más rápido?

Actualización: He probado todas las sugerencias hasta la fecha, incluso scipy. Aquí está la versión para vencer, incorporando sugerencias de Mike y 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

Si puede usar SciPy, puede usar coseno de spatial.distance :

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

Si no puede usar SciPy, puede intentar obtener una pequeña aceleración reescribiendo su Python (EDITAR: pero no funcionó como pensaba, ver más abajo).

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

Es mejor generar una excepción cuando las longitudes de ayb no coinciden.

Al usar expresiones generadoras dentro de las llamadas a sum () puede calcular sus valores con la mayor parte del trabajo realizado por el código C dentro de Python. Esto debería ser más rápido que usar un bucle for .

No he cronometrado esto, así que no puedo adivinar qué tan rápido podría ser. Pero es casi seguro que el código SciPy está escrito en C o C ++ y debería ser lo más rápido posible.

Si está haciendo bioinformática en Python, debería utilizar SciPy de todos modos.

EDITAR: Darius Bacon cronometró mi código y lo encontró más lento. Así que cronometré mi código y ... sí, es más lento. La lección para todos: cuando intente acelerar las cosas, no adivine, mida.

Estoy desconcertado sobre por qué mi intento de poner más trabajo en los componentes internos de Python es más lento. Lo probé para listas de longitud 1000 y todavía era más lento.

No puedo pasar más tiempo intentando piratear Python de manera inteligente. Si necesita más velocidad, le sugiero que pruebe SciPy.

EDITAR: acabo de probar a mano, sin tiempo. Encuentro que para abreviar ayb, el código antiguo es más rápido; para largo ayb, el nuevo código es más rápido; en ambos casos la diferencia no es grande. (Ahora me pregunto si puedo confiar en timeit en mi computadora con Windows; quiero probar esta prueba nuevamente en Linux). No cambiaría el código de trabajo para intentar hacerlo más rápido. Y una vez más te insto a que pruebes SciPy. :-)

Answer 2

(originalmente pensé) no vas a acelerarlo mucho sin romper a C (como numpy o scipy) o cambiar lo que calculas. Pero así es como lo intentaría, de todos modos:

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))))

Es aproximadamente el doble de rápido en Python 2.6 con matrices de 500k elementos. (Después de cambiar el mapa a imap, siguiendo a Jarret Hardie).

Aquí hay una versión modificada del código revisado del póster original:

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)

Es feo, pero sale más rápido. . .

Editar: ¡Y pruebe Psyco ! Acelera la versión final por otro factor de 4. ¿Cómo podría olvidarlo?

Answer 3

No es necesario tomar abs () de a [i] y b [i] si lo cuadra.

Almacene a [i] y b [i] en variables temporales, para evitar la indexación más de una vez. Tal vez el compilador pueda optimizar esto, pero tal vez no.

Regístrese en el operador ** 2 . ¿Lo está simplificando en una multiplicación, o está utilizando una función de potencia general (log - multiplicar por 2 - antilog)?

No hagas sqrt dos veces (aunque el costo es pequeño). Haga sqrt (denoma * denomb) .

Answer 4

Similar a la respuesta de Darius Bacon, he estado jugando con el operador y las herramientas de iterto para producir una respuesta más rápida. Lo siguiente parece ser 1/3 más rápido en una matriz de 500 elementos según 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

Esto es más rápido para conjuntos de más de 1000 elementos.

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

El uso del código C dentro de SciPy gana mucho para las matrices de entrada largas. Usando Python simple y directo gana para arreglos de entrada cortos; El código basado en izip () de Darius Bacon se comparó mejor. Por lo tanto, la solución final es decidir cuál usar en tiempo de ejecución, en función de la longitud de las matrices de entrada:

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)

Hice un arnés de prueba que probó las funciones con diferentes entradas de longitud, y descubrí que alrededor de la longitud 200 la función SciPy comenzó a ganar. Cuanto más grandes son las matrices de entrada, más grande gana. Para matrices de longitud muy corta, digamos longitud 3, gana el código más simple. Esta función agrega una pequeña cantidad de sobrecarga para decidir de qué manera hacerlo, luego lo hace de la mejor manera.

En caso de estar interesado, aquí está el arnés de prueba:

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

Su solución actualizada todavía tiene dos raíces cuadradas. Puede reducir esto a uno reemplazando la línea sqrt con:

  

resultado = 1 - numerador /   (sqrt (denoma * denomb))

Una multiplicación suele ser bastante más rápida que un sqrt. Puede que no parezca mucho, ya que solo se llama una vez en la función, pero parece que está calculando muchas distancias de coseno, por lo que la mejora se sumará.

Parece que su código debería estar maduro para optimizaciones vectoriales. Entonces, si el soporte multiplataforma no es un problema y desea acelerarlo aún más, puede codificar el código de distancia cosenoidal en C y asegurarse de que su compilador vectorice agresivamente el código resultante (incluso Pentium II es capaz de alguna vectorización de punto flotante )