O script do documento Quarto MarkDown começa com o YAML, que não aparece no documento renderizado. Não apague, nem altere o YAML. Ele contem instruções sobre a renderização do script para o formato html, entre outras funções.
Prefira uso de projeto do R, um avanço sobre o uso de setwd()
Usar um projeto R torna seu código portátil, pois usa caminhos relativos e não absolutos. A pasta do projeto pode ficar em qualquer local no seu disco, no disco de colaboradores, de usuários ou alunos. Para compartilhar uma pasta de projeto, compactar a pasta com zip e envie a pasta por email ou guarde a pasta no GDrive e compartilhe um link.
Se você recebeu ou baixou uma pasta contendo o arquivo “mapas_precip_brasil.Rproj”, pular os três passos numerados abaixo.
Passos somente se estiver criando um novo projeto:
No RStudio: File - New Project - nova pasta ou selecione uma pasta existente vazia. Será criado nesta pasta um arquivo “nome_da_pasta.Rproj” e uma pasta escondida com nome .Rproj.user. Não apague o arquivo .Rproj nem a pasta escondida.
O caminho entre o diretório raiz até a pasta do projeto não deve conter espaços nem caracteres especiais.
Crie ou coloque seu script qmd (Quarto markdown) que contem o código deste exercício na pasta do projeto. Criar duas subpastas, com os nomes data e docs.
Se você recebeu ou baixou uma pasta contendo o arquivo “mapas_precip_brasil.Rproj”, pular os três passos numerados abaixo.
Passos somente se estiver criando um novo projeto:
No RStudio: File - New Project - nova pasta ou selecione uma pasta existente vazia. Será criado nesta pasta um arquivo “nome_da_pasta.Rproj” e uma pasta escondida com nome .Rproj.user. Não apague o arquivo .Rproj nem a pasta escondida.
O caminho entre o diretório raiz até a pasta do projeto não deve conter espaços nem caracteres especiais.
Crie ou coloque seu script qmd (Quarto markdown) que contem o código deste exercício na pasta do projeto. Criar duas subpastas, com os nomes data e docs.
Como baixar a pasta do projeto e colocar em qualquer lugar no seu micro
Ao baixar a pasta do GDrive, certifique que baixou a pasta inteira e não o conteúdo da pasta.
Basta colocar a pasta baixada em qualquer lugar no seu disco. No entanto, é bom certificar que o caminho até a pasta de projeto não contem espaços, acentos ou caracteres especiais, e que o caminho não é muito comprido.
Ativar o projeto e abra o script com extensão qmd
Ativar o projeto (File - Open Project - navegue até a pasta do projeto (pasta ‘worldclim_map’) – clicar duas vezes no arquivo com extensão Rproj
Abrir o arquivo qmd contendo o script (File - Open File - navegue até a pasta do projeto – clicar duas vezes no arquivo com a extensão qmd).
O arquivo com extensão qmd é um arquivo texto contendo regras de formatação quando renderizado ou pre-renderizado. Ele inclui YAML no início, cabeçalhos de diferentes níveis, figuras inseridas, listas e o mais importante – blocos de código intercalados com texto (como este texto que você está lendo).
Pode-se visualizar uma prévia da renderização do qmd clicando no botão “Visual”, no painel Source de RStudio. Você pode editar o arquivo qmd usando ou a aba “Source”, ou a aba “Visual”. Lembre de salvar o arquivo qmd periodicamente, clicando no ícono de disquete.
Para inserir um novo bloco de código no script qmd, clicar no ícone verde +C, no cabeçalho do painel Source
Passo a Passo
1. Instale e carregue pacotes do R
Rodar este priméiro bloco de código, clicando na seta verde “Run Current Chunk”. Pode ignorar as advertências. Este bloco lista os pacotes a serem usados, avalia se já estão instalados (baixados), e instala aqueles pacotes não instalados. Finalmente, todos os pacotes da lista são carregados na memória de seu computador. Este bloco tem a mesma função de rodar manualmente install.packages() e library().
Um bloco de código inicia com uma linha contendo tres crazes, seguidos pela letra r entre chaves. Termina com outra linha contendo três crazes. As crazes podem ser vistas no seu código usando a aba “Source” do painel Source. Opcionalmente, depois da letra r e dentro das chaves, pode se dar um espaço e inserir um nome qualquer para o bloco de código. Não é permitido usar o mesmo nome em dois blocos de código diferentes.
# Este é um comentário, dentro do bloco de código# Defina os nomes dos pacotesnomes_pacotes <-c("rnaturalearth", "terra", "geodata", "tidyverse")# Verifica se cada pacote está instalado e instala se necessáriofor (pkg in nomes_pacotes) {if (!pkg %in%installed.packages()[,"Package"]) {install.packages(pkg)}# Carregue os pacoteslibrary(pkg, character.only =TRUE)}
Warning: package 'rnaturalearth' was built under R version 4.3.3
Warning: package 'terra' was built under R version 4.3.2
terra 1.7.65
Warning: package 'geodata' was built under R version 4.3.3
Warning: package 'ggplot2' was built under R version 4.3.3
Warning: package 'readr' was built under R version 4.3.2
Warning: package 'purrr' was built under R version 4.3.2
Warning: package 'dplyr' was built under R version 4.3.2
Warning: package 'stringr' was built under R version 4.3.2
Warning: package 'lubridate' was built under R version 4.3.2
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ tidyr::extract() masks terra::extract()
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag() masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
2. Usar pacote geodata para baixar do site WorldClim a precipitação mensal do Brasil
# confira que o objeto prec_meses é um arquivo TIFF da classe Spatial Raster, com 12 bandas. Cada banda contem a chuva mensal. A resolução espacial é 0.00833 graus, o que corresponde a 925 metros perto do equador. O arquivo raster é grande, com 4740 linhas por 5460 colunas em cada uma das 12 bandas.
O arquivo tem 137 mb. Tenha paciência se sua conexão é lenta. Acompanhe o progresso clicando no ícone de nova janela informativa que aparece na sua barra de tarefas. Usando Windows File Explorer, verifique se foi criado este caminho e esta imagem tif. O ponto inicial representa o caminho até a pasta de seu Projeto R: “./Output/wc2.1_country/BRA_wc2.1_30s_prec.tif”.
O próximo bloco de código mostra as 12 “bandas”, uma para cada mês. A função plot(), que é do Base R, já produz um mapa rudimentar para cada banda contendo título, moldura, coordenadas de latitude e longitude e uma rampa de cores adequada para dados contínuos.
plot(prec_meses) # serão mostradas 12 imagens de chuva mensal
Para poupar memória nos processamentos mais a frente, vamos reamostrar o spatial rastar prec_meses para um pixel 5 vezes maior no sentido X e no sentido y, o que reduz em 25x o tamanho de cada banda. O pixel novo terá como valor a média dos 25 pixels originais nele contidos. Pixels com valor NA serão desconsiderados.
A expressão “terra::aggregate” significa “rodar a função aggregate, do pacote terra”. Isto evita confusão se outro pacote instalado oferecer uma função com o mesmo nome.
rm(prec_meses) # Remove da memoria o Spatial Raster com resolução finagc() # "garbace collection", efetue a liberação de memória
used (Mb) gc trigger (Mb) max used (Mb)
Ncells 1762921 94.2 3428158 183.1 2377697 127.0
Vcells 2597080 19.9 8388608 64.0 8352911 63.8
3. Usar o argumento “fun=sum” dentro da função app, para obter a precipitação anual.
class : SpatRaster
dimensions : 948, 1092, 1 (nrow, ncol, nlyr)
resolution : 0.04166667, 0.04166667 (x, y)
extent : -74, -28.5, -34, 5.5 (xmin, xmax, ymin, ymax)
coord. ref. : lon/lat WGS 84 (EPSG:4326)
source(s) : memory
name : sum
min value : 0.00
max value : 6090.48
plot(prec_anual) # será mostrado um mapa preliminar de chuva anual
O argumento “fun=sum” de terra::app faz a soma das 12 “bandas” contidas no objeto prec_meses, que é um Spatial Raster. Será produzido um novo Spatial Raster (prec_anual) contendo uma única banda, que é o resultado da soma em cada posição de pixel. Veja esta explicação de álgebra de mapas:
Caso queira salvar o raster de prec_anual como imagem TIFF de uma única banda, rodar este comando: writeRaster(prec_anual, “prec_anual.tif”). Mas isto não é necessário para noso fluxo de trabalho.
4. Spatial Raster e mapa preliminar do mes mais chuvoso: 1=jan,2=fev,3=mar…. 12=dez
Faremos outro uso da função terra::app, para transformar o Spatial Raster prec_meses, contendo 12 bandas, em um outro Spatial Raster com uma única camada, que reporta para cada célula, o índice da banda com o maior valor (ou seja, o mês mais chuvoso daquela célula). Os índices das bandas vão de 1 a 12, em órdem cronológica no Spatial Raster. Observe que usanos o argumento fun=which.max para identificar o index do mês maix chuvos em cada posição no mapa.
class : SpatRaster
dimensions : 948, 1092, 1 (nrow, ncol, nlyr)
resolution : 0.04166667, 0.04166667 (x, y)
extent : -74, -28.5, -34, 5.5 (xmin, xmax, ymin, ymax)
coord. ref. : lon/lat WGS 84 (EPSG:4326)
source(s) : memory
name : which.max
min value : 1
max value : 12
plot(idx_mes_max) # 1=jan,2=fev,3=mar.... 12=dez
5. Obter polígono vetorial do tipo sf para Brasil
Usar Natural Earth para extrair países da America do Sul como “simple features” vetoriais, filtrar para Brasil e plotar o polígono do Brasil. Na tabela de atributos, Natural Earth fornece muitos atributos. Serão plotados os primeiros nove mapas do Brasil entre 168 possíveis, pois tem 168 colunas de atributos na tabela do polígono. Nenhum destes atributos nos interessam no momento, mas também não é necessrio remover eles, pois o objeto sf não usa muita memória.
Warning: plotting the first 9 out of 168 attributes; use max.plot = 168 to plot
all
# Exportando o shapefile do Brasilsf::st_write(brasil, "./Output/brasil_shapefile.shp", overwrite =TRUE)
Writing layer `brasil_shapefile' to data source
`./Output/brasil_shapefile.shp' using driver `ESRI Shapefile'
Writing 1 features with 168 fields and geometry type Multi Polygon.
Warning in CPL_write_ogr(obj, dsn, layer, driver,
as.character(dataset_options), : GDAL Message 1: Value 211049527 of field
pop_est of feature 0 not successfully written. Possibly due to too larger
number with respect to field width
Warning in CPL_write_ogr(obj, dsn, layer, driver,
as.character(dataset_options), : GDAL Message 1: Value 1159320441 of field
ne_id of feature 0 not successfully written. Possibly due to too larger number
with respect to field width
Warning in CPL_write_ogr(obj, dsn, layer, driver,
as.character(dataset_options), : GDAL Message 1: One or several characters
couldn't be converted correctly from UTF-8 to ISO-8859-1. This warning will
not be emitted anymore.
6. Transformar os dois Spatial Rasters em dataframes
Faremos os dois mapas finais com a função geom_tile do ggplot2. Daqui em diante, não iremos usar imagens raster no R. Usaremos o polígono do Brasil e dois dataframes derivados de nossos dois mapas de chuva (annual e mes mais chuvoso). Estes dois dataframes são tabelas com três colunas: latitude, longitude e o atributo a ser mapeado, herdado dos dous respectivos Spatial Rasters que criamos acima com argumentos do tipo fun=, dentro da função terra::app.
A razão para transformar os rasters em dataframes é que as funções geom_raster e geom_tile – que fazem parte do pacote ggplot2 e que usaremos para preencher o polígono do Brasil com dados raster – somente aceitam os valores das células em formato de um dataframe, contendo uma coluna com o valor de cada célula e outras duas colunas para as coordenadas x,y representando a posição de cada célula. O número de linhas em cada dataframe é equivalente ao número de células contendo dados no respectivo Spatial Raster
prec_anual_df <-as.data.frame(prec_anual, xy=TRUE)%>%drop_na()idx_mes_max_df <-as.data.frame(idx_mes_max, xy=TRUE)%>%drop_na()# Examine as primeiras linhas de cada dataframe criadahead(prec_anual_df)
x y sum
1 -73.97917 5.479167 1917.32
2 -73.93750 5.479167 1599.96
3 -73.89583 5.479167 1099.60
4 -73.85417 5.479167 1112.96
5 -73.81250 5.479167 1261.64
6 -73.77083 5.479167 1337.16
A coluna ‘sum’ do dataframe prec_anual_df, tem a precipitação anual. A coluna ‘which.idx’ do dataframe idx_mes_max_df tem os números de 1 a 12 (janeiro a dezembro), que neste caso é o mês mais chuvoso do ano para cada célula.
#Exportar o raster de precipitacao anualterra::writeRaster(prec_anual,filename="./Output/prec_anual_df.tif", overwrite=TRUE) #Exportar o raster dos meses mais chuvosos do anoterra::writeRaster(idx_mes_max,filename="./Output/idx_mes_max.tif", overwrite=TRUE, datatype ="INT2S")
Cheque se o raster foi devidamente exportado, abra no Qgis para confirmar e cheque as propriedades.
rm(prec_anual, idx_mes_max) #remover Spatial Rasters. Não serão usados maisgc() # liberação de memória
used (Mb) gc trigger (Mb) max used (Mb)
Ncells 1814780 97 3428158 183.1 3157826 168.7
Vcells 6805876 52 18453968 140.8 18453968 140.8
7. Criar mapa de precipitação anual usando rampa contínua viridis_c
Na primeira linha, atribuimos o mapa a um objeto com nome. Isto impede que o mapa seja mostrado auomaticamente. Isto é feito propositalmente, pois queremos ver como fica o mapa salvo em disco. Gravar o mapa em disco é feito com a função ‘ggsave’, que requer que o output de ggplot seja um objeto com nome.
Na terceira linha do bloco de código abaixo, a função geom_sf usa o objeto sf (um polígono vetorial do Brasil) obtido do site Natural Earth e inclui o símbolo de ‘graus’ e outra formatação dos rótulos dos eixos.
mapa_precip_annual <-ggplot() +geom_tile(aes(x = x, y = y, fill = sum), data = prec_anual_df) +# sum é campo no dataframegeom_sf(data = brasil, fill ='transparent', color ='white', size =3) +# borda Brasil scale_fill_viridis_c(name ='mm/ano\n(cinza, >3500)', direction =-1,limits =c(0, 3500)) +ggtitle("Chuva anual no Brasil") +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ítulolegend.key.size =unit(0.6, "cm"), # tamanho da legendalegend.title =element_text(size =8), # tamanho da fonte na legendalegend.position.inside =c(0.86, 0.20), # Posiciona legenda dentro da figuralegend.background =element_rect(fill ="white")) # fundo brancoplot(mapa_precip_annual)
Salvar em disco o mapa de precip anuial com a função ggsave, especificando tamanho, resolução e fundo branco. Usar este produto para ajustes finos de posição, tamanho fonte etc dos elementos no mapa dentro do bloco de código logo acima.
A maior parte do código será igual ao usado para chuva anual. Os valores de cada célula vão de 1 a 12 (janeiro a dezembro). Para a rampa de cores, usamos a opção “turbo” onde 1 e 12 (janeiro e dezembro) tem cores similares, pois estes dois meses são vizinhos no tempo.
meses <-c("Jan", "Fev", "Mar", "Abr", "Mai", "Jun", "Jul", "Ago", "Set", "Out", "Nov", "Dez")mes_mais_chuvoso <-ggplot() +geom_tile(aes(x = x, y = y, fill =as.factor(which.max)), data = idx_mes_max_df) +geom_sf(data = brasil, fill ='transparent', color ='white', size =3) +scale_fill_manual(values = viridis::viridis(12, option ="turbo", direction =-1),name ='Mês mais chuvoso', labels = meses) +ggtitle("Mês Mais Chuvoso no Brasil") +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ítulolegend.key.size =unit(0.6, "cm"), # tamanho da legendalegend.title =element_text(size =8), # tamanho da fonte na legendalegend.position.inside =c(0.86, 0.20), # Posiciona legenda dentro da figuralegend.background =element_rect(fill ="white")) # fundo brancoplot(mes_mais_chuvoso)
Com as imagens raster baixadas, abra no Qgis e faca um mapa com elementos necessários (i.e. legenda, escala, orientacao, fonte de dados, projecao cartográfica, basemap).
Assista o vídeo ensinando a fazer Layout do QGis na prática aqui
-Layout de chuva anual:
Você pode aplicar coloração por quartis ou em uma escala contínua diretamente no QGIS ao configurar a simbologia de um raster de precipitação anual. Abaixo, explico ambos os métodos:
1. Cores em Quartis
Abrir Propriedades do Raster:
Clique com o botão direito no raster e escolha Propriedades.
Configurar Simbologia:
Na aba Simbologia, escolha o tipo Renderizador Contínuo (Singleband pseudocolor).
No campo Modo, selecione Intervalos definidos por estatísticas.
Em Divisão, escolha Quantile (quartis).
Defina número de classes como 4 (para quartis).
Selecionar Paleta de Cores:
Escolha uma paleta em Cor, como Viridis, Spectral, ou outra.
Ajuste se necessário (inverter direção, cores personalizadas etc.).
Aplicar e Visualizar:
-Layout de meses mais chuvosos:
Em simbologia, modifique para paletizado/valores únicos e defina o nome de cada mes em funcao do número do mes.
Abrir Propriedades do Raster:
Clique com o botão direito no raster e selecione Propriedades.
Configurar Simbologia Categórica:
Vá até a aba Simbologia.
Escolha o tipo Valores Únicos (Categórico).
No campo de valores, use os números de 1 a 12 (ou os valores existentes no raster).
Para cada valor, adicione os nomes dos meses como rótulos: