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.
Aviso importante!
Este tutorial necessita de conhecimentos básicos da linguagem por parte do usuário. Desta forma, caso você queira uma boa fonte para se iniciar na linguagem, recomendo visitar o sítio web da disciplina Preparação de dados para Análise Estatística do Programa de Pós-graduação em Ciências Biológicas (Botânica) do INPA. Lá você encontrará rica oferta de materiais para se familiarizar com o R e ser capaz de seguir este tutorial.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
Boa!! Obrigado Labotan!! Passo a passo bem organizado!
ResponderExcluirValeu pelo elogio. Fique ligado que a intenção é postar frequentemente ações de rotina ou coisas interessantes com o R, então, dê sempre sua visita e opinião, será sempre bem vinda.
ResponderExcluirNossa que legal! Vou tentar fazer um mapa, ajuda bem vinda na hora certa. Muito obrigada!Qq coisa heip me rsrs
ResponderExcluirmuito obrigada, vamos testar a ferramenta!
ResponderExcluirTeste mesmo. E me mantenha informado para saber como foi. Alguns colegas não obtiveram sucesso com a importação dos dados. Eu testei em meu MAC e já testei em um Ubuntu 14.04, e deu certo. Abraços
ResponderExcluirR.
Saciiii eu te amo!!
ResponderExcluirMuito Obrigada!!!
ResponderExcluir