Seção 1: Introdução
Em março de 2019, o Projeto de Mapeamento Anual da Cobertura e Uso do Solo do Brasil (MapBiomas) (http://mapbiomas.org/) lançou o banco de dados Mapbiomas v.3.1 sobre uso e cobertura da terra para o período de 1985 a 2017.
Este conjunto de dados contém informações sobre a área total (em hectares) de formações florestais, formações naturais não-florestais, agricultura, áreas sem vegetação e corpos d’água - informações detalhadas sobre estas formações e suas subdivisões podem ser encontradas aqui.
E ainda que o site forneça ferramentas para qualquer um criar seus próprios mapas, decidi criar meu próprio mapa animado usando R e os pacotes ggplot2 e gifski porque queria ver a porcentagem de uma área da cidade que é coberta por vegetação nativa.
Eu também estou usando o pacote brazilmaps para traçar o mapa do Brasil.
Eu queria me concentrar apenas nos dados da Formação Florestal Natural e ver como sua área total evoluiu desde 1985.
Este tutorial tem 3 seções, além desta introdução. A segunda seção é apenas sobre os pacotes e fontes que eu usei e como baixar essas fontes. Na Seção 3, examino o conjunto de dados e falo sobre suas variáveis e suas definições. Eu também passo por etapas de importar o conjunto de dados e fazer algumas manipulações de dados para obter as informações que desejo. A seção 4 fornece o código para o mapa final.
Você pode baixar o arquivo .R com o código completo do meu GitHub.
Seção 2: Pacotes
Começarei examinando as fontes que usei e como instalá-las. Eu uso o Windows, então não sei se os mesmos passos se aplicam a outros sistemas (desculpe usuários Mac e Unix).
Estou usando Lato, uma família de fontes de código aberto. Então, certifique-se de baixar e instalar a fonte Lato antes de começar ou você pode deixar R usar a fonte padrão do sistema - você não terá que alterar o código, você receberá algumas mensagens de aviso, mas R ainda irá traçar o mapa .
Você pode baixá-la aqui.
A primeira coisa que faço é ter certeza de que R pode usar essa fonte, então vou usar o pacote extrafont para isso. Se é a primeira vez que você usa essas fontes, você terá que executar o comando font_import () (isso pode demorar um pouco dependendo do número de fontes instaladas que você possui, o bom é que você terá que executar este comando apenas uma vez).
# ************************************************************************* ----
# Pacote extrafont ----
# ************************************************************************* ----
#install.packages("extrafont")
library("extrafont")
# ************************************************************************* ----
# Importar Fonte ----
# ************************************************************************* ----
# Instale a fonte Lato
# Rode font_import apenas uma vez:
# font_import()
# Verifique se Lato foi incluida
fonts()[grep("Lato", fonts())]
# Registre as fontes no R
loadfonts(device = "win")
# ************************************************************************* ----
# Outros Pacotes ----
# ************************************************************************* ----
#install.packages("openxlsx")
#install.packages("tidyverse")
#install.packages("brazilmaps")
#install.packages("gifski")
library("openxlsx")
library("tidyverse")
library("brazilmaps")
library("gifski")
Section 3: Os Dados
Dando os devidos créditos
Os dados do Projeto MapBiomas são públicos, abertos e gratuitos através da simples referência da fonte:
“Projeto MapBiomas – Coleção 3.1 da Série Anual de Mapas de Cobertura e Uso de Solo do Brasil, acessado em 24/04/2019 através do link:(http://mapbiomas.org/pages/estatisticas)"
“Projeto MapBiomas - é uma iniciativa multi-institucional para gerar mapas anuais de cobertura e uso do solo a partir de processos de classificação automática aplicada a imagens de satélite. A descrição completa do projeto encontra-se em http://mapbiomas.org“
Também estou usando dados sobre a área das cidades no Brasil do IBGE (Instituto Brasileiro de Geografia e Estatística). Devido às mudanças no tamanho das cidades ao longo dos anos, decidi usar os dados de 2016. Eu estava obtendo alguns resultados incorretos usando dados de 2017 ou 2018, como cidades com mais de 100% de sua área coberta por vegetação nativa. Claramente, essas cidades foram divididas, mas isso não foi contabilizado no conjunto de dados do MapBiomas. Não foram observadas cidades com mais de 100% de sua área coberta por vegetação nativa ao utilizar dados de 2016.
O código abaixo importa os conjuntos de dados para R, mas se você quiser, pode baixar os arquivos xlsx dos sites oficiais clicando em aqui para o banco de dados do MapBiomas e aqui para os dados do IBGE.
Para manter meu código, também estou fazendo o upload dos bancos de dados originais para o meu repositório GitHub e esses são os dados que vou usar no meu código abaixo.
Dicionário de Dados
O arquivo xlsx do MapBiomas possui as seguintes planilhas:
- INSTRUCTIONS;
- LAND COVER - BIOMAS e UF;
- LAND COVER - MUN_UF;
- CONSULTA BIOMA-UF; e
- CONSULTA MUN-UF.
Usarei a planilha LAND COVER - MUN_UF. As variáveis dessa planilha são:
| Variável | Classe | Descrição |
|---|---|---|
| COD_ESTADO | numeric | identificador do Estado |
| estado | character | nome do Estado |
| CODIBGE | numeric | identificador do municipio |
| municipio | character | nome do municipio |
| cod.classe | numeric | identificador do tipo de classe de cobertura terrestre |
| nivel1 | character | classe nível 1 |
| nivel2 | character | classe nível 2 |
| nivel3 | character | classe nível 3 |
| 1985 | numeric | área coberta em hectares(ha) em 1985 |
| … | … | … |
| 2017 | numeric | área coberta em hectares(ha) em 2017 |
O arquivo xlsx do IBGE possui as seguintes planilhas:
- AR_BR_MUN_2016;
- AR_BR_UF_2016;
- AR_BR_RG_2016; e
- AR_BR_2016.
Usarei a planilha AR_BR_MUN_2016. As variáveis dessa planilha são:
The variables in IBGE’s AR_BR_MUN_2016 data set are:
| Variável | Classe | Descrição |
|---|---|---|
| ID | numeric | ID |
| CD_GCUF | character | identificador do Estado |
| NM_UF | character | nome do Estado |
| NM_UF_SIGLA | character | UF |
| CD_GCMUN | character | identificador do municipio |
| NM_MUN_2016 | character | nome do municipio |
| AR_MUN_2016 | numeric | área em km quadrados |
# ************************************************************************* ----
# Importar os Dados ----
# ************************************************************************* ----
# Mudar formato dos numeros
options(scipen=10)
# Caminho para o Arquivo no GitHub
file_path_forest_area = paste0(
"https://github.com/gu-stat/forest-area/blob/master/",
"data/Dados_Cobertura_MapBiomas_3.1_Biomas-UF-Municipios_%20SITE_UNPROTECTED.xlsx",
"?raw=true")
file_path_city_area = paste0(
"https://github.com/gu-stat/forest-area/blob/master/",
"data/AR_BR_RG_UF_MUN_2016.xlsx",
"?raw=true")
# Importacao
df.original <- openxlsx::read.xlsx(file_path_forest_area,
sheet = "LAND COVER - MUN_UF")
df.city.area <- openxlsx::read.xlsx(file_path_city_area,
sheet = "AR_BR_MUN_2016")
## \__ Dados para Criar Mapa por Cidade ----
brazil.map <- get_brmap(geo = "City", class = "sf")
Manipulação de dados
# ************************************************************************* ----
# Manipulacao dos Dados ----
# ************************************************************************* ----
## \__ Transforma Codigo da Cidade em Numerico ----
df.city.area2 <-
df.city.area %>%
dplyr::mutate(CODIBGE = as.numeric(CD_GCMUN)) %>%
dplyr::select(CODIBGE, AR_MUN_2016)
## \__ Dados de Area de Cobertura de Florestas Natuais ----
df.nat.forest <-
df.original %>%
filter(nivel2 == "Floresta Natural") %>%
select(COD_ESTADO, estado, CODIBGE, nivel2, nivel3,
paste0(seq(1985, 2017))
) %>%
# SOMA OS DIFERENTES SUB-GRUPOS DE Floresta Natural
group_by(COD_ESTADO, estado, CODIBGE) %>%
summarize_at(vars(-nivel2,-nivel3),sum) %>%
ungroup()
## \__ Juntar com Dados de Area da Cidade ----
df.nat.forest.AREA <-
df.nat.forest %>%
# ADICIONA AREA DE CADA CIDADE NA BASE COM DADOS DE VEGETACAO
dplyr::left_join(x = ., y = df.city.area2, by = "CODIBGE") %>%
# CRIA PERCENTUAIS
## FAZ x/100 PARA TRANSFORMAR DE HECTARES PARA KM2
dplyr::mutate_at(.vars=vars(paste0(seq(1985,2017))),
.fun = function(x) round(100*((x/100)/.$AR_MUN_2016),3)
)
# VERIFICA SE ALGUM MUNICIPIO TEM MAIS DE 100%
greater.1 <- NULL
for (i in 1985:2017) {
greater.1 <- append(greater.1, which(df.nat.forest.AREA[,paste0(i)] >= 100))
}
greater.1 <- greater.1 %>% unlist() %>% unique()
## \__ Cria os limites de percentuais ----
df.nat.forest.AREA <-
df.nat.forest.AREA %>%
dplyr::mutate_at(.vars=vars(paste0(seq(1985,2017))),
.fun = function(x) cut(x,breaks=c(-Inf, 0, 1, seq(0, 100, 10)[-1]))
)
Seção 4: O Mapa
# ************************************************************************* ----
# Mapa ----
# ************************************************************************* ----
## \__ Funcao que junta os dados de vegetacao com os de posicao das cidades ----
create.map.df <- function(DATA){
map <- join_data(map = brazil.map, data = DATA, by = c("City" = "CODIBGE"))
map.ggplot <- as(map, "Spatial")
map.ggplot@data$id = row.names(map.ggplot@data)
df.map.ggplot <- suppressMessages(ggplot2::fortify(map.ggplot))
df.map.ggplot <- dplyr::left_join(df.map.ggplot, map.ggplot@data, by = "id")
return(df.map.ggplot)
}
## \____ Base de dados Para criar o Mapa ----
df.map.ggplot <- create.map.df(df.nat.forest.AREA)
## \__ Map Colors ----
# USEI https://gka.github.io/palettes PARA OBTER UMA ESCALA DE CORES
# COM 12 INTERVALOS (MESMO NUMERO DE INTERVALOS DE
# cut(x,breaks=c(-Inf, 0, 1, seq(0, 100, 10)[-1])) )
# QUE VAI DE white ATE darkgreen QUE EH APROPRIADA PARA PESSOAS COM
# VARIOS TIPOS DE DALTONISMO
map.colors <- c('#ffffff','#eaf1e7','#d4e2cf','#bfd4b8','#abc5a1','#96b78a',
'#81a974','#6c9b5f','#588d4a','#428034','#28711d','#006400')
### \__ Rotulos para a Legenda ----
map.labels <- c("0", "", "10%", rep("",3), "50%", rep("",3), "90%", "")
## \__ Criar Tema para o Mapa ----
map.theme <-
theme(
text = element_text(size = 12, family = "Lato", color = "#f2efef"),
axis.text = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
plot.background = element_rect(fill = "#3a3935", color = NA),
panel.background = element_rect(fill = "#3a3935", color = NA),
panel.grid = element_blank(),
plot.margin = margin(t = 0.25,r = 0.25,b = 0.25,l = 0.25,unit = "cm"),
plot.title = element_text(margin = margin(t = 0.25, unit = "cm")),
plot.caption = element_text(size = 10),
legend.background = element_rect(fill = "#3a3935", color = "#3a3935"),
legend.text = element_text(size = 11),
legend.title = element_text(size = 12),
legend.position = "left",
legend.key = element_rect(color = "transparent"),
legend.key.size = unit(0, units = "mm"),
legend.key.width = unit(2, units = "mm")
)
## \__ Plotar Mapa ----
### \____ Funcao para o GIF ----
plot.map.forest <- function(){
for (i in 1985:2017) {
lost.area.year.km2 <-
round(
sum(-(df.nat.forest %>% summarize_at(paste0(i),sum)/100),
df.nat.forest %>% summarize_at(paste0(1985),sum)/100)
,2
)
lost.area.year.mi2 <- round(lost.area.year.km2/2.59, 2)
map.forest <-
# MAPA
ggplot(mapping = aes(x = long, y = lat, group = group)) +
# MAPA COM LINHAS DA CIDADE TRANSPARENTES
geom_polygon(
data = df.map.ggplot,
aes_string(fill = paste0("X",i)),
color = "transparent"
) +
# LEGENDA
scale_fill_manual(
values = map.colors,
drop = FALSE,
na.value = "black",
na.translate = FALSE,
name = "",
label = map.labels,
guide = guide_legend(
direction = "vertical",
ncol = 1,
reverse = TRUE,
label.vjust = 2.5
)
) +
# LIMITES DO MAPA
xlim(range(df.map.ggplot$long)) +
ylim(range(df.map.ggplot$lat)) +
coord_map() +
# ADICIONAR TITULO E CAPTIONS
labs(
title = "Percentual da Área da Cidade coberta por Vegetação Nativa, Brasil.",
caption = paste0("Mapa: Gustavo Varela-Alvarenga - ogustavo.com/pt \n \n ",
"Dados: MapBiomas v.3.1 - mapbiomas.org & ",
"IBGE - www.ibge.gov.br ")
) +
# ADICIONAR TEMA
map.theme +
# ADICIONAR ANNOTATIONS
## PALAVRA ANO
annotate(
geom = "text", x = -44.1, y = 2.4, size = 4, family = "Lato",
color = "#f2efef", label = "ano", angle = 90
) +
## VALOR DO ANO
annotate(
geom = "text", x = -40, y = 2.5, size = 10, family = "Lato",
color = "#f2efef", label = paste0(i)
) +
## TITULO - TABELA A ESQUERDA
annotate(
geom = "text", x = -68.5, y = -26, size = 4, family = "Lato",
color = "#f2efef", hjust = 0.5,
label =
paste0("Área Total de Vegetação \n", "Nativa Perdida Desde 1985")
) +
## LINHA ABAIXO TITULO
annotate(
geom = "segment", x = -64.5, xend = -72.5, y = -28, yend = -28,
linetype = "solid", color = "#f2efef"
) +
## VALORES EM KMs QUADRADOS COM . PARA MARCAR OS MILHARES
annotate(
geom = "text", x = -68.5, y = -29, size = 4.5, family = "Lato",
color = "#f2efef", hjust = 0.5,
label =
ifelse(lost.area.year.km2 == 0, "0",
paste0(prettyNum(round(lost.area.year.km2, 0),
big.mark=".", decimal.mark = ","),
" km2")
)
) +
## VALORES EM MILHAS QUADRADAS COM , PARA MARCAR OS MILHARES
annotate(
geom = "text", x = -68.5, y = -31, size = 4.5, family = "Lato",
color = "#f2efef", hjust = 0.5,
label =
ifelse(lost.area.year.mi2 == 0, "0",
paste0(prettyNum(round(lost.area.year.mi2, 0), big.mark=","),
" milhas2")
)
) +
## LINHA ABAIXO TABELA
annotate(
geom = "segment", x = -64.5, xend = -72.5, y = -32, yend = -32,
linetype = "solid", color = "#f2efef"
) +
## TITULO - TABELA A DIREITA
annotate(
geom = "text", x = -39, y = -26, size = 4, family = "Lato",
color = "#f2efef", hjust = 0.5,
label = paste0("Área Perdida em \n", "% da Área do(a)")
) +
## LINHA ABAIXO TITULO
annotate(
geom = "segment", x = -35, xend = -43, y = -28, yend = -28,
linetype = "solid", color = "#f2efef"
) +
## ADICIONAR PALAVRA TEXAS
annotate(
geom = "text", x = -41, y = -29, size = 4.5, family = "Lato",
color = "#f2efef", label = "Texas : "
) +
## ADICIONAR PERCENTUAL PARA TEXAS
annotate(
geom = "text", x = -34, y = -29, size = 4.5, family = "Lato",
color = "#f2efef", hjust = 1,
label =
ifelse(
lost.area.year.km2 == 0, "0%",
paste0(round(lost.area.year.km2/696241, 2) * 100, "%")
)
) +
## ADICIONAR PALAVRA ALEMANHA
annotate(
geom = "text", x = -41, y = -31, size = 4.5, family = "Lato",
color = "#f2efef", label = "Alemanha : "
) +
## ADICIONAR PERCENTUAL PARA ALEMANHA
annotate(
geom = "text", x = -34, y = -31, size = 4.5, family = "Lato",
color = "#f2efef", hjust = 1,
label =
ifelse(
lost.area.year.km2 == 0, "0%",
paste0(round(lost.area.year.km2/357386, 2) * 100, "%")
)
) +
## LINHA ABAIXO DA TABELA
annotate(
geom = "segment", x = -35, xend = -43, y = -32, yend = -32,
linetype = "solid", color = "#f2efef"
)
print(map.forest)
}
}
### \____ Criar GIF ----
### O comando abaixo salvara um arquivo gif de nome forest_animation_pt-br.gif
### no seu working directory
gif_file <-
gifski::save_gif(
expr = plot.map.forest(),
gif_file = "forest_animation_pt-br.gif",
delay = 0.75, width = 738, height = 788, res = 100
)
utils::browseURL(gif_file)
Se você quiser, você pode usar o código abaixo para obter imagens estáticas de alta qualidade dos mapas. Esteja ciente de que isso salvará os gráficos em seu diretório de trabalho.
### \____ Criar Mapas Estaticos ----
#install.packages("Cairo")
library("Cairo")
for (i in 1985:2017) {
lost.area.year.km2 <-
round(
sum(-(df.nat.forest %>% summarize_at(paste0(i),sum)/100),
df.nat.forest %>% summarize_at(paste0(1985),sum)/100)
,2
)
lost.area.year.mi2 <- round(lost.area.year.km2/2.59, 2)
map.forest <-
# MAPA
ggplot(mapping = aes(x = long, y = lat, group = group)) +
# MAPA COM LINHAS DA CIDADE TRANSPARENTES
geom_polygon(
data = df.map.ggplot,
aes_string(fill = paste0("X",i)),
color = "transparent"
) +
# LEGENDA
scale_fill_manual(
values = map.colors,
drop = FALSE,
na.value = "black",
na.translate = FALSE,
name = "",
label = map.labels,
guide = guide_legend(
direction = "vertical",
ncol = 1,
reverse = TRUE,
label.vjust = 2.5
)
) +
# LIMITES DO MAPA
xlim(range(df.map.ggplot$long)) +
ylim(range(df.map.ggplot$lat)) +
coord_map() +
# ADICIONAR TITULO E CAPTIONS
labs(
title = "Percentual da Área da Cidade coberta por Vegetação Nativa, Brasil.",
caption = paste0("Mapa: Gustavo Varela-Alvarenga - ogustavo.com/pt \n \n ",
"Dados: MapBiomas v.3.1 - mapbiomas.org & ",
"IBGE - www.ibge.gov.br ")
) +
# ADICIONAR TEMA
map.theme +
# ADICIONAR ANNOTATIONS
## PALAVRA ANO
annotate(
geom = "text", x = -44.1, y = 2.4, size = 4, family = "Lato",
color = "#f2efef", label = "ano", angle = 90
) +
## VALOR DO ANO
annotate(
geom = "text", x = -40, y = 2.5, size = 10, family = "Lato",
color = "#f2efef", label = paste0(i)
) +
## TITULO - TABELA A ESQUERDA
annotate(
geom = "text", x = -68.5, y = -26, size = 4, family = "Lato",
color = "#f2efef", hjust = 0.5,
label =
paste0("Área Total de Vegetação \n", "Nativa Perdida Desde 1985")
) +
## LINHA ABAIXO TITULO
annotate(
geom = "segment", x = -64.5, xend = -72.5, y = -28, yend = -28,
linetype = "solid", color = "#f2efef"
) +
## VALORES EM KMs QUADRADOS COM . PARA MARCAR OS MILHARES
annotate(
geom = "text", x = -68.5, y = -29, size = 4.5, family = "Lato",
color = "#f2efef", hjust = 0.5,
label =
ifelse(lost.area.year.km2 == 0, "0",
paste0(prettyNum(round(lost.area.year.km2, 0),
big.mark=".", decimal.mark = ","),
" km2")
)
) +
## VALORES EM MILHAS QUADRADAS COM , PARA MARCAR OS MILHARES
annotate(
geom = "text", x = -68.5, y = -31, size = 4.5, family = "Lato",
color = "#f2efef", hjust = 0.5,
label =
ifelse(lost.area.year.mi2 == 0, "0",
paste0(prettyNum(round(lost.area.year.mi2, 0), big.mark=","),
" milhas2")
)
) +
## LINHA ABAIXO TABELA
annotate(
geom = "segment", x = -64.5, xend = -72.5, y = -32, yend = -32,
linetype = "solid", color = "#f2efef"
) +
## TITULO - TABELA A DIREITA
annotate(
geom = "text", x = -39, y = -26, size = 4, family = "Lato",
color = "#f2efef", hjust = 0.5,
label = paste0("Área Perdida em \n", "% da Área do(a)")
) +
## LINHA ABAIXO TITULO
annotate(
geom = "segment", x = -35, xend = -43, y = -28, yend = -28,
linetype = "solid", color = "#f2efef"
) +
## ADICIONAR PALAVRA TEXAS
annotate(
geom = "text", x = -41, y = -29, size = 4.5, family = "Lato",
color = "#f2efef", label = "Texas : "
) +
## ADICIONAR PERCENTUAL PARA TEXAS
annotate(
geom = "text", x = -34, y = -29, size = 4.5, family = "Lato",
color = "#f2efef", hjust = 1,
label =
ifelse(
lost.area.year.km2 == 0, "0%",
paste0(round(lost.area.year.km2/696241, 2) * 100, "%")
)
) +
## ADICIONAR PALAVRA ALEMANHA
annotate(
geom = "text", x = -41, y = -31, size = 4.5, family = "Lato",
color = "#f2efef", label = "Alemanha : "
) +
## ADICIONAR PERCENTUAL PARA ALEMANHA
annotate(
geom = "text", x = -34, y = -31, size = 4.5, family = "Lato",
color = "#f2efef", hjust = 1,
label =
ifelse(
lost.area.year.km2 == 0, "0%",
paste0(round(lost.area.year.km2/357386, 2) * 100, "%")
)
) +
## LINHA ABAIXO DA TABELA
annotate(
geom = "segment", x = -35, xend = -43, y = -32, yend = -32,
linetype = "solid", color = "#f2efef"
)
### Observe que o comando abaixo ira salvar varios arquivos
### png com nome Y*_pt-br.png, onde * eh o ano, no seu
### working directory
ggsave(
paste0("Y",i,"_pt-br.png"), device="png", type="cairo",
width = 7, height = 7
)
}
