Un amigo y lector del blog me ha pedido un mapa de códigos postales donde poder identificar los centroides para andar calculando distancias a otros puntos. Yo no tengo un mapa de España por códigos postales para poder usar con fines comerciales, pero sí cuento en el blog cómo poder obtenerlo bajo ciertas condiciones. Lo que sí puedo contar a Juan es cómo hacer un mapa por municipios con SAS; aunque ya he hablado de ello, hay ciertos aspectos que pueden ser interesantes. Y todo empieza donde siempre: http://www.gadm.org/country, la web donde tenemos los mapas «libres» por países. Seleccionáis Spain y el formato shapefile; una vez descargados los mapas en vuestros equipos, empezamos con el trabajo en SAS:
proc mapimport datafile = "\\directorio\\mapa\\ESP_adm_shp.shp"
out = work.espania;
run;
proc contents;
run;
quit;
El procedimiento MAPIMPORT ha creado un conjunto de datos SAS donde tenemos caracterizados todos los polígonos que componen el shapefile. Entonces, si tenemos que calcular el centroide de un municipio con SAS, sugiero realizar un PROC SQL de la siguiente forma:
proc sql;
create table centroides as
select id_4,
avg(x) as centroide_x,
avg(y) as centroide_y
from espania
group by 1;
quit;
¡A que molo! El centroide de cada municipio será la media de las coordenadas que lo definen y ahora viene el motivo por el cual he creado una entrada en el blog en vez de llamar por teléfono a Juanito, que se ha ido muy lejos a trabajar y no me invita. Vamos a realizar el mapa de la provincia de Barcelona con SAS con el procedimiento SGPLOT, que tiene una serie de matices que me hacen sugerir que empleéis R para la realización de este tipo de mapas:
* Creación de un mapa de Barcelona;
data barna;
set espania;
if name_2 eq "Barcelona";
orden = _n_;
run;
proc sql;
create table barna as
select distinct a.*, b.*
from barna a
left join centroides b on a.id_4 eq b.id_4
order by orden;
quit;
proc sgplot data = barna;
polygon x = x y = y ID = id_4 / fill outline;
xaxis display = none;
yaxis display = none;
scatter x = centroide_x y = centroide_y;
run;
quit;
Matices y rarezas que tiene este código SAS: necesita el orden del shapefile; no puedo explicar por qué, lo investigaré, pero los primeros mapas no me salieron e intuí que por ahí iban los problemas. Además, el left join ha necesitado un distinct, algo que no entiendo muy bien ya que se puede realizar el mapa con los mismos datos sin el left join y este no genera duplicados; esto sí lo he investigado y el resultado ha sido PROBLEMA 1 – RAÚL 0, volveré con ello. Por otro lado, tenemos el PROC SGPLOT, que emplea POLYGON para la representación de los shapefile y, como deseamos que nos pinte un punto para identificar los centroides, usamos la sentencia SCATTER. Bueno, espero que a Juan y a otros os sirva esta entrada; imagino que sí porque ahora será capaz de representar en mapas la distancia de los centroides de un municipio a un punto determinado o incluso la adyacencia entre municipios. Saludos.