Mapa de Vegetação de Areia Branca na Amazônia Brasileira

Author

dados do IBGE, código R com apoio de GPT4

Objetivos de aprendizagem

Usar pacotes R (sf, ggplot2, ggspatial, cowplot) para produzir um mapa de vegetação.

Manipular tabela de atributos
Manipular simbologia e formatação, visando uma figura para artigo ou tese
Preparar e inserir elementos cartográficos: mapa de localização, barra de escala, seta de norte
Salvar em disco o mapa final em fornato png

Criar fluxo de trabalho reproduzível e didático com Quarto MarkDown

Acessar pasta de projeto compartilhada
Formatação do corpo de texto e dos blocos de código do script qmd
Conceitos de pre-visualização e renderização do script qmd
Entender alguns comandos do YAML (parte inicial do script qmd)
Inserir figuras didáticas no script qmd (criadas no Paint, Snagit, Photoshop etc.)
Solicitar apoio de GPT4 para adaptar o script qmd e para comentar códigos recebidos
Como postar script renderizado no site Rpubs

Passos resumidos

Baixar e ativar a pasta deste projeto, disponibilizada no GDrive
Localizar e acessar shapefile de fonte oficial (IBGE)
Transformar o shapefile em Simple Features Collection
Filtrar e manipular a tabela de atributos
Usar ggplot2 para a construção incremental do mapa – cerne do exercício!
Criar e incluir elementos de mapa
Salvar o mapa completo em formato png e inserir como figura no QMD
Renderizar o documento Markdown em formato html
Publicar o documento Markdown renderizado no site rpubs.com
Veja mais informações sobre ggplot2 e mapas vetoriais

Preparativos

Instalar R e RStudio mais recentes

Obter e instalar R Uma vez instalado, faça o update automático de seu RStudio (Menu de Ajuda - Procurar Atualizações), ou baixar RStudio pela primeira vez no site da posit

Guia de Quarto Markdown

Consultar aqui

Ajustar aparência de RStudio

View - Panes - ativar “View All Panes”, ativar “Console on the Left”.

Tools - Global Options - Appearance - IPlastic. Os blocos de código ficarão mais visíveis.

Organizar seus dados em uma pasta de projeto de RStudio

Uma pasta definida como projeto pelo RStudio pode ficar em qualquer local no seu disco, no disco de colaboradores ou de alunos. Um projeto de RStudio torna seu código mais portátil. Estimula o usuário a manter juntos todos os dados de entrada e de saída de um script. Aqui usaremos um pasta já definida como projeto. (Para começar um novo projeto, os passos seriam File - New Project - New Directory - New Project - Dê um nome e definir o local)

Acessar a pasta de projeto RStudio compartilhada, contendo este script .Qmd

A pasta de projeto com nome “veg_areia_amazonia” foi salva em nuvem. O acesso é aberto pelo link abaixo.

https://drive.google.com/drive/folders/1npFt2UfGxGuCHZ4QKDX9IIvu2_M4yvL2?usp=sharing

Colar o link em seu navegador de internet e baixar a pasta inteira (não apenas os conteudos, mas sim a pasta, que tem o nome veg_areia_amazonia). Siga as respectivas figuras abaixo para Chrome e para Firefox. (Estas figuras didáticas vão aparecer no seu script qmd apenas depois baixar de baixar a pasta e designar como projeto. Elas estão salvas n a subpasta ./images, da pasta de projeto. O comando para inserir estas figuras pode ser visualizado clicando na aba “Source’, para ver a versão texto do script QMD.)

Chrome:

Firefox:

Uma pasta zip cujo nome raiz termina com um número comprido estará nos downloads. A pasta do projeto estará dentro desta; Veja figura abaixo

Depois de colocar a pasta do projeto no lugar desejado, basta reconhecer como projeto ativo. Usando RStudio, vá em File - Open file - navegue até a pasta veg_areia_amazonia e clicar no arquivo veg_areia_amazonia.Rproj.

Agora abra o arquivo .qmd usando RStudio (File - open file - veg_areia_amazonia.qmd). O script deve aparecer no painel Source. Há duas opções de visualizar o script, usando os botões “Source” ou “Visual” , respectivamente sem ou com a renderização prévia (figura abaixo). O arquivo qmd em si é um arquivo texto contendo comentários, script de código organizado em blocos, e as regras de formatação quando renderizado para os formatos html, pdf, ou word.doc. O arquivo qmd inicia com um YAML (linhas 1 a 11 na figura), que não deve ser alterado.

Blocos de Código

Para este exercício, os blocos de código já foram criados (com apoio de GPT4). Escolhendo a opção de pre-renderizado do script, pode se ver os blocos de código com coloração mais escura. Pode rodar o código do primeiro bloco abaixo, clicando na seta verde “Run Current Chunk”.

No primeiro bloco, definimos os pacote do R necessários para o restante deste projeto; R verifica se já foram instalados (baixados); instala apenas aqueles que faltam, com o comando “install.packages”; e finaliza carregando os pacotes com o comando “library”. Este bloco lhe servirá de modelo para constriur o primeiro bloco de todos seus códigos futuros. Não é necessário aprender ou lembrar o síntax exato de seus blocos de código. Para isso, existem chat.openai.com e stackexchange.com !

# Defina os nomes dos pacotes
package_names <- c("sf", "tidyverse", "ggspatial", "mapdata", "cowplot")

# Verifica se cada pacote está instalado e instala se necessário
for (pkg in package_names) {
  if (!pkg %in% installed.packages()[,"Package"]) {
    install.packages(pkg)}
  
  # Carregue os pacotes
  library(pkg, character.only = TRUE)}

Baixar shapefiles do IBGE (vegetação do Brasil) e de Natural Earth (paises, para o mapa de localização)

Cada shapefile é, na verdade, uma coletânea de pelo menos quatro arquivos com o mesmo nome raiz e com diferentes extensões depois do ponto no nome. Pelo site do IBGE o conjunto de arquivos que compoem um shape é oferecido já consolidado em único arquivo zip.

url <- "https://geoftp.ibge.gov.br/informacoes_ambientais/vegetacao/vetores/brasil_5000_mil/Vegetacao_5000mil.zip" #mapa de vegetação do Brasil, escala 1:5 milhões

# repita para o shape de paises do mundo, para o mapa de localização
url_2 <- "https://www.naturalearthdata.com/http//www.naturalearthdata.com/download/110m/cultural/ne_110m_admin_0_countries.zip"

# indique novo nome de arquivo zipado e a subpasta dentro da pasta do projeto onde será guardado
zipped_shape_veg <- "./data/veget_br_5milhoes.zip"
zipped_shape_paises <- "./data/ne_50m_admin_0_countries.zip"

# baixar os shapes compactados e salvar 
download.file(url, zipped_shape_veg)
download.file(url_2, zipped_shape_paises)

# Descompactar no mesmo diretório. Qual é o signficado de "./data/" ?
unzip(zipped_shape_veg, exdir = "./data/")
unzip(zipped_shape_paises, exdir = "./data/")

Gerenciar aquivos shapefile baixados

Usando windows field explorer, examine os shapefiles extraidos do zip. Um dos componentes do shapefile de vegetação tem o sufixo .prj.txt e precisa ser apenas .prj Corriga e excluia os arquivos zip baixados.

file.rename("./data/Vegetacao_5000.prj.txt", "./data/Vegetacao_5000.prj")

[1] TRUE

file.remove("./data/veget_br_5milhoes.zip")

[1] TRUE

file.remove("./data/ne_50m_admin_0_countries.zip")

[1] TRUE

Acessar o arquivo com extensão proj de cada shapefile. Abre cada um como arquivo texto, usando bloco de notas, para saber qual o Sistema de Coordenadas (CRS) de cada um. Para o mapa de vegetação, temos “GCS_SIRGAS_2000”. Para o mapa de paises fornecido pelo Natural Earth, temos “GCS_WGS_1984”. A abreviação GCS significa que ambos usam um Geographic Coordinate System, ou seja, graus de longitude e latitude representam as posiçãoes X e Y em um sistema cartesiano de coordenadas.

Uso de latitude e longitide como posição em uma tela ou um papel representa bem as distâncias e as formas perto do equador. Longe do equador, as distâncias no sentido Norte-Sul são corretas, mas as distâncias no sentido Leste-Oeste são exageradas, alargando e ampliando a area dos polígonos de países. “SIRGAS_2000” e “WGS_1984” indicam o respectivo datum – o modelo usado para a forma do esferóide da Terra. Estes dois “data” são quase idênticos, com diferença menor que um metro na posição.

Feche o Bloco de Notas, sem salvar ou alterar os arquivos prj dos dois shapefiles.

Converta o dois shapefiles em dois sf (simple features collections)

Especifique o caminho até o shapefile; leia e converta para formato sf. Este é o novo formato padrão no R para dados vetoriais (coleções de pontos, de linhas ou de polígonos) com respectivas tabelas contendo os atributos dos pontos, das linhas ou dos polígonos.

# converte o shapefile de tipos de vegetação em um objeto do tipo "simple features collection"
shape_veg_path <- "./data/Vegetacao_5000.shp" # especifica o caminho
sf_veget_br_5milhoes <- st_read(shape_veg_path)

Reading layer `Vegetacao_5000' from data source 
  `C:\R_projects\veg_areia_amazonia\data\Vegetacao_5000.shp' 
  using driver `ESRI Shapefile'
Simple feature collection with 2170 features and 9 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: -73.9806 ymin: -35.8075 xmax: -25.7463 ymax: 7.0518
Geodetic CRS:  SIRGAS 2000

# Repita para o shapefile de paises do mundo
shape_paises_path <- "./data/ne_110m_admin_0_countries.shp" 
sf_paises <- st_read(shape_paises_path)

Reading layer `ne_110m_admin_0_countries' from data source 
  `C:\R_projects\veg_areia_amazonia\data\ne_110m_admin_0_countries.shp' 
  using driver `ESRI Shapefile'
Simple feature collection with 177 features and 168 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -180 ymin: -90 xmax: 180 ymax: 83.64513
Geodetic CRS:  WGS 84

# Na aba Environment, clicar duas vezes no objeto sf_paises e examine a coluna CONTINENT.  Extrair apenas os polígonos de paises da América do Sul
sf_south_america <- sf_paises[sf_paises$CONTINENT == 'South America', ]

O próximo bloco de código cria o mapa de localização

# Ajuste a coordenada do limite sul perto do extremo sul do Brasil
bbox <- st_bbox(sf_south_america)
bbox["ymin"] <- -34.7

# Converta o bounding box em um objeto sf
bbox_sf <- st_as_sfc(bbox, crs = st_crs(sf_south_america))

# Recorte sf_south_america usando o bounding box
sf_south_america_clipped <- st_crop(sf_south_america, bbox_sf)

# Imprimir o bounding box do resultado
print(st_bbox(sf_south_america_clipped))

     xmin      ymin      xmax      ymax 
-81.41094 -37.01564 -34.72999  12.43730

# Depois de clipado, o limite sul do mapa inset está aprox três graus ao sul do retângulo usado para clipar o mapa. Não sabemos porque.

Criar nomes e cores para os tipos de cobertura de interesse

Prosseguimos com a manipulação da tabela de atributos do mapa principal.

Cada linha na tabela de atributos do objeto ‘sf_veget_br_5milhoes’ representa um polígono, que é uma mancha de determinado tipo de vegetação. Na aba Environment, clicar duas vezes no objeto ‘sf_veg_br_5milhoes’. Sua tabela de atributos deve aparecer como nova aba, no painel Source. Na tabela de atributos, ordene alfabeticamente a coluna COD_LEG_CA, clicando no cabeçalho. Identifique as cinco categorias que começam com a letra “L”. Estas representam os polígonos dos cinco tipos de vegetação de areia branca. Identifique também todas as categorias para água aberta.

No sistema de classificação da vegetação brasileira adotado pelo IBGE, toda a vegetação em solo podzólico arenoso é chamada de “Campinarana”. Há cinco sub-classes em ordem crescente de sua biomassa lenhosa. São elas: Campinarana gramíneo-lenhosa (código Lg), Campinarana arbustiva (código Lb), Campinarana arbórea (código La), Campinarana florestada (código Ld) e Ecótono (código LO). O ecótono é uma transição entre a floresta de areia branca e a floresta de solo argiloso. O ecótono ou é gradual ou tem encraves pequenos de um dos dois tipos em uma matriz do outro tipo.

# Na tabela de atributos do objeto tipo sf de vegetação, crie uma nova coluna chamada tipo_veg_areia contendo nomes para os cinco tipos de vegetação de solo arenoso; o nome "água" para todos os tipos de água aberta; e o nome "Outro" no restante das linhas da tabela.
sf_veget_br_5milhoes$tipo_veg_areia <- case_when(
  sf_veget_br_5milhoes$COD_LEG_CA %in% c("Lg") ~ "herbácea",
  sf_veget_br_5milhoes$COD_LEG_CA %in% c("Lb") ~ "arbustiva",
  sf_veget_br_5milhoes$COD_LEG_CA %in% c("La") ~ "árvores esparsas",
  sf_veget_br_5milhoes$COD_LEG_CA %in% c("Ld") ~ "árvores densas",
  sf_veget_br_5milhoes$COD_LEG_CA %in% c("LO") ~ "ecótono",
  sf_veget_br_5milhoes$COD_LEG_CA %in% c("00Magua", "01Magua", "02Magua", "03Magua") ~ "água",
  TRUE ~ "Outro"
)

# Forçar a nova coluna a ser do tipo "factor", tendo "níveis" na mesma ordem desejada para a legenda do mapa (ordem crescente de biomassa lenhosa)
sf_veget_br_5milhoes$tipo_veg_areia <- factor(sf_veget_br_5milhoes$tipo_veg_areia,    levels = c("herbácea", "arbustiva", "árvores esparsas", "árvores densas", "ecótono", "água", "Outro"))

# Curioso para saber porque usamos factor() e não as.factor()?  Consulte GPT4 para saber a diferença.

Criar o mapa de localização de maneira incremental usando ggplot

Primeiro, criar mapa inset. O primeiro geom_sf() apenas desenha a borda para o Brasil (e preenche com NA para tornar o preenchimento transparente), enquanto o segundo geom_sf() preenche todos os países (incluindo o Brasil) com base no nome do país. A função ifelse() é usada para criar um novo fator que só tem dois níveis: “Brazil” e “Other”.

# Faça um mapa de inset
x_min <- -81.41094
x_max <- -34.72999
y_min <- -37.01564
y_max <- 12.43730

south_america <- ggplot() +
  # Adicione um retângulo branco como fundo
  annotate("rect", xmin = x_min, xmax = x_max, ymin = y_min, ymax = y_max, fill = "white") +
  geom_sf(data = sf_south_america_clipped[sf_south_america_clipped$SOVEREIGNT == "Brazil", ], fill = NA, color = "black") +
  geom_sf(data = sf_south_america_clipped, aes(fill = ifelse(SOVEREIGNT == "Brazil", "Brazil", "Other"))) +
  # Adicione um retângulo preto como borda
  annotate("rect", xmin = x_min, xmax = x_max, ymin = y_min, ymax = y_max, fill = NA, color = "black", lwd = 0.5) +
  scale_fill_manual(values = c("Brazil" = "lightgrey", "Other" = "white")) +
  theme_void() +  # Remove os eixos e texto
  theme(legend.position = "none")  # Remove a legenda

# Converte o mapa do inset para um grob (grid object)
south_america_grob <- ggplotGrob(south_america)

Criar o mapa principal de maneira incremental usando ggplot

Ajustar os limites de longitude xlim e latitude ylim para alterar o zoom para focar na Amazônia brasileira.

# Faça o mapa principal com passos incrementais
mapa_principal <- ggplot() +
  geom_sf(data = sf_veget_br_5milhoes, aes(fill = tipo_veg_areia), color = NA) +
  ggtitle("Mapa de Vegetação do Brasil, escala 1:5,000,000") +
  labs(x = "Longitude", y = "Latitude") + 
  theme_minimal() +
  theme(panel.grid.major = element_line(colour = "black", linetype = "dashed"),
        panel.grid.major.x = element_line(color = "black", linewidth = 0.5),
        panel.grid.major.y = element_line(color = "black", linewidth = 0.5),
        axis.text.x = element_text(angle = 0, vjust = 0.5, color = "black"),
        axis.text.y = element_text(angle = 0, hjust = 1, color = "black"),
        plot.title = element_text(hjust = 0.5), # Centraliza o título
        legend.position = c(0.9, 0.18), # Posiciona legenda dentro da figura: 
        # ((0,0) é o canto inferior esquerdo; (1,1) é o canto superior direito)
        legend.background = element_rect(fill = "white")) + # fundo branco
  scale_x_continuous(breaks = seq(-75, -45, 5)) + # Extremos e intervalos de longitude
  scale_y_continuous(breaks = seq(-15, 10, 5)) +  # Extremos e intervalos de latitude
  scale_fill_manual(values = c("herbácea" = "#ccffcc",
                             "arbustiva" = "#99ff99",
                             "árvores esparsas" = "#33cc33",
                             "árvores densas" = "#009900",
                             "ecótono" = "#004d00",
                             "água" = "blue",
                             "Outro" = "lightgrey"),
        name = "Tipos de Vegetação\nem Solo Arenoso", # "\n" quebra o título
        guide = guide_legend(override.aes = list(colour = "darkgrey"))) +
  coord_sf(xlim = c(-75, -45), ylim = c(-15, 10), expand = FALSE) + # enquadramento da legenda
  annotate("text", x = -55, y = -14, label = "Fonte: https://geoftp.ibge.gov.br", size = 3) + # cria e posiciona o texto informando fonte é IBGE
  annotation_scale(location = "bl") + # barra de escala no bottom left 
  annotation_north_arrow(location = "tr", which_north = "true") # seta de norte no top right

Adicionar o mapa de localização como um grob e plotar os dois mapas com ggdraw, para obter o objeto ‘mapa_final’

# Adicione o mapa da América do Sul como um inset. 
mapa_final <- ggdraw() + 
  draw_plot(mapa_principal) + 
  draw_grob(south_america_grob, width = 0.3, height = 0.3, x = 0.063, y = 0.665)
# Experimente diferentes valores x,y para posicionar o mapa de localização onde desejar. Valores de x,y vão de (0,0) que é o canto inferior esquerdo no frame do mapa principal, até (1,1) que é o canto superior direito. Width e height indicam o grau de encolhimento do mapa de localização

Salvar o mapa final em disco como um arquivo de imagem. Ajustar os valores de largura e altura em pólegadas, para obter a proporção desejada e o nível de detalhe. Ajuste o valor de dpi para a resolução desejada.

ggsave(filename = "mapa_veg_areia_amazonia.png", plot = mapa_final, width = 9, height = 7, dpi = 300, bg = "white")

Mapa final está pronto e salvo!

Clicar duas vezes sobre o nome do arquivo salvo pelo ggsave. Este é seu mapa final, que pode ser inserido em documento Word, powerpoint etc

Inserir o mapa final neste QMD

Uma vez criada a imagem png do mapa, ativar a aba “Visual” no painel “Source”. Clicar no ícone de Figure/image (na barra de ferramentas na parte superior do painel Source). Navegue até o arquivo png do mapa, e insere. Observe que este comando será no texto do QMD e não dentro de um bloco de código

No entanto, o comando já existe logo abaixo para inserir o mapa – portanto não precisa fazer.

Renderizar

Clicar na seta Render, na barra de ferramentas do painel Source, com a aba “Visual” selecionada. Será criado um arquivo html e uma pasta com o mesmo nome deste arquivo, ambos dentro da pasta do projeto. Mantenham esta pasta e o arquivo html juntos. Clicar duas vezes sobre o arquivo html para ver o script renderizado no seu navagador. Pode também visualizar o script renderizado na aba Viewer do painel de Files, Plots, Packages, Help etc.

Publicar seu documento Quarto Markdown no rpubs

Na aba view, do painel onde aparece seu documento renderizado, clicar em “publish”. O documento será salvo gratuitamente na sua conta de rpubs.com É necssário primeiro criar uma conta com senha neste site.

FIM DA AULA

Exercício Prático de mapa feito com R

Usando o código fornecido para vegetação de areia branca, produzir um mapa similar, que mostra no país inteiro os tipos de vegetação de savana
Pode eliminar o mapa de localização, se quiser
Pode deixar a legenda na sua posição padrão, se quiser, que é fora do mapa, no lado direito
Criar uma nova pasta de projeto, para guardar este código novo: “veg_savana_brasil”
Cria as subpastas: “data” , “docs”, “images”;
copiar quatro figuras didática para a subpasta “images” do novo projeto (estão na subpasta “images” do projeto veg_areia_amazonia);
Copiara o arquivo qmd da pasta de projeto “veg_areia_amazonia”, para a pasta “veg_savana_brasil”. Renomeia o arquivo qmd para “veg_savana_brasil.qmd”.
Faça as alterações necessárias no código para obter um mapa de todos os tipos de savana, enquadrando o pais inteiro. Isto inclui os seguintes passos
- Rodar os blocos de código até o ponto de ter o shapefile do IBGE “veget_br_5milhoes.shp” e todos seus arquivos ancilares, descompactados na pasta “data’.
- Abrir o shapefile no QGis, visualize sua tabela de atributos, e ordene alfabeticamente pelo campo “COD_LEG_CA” (basta clicar no cabeçalho) e localize os tipos de vegetação listados na tabela abaixo, onde estão re-ordenados para representar a progressão de biomassa lenhosa.

COD_LEG_CA	novo nome em novo campo criado	cor
Sg	“herbácea”
Sp	“savana parque”
Sa	“árvores esparsas”
Sd	“árvores densas”
ST	“ecótono savana - savana estépica”
STN	“ecótono savana - sav estép - floresta estacional”
SN	“ecótono savana - floresta estacional”
SO	“ecótono savana - floresta ombrófila”

Usando um site com color picker, identifique uma progressão de vermelho claro para vermelho saturado para as quatro classes que não são ecótonos. Faça a mesma coisa para as quatro classes de ecótionos, passando de laranja clara para laranja saturada (ou escolhe uma outra cor base para representar a porgressão de quatro ecótonos). Anote os oito códigos hexadecimais correspondentes de cada cor e colar estes oito códigos na última coluna da tabela acima.

O enquadramento será o pais de Brasil inteiro. Não é obrigatório produzir um mapa de localização do Brasil. Se não incluir, pode excluir todos as linhas de código que tratam do mapa_inset de localização.