Cómo calcular la correlación de cada fila en una matriz 2D con una matriz 1D de la misma longitud

He verificado esto con la función corrcoef de Numpy, pero usaré esto como una oportunidad para comprender y practicar las funciones de vectorización usando numpy. Estoy seguro de que existen algos e implementaciones más eficientes

Deje que esta matriz de 1D de interés sea la siguiente:

x = np.random.randn (10)

La matriz 2D de interés será la siguiente:

y = np.random.randn (100, 10)

Ahora queremos usar la fórmula del coeficiente de Pearson repetidamente para la fila de la matriz y , con x. No seas flojo, búscalo o el resto de esta respuesta no tiene sentido (accesorios si lo sabes de memoria, no lo hice). Correlación y dependencia – Wikipedia

Ahora calculemos la parte [math] x_i – \ bar {x} [/ math] de la fórmula. Esto es facil. Numpy permite la transmisión, lo que básicamente significa que suma y resta a lo largo de una dimensión que “tiene sentido” dadas dos matrices de dimensiones desiguales. Se puede encontrar una definición más precisa aquí Broadcasting – NumPy v1.10 Manual

x_bar = np.mean (x)
x_intermediate = x – x_bar

Hacer lo mismo para y no es tan fácil.

y_bar = np.mean (y, axis = 1) # esto aplana y para que sea (100,) que es una matriz 1D. El problema es que y es 100, 10, por lo que la transmisión de numpy no sabe qué eje tratar como el que se va a transmitir.
y_bar = y_bar [:, np.newaxis] # Al agregar esta dimensión adicional, estamos obligando a numpy a tratar el eje 0 como el que se transmite, lo que hace posible el siguiente paso. y_bar ahora es 100, 1
y_intermediate = y – y_bar

La función np.mean toma la media sobre el eje especificado pero aplana la matriz y borra la segunda dimensión. Borrar la segunda dimensión no proporciona suficiente información para que funcione el mecanismo de transmisión de numpy, que puede solucionarse como se muestra.

Ahora para encontrar los numeradores del coeficiente de Pearson. Para esto, queremos un producto punto de cada fila de y_intermediate con x_intermediate.

numeradores = y_intermediate.dot (x_intermediate) # o x_intermediate.dot (y_intermediate.T)

y para los denominadores queremos la suma de [math] \ sum_i [/ ​​math] [math] (x_intermediate) ^ 2 [/ math].

x_sq = np.sum (np.square (x_intermediate))

Lo mismo para y_intermediate, excepto que esta es una matriz. Entonces, cuando hacemos la suma, necesitamos una suma sabia de columnas.

y_sqs = np.sum (np.square (y_intermediate), axis = 1)

El denominador es el producto de estos dos.

denoms = np.sqrt (x_sq * y_sqs) # vector de tiempos escalares

Finalmente:

pearsons = (numeradores / denominadores) # numeradores es forma (100,) y denominadores es forma (100,)

Editar: Tal vez una palabra sobre por qué esto es eficiente está en orden. Las funciones de Numpy como mean y sum están vectorizadas. También lo son las operaciones de álgebra lineal como los productos de punto. Esto está disponible mediante bibliotecas eficientes que Numpy utiliza bajo el capó como BLAS y LAPACK.

Si ‘a’ es su matriz 2d y ‘b’ es su matriz 1d:

resultado = mapa (lambda x: numpy.correlate (x, b))

El resultado contendrá una matriz 1d donde el enésimo elemento es la correlación de la enésima fila de ‘a’ con ‘b’. Si solo desea el índice de qué fila tiene la mayor correlación, puede encadenarlo así:

nunmpy.argmax (mapa (lambda x: numpy.correlate (x, b))

El mapa es mucho más rápido que un bucle for.

Si tiene muchos datos, puede echar un vistazo a ufuncs vectorizados o numpy como se menciona en esta publicación: La forma más eficiente de mapear la función sobre una matriz numpy