Análisis Cluster (II): Clasificación no supervisada mediante clasificación jerárquica aglomerativa

(Este es el segundo capítulo de la mini-serie de artículos sobre análisis cluster que estamos haciendo en pybonacci, si todavía no has leído el artículo inicial échale un ojo ahora).

Como vimos anteriormente, existen diferentes formas de hacer clustering y, como también comentamos anteriormente, una de las más habituales es el clustering jerárquico.

El clustering jerárquico asociativo pretende, partiendo de m observaciones, ir encontrado agrupaciones de forma jerarquizada. Para ello, primero busca similitudes entre las observaciones y después procura asociar en grupos o 'clusters' las observaciones que se encuentran 'más cercanas' o presentan mayor similitud.

Si os acordáis, en el primer capítulo de esta mini-serie, entrecomillamos la palabra ‘similitud’. Vamos a ver qué significa esto de similitud en general y en nuestro ejemplo concreto (ver capítulo anterior de la serie para conocer el ejemplo). Dependiendo del problema concreto necesitaremos asociar las variables para poder medir como son de similares o a qué distancia se encuentran entre sí. Estas medidas de asociación (similitud o distancia) dependerán del problema concreto con el que nos encontremos y no se usaría lo mismo para variables booleanas, binarias, reales,... Gracias al módulo scipy.spatial.distance podemos ver muchas de estas medidas de asociación en las que no voy a entrar en detalle. Un resumen de la mayoría de ellas lo podéis ver en el siguiente enlace del profesor José Ángel Gallardo San Salvador, nuevamente. En el ejemplo propuesto usaremos la correlación como medida de asociación ya que es una medida invariante aunque las variables se escalen y/o se les sumen/resten parámetros (algo muy útil para el caso que nos ocupa).

Una vez que hemos calculado todas las medidas de asociación entre las m observaciones hemos de establecer un método para las uniones de las observaciones para así ir creando los clusters. Nuevamente, scipy acude al rescate para poder establecer esto mediante el módulo scipy.cluster y la función scipy.cluster.hierarchy.linkage. El 'linkage', enlazado o amalgamiento se puede realizar mediante diferentes métodos. En la función anterior podemos elegir entre algunos de ellos como el método simple, completo, promedio, pesado, de Ward,... ¿Dónde los podéis ver explicados de forma detallada y excelente? Pues sí, en los apuntes del profesor José Ángel Gallardo San Salvador, de la escuela de estadística de la Universidad de Granada (perded 10 minutos en leer el enlace anterior y lo que sigue a continuación tendrá más sentido). El enlazado se realiza en el momento inicial, cuando se dispone de las m observaciones, y después de hacer un nuevo grupo o cluster usando el método de enlazado (o linkage o amalgamiento) que se haya seleccionado. En las siguientes transparencias se intenta detallar un ejemplo, paso a paso, de un enlazado simple para una distancia mínima:

[slideshare id=15253227&w=427&h=356&sc=no]

Bueno, bueno,..., mucha palabrería, mucho enlace a apuntes,... Falta una cosa, ¡¡!SHOW ME THE CODE!! A eso vamos. Para ello vamos a hacer uso de unos datos de temperatura que descargaremos de la siguiente url (pinchad sobre 'FTP a copy of the file'). Una vez que tengáis el fichero descargado vamos a importar todo lo que vamos a necesitar:
from scipy import cluster
import numpy as np
from matplotlib import pyplot as plt
import netCDF4 as nc
from mpl_toolkits import basemap as bm

Ahora vamos a preparar los datos:
print('leyendo datos nc')
tsfc = nc.Dataset('X83.42.0.38.323.14.41.15.nc') #tsfc 'por temperature at surface'
lat = tsfc.variables['lat'][:]
lon = tsfc.variables['lon'][:]
tmp = tsfc.variables['air'][:]

Tenemos las series de temperatura (tmp) y su posición (lon, lat). Vamos a ver los nodos (series) que tenemos sobre un mapa.

m = bm.Basemap(llcrnrlon = np.min(lon) - 1, llcrnrlat = np.min(lat) - 1,
               urcrnrlon = np.max(lon) + 1, urcrnrlat = np.max(lat) + 1,
projection = 'mill')
lon, lat = np.meshgrid(lon, lat)
x, y = m(lon, lat)
m.scatter(x, y, color= 'y')
m.drawcountries()
m.drawlsmask(land_color = 'g', ocean_color = 'c')
plt.show()

Bueno, quizá son demasiadas series para ver el ejemplo pero luego podéis toquetear el código para usar menos series o vuestras propias series. Ahora es cuando hacemos los cálculos propios del análisis cluster:
tmp = tmp.reshape(tmp.shape[0], tmp.shape[1] * tmp.shape[2])
print('calculando grupos')
enlaces = cluster.hierarchy.linkage(tmp.T, method = 'single', metric = 'correlation')
print('dibujando dendrograma')
cluster.hierarchy.dendrogram(enlaces, color_threshold=0)
plt.show()

Y veremos algo como lo siguiente, que se conoce como dendrograma:

Desgraciadamente no se ve muy bien puesto que son muchas observaciones (682), lo dicho, toquetead para hacerlo con menos series y, si queréis, cambiando el método de 'linkage' y de distancias (en el ejemplo se usa el método simple y la correlación, respectivamente). Vemos que en el eje y hay unos valores, Pues bien, si cortamos el dendrograma horizontalmente por un valor nos quedaremos con tantos grupos como 'barritas verticales' o ramas cortemos. Vamos a hacer una prueba usando el valor 0.15 como valor de corte y vamos a representar a qué grupo pertenece cada una de las observaciones:
clusters = cluster.hierarchy.fcluster(enlaces, 0.15, criterion = 'distance')
clusters = clusters.reshape(lat.shape)
m = bm.Basemap(llcrnrlon = np.min(lon) - 1, llcrnrlat = np.min(lat) - 1,
               urcrnrlon = np.max(lon) + 1, urcrnrlat = np.max(lat) + 1,
               projection = 'mill')
colores = np.linspace(0,1,np.max(clusters))
for j in range(clusters.shape[0]):
    for i in range(clusters.shape[1]):
        plt.text(x[j,i],y[j,i], str(clusters[j,i]),
                 color = str(colores[clusters[j,i] - 1]),
                 horizontalalignment='center',
                 verticalalignment='center')
plt.title(u'Número de grupos: %03d' % np.max(clusters))
m.drawcountries()
m.drawlsmask(land_color = 'g', ocean_color = 'c')
plt.show()

Y veremos un número sobre cada nodo que es el grupo al que pertenece cada observación si cortamos en 0.15 (los colores no indican nada, solo sirve para poder visualizar e identificar un poco más fácilmente los grupos).

Y ahora es donde vendría el trabajo del experto para interpretar el resultado, para conocer si con 19 grupos es suficiente, es demasiado,... Haciendo una interpretación (no os creáis nada de lo que viene a continuación) de estos resultados podríamos decir que hay una clara separación entre temperaturas en el ecuador y el Caribe (grupo 11) y temperaturas por debajo del subtrópico (grupo 8), parece que hay varios grupos (17, 18, 6, 5, 12) que podrían estar diciendo que parece que estamos identificando un fenómeno de ¿El Niño?,...

He intentado sintetizar excesivamente y no sé si habréis entendido algo, pero bueno, si queréis profundizar en el tema, tenéis enlaces para leer información más en profundidad, hemos localizado bibliotecas python que nos ayudan a hacer este tipo de análisis y tenéis los comentarios para corregir, discutir, criticar, preguntar..., lo que consideréis oportuno.

Hasta la próxima!!

Kiko Correoso

Licenciado y PhD en Ciencias Físicas, especializado en temas de física, meteorología, climatología, energías renovables, estadística, aprendizaje automático, análisis y visualización de datos. Apasionado de Python y su comunidad. Fundador de pybonacci y editor del sitio en el que se divulga Python, Ciencia y el conocimiento libre en español.

More Posts

Follow Me:
TwitterLinkedIn

2 thoughts on “Análisis Cluster (II): Clasificación no supervisada mediante clasificación jerárquica aglomerativa

Leave a Reply

Your email address will not be published. Required fields are marked *