4) Operaciones con vectores en arrays o vectorización.

Anteriormente. Matématicas con Numpy

Operaciones con vectores en arrays o vectorización.

El trabajar con bucles "for" en python resulta muy lento. Si tienes que repetir una operación matemática en muchos elementos consecutivos de una matriz, siempre resulta mejor utilizar una operación vectorizada si es posible, es hasta un 80% más rápido. 

En la práctica una operación vectorizada, significa rehacer nuestro código para evitar el uso de bucles y utilizar en su lugar slicing, rebanadas o vectores de arrays de Numpy para aplicar la operación a todo o parte del array.

Vamos a ver la diferencia que existe con un ejemplo. Supongamos que tienes que calcular la diferencia existente entre los elementos consecutivos de una matriz. 

Como en este ejemplo solo queremos ver la diferencia de tiempos de ejecución que hay con ambos sistemas vamos a crear una matriz "a" que va a contener todos los números desde el cero al 999. 

La matriz "a" seria algo así:

a = [0,1,2,3,4,5,6,..........,998,999], tiene 1000 elementos

Lo que queremos hacer es obtener otra matriz de resultados, a la que llamaremos "b", con la diferencia entre los elementos consecutivos de la matriz. Es decir, la matriz "b" seria algo como esto:

b = [(1-0),(2-1),(3-2),(4-3).....................(998-997),(999-998)] que tiene 999 elementos.

Como ves en esta matriz, todos sus elementos van a ser unos. Pero esto no nos interesa, no nos interesa el resultado, sino el proceso que hay que seguir en python para hacerlo.

Para ello vamos a crear un archivo en python, yo lo llamaré consecutivos.py, con el siguiente código:

# importamos la libreria numpy para crear los array
import numpy as np

# Creamos la matriz a con 1000 números consecutivos del 0 al 999
a = np.arange(1000)

# Creamos la matriz b rellena con ceros para luego almacenar los resultados.
# Al estar formada con la diferencia entre 2 elementos consecutivos, tendrá
# 1000 - 1 elementos, es decir 999
b = np.zeros(999, int)

# Bucle for que calcula la diferencia entre un elemento (i) y el anterior
# (i-1), recorre haciendo el calculo con un bucle for la matriz a
# y asigna el resultado en la matriz b.
# e.j b[0]= a[1] que es 1 menos a[0] que es el cero 
for i in range(1,len(a)):
    b[i-1]=a[i]-a[i-1]

Pues bien, esto mismo se puede conseguir usando una operación vectorizada, con el siguiente código:

c = a[1:] - a[:-1]

Con esto estaríamos restando, elemento a elemento, la siguiente rebanada de la matriz a:

a[1:] = [1,2,3,4,...,998,999] desde el elemento segundo, cuyo valor es 1, hasta el final. (se empieza a contar en el cero)

a[: -1] = [0,1,2,3,....,997,998] desde el primer elemento que es el 0 hasta el penúltimo, el 998.

Si no lo acabas de ver bien, imagínatelo en pequeñito, con este ejemplo.

[1,2,3,4,5]

arr[1:]   → [2,3,4,5]

arr[:-1] →  [1,2,3,4]

dif =           [1,1,1,1]

Vamos ahora a juntarlo todo y ver lo que se tarda en ejecutar de las dos formas. Pero para que se note un poco más y se vea claramente la diferencia vamos a realizarlo, no con una matriz no de 1.000 elementos sino una de 1.000.000 de números.

# importamos la libreria numpy para crear los array
import numpy as np
# importamos time para medir el tiempo de ejecución
from time import time

# Creamos la matriz a con 1.000.000 de números consecutivos del 0 al 999.999
# en python el guion bajo _ nos permite separar los miles para verlo mejor.
a = np.arange(1_000_000)

# Creamos la matriz b rellena con ceros para luego almacenar los resultados.
# Al estar formada con la diferencia entre 2 elementos consecutivos, tendrá
# 1.000.000 - 1 elementos, es decir 999.999
b = np.zeros(len(a-1), int)

# Bucle for que calcula la diferencia entre un elemento (i) y el anterior
# (i-1), recorre haciendo el calculo con un bucle for la matriz a
# y asigna el resultado en la matriz b.
# e.j b[0]= a[1] que es 1 menos a[0] que es el cero 
t0 = time()
for i in range(1,len(a)):
    b[i-1]=a[i]-a[i-1]
t1 = time()

# Utilizando vectores de la matriz o slicing en el array.
t2 = time()
c = a[1:]-a[:-1]
t3 = time()

print("tiempo del bucle for ", (t1-t0))
print("tiempo con vectorización", (t3-t2))
print((t1-t0)>(t3-t2))

Yo, en mi equipo obtengo los siguientes resultados:

tiempo del bucle for  0.6465959548950195
tiempo con vectorización 0.0030362606048583984
True

Como puedes ver, la diferencia de tiempo es MUY significativa.

Seleccionar elementos que cumplan una determinada condición en un vector.

Una tarea que se realiza con bastante frecuencia es seleccionar elementos de un vector. Esta selección se pude llevar a cabo en base a índices o en base a un determinado criterio (por ejemplo que los elementos a seleccionar sean mayores que un determinado valor o que se encuentren en un determinado rango)

Vamos a realizar un ejemplo. Creamos una matriz de ejemplo en Numpy:

>>> import numpy as np

>>> matriz = np.arange(10)

>>> matriz
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])

Ahora vamos a seleccionar todos los elementos del array que sean mayores que 5. Al utilizar operadores de igualdad o comparación vamos a obtener un nuevo array con valores booleanos.

>>> matriz > 5
array([False, False, False, False, False, False,  True,  True,  True,
        True])

Como ves obtenemos un nuevo vector en el que los 6 primeros elementos son falsos y los 4 últimos verdaderos. Este resultado se puede guardar en una variable y utilizarlo para seleccionar los elementos que cumplen la condición.

Si a un array de Numpy se le pasa un vector de valores booleanos el resultado es un nuevo vector en el que solo estarán los elementos cuyas condiciones sean ciertas. Vamos a verlo en nuestro ejemplo:

>>> mayor_5 = matriz > 5

# seleccionamos los elementos de matriz en base a matriz > 5
>>> matriz[mayor_5]
array([6, 7, 8, 9])

Y como puedes ver hemos seleccionado los elementos de nuestra matriz original que cumplen la condición de ser mayores de 5.

Seleccionar elementos en base a múltiples condiciones en Numpy

Una vez que hemos visto como seleccionar elementos de un vector en base a una condición, podemos hacer las selecciones complejas que queramos utilizando operadores lógicos igual que en Python. Por ejemplo vamos como obtendríamos todos los números mayores que 5 pero a la vez que sean pares. Recordemos que para que un número sea para el resto de dividir ese número entre dos tiene que ser cero. Pues bien vamos a utilizar la función np.mod(matriz, divisor) que nos va a proporcionar el resto de dividir cada elemento de la matriz entre el divisor que elijamos.

Con todo lo anterior para elegir cual de los elementos son pares y mayores que 5 el código sería el siguiente:

>>> matriz[np.mod(matriz,2) & matriz > 5]
array([6, 8])

Le estamos diciendo Numpy que selecciones los elementos de la matriz que cumplan 2 criterios:

- Que el resto de su división entre 2 sea cero o lo que es lo  mismo que el número es par.

- Que sea mayor que 5.

Como hemos visto, con está técnica podemos filtrar grandes cantidades de datos en python de forma muy eficiente ya podemos implementar fácilmente filtros complejos. 

EJERCICIOS.

Calculo Diferencial

Con este ejercicio vamos a practicar también la vectorización, lo cual es crucial para obtener buenos resultados en términos de velocidad de cálculo y rendimiento con Numpy.

Para el que lo haya visto recordar que las derivadas, en matemáticas, se pueden calcular numéricamente con el método de las diferencias finitas como:

f′(xi)=(f(xi+Δx)−f(xi−Δx))/2Δx

Vamos a construir una matriz unidimensional que contenga los valores de xi en el intervalo 0 y π/2 con incrementos de 0.10.

Luego calcularemos numéricamente la derivada de la función seno y coseno en ese intervalo (excluyendo los puntos que constituyen los extremos) y para terminar las representaremos gráficamente.

Primeramente importamos las librerías necesarias:

>>> import numpy as np
>>> import matplnotlib.pyplot as plt

Con un array de Numpy realizamos los cálculos y construimos las funciones seno y coseno entre 0 y pi/2.

>>> dx=0.10
>>> x=np.arange(0,np.pi/2,dx) # Intervalo entre 0 y pi/2 con incrementos de 0.10
>>> f=np.sin(x) # función seno en ese intervalo
>>> g=np.cos(x) # función coseno en ese intervalo.

>>> print(x)
[0.  0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.  1.1 1.2 1.3 1.4 1.5]

>>> print(f) # Seno
[0.         0.09983342 0.19866933 0.29552021 0.38941834 0.47942554
 0.56464247 0.64421769 0.71735609 0.78332691 0.84147098 0.89120736
 0.93203909 0.96355819 0.98544973 0.99749499]

>>> print(g) # Coseno
[1.         0.99500417 0.98006658 0.95533649 0.92106099 0.87758256
 0.82533561 0.76484219 0.69670671 0.62160997 0.54030231 0.45359612
 0.36235775 0.26749883 0.16996714 0.0707372 ]

Ahora calcularemos las derivadas de ambas funciones por el método de las diferencias finitas.

>>> df_sen=(f[2:]-f[:-2])/2*dx
>>> print(df_sen)
[0.00993347 0.00978434 0.00953745 0.00919527 0.00876121 0.00823961
 0.00763568 0.00695546 0.00620574 0.00539402 0.00452841 0.00361754
 0.00267053 0.00169684]

>>> df_cos=(g[2:]-g[:-2])/2*dx
>>> print(df_cos)
[-0.00099667 -0.00198338 -0.00295028 -0.0038877  -0.00478627 -0.00563702
 -0.00643145 -0.00716161 -0.00782022 -0.00840069 -0.00889723 -0.00930486
 -0.00961953 -0.00983808]

Para entenderlo mejor vamos a desarrollar los arreglos de la función seno para que veas como funcionan:

# función f(seno)
>>> print(f)
[0.         0.09983342 0.19866933 0.29552021 0.38941834 0.47942554
 0.56464247 0.64421769 0.71735609 0.78332691 0.84147098 0.89120736
 0.93203909 0.96355819 0.98544973 0.99749499]

>>> print(f[2:])
[0.19866933 0.29552021 0.38941834 0.47942554 0.56464247 0.64421769
 0.71735609 0.78332691 0.84147098 0.89120736 0.93203909 0.96355819
 0.98544973 0.99749499]

>>> print(f[:-2])
[0.         0.09983342 0.19866933 0.29552021 0.38941834 0.47942554
 0.56464247 0.64421769 0.71735609 0.78332691 0.84147098 0.89120736
 0.93203909 0.96355819]

# Como ves ambos arreglos tienen el mismo número de elementos para poder
# luego restarlos y acabar de aplicar la fórmula.

Para dibujar las funciones y su relación con las derivadas:

plt.plot(x[1:-1], df_sen)
plt.plot(x[1:-1], f[1:-1])
plt.plot(x[1:-1], df_cos)
plt.plot(x[1:-1], g[1:-1])
plt.show()

Próximo Post. Manipulación de Arrays y Transmisión en Numpy.