O R é uma
linguagem ou software que tem se tornado importante ferramenta para análise de dados espaciais e produção de mapas. Como exemplo, sugiro visitar o sítio web
R-bloggers e fazer uma pequena busca utilizando as palavras “plot maps”. Vejam quantos resultados aparecerão. Para mais usos e maneiras de produzir mapas no R, veja também o sítio web
Spatial data in R: Using R as a GIS.
Nesta postagem, quero demonstrar o quanto podemos produzir mapas informativos de uma maneira bem simples, com poucas linhas de comando. Como exemplo, vejam que com duas linhas apenas, podemos plotar o mapa mundi.
library(maps)
map()
O tutorial visa a produção de um mapa de distribuição geográfica utilizando dados que abrangem a distribuição espacial de duas espécies de
Burseraceae,
Protium aracouchini (Aubl.) Marchand e
P. heptaphyllum (Aubl.) Marchand. A primeira planta faz parte de um complexo de espécies informalmente denominado completo
Protium aracouchini que é o modelo de estudo de
minha tese de doutorado. A segunda espécie é o modelo de estudo de um amigo,
Gabriel Damasco, colaborador do LABOTAM e sob orientação do Dr.
Paul Fine, também meu orientador.
Os dados foram obtidos através da base de dado do
herbário virtual do Jardim Botânico de Nova Iorque.
Os tópicos a serem explorados neste tutorial são:
- Importar os dados para o R;
- Verificar e limpar os dados;
- Preparação de variáveis a serem utilizadas no mapa;
- Plotar o mapa.
Passo 1 - Instalar o R
O passo essencial deste tutorial é ter a última versão do R instalado em seu computador. Para obtê-lo, visite a página do projeto
R-project.
Caso haja dúvidas sobre o processo de instalação, favor visitar este
link.
Passo 2 - Chamar os pacotes
Faremos uso dos pacotes
maps (para plotar o mapa),
RColorBrewer (paleta de cores),
magrittr (pipeline) e
dplyr (manipulação de ‘dataframes’) . É imprescindível que todos estejam atualizados e com as dependências instaladas.
Iniciando a sessão, execute os comandos abaixo para limpar a área de trabalho E chamar os pacotes para dentro da sessão.
rm(list=ls()) #limpa COMPLETAMENTE o ambiente de trabalho
#chama os pacotes
libs <- c('maps','RColorBrewer','magrittr','dplyr')
sapply(libs,library,character.only=T,logical.return=T)
Passo 3 - Importar os dados para a sessão
Os dados estão reunidos em dois arquivos ‘.csv’ separados por tabulação e codificados em UTF-8. Baixe-os para seu computador:
Cabe ressaltar que vários formatos podem ser utilizados para importar dados ao R. Utilizo o ‘
CSV’ porque é mais simples.
Seguindo o tutorial, geralmente estabeleço inicialmente o diretório de trabalho usando a função
setwd(). Desta forma, posso executar o script por inteiro, sem me preocupar em escolher o diretório utilizando botões, enfim, é mais prático. Ressalto que isso deve alterado por você para o seu diretório de trabalho!
#estabelece o diretorio de trabalho
setwd('/Users/ricoperdiz/Documents/DOC/PROJETO_DOC/LABOTAM/data')
Chamamos o arquivo para dentro da sessão e utilzamos o operador %>% para passar o resultado da importação para a próxima ação, que é a de ordenar esses dados segundo o epíteto específico, nome do coletor e número de coleta. Ao resultado dessa ação é atribuído um objeto,
‘cpa’. Veja como estão as primeiras linhas desse objeto após a ordenação.
#importa os dados - importante que haja uma coluna com dados de lat e outra de long
#protium aracouchini
read.table("nybg_paracouchini.csv", header = T, as.is = T, sep = '\t', dec = '.') %>%
select(recordedBy, recordNumber, decimalLatitude, decimalLongitude, identifiedBy, specificEpithet) %>%
arrange(recordedBy, recordNumber) -> aracouch
#protium heptaphyllum
read.table("nybg_pheptaphyllum.csv", header = T, as.is = T, sep = '\t', dec = '.') %>%
select(recordedBy, recordNumber, decimalLatitude, decimalLongitude, identifiedBy, specificEpithet) %>%
arrange(recordedBy, recordNumber) -> heptaphy
Passo 4 - Verifica os dados
Em nosso caso, devemos verificar: * existência valores vazios; em caso positivo, devemos eliminá-los; * confiabilidade dos valores de latitude e longitude, às vezes, por diversos fatores, há troca de sinais (negativos e positivos) ocasionando equívocos quanto à ocorrência exata da amostra. Se for percebido algo assim, é bom checar os dados e buscar corrigí-los.
############################################################
############################################################
### Limpeza e verificação de dados
#une os dados e passa para a proxima acao
rbind(aracouch,heptaphy) %>%
#elimina os registros vazios de coletor
filter(recordedBy != '') %>%
#elimina os registros sem lat ou long
filter(decimalLatitude != '' | decimalLongitude != '') %>%
#filtra apenas os especimes identificados pelo especialista
filter(identifiedBy == 'D. C. Daly') -> dados
#ha registros duplicados
#busca-se entao apenas os dados unicos
#limpa os dados de coletor e numero para poder criar um identificador
#limpeza consiste em eliminar '.', espacos vazios, apostrofe e '_' duplos
gsub('\\.', '_', dados$recordedBy) %>%
gsub(' ', '_', .) %>%
gsub("'", '_', .) %>%
gsub('__', '_', .) -> dados$recordedBy
dados$recordedBy
#faz-se o mesmo para os numeros de coleta
gsub('/','_', dados$recordNumber) %>%
gsub(' ', '_', .) %>%
gsub('\\.', '_', .) -> dados$recordNumber
dados$recordNumber
#cria o identificador de coleta e especie
dados %>%
mutate( ID = paste(recordedBy,recordNumber, sep = '_'), Species = paste('Protium', specificEpithet, sep = ' ')) -> dados
dados$ID <- gsub('__','_',dados$ID)
#quem sao os dados unicos
unicos <- unique(dados$ID)
unicos
#agora filtra os dados unicos no dataframe, eliminando os duplicados
#faz-se uso da funcao match para obter esse resultado
match(unicos, dados$ID) %>%
dados[.,] -> prot
#verifica a cobertura de lat e long para ver se estao dentro
# dos limites da America do Sul
range(prot$decimalLatitude)
range(prot$decimalLongitude) #aqui tem algo estranho
#percebe-se aqui que ha valores que caem fora da Am Sul
#limite e pouco mais de -80
head(sort(prot$decimalLongitude))
#devemos eliminar
prot %>%
filter(decimalLongitude > -80) -> prot
Passo 5 - Prepara variáveis para plotar o mapa
#vetor para lat e long
lat <- prot$decimalLatitude
long <- prot$decimalLongitude
#vetor com epitetos especificos
spp <- unique(prot$Species)
#vetor de cores para utilizar no mapa
#para cada especie, uma cor
prot$cores.map <- ifelse(prot$Species == spp[1], 'red','black')
cores.map <- prot$cores.map
#cria um vetor de tamanho para cada especie
#como uma coluna de prot
#P heptaphyllum possui uma distribuicao mais ampla
#por isso atribuo um tamanho menor pra ela
prot$cex.p <- ifelse(prot$Species == spp[1], 1, 0.8)
#vetor para lat e long
y1 <- range(lat) + c(-1,1)
x1 <- range(long) + c(-1,1)
#vetor para os simbolos pch
prot$pontos <- ifelse(prot$Species == spp[1], 21, 24)
#vetor para paises
paises <- c('Brazil','Argentina','Peru','Paraguay','Ecuador','Chile','Uruguay','French Guiana','Suriname','Venezuela','Colombia','Guyana','Bolivia','Panama','Costa Rica')
Passo 6 - Plota o mapa
#plota o mapa
map(regions = paises, fill = F, xlim = x1, ylim = y1)
#plota os pontos de ocorrencia de cada especie
points(long,lat, pch = prot$pontos, col = cores.map, bg = cores.map, cex = prot$cex.p)
#coloca eixos das coordenadas
map.axes()
axis(side=4,las=1)
axis(side=3,las=1)
#coloca escala do mapa
par(cex=1, las=1)
map.scale(max(long) - 12, ratio = F, cex = 1, metric = T)
#plota a linha do equador
abline(h=0,lwd=0.5,lty="dotted")
#nomeia a linha do equador
text(x = max(long) - 1, y = 1,labels="Equador", cex=1, adj=c(1,0.5))
#plota uma legenda
legend(max(long) - 16 , max(lat),legend = spp, pch = unique(prot$pontos), pt.bg = unique(cores.map), cex = 0.8, x.intersp = 0.8, text.font = 3)
Pode-se salvar também este mapa em um pdf, utilizando a função
pdf(). Veja abaixo.
### Cria um pdf e Plota o mapa
#inicia o pdf
pdf('meu_mapa_complexo_protium_aracouchini.pdf',height=8,width=10)
map(regions = paises, fill = F, xlim = x1, ylim = y1)
points(long,lat, pch = prot$pontos, col = cores.map, bg = cores.map, cex = prot$cex.p)
map.axes()
axis(side=4,las=1)
axis(side=3,las=1)
par(cex=1, las=1)
map.scale(max(long) - 12, ratio = F, cex = 1, metric = T)
abline(h=0,lwd=0.5,lty="dotted")
text(x = max(long) - 1, y = 1,labels="Equador", cex=1, adj=c(1,0.5))
legend(max(long) - 16 , max(lat),legend = spp, pch = unique(prot$pontos), pt.bg = unique(cores.map), cex = 0.8, x.intersp = 0.8, text.font = 3)
dev.off()
Conclusão
Simples demais… Brincadeira! Com o tempo fica fácil.
Todavia, há muito outras possibilidades, como, por exemplo, pode-se importar
shapefiles, plotar rios, colocar cores etc. Em próximo post, pretendo demonstrar como importar
shapefiles para dentro da sessão e utilizá-los para incrementar o mapa.
Dúvidas, sugestões, favor entrar em contato via comentários.
O script completo deste tutorial pode ser encontrado
aqui.
Esta postagem foi escrita em ambiente R utilizando os pacotes
R Markdown e
knitr. Você pode reproduzir a postagem baixando o arquivo ‘Rmd’ em seu computador (baixe-o
aqui, em arquivo zipado; faça a descompressão antes de fazer uso do mesmo) e processando-a com o pacote ‘knitr’. Por questões de portabilidade entre R e a plataforma
Blogger, este arquivo ‘Rmd’ não contem as imagens desta postagem.
Divirta-se!!!
Para saber mais
