Rggplot2 mescla com arquivos shapefile e csv para preencher polígonos

Produzimos diariamente mapas que mostram um nível calculado para a temperatura em 30 áreas distintas da nossa região, cada área é preenchida com uma cor diferente, dependendo do nível. Esses mapas parecem

insira a descrição da imagem aqui

Agora eu quero mudar a geração de mapas para R. Eu baixei limites provinciais e municipais (você pode encontrar limites para toda a Espanha ou aqui o subconjunto para minha região ) e consegui plotar eles com ggplot2 seguindo o exemplo de Hadley.

Eu também posso produzir um arquivo ascii que contém duas colunas: identificador (CODINE) e nível diário. Você pode baixar aqui .

Este é meu primeiro script tentando plotar shapefiles com R e ggplot2, então pode haver erros e com certeza ele pode ser melhorado, sugestões bem-vindas. O código a seguir (baseado no Hadley mencionado anteriormente) funciona para mim:

> require("rgdal") > require("maptools") > require("ggplot2") > require("plyr") # Reading municipal boundaries esp = readOGR(dsn=".", layer="lineas_limite_municipales_etrs89") muni=subset(esp, esp$PROV1 == "46" | esp$PROV1 == "12" | esp$PROV1 == "3") muni@data$id = rownames(muni@data) muni.points = fortify(muni, region="id") muni.df = join(muni.points, muni@data, by="id") # Reading province boundaries prov = readOGR(dsn=".", layer="poligonos_provincia_etrs89") pr=subset(prov, prov$CODINE == "46" | prov$CODINE == "12" | prov$CODINE == "03" ) pr@data$id = rownames(pr@data) pr.points = fortify(pr, region="id") pr.df = join(pr.points, pr@data, by="id") ggplot(muni.df) + aes(long,lat,group=group) + geom_path(color="blue") + + coord_equal()+ geom_path(data=pr.df, + aes(x=long, y=lat, group=group),color="red", size=0.5) 

Este código traça um bom mapa com todos os limites insira a descrição da imagem aqui

Para o preenchimento de polígonos por nível, tentei ler e depois mesclar como sugerido em http://tormodboe.wordpress.com/2011/02/22/g%C3%B8y-med-kart-2/

level = read.csv (“levels.dat”, header = T, sep = “”)
munlevel = mesclar (muni.df, nível, por = “CODINE”)

mas dá um erro

Erro en fix.by (by.x, x): ‘por’ deve especificar uma coluna unicamente válida

Eu não estou familiarizado com shapefiles, talvez eu precise aprender mais sobre os atributos de dados shp para encontrar a escolha certa para mesclar os dois conjuntos de dados. Como posso mesclar dados para poder traçar as linhas (limites municipais) e depois preenchê-los com níveis?

[NB: Esta pergunta foi feita há mais de um mês, então o OP provavelmente encontrou uma maneira diferente de resolver seu problema. Eu tropecei nisso enquanto trabalhava nessa questão relacionada . Esta resposta está incluída na esperança de que isso beneficiará alguém.]

Isso parece ser o que o OP está pedindo …

… e foi produzido com o seguinte código:

 require("rgdal") require("maptools") require("ggplot2") require("plyr") # read temperature data setwd("") temp.data <- read.csv(file = "levels.dat", header=TRUE, sep=" ", na.string="NA", dec=".", strip.white=TRUE) temp.data$CODINE <- str_pad(temp.data$CODINE, width = 5, side = 'left', pad = '0') # read municipality polygons setwd(" 

Explicação:

Shapefiles importados para R com readOGR(...) são do tipo SpacialDataFrame e têm duas seções principais: uma seção de polígonos que contém as coordenadas de todos os pontos em cada polígono e uma seção de dados que contém informações sobre cada polígono (portanto, um linha por polígono). Estes podem ser referenciados, por exemplo, usando muni@polygons e muni@data . A function de utilidade fortify(...) converte a seção de polígonos em uma estrutura de dados organizada para plotagem com ggplot . Portanto, o stream de trabalho básico é:

 [1] Import temperature data file (temp.data) [2] Import polygon shapefile of municipalities (muni) [3] Convert muni polygons to a data frame for plotting (muni.df <- fortify(...)) [4] Join columns from muni@data to muni.df [5] Join columns from temp.data to muni.df [6] Make the plot 

As junções devem ser feitas em campos comuns, e é aí que a maioria dos problemas aparece. Cada polígono no arquivo de forma original possui um atributo de ID exclusivo. Running fortify(...) no shapefile cria uma coluna, id , que é baseada nisso. Mas não há coluna de ID na seção de dados. Em vez disso, os IDs de polígono são armazenados como nomes de linha. Então, primeiro devemos adicionar uma coluna id para muni@data seguinte forma:

 muni@data$id <- rownames(muni@data) 

Agora temos um campo de id em muni@data e um campo de id correspondente em muni.df , para que possamos fazer a junit:

 muni.df <- join(muni.df, muni@data, by="id") 

Para criar o mapa, precisamos definir colors de preenchimento com base no nível de temperatura. Para fazer isso, precisamos nos juntar à coluna LEVEL de temp.data para muni.df Em temp.data existe um campo CODINE que identifica o município. Há também, agora, um campo correspondente CODIGOINE em muni.df Mas há um problema: CODIGOINE é char(5) , com zeros à esquerda, enquanto CODINE é inteiro, o que significa falta de zeros à esquerda (importados do Excel, talvez?). Então apenas juntar-se a esses dois campos não produz nenhum resultado. Devemos primeiro converter CODINE em char(5) com zeros à esquerda:

 temp.data$CODINE <- str_pad(temp.data$CODINE, width = 5, side = 'left', pad = '0') 

Agora podemos unir o temp.dat ao muni.df base nos campos correspondentes.

 muni.df <- merge(muni.df, temp.data, by.x="CODIGOINE", by.y="CODINE", all.x=T, a..ly=F) 

Usamos merge(...) invés de join(...) porque os campos de junit possuem nomes diferentes e join(...) requer que eles tenham o mesmo nome. (Note, no entanto, que join(...) é mais rápido e deve ser usado se possível). Então, finalmente, temos um quadro de dados que contém todas as informações para traçar os polígonos e a temperatura LEVEL que pode ser usada para estabelecer a cor de preenchimento para cada polígono.

Algumas notas no código original do OP:

  1. O primeiro mapa do OP (o verde no topo) identifica "30 áreas distintas para a nossa região ...". Eu não encontrei shapefile identificando essas áreas. O arquivo municipal identifica 543 municípios, e eu não via como agrupá-los em 30 áreas. Além disso, o arquivo de nível de temperatura tem 542 linhas, uma para cada município (mais ou menos).

  2. OP estava importando shapefiles de linha para o município para desenhar os limites. Você não precisa disso porque geom_polygon(...) irá desenhar (e preencher) os polígonos e geom_path(...) irá desenhar os limites.