xander_bakker

Usar gráficos de radar para evaluar el impacto

Blog Post created by xander_bakker on Sep 28, 2015

En la presentación sobre Python en ArcGIS Pro de la conferencia CCU2015 mostramos un gráfico de radar para evaluar el impacto de cada componente (ingeniería, entorno construido, socio-ambiental y la suma de los tres) sobre cada uno de ellos.

Para este fin utilizamos el ejemplo de código de Matplotlib de la siguiente página:

api example code: radar_chart.py — Matplotlib 1.4.3 documentation

 

Los datos raster y vectoriales son basado sobre el cálculo de corredores (según el método de EPRI). Los rasters representen raster de costos por celda que indican si es fácil (barato) o difícil (costoso) para viajar a través de esta celda. Los polígonos corresponden al corredor de mejor aptitud para la construcción de la linea de transmisión (en este caso 10% del superficie) conectando el punto de inicio con el punto final para cada componente (ingeniería, entorno construido, socio-ambiental y la suma de los tres):

Abajo el código Python que se utilizó en ArcGIS Pro para lograr el resultado:

 

import numpy as np import matplotlib.pyplot as plt from matplotlib.path import Path from matplotlib.spines import Spine from matplotlib.projections.polar import PolarAxes from matplotlib.projections import register_projection  def radar_factory(num_vars, frame='circle'):     """Create a radar chart with 'num_vars' axes.      This function creates a RadarAxes projection and registers it.      Parameters     ----------     num_vars : int         Number of variables for radar chart.     frame : {'circle' | 'polygon'}         Shape of frame surrounding axes.      """     # calculate evenly-spaced axis angles     theta = 2*np.pi * np.linspace(0, 1-1./num_vars, num_vars)     # rotate theta such that the first axis is at the top     theta += np.pi/2      def draw_poly_patch(self):         verts = unit_poly_verts(theta)         return plt.Polygon(verts, closed=True, edgecolor='k')      def draw_circle_patch(self):         # unit circle centered on (0.5, 0.5)         return plt.Circle((0.5, 0.5), 0.5)      patch_dict = {'polygon': draw_poly_patch, 'circle': draw_circle_patch}     if frame not in patch_dict:         raise ValueError('unknown value for `frame`: %s' % frame)      class RadarAxes(PolarAxes):          name = 'radar'         # use 1 line segment to connect specified points         RESOLUTION = 1         # define draw_frame method         draw_patch = patch_dict[frame]          def fill(self, *args, **kwargs):             """Override fill so that line is closed by default"""             closed = kwargs.pop('closed', True)             return super(RadarAxes, self).fill(closed=closed, *args, **kwargs)          def plot(self, *args, **kwargs):             """Override plot so that line is closed by default"""             lines = super(RadarAxes, self).plot(*args, **kwargs)             for line in lines:                 self._close_line(line)          def _close_line(self, line):             x, y = line.get_data()             # FIXME: markers at x[0], y[0] get doubled-up             if x[0] != x[-1]:                 x = np.concatenate((x, [x[0]]))                 y = np.concatenate((y, [y[0]]))                 line.set_data(x, y)          def set_varlabels(self, labels):             self.set_thetagrids(theta * 180/np.pi, labels)          def _gen_axes_patch(self):             return self.draw_patch()          def _gen_axes_spines(self):             if frame == 'circle':                 return PolarAxes._gen_axes_spines(self)             # The following is a hack to get the spines (i.e. the axes frame)             # to draw correctly for a polygon frame.              # spine_type must be 'left', 'right', 'top', 'bottom', or `circle`.             spine_type = 'circle'             verts = unit_poly_verts(theta)             # close off polygon by repeating first vertex             verts.append(verts[0])             path = Path(verts)              spine = Spine(self, spine_type, path)             spine.set_transform(self.transAxes)             return {'polar': spine}      register_projection(RadarAxes)     return theta  def unit_poly_verts(theta):     """Return vertices of polygon for subplot axes.      This polygon is circumscribed by a unit circle centered at (0.5, 0.5)     """     x0, y0, r = [0.5] * 3     verts = [(r*np.cos(t) + x0, r*np.sin(t) + y0) for t in theta]     return verts   def example_data():     import arcpy, os      fcs = [r"D:\Xander\LineasTransmision\SanLorenzo\gdb\MacroCorridor.gdb\MacroCorridor_total",            r"D:\Xander\LineasTransmision\SanLorenzo\gdb\MacroCorridor.gdb\MacroCorridor_ingenieria",            r"D:\Xander\LineasTransmision\SanLorenzo\gdb\MacroCorridor.gdb\MacroCorridor_construido",            r"D:\Xander\LineasTransmision\SanLorenzo\gdb\MacroCorridor.gdb\MacroCorridor_socioamb"]      lst_names = ['Corredor total', 'Corredor ingenieria',                  'Corredor entorno construido', 'Corredor socio ambiental']      rasters = [r"D:\Xander\LineasTransmision\SanLorenzo\gdb\AnalysisRaster.gdb\cmp_ingenieria",                r"D:\Xander\LineasTransmision\SanLorenzo\gdb\AnalysisRaster.gdb\cmp_construido",                r"D:\Xander\LineasTransmision\SanLorenzo\gdb\AnalysisRaster.gdb\cmp_socioamb",                r"D:\Xander\LineasTransmision\SanLorenzo\gdb\AnalysisRaster.gdb\cmp_total"]      fld_val = "Value"     fld_cnt = "Count"      all_lists = []     for fc in fcs:         ws, fc_name = os.path.split(fc)         arcpy.env.mask = fc          lsts = []         for ras in rasters:             rasi = arcpy.sa.Int(ras)             dct = {r.getValue(fld_val): r.getValue(fld_cnt) for r in arcpy.SearchCursor(rasi)}             tot =  sum(dct.values())             lst = []             for i in range(1, 10):                 if i in dct:                     val = float(dct[i]) * 100 / float(tot)                     lst.append(val)                 else:                     lst.append(0.0)              lsts.append(lst)         all_lists.append(lsts)      data = {         'column names':             ['1 (mas economico)', '2', '3', '4', '5 (costos elevados)', '6',              '7', '8', '9 (mas costoso)'],             lst_names[0] : all_lists[0],             lst_names[1] : all_lists[1],             lst_names[2] : all_lists[2],             lst_names[3] : all_lists[3]             }     return data   if __name__ == '__main__':     N = 9     theta = radar_factory(N, frame='polygon')      data = example_data()     spoke_labels = data.pop('column names')      fig = plt.figure(figsize=(15, 10))     fig.subplots_adjust(wspace=0.5, hspace=0.50, top=0.85, bottom=0.05)      colors = ['b', 'r', 'g', 'm', 'y']      # Plot the four cases from the example data on separate axes     for n, title in enumerate(data.keys()):         ax = fig.add_subplot(2, 2, n+1, projection='radar')         plt.rgrids([0.2, 0.4, 0.6, 0.8])         ax.set_title(title, weight='bold', size='medium', position=(0.5, 1.1),                      horizontalalignment='center', verticalalignment='center')         for d, color in zip(data[title], colors):             ax.plot(theta, d, color=color)             ax.fill(theta, d, facecolor=color, alpha=0.25)         ax.set_varlabels(spoke_labels)      # add legend relative to top-left plot     plt.subplot(2, 2, 1)     labels = ('Ingenieria', 'Entorno construido', 'Socio Ambiental', 'Suma 3 componentes')      legend = plt.legend(labels, loc=(0.9, .95), labelspacing=0.1)     plt.setp(legend.get_texts(), fontsize='small')      plt.figtext(0.5, 0.965, 'macro corredores vs costos de componentes',                 ha='center', color='black', weight='bold', size='large')     plt.show()

 

Para los primeros 100 lineas de código, no apliqué ningún cambio en el código Python. La función "example_data" (empezando en la linea 101) es donde puse el código para conectarse a mis datos (raster y vectorial), con el fin de:

  • crear un raster por componente que contiene valores enteros (1-9) de los costos
  • generar un diccionario por componente con los datos de valor de celda vs conteo de celdas con este valor
  • generar una lista del porcentaje de cada valor (1 a 9) por combinación de polígono de corredor y raster (4 x 4 = 16)
  • en la linea 141, se crea un diccionario, siguiendo el formato de los datos que se está utilizando en este ejemplo de radar.

 

En la función main, se cambiaron unos rótulos y el tamaño del gráfico y en verdad eso fue todo.

 

Si miramos con más detalla a un ejemplo de radar:

... podemos observar que el escenario del corredor con costos solamente basados en el entorno construido, tiene un impacto grande sobre el componente socio ambiental (ver el pique verde en el valor 7).

Outcomes