Este breve tutorial tem como objetivo apresentar o processo de aquisição e visualização dos dados disponibilizados pelo CPTEC/INPE e UPF, que visa o desenvolvimento de uma ferramenta para acesso aos cenários de mudanças climáticas regionalizados no Brasil PROJETA.
O material mostrado neste documento serviu também de base para as publicações:
“The impact of climate change on renewable energy in the brazilian northeast” -> Pieter de Jong; Clemente A S Tanajura; Tarssio B Barreto; Asher Kipersok; Ednildo A Torres. Apresentado no XI Congresso Brasileiro de Planejamento Energético.
“The impact of climate change on hydroelectricity in Brazil” -> Pieter de Jong; Tarssio B Barreto; ASher Kipersok; Ednildo A Torres. Apresentado em: Advanced School on Water-Energy-Food Nexus.
“Estimating the Impact of Climate Change on Wind and Solar Energy in Brazil using a South American regional climate model” -> Pieter de Jong; Clemente A S Tanajura; Tarssio B Barreto; ASher Kipersok; Ednildo A Torres; Karla Esquerre. Aprovado na revista Renewable Energy.
Uma das questões levantadas sobre a mudança climática diz respeito a necessidade de maior geração de energia renovável para mitigação destes impactos. Porém, estas fontes de energias também serão impactadas por estas mudanças, como a intensidade da radiação solar.
Apesar disto, há poucos estudos que examinaram o impacto das mudanças climáticas na radiação solar e não há nenhum que realizou este estudo para a América do Sul. Esta série de artigos publicados visam contribuir, justamente, nesta lacuna do conhecimento.
Os dados podem ser adquiridos no PROJETA. Neste breve tutorial trabalharemos com dois períodos a fim de comparar duas normas climáticas, sendo a primeira entre 1960 e 1990 e a segunda entre 2070 e 2099. Utilizaremos os dados referentes a radiação solar dado o modelo HADGEM2 ES através do cenário RCP 8.5.
Para facilitar a comunicação entre os membros que trabalharam nestas publicações, os dados foram adquiridos em xls. Porém, para inclusão no R foi mais pertinente os transformar em txt.
data_f <- read.delim("dataf.txt", sep = ";") Com os pacotes do tidyverse iremos fazer pequenos ajustes. O primeiro deles diz respeito aos dias referentes as observações, uma breve análise no banco de dados apontou que os dados válidos estão indexados no dia 01/01. Dado isto, transformaremos nossa coluna de data para o formato DATE, retornaremos apenas o dia da observação e filtraremos as observações referentes ao dia 01.
É necessário também que utilizemos a função melt do pacote reshape2 ou outra similar como gather (tidyr). Desejamos criar, enfim, um banco de dados onde as duas primeiras colunas sejam Latitude e Longitude. A seguir, temos uma coluna com o valor da intensidade da radiação solar (OCIS) e o ano da observação está na coluna 5.
Este procedimento é visto em 1 (para os dados do futuro) e 2 (para os dados históricos):
dat_f <- data_f %>%
mutate(date = lubridate::ymd(Data)) %>%
mutate(day = lubridate::day(date)) %>%
filter(day == 1) %>%
mutate(year = lubridate::year(date)) %>%
dplyr::select(-c(day, Data, Hora, date)) %>%
mutate(Latitude = as.numeric(as.character(Latitude))) %>%
melt(id.vars = c("Latitude", "Longitude", "OCIS"),
measure.variable = c("year")) %>%
na.omit()## Warning: 2119947 failed to parse.
data_h <- read.delim("datah.txt", sep = ";")
dat_h <- data_h %>%
mutate(date = lubridate::ymd(Data)) %>%
mutate(day = lubridate::day(date)) %>%
filter(day == 1) %>%
mutate(year = lubridate::year(date)) %>%
dplyr::select(-c(day, Data, Hora,date)) %>%
mutate(Latitude = as.numeric(as.character(Latitude))) %>%
melt(id.vars = c("Latitude", "Longitude", "OCIS"),
measure.variable = c("year")) %>%
na.omit()## Warning: 2088306 failed to parse.
Uniremos as linhas destes bancos de dados:
dat_final <- rbind(dat_f, dat_h)Partimos agora para encontrar as normas climáticas, através da média, nos dois intervalos de tempo que foram trabalhos nestas publicações. Uma forma de fazermos este corte é utilizando a função cut e determinando breaks, outras alternativas como o case_when do pacote dplyr me parecem, hoje, mais racional.
dat_a <- dat_final %>%
mutate(factor = cut(value, breaks = c(-Inf, 1990, 2070, 2099),
labels = c("1", "2", "3"))) %>%
filter(factor == 1) %>%
group_by(Latitude, Longitude) %>%
summarise(mean_1 = mean(OCIS))dat_b <- dat_final %>%
mutate(factor = cut(value, breaks = c(-Inf, 1990, 2070, 2099),
labels = c("1", "2", "3"))) %>%
filter(factor == 3) %>%
group_by(Latitude, Longitude) %>%
summarise(mean_2 = mean(OCIS)) %>%
ungroup() Uniremos os dois bancos de dados (inner_join) e criaremos uma variável que representa a variação em porcentagem da norma climática ponto a ponto.
var <- dat_a %>%
inner_join(dat_b, by = c("Latitude" = "Latitude", "Longitude" = "Longitude")) %>%
mutate(var = (mean_2 - mean_1)*100/ mean_1) %>%
mutate(group = 0)Chegamos ao fim com um banco de dados com cerca de 40 mil linhas, perceba que nos primeiros bancos de dados tinha cerca de 1 milhão de linhas. Este processo foi fundamental para que conseguíssemos criar uma visualização interpretável dos dados. Vamos agora a esta etapa
O primeiro passo é encontrar um shape adequado para este trabalho, a princípio, como queremos trabalhar com estados, sugiro o uso deste.
shape <- readOGR("regioes_2010.shp")## OGR data source with driver: ESRI Shapefile
## Source: "C:\Users\barreto\Desktop\Nova pasta\Pieter\regioes_2010.shp", layer: "regioes_2010"
## with 5 features
## It has 3 fields
## Integer64 fields read as strings: id
aux <- fortify(shape) # Transformando nosso mapa em um dataframe## Regions defined for each Polygons
O primeiro passo é criar breakpoints para determinar nossa escala de cores. Criaremos para isto um histograma da variável OCIS.
p1 <- ggplot(data = var, aes(x = var)) +
geom_histogram(fill = "darkblue", col = "white") +
labs(x = "Variação (%)", y = "") +
geom_density() +
theme_bw()
plotly::ggplotly(p1)## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
Optamos, depois de algumas tentativas, pelo seguinte formato:
pretty_breaks <- c(-10, -8, -6, -4, -3, -2, -1, 0,
1, 2, 3, 4, 6, 8, 10, 15, 20, 30)Para montarmos a nossa escala, vamos utilizar uma sequência de funções. Nesta sequência determinaremos os extremos, as labels e o arredondamento destas:
# Encontrando Extremos
minVal <- round(min(var$var, na.rm = T),0)-1
# Labels
labels <- c()
brks <- c(minVal, pretty_breaks)
# Arrendondamento
for(idx in 1:length(brks)){
labels <- c(labels,round(brks[idx + 1], 2))
}
labels <- labels[1:length(labels)-1]
var$brks <- cut(var$var,
breaks = brks,
include.lowest = TRUE,
labels = labels)
brks_scale <- levels(var$brks)
labels_scale <- rev(brks_scale)Então, agora uma breve colaboração de cores dada por Isabela Almeida:
colb <- c("red4", "#800026","#B10026", "#E31A1C", "#FC4E2A", "#FD8D3C", "#FEB24C",
"#FED976", "khaki1", "#FCFCBD","#F7FCF0","#B7FCE3","#CCEBC5","#A8DDB5",
"#7BCCC4" ,"#4EB3D3", "#2B8CBE","#0868AC","#084081")ggplot() +
geom_raster(data = var, aes(x = Longitude, y = Latitude, fill = brks)) +
geom_polygon(data = aux, aes(x = long, y = lat, group = group),
fill = NA, color = "black", size = 1) +
labs(fill = "Variation of OCIS (%)") +
coord_equal() +
theme(legend.position = "bottom") +
theme_bw()Nos artigos em questão, trabalhamos a questão das usinas geradoras de energia solar. É de destaque, por exemplo, a cidade de Salvador, com acréscimo de 4,7%, atingindo um índice médio anual de 261 W/m².Em Pirapora, em Minas Gerais, onde este acréscimo é de 6.3%, alcançando cerca de 300 W/m².
A seguir vamos encontrar o ponto de maior acréscimo e aquele que possuí, ao final, a maior intensidade de radiação solar. Esta parte não foi de interesse da publicação, mas exploraremos mais algumas possibilidades, já que parte do que foi feito ainda não está disponível ou publicado.
maximo <- max(var$var)
var %>%
filter(var == maximo)## # A tibble: 1 x 7
## # Groups: Latitude [1]
## Latitude Longitude mean_1 mean_2 var group brks
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <fct>
## 1 -1.50 -58.5 211. 258. 22.5 0 30
Se pesquisarmos esta coordenada, chegamos a um local próximo a Manaus, mostrando o impacto que as mudanças climáticas teriam na Região Amazônica, pode-se conferir no mapa esta afirmação.
maximo <- max(var$mean_2)
var %>%
filter(mean_2 == maximo)## # A tibble: 1 x 7
## # Groups: Latitude [1]
## Latitude Longitude mean_1 mean_2 var group brks
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <fct>
## 1 -13.1 -73.1 339. 336. -0.705 0 0
No período de 2070-2090, o local com maior intensidade de radiação solar será no Peru, próximo a cidade histórica de Machu Pichu. É interessante, então, observar qual seria, no Brasil, a localidade que assumiria este posto.
Para isto, faremos um “merge” entre o banco de dados que representa nosso mapa do Brasil e o que representa os nossos dados. Desta forma, teremos um mapa que será limitado as fronteiras nacional:
aux2 <- aux %>%
mutate(lat = round(lat, digits = 2)) %>%
mutate(long = round(long, digits = 2)) %>%
inner_join(var, by = c("lat" = "Latitude", "long" = "Longitude"))Agora, repetiremos o processo para determinar o local no Brasil que terá a maior média de incidência solar no período de 2070 a 2090.
maximo <- max(aux2$mean_2)
aux2 %>%
filter(mean_2 == maximo)## long lat order hole piece id group.x mean_1 mean_2 var
## 1 -43.9 -14.3 5383 FALSE 1 1 1.1 300.1925 310.1696 3.323563
## 2 -43.9 -14.3 5385 FALSE 1 1 1.1 300.1925 310.1696 3.323563
## 3 -43.9 -14.3 2336 FALSE 1 3 3.1 300.1925 310.1696 3.323563
## 4 -43.9 -14.3 2338 FALSE 1 3 3.1 300.1925 310.1696 3.323563
## group.y brks
## 1 0 4
## 2 0 4
## 3 0 4
## 4 0 4
O local que apresenta maior OCIS fica próximo a bacia do Rio São Francisco e ao município de Carinhanha.
Percebam que pela imagem indicado e pela localização temos uma zona de semiárido que já é bastante assolada pela questão do balanço hídrico que ficará ainda mais comprometido com o aumento da radiação solar.
Para fechar, vamos comparar a série temporal de radiação solar antes e depois para o município de Salvador. Para filtrar apenas este município usaremos o “spdplyr”.
mun_shapes <- readOGR("DPA_A_100K_2017_06_14_GCS_SIR_SEI.shp") %>%
filter(MUNICIPIO == "Salvador")## OGR data source with driver: ESRI Shapefile
## Source: "C:\Users\barreto\Desktop\Nova pasta\Pieter\DPA_A_100K_2017_06_14_GCS_SIR_SEI.shp", layer: "DPA_A_100K_2017_06_14_GCS_SIR_SEI"
## with 417 features
## It has 2 fields
ssa <- fortify(mun_shapes) # Transformando nosso mapa em um dataframe## Regions defined for each Polygons
# Filtrando os valores referentes ao município de Salvador
ssa_value <- ssa %>%
mutate(lat = round(lat,digits = 1)) %>%
mutate(long = round(long, digits = 1)) %>%
inner_join(dat_final, by = c("lat" = "Latitude", "long" = "Longitude"))
# Encontrando a média anual do OCIS para o município de Salvador
ssa_hist <- ssa_value %>%
group_by(value) %>%
summarise(mean = mean(OCIS))
# Plotting
p3 <- ssa_hist %>% mutate(factor = cut(value, breaks = c(-Inf, 1990, 2070, 2099),
labels = c("1", "2", "3"))) %>%
filter(factor == 1 | factor == 3) %>%
ggplot(aes(x = value, y = mean)) +
geom_line(aes(fill = factor)) +
facet_wrap(~factor, scales = "free_x", ncol = 1) +
theme_bw()## Warning: Ignoring unknown aesthetics: fill
ggplotly(p3)Chegamos, então, ao fim deste breve tutorial de como adquirir e visualizar dados referentes ao PROJETA (CPTEC/INPE). Qualquer coisa, segue meu email: tarssioesa at gmail dot com
Referência:
-> CPTEC/INPE.- Centro de Previsão de Tempo e Estudos Climáticos / Instituto Nacional de Pesquisas Espaciais. Plataforma PROJETA - Projeções de mudança do clima para a América do Sul regionalizadas pelo Modelo Eta, 2018. Available at: https://projeta.cptec.inpe.br/ Accessed on: 25/04/2018.