Los cuadernos de Jupyter II Parte

Comenzamos con la II parte de los Cuadernos de Jupyter aplicado a las Geo-Ciencias, si necesitas mayor información no dudes en escribir.

DATOS GEOGRÁFICOS CON BASEMAP

Un tipo común de visualización en la ciencia de datos es el de los datos geográficos. La herramienta principal de Matplotlib para este tipo de visualización es el kit de herramientas Basemap, que es uno de varios kits de herramientas Matplotlib que vive bajo el espacio de nombres mpl_toolkits. Es cierto que Basemap se siente un poco torpe de usar, y muchas veces incluso las visualizaciones simples tardan mucho más tiempo en renderizarse de lo que esperas. Las soluciones más modernas como el leaflet o API de Google Maps pueden ser una mejor opción para visualizaciones de mapas más intensivas. Aún así, Basemap es una herramienta útil para los usuarios de Python en sus herramientas virtuales. En esta sección, mostraremos varios ejemplos del tipo de visualización de mapas que es posible con este kit de herramientas.

La instalación de Basemap es simple; Si está utilizando conda, puede escribir esto y el paquete se descargará:

$ conda install basemap

Agregamos solo una nueva importación a nuestro estándar repetitivo:

In [1]:

%matplotlib inline

import numpy as np

import matplotlib.pyplot as plt

from mpl_toolkits.basemap import Basemap

Una vez que haya instalado e importado el juego de herramientas Basemap, las parcelas geográficas están a solo unas líneas de distancia (los gráficos a continuación también requieren el paquete PIL en Python 2, o el paquete pillow en Python 3):

In [2]:

plt.figure(figsize=(8, 8))

m = Basemap(projection='ortho', resolution=None, lat_0=0, lon_0=-73)

m.bluemarble(scale=0.6);

Lo útil es que el globo que se muestra aquí no es una mera imagen; ¡es un eje Matplotlib completamente funcional que comprende coordenadas esféricas y que nos permite sobrerreglar fácilmente datos en el mapa! Por ejemplo, podemos usar una proyección de mapa diferente, acercarnos a Suramérica y trazar la ubicación de Quito. Usaremos una imagen etopo (que muestra las características topográficas tanto en tierra como debajo del océano) como fondo del mapa:

In [13]:

fig = plt.figure(figsize=(8, 8))

m = Basemap(height=16700000,width=12000000,

            resolution='l',area_thresh=1000.,projection='omerc',\

            lon_0=-100,lat_0=15,lon_2=-120,lat_2=65,lon_1=-50,lat_1=-55)

m.drawcoastlines()

m.fillcontinents(color='coral',lake_color='aqua')

# draw parallels and meridians.

m.drawparallels(np.arange(-80.,81.,20.))

m.drawmeridians(np.arange(-180.,181.,20.))

m.drawmapboundary(fill_color='aqua')

# Map (long, lat) to (x, y) for plotting

x, y = m(-78.503722, -0.216471)

plt.plot(x, y, 'ok', markersize=5)

plt.text(x, y, ' Quito', fontsize=12);

Esto le da una breve visión del tipo de visualizaciones geográficas que son posibles con solo unas pocas líneas de Python. Ahora analizaremos las características del mapa base con más profundidad y proporcionaremos varios ejemplos de visualización de datos de mapas. Usando estos ejemplos breves como bloques de construcción, debería poder crear casi cualquier visualización de mapa que desee.

Proyecciones de mapas

Lo primero que debe decidir al usar mapas es qué proyección usar. Probablemente esté familiarizado con el hecho de que es imposible proyectar un mapa esférico, como el de la Tierra, sobre una superficie plana sin distorsionarlo o romper su continuidad. Estas proyecciones se han desarrollado en el transcurso de la historia humana, ¡y hay muchas opciones! Dependiendo del uso previsto de la proyección del mapa, existen ciertas características del mapa (por ejemplo, dirección, área, distancia, forma u otras consideraciones) que son útiles para mantener.

El paquete de mapa base implementa varias docenas de tales proyecciones, todas referenciadas por un código de formato corto. Aquí mostraremos brevemente algunos de los más comunes.

Comenzaremos por definir una rutina de conveniencia para dibujar nuestro mapa del mundo junto con las líneas de longitud y latitud:

In [30]:

from itertools import chain

​

def draw_map(m, scale=0.2):

    # draw a shaded-relief image

    m.shadedrelief(scale=scale)

   

    # lats and longs are returned as a dictionary

    lats = m.drawparallels(np.linspace(-90, 90, 13))

    lons = m.drawmeridians(np.linspace(-180, 180, 13))

​

    # keys contain the plt.Line2D instances

    lat_lines = chain(*(tup[1][0] for tup in lats.items()))

    lon_lines = chain(*(tup[1][0] for tup in lons.items()))

    all_lines = chain(lat_lines, lon_lines)

   

    # cycle through these lines and set the desired style

    for line in all_lines:

        line.set(linestyle='-', alpha=0.3, color='w')

Proyecciones cilíndricas

Las proyecciones de mapa más simples son proyecciones cilíndricas, en las que las líneas de latitud y longitud constantes se asignan a líneas horizontales y verticales, respectivamente. Este tipo de mapeo representa bastante bien las regiones ecuatoriales, pero produce distorsiones extremas cerca de los polos. El espaciado de las líneas de latitud varía entre las diferentes proyecciones cilíndricas, que conducen a diferentes propiedades de conservación y a una distorsión diferente cerca de los polos. En la siguiente figura mostramos un ejemplo de la proyección cilíndrica equidistante, que elige una escala de latitud que preserva las distancias a lo largo de los meridianos. Otras proyecciones cilíndricas son las proyecciones de Mercator (proyección = 'merc') y área cilíndrica igual (proyección = 'cea') .:

In [31]:

fig = plt.figure(figsize=(8, 6), edgecolor='w')

m = Basemap(projection='cyl', resolution=None,

            llcrnrlat=-90, urcrnrlat=90,

            llcrnrlon=-180, urcrnrlon=180, )

draw_map(m)

Los argumentos adicionales al mapa base para esta vista especifican la latitud (lat) y la longitud (lon) de la esquina inferior izquierda (llcrnr) y la esquina superior derecha (urcrnr) para el mapa deseado, en unidades de grados.

Proyecciones pseudo-cilíndricas

Las proyecciones pseudo-cilíndricas relajan el requisito de que los meridianos (líneas de longitud constante) permanezcan verticales; esto puede dar mejores propiedades cerca de los polos de la proyección. La proyección de Mollweide (proyección = 'moll') es un ejemplo común de esto, en el que todos los meridianos son arcos elípticos. Está construido para preservar el área a lo largo del mapa: aunque hay distorsiones cerca de los polos, el área de pequeños parches refleja el área verdadera. Otras proyecciones pseudo-cilíndricas son las proyecciones sinusoidales (proyección = 'sinu') y Robinson (proyección = 'robin').

In [32]:

fig = plt.figure(figsize=(8, 6), edgecolor='w')

m = Basemap(projection='moll', resolution=None,

            lat_0=0, lon_0=0)

draw_map(m)

Los argumentos adicionales al mapa base aquí se refieren a la latitud central (lat_0) y la longitud (lon_0) para el mapa deseado.

Proyecciones de perspectiva

Las proyecciones de perspectiva se construyen usando una opción particular de punto de perspectiva, similar a la Tierra desde un punto particular en el espacio (un punto que, para algunas proyecciones, se encuentra técnicamente dentro de la Tierra). Un ejemplo común es la proyección ortográfica (proyección = 'orto'), que muestra un lado del globo visto desde un espectador a una distancia muy larga. Como tal, puede mostrar solo la mitad del globo a la vez. Otras proyecciones basadas en perspectiva incluyen la proyección gnomónica (proyección = 'gnom') y la proyección estereográfica (proyección = 'estere'). Estos suelen ser los más útiles para mostrar pequeñas porciones del mapa.

Aquí hay un ejemplo de la proyección ortográfica:

In [33]:

fig = plt.figure(figsize=(8, 8))

m = Basemap(projection='ortho', resolution=None,

            lat_0=50, lon_0=0)

draw_map(m);

Proyecciones cónicas

Una proyección Cónica proyecta el mapa en un solo cono, que luego se desenrolla. Esto puede conducir a propiedades locales muy buenas, pero las regiones alejadas del punto de enfoque del cono pueden distorsionarse mucho. Un ejemplo de esto es la proyección Cónica Conforme de Lambert (proyección = 'lcc'). Proyecta el mapa en un cono dispuesto de tal manera que dos paralelos estándar (especificados en el mapa base por lat_1 y lat_2) tienen distancias bien representadas, con una escala decreciente entre ellos y un aumento fuera de ellos. Otras proyecciones cónicas útiles son la proyección cónica equidistante (proyección = 'eqdc') y la proyección de área igual de Albers (proyección = 'aea'). Las proyecciones cónicas, como las proyecciones de perspectiva, tienden a ser buenas opciones para representar parcelas pequeñas o medianas del globo.

In [34]:

fig = plt.figure(figsize=(8, 8))

m = Basemap(projection='lcc', resolution=None,

            lon_0=0, lat_0=50, lat_1=45, lat_2=55,

            width=1.6E7, height=1.2E7)

draw_map(m)

Otras proyecciones

Si vas a hacer mucho con visualizaciones basadas en mapas, te recomiendo que leas otras proyecciones disponibles, junto con sus propiedades, ventajas y desventajas. Lo más probable es que estén disponibles en el paquete de mapa base http://matplotlib.org/basemap/users/mapsetup.html. Si profundizas lo suficiente en este tema, encontrarás una increíble subcultura de geeks de geo-viz que estarán listos para argumentar fervientemente en apoyo de su proyección favorita para cualquier aplicación dada.

Dibujando un fondo de mapa

Anteriormente vimos los métodos bluemarble () y shadedrelief () para proyectar imágenes globales en el mapa, así como los métodos drawparallels () y drawmeridians () para dibujar líneas de latitud y longitud constantes. El paquete de mapa base contiene una gama de funciones útiles para dibujar bordes de características físicas como continentes, océanos, lagos y ríos, así como límites políticos, como países y estados y condados de EE. UU. Las siguientes son algunas de las funciones de dibujo disponibles que puede explorar con las funciones de ayuda de IPython:

Límites físicos y cuerpos de agua

drawcoastlines (): dibuja líneas de costa continental drawlsmask (): dibuja una máscara entre la tierra y el mar, para utilizar con la proyección de imágenes en uno u otro drawmapboundary (): dibuja el límite del mapa, incluido el color de relleno para los océanos. drawrivers (): Dibuja ríos en el mapa fillcontinents (): llena los continentes con un color dado; opcionalmente llenar lagos con otro color Límites políticos

drawcountries (): dibujar los límites del país drawstates (): dibuja los límites del estado de EE. UU. drawcounties (): Dibuja los límites del condado de EE. UU. Características del mapa

drawgreatcircle (): Dibuja un gran círculo entre dos puntos drawparallels (): Dibuja líneas de latitud constante drawmeridians (): Dibuja líneas de longitud constante drawmapscale (): Dibuja una escala lineal en el mapa

Imágenes de globo entero

bluemarble (): Proyecte la imagen de mármol azul de la NASA en el mapa shadedrelief (): Proyecte una imagen sombreada en relieve en el mapa etopo (): Dibuje una imagen de relieve etopo en el mapa warpimage (): proyecta una imagen proporcionada por el usuario en el mapa

Para las características basadas en límites, debe establecer la resolución deseada al crear una imagen de mapa base. El argumento de resolución de la clase Basemap establece el nivel de detalle en los límites, ya sea 'c' (crudo), 'l' (bajo), 'i' (intermedio), 'h' (alto), 'f' (completo) , o Ninguno si no se usarán límites. Esta elección es importante: establecer límites de alta resolución en un mapa global, por ejemplo, puede ser muy lento.

Aquí hay un ejemplo de cómo dibujar los límites tierra / mar, y el efecto del parámetro de resolución. Crearemos un mapa de baja y alta resolución de la bella isla de Skye de Escocia. Está ubicado a 57.3 ° N, 6.2 ° W, y un mapa de 90,000 × 120,000 kilómetros:

In [35]:

fig, ax = plt.subplots(1, 2, figsize=(12, 8))

​

for i, res in enumerate(['l', 'h']):

    m = Basemap(projection='gnom', lat_0=57.3, lon_0=-6.2,

                width=90000, height=120000, resolution=res, ax=ax[i])

    m.fillcontinents(color="#FFDDCC", lake_color='#DDEEFF')

    m.drawmapboundary(fill_color="#DDEEFF")

    m.drawcoastlines()

    ax[i].set_title("resolution='{0}'".format(res));

Tenga en cuenta que las líneas costeras de baja resolución no son adecuadas para este nivel de zoom, mientras que las de alta resolución funcionan bien. Sin embargo, el nivel bajo funcionaría bien para una vista global, ¡y sería mucho más rápido que cargar los datos de borde de alta resolución para todo el mundo! Puede requerir cierta experimentación encontrar el parámetro de resolución correcto para una vista determinada: la mejor ruta es comenzar con una gráfica rápida de baja resolución y aumentar la resolución según sea necesario.

Trazar datos en mapas

Tal vez la pieza más útil del conjunto de herramientas Basemap es la capacidad de sobreponer una variedad de datos en un fondo de mapa. Para gráficos y texto simples, cualquier función plt funciona en el mapa; puede usar la instancia del mapa base para proyectar coordenadas de latitud y longitud a coordenadas (x, y) para trazar con plt.

Además de esto, hay muchas funciones específicas del mapa disponibles como métodos de la instancia del mapa base. Estos funcionan de forma muy similar a sus contrapartidas Matplotlib estándar, pero tienen un argumento booleano latlon adicional, que si se establece en True le permite pasar latitudes y longitudes sin formato al método, en lugar de las coordenadas proyectadas (x, y).

Algunos de estos métodos específicos del mapa son:

contour () / contourf (): Dibuja líneas de contorno o contornos rellenos imshow (): Dibuja una imagen pcolor () / pcolormesh (): dibuja un gráfico de pseudocolor para mallas irregulares / regulares plot (): Dibuja líneas y / o marcadores. scatter (): Dibuja puntos con marcadores. carcaj (): dibujar vectores. púas (): dibuja barbas de viento. drawgreatcircle (): Dibuja un gran círculo.

Veremos algunos ejemplos de algunos de ellos mientras continuamos. Para obtener más información sobre estas funciones, incluidos varios diagramas de ejemplo, consulte la documentación en línea del mapa base http://matplotlib.org/basemap/

Ejemplo: Ciudades de California

Recordemos que en Customizing Plot Legends, demostramos el uso del tamaño y el color en un diagrama de dispersión para transmitir información sobre la ubicación, el tamaño y la población de las ciudades de California. Aquí, crearemos este gráfico de nuevo, pero usando Basemap para poner los datos en contexto.

Comenzamos cargando los datos, como lo hicimos antes:

In [48]:

import pandas as pd

cities = pd.read_csv('data/california_cities.csv')

​

# Extract the data we're interested in

lat = cities['latd'].values

lon = cities['longd'].values

population = cities['population_total'].values

area = cities['area_total_km2'].values

A continuación, configuramos la proyección del mapa, dispersamos los datos y luego creamos una barra de colores y una leyenda:

In [49]:

# 1. Draw the map background

fig = plt.figure(figsize=(8, 8))

m = Basemap(projection='lcc', resolution='h',

            lat_0=37.5, lon_0=-119,

            width=1E6, height=1.2E6)

m.shadedrelief()

m.drawcoastlines(color='gray')

m.drawcountries(color='gray')

m.drawstates(color='gray')

​

# 2. scatter city data, with color reflecting population

# and size reflecting area

m.scatter(lon, lat, latlon=True,

          c=np.log10(population), s=area,

          cmap='Reds', alpha=0.5)

​

# 3. create colorbar and legend

plt.colorbar(label=r'$\log_{10}({\rm population})$')

plt.clim(3, 7)

​

# make legend with dummy points

for a in [100, 300, 500]:

    plt.scatter([], [], c='k', alpha=0.5, s=a,

                label=str(a) + ' km$^2$')

plt.legend(scatterpoints=1, frameon=False,

           labelspacing=1, loc='lower left');

Esto nos muestra aproximadamente donde poblaciones más grandes de personas se han asentado en California: se agrupan cerca de la costa en las áreas de Los Ángeles y San Francisco, se extienden a lo largo de las carreteras en el plano valle central y evitan casi por completo las regiones montañosas a lo largo de las fronteras el estado.

Ejemplo: datos de temperatura superficial

Como ejemplo de visualizar algunos datos geográficos más continuos, consideremos el "vórtice polar" que golpeó la mitad oriental de los Estados Unidos en enero de 2014. Una gran fuente para cualquier tipo de datos climáticos es el Instituto Goddard de Estudios Espaciales de la NASA. Aquí usaremos los datos de temperatura de GIS 250, que podemos descargar usando comandos de shell (estos comandos pueden tener que ser modificados en máquinas con Windows). Los datos utilizados aquí se descargaron el 6/12/2016 y el tamaño del archivo es de aproximadamente 9 MB:

In [ ]:

# !curl -O http://data.giss.nasa.gov/pub/gistemp/gistemp250.nc.gz

# !gunzip gistemp250.nc.gz

Los datos vienen en formato NetCDF, que la biblioteca netCDF4 puede leer en Python. Puede instalar esta biblioteca como se muestra aquí

$ conda install netcdf4

Leemos los datos de la siguiente manera :

In [50]:

from netCDF4 import Dataset

data = Dataset('data/gistemp250.nc')

El archivo contiene muchas lecturas de temperatura global en una variedad de fechas; tenemos que seleccionar el índice de la fecha que nos interesa, en este caso, 15 de enero de 2014:

In [51]:

from netCDF4 import date2index

from datetime import datetime

timeindex = date2index(datetime(2014, 1, 15),

                       data.variables['time'])

Ahora podemos cargar los datos de latitud y longitud, así como la anomalía de temperatura para este índice:

In [52]:

lat = data.variables['lat'][:]

lon = data.variables['lon'][:]

lon, lat = np.meshgrid(lon, lat)

temp_anomaly = data.variables['tempanomaly'][timeindex]

Finalmente, usaremos el método pcolormesh () para dibujar una malla de color de los datos. Veremos América del Norte y utilizaremos un mapa en relieve sombreado en el fondo. Tenga en cuenta que para estos datos elegimos específicamente un mapa de color divergente, que tiene un color neutral en cero y dos colores que contrastan en valores negativos y positivos. También dibujaremos ligeramente las costas sobre los colores para referencia:

In [53]:

fig = plt.figure(figsize=(10, 8))

m = Basemap(projection='lcc', resolution='c',

            width=8E6, height=8E6,

            lat_0=45, lon_0=-100,)

m.shadedrelief(scale=0.5)

m.pcolormesh(lon, lat, temp_anomaly,

             latlon=True, cmap='RdBu_r')

plt.clim(-8, 8)

m.drawcoastlines(color='lightgray')

​

plt.title('January 2014 Temperature Anomaly')

plt.colorbar(label='temperature anomaly (°C)');

Los datos muestran una imagen de las anomalías localizadas de temperatura extrema que ocurrieron durante ese mes. La mitad oriental de los Estados Unidos era mucho más fría de lo normal, mientras que la mitad occidental y Alaska eran mucho más cálidas. Las regiones sin temperatura registrada muestran el fondo del mapa.