Dados do comércio internacional no R

By | Data Science, Economia

A base de dados UN Comtrade fornece acesso gratuito às informações de comércio global. É possível obter dados usando sua interface de extração ou API. Por sorte, existe um pacote chamado comtradr que facilita o uso da API no R. Neste post, iremos mostrar um pouco da funcionalidade dele.

 

 
library(comtradr)
library(tidyverse)
library(splitstackshape)
library(inlmisc)

Primeiramente, iremos selecionar os dados de exportações brasileiras com quatro países: China, Estados Unidos, Argentina e Holanda. Como não selecionamos o período, todos os dados desde 1990 são retornados. Além disso, é possível filtrar por tipo de produto. Neste caso, estamos selecionando o total exportado.

 

 
braziltop <- ct_search(reporters = "Brazil", 
partners = c("USA", "China", "Argentina", "Netherlands"), 
trade_direction = "exports") %>%
ct_use_pretty_cols()

Veja que o formato de saída torna muito simples fazer um gráfico com a evolução das exportações ao longo dos anos. Fica evidente a rápida expensão chinesa para se tornar principal destino dos produtos brasileiros.

 

 
ggplot(braziltop, aes(Year, `Trade Value usd`, color = factor(`Partner Country`))) +
geom_point(size = 2) +
geom_line(size = 1) +
scale_color_manual(values = c("#011627", "#FF3366", "#2EC4B6", "#20A4F3"), 
name = "País") +
scale_shape_discrete(name = "País") +
labs(title = "Destino das exportações brasileiras",
y = "Total em dólares",
x = "Ano") +
theme_minimal()

Agora, ao invés de selecionarmos todos os produtos, iremos escolher apenas as exportações de peixe. Para fazer a filtragem, é preciso utilizar os códigos do sistema harmonizado (SH), que cataloga os produtos em categorias gerais e específicas. Essa filtragem é feita no parâmetro "commod_codes" da função. Iremos extrair apenas os dados referentes a 2020.

Como nós selecionamos os produtos, a API retorna o dado individual de cada categoria. Assim, iremos somar o valor exportado por país de destino.

export_peixe <- ct_search(reporters = "Brazil", 
partners = "All", 
trade_direction = "exports",
start_date = 2020, 
end_date = 2020,
commod_codes = c("0301","0302","0303","0304","0305")
) %>%
ct_use_pretty_cols() %>%
group_by(`Partner Country`) %>%
summarise("Total" = sum(`Trade Value usd`)) %>%
filter(`Partner Country` != "World")

A partir destes dados, iremos criar um gráfico que mostre a composição relativa de cada país. Para isso, é preciso antes fazer algumas alterações na formatação dos dados, de modo que o dataframe final tenha a o número de quadrantes ocupado por cada país e suas posições.

</pre>
export_peixe$prop = 100*export_peixe$Total/sum(export_peixe$Total)
export_peixe = expandRows(export_peixe, "prop")
export_peixe <- rbind(export_peixe,
export_peixe[sample(nrow(export_peixe), 100- nrow(export_peixe)), ])

var = export_peixe$`Partner Country`
nrows <- 10
categ_table <- round(table(var) * ((nrows*nrows)/(length(var))))

base <- expand.grid(y = 1:nrows, x = 1:nrows)
base$category <- factor(rep(names(categ_table), categ_table))

ggplot(base, aes(x = x, y = y, fill = category)) +
geom_tile(color = "black", size = 0.5) +
scale_fill_manual(values = as.character(GetColors(n = 9,
scheme = "roma",
alpha = 0.9))) +
theme_void() +
labs(title = "Destino das exportações de Peixe pelo Brasil (2020)") +
scale_x_continuous(expand = c(0, 0)) +
scale_y_continuous(expand = c(0, 0), trans = 'reverse') +
theme(panel.border = element_rect(size = 2, fill = NA))
<pre>
 

 

 

Com isso, vemos uma predominância dos Estados Unidos como destino das exportações brasileiras de peixe.

Hackeando o R: calculando o carry-over estatístico de uma variável

By | Hackeando o R

No Hackeando o R de hoje, vamos mostrar como fazer a visualização do carry-over estatístico de uma série temporal. Essa estatística pode ser uma importante ferramenta para a análise de agregações de dados ao longo do tempo pois facilita identificar a variação que ocorreu apenas no período corrente, extraindo a variação que é apenas resíduo do período anterior, chamada de carry-over. Abaixo, visualizamos os dois efeitos teóricos com uma simulação de dados

library(RcppRoll)
library(ggplot2)

set.seed(1)
x = data.frame(valor = (1:100) + rnorm(100), t = 1:100)

ggplot(x[49:72,], aes(x = t, y = valor))+geom_bar(stat = 'identity') +
coord_cartesian(ylim = c(45, 73))+
geom_segment(aes(x=49, xend=60, y = 54.68361, yend=54.68361), size = 1.2)+
geom_label(aes(x=49, y=54.68361, label = 'A'))+
geom_segment(aes(x=61, xend=72, y = 59.86495, yend=59.86495), size = 1.2)+
geom_label(aes(x=61, y=59.86495, label = 'B'))+
geom_segment(aes(x=61, xend=72, y = 66.85647, yend=66.85647), size = 1.2)+
geom_label(aes(x=61, y=66.85647, label = 'C'))+
labs(x='', y = '')+
theme_bw()

No exemplo acima, A é a média do ano anterior, C a média do ano corrente, e B é o valor da última observação do ano anterior, repetido para o ano corrente, ou seja, a média do ano corrente caso não houvesse crescimento. Ao compararmos a variação interanual dos dois períodos, podemos decompor esse valor em duas partes: a variação percentual de A a B, chamada de carry-over, e a variação percentual de B a C (mensurada no nível de A), que é o crescimento que ocorreu apenas a partir da última observação do ano anterior. A função abaixo calcula tais valores para uma variável mensal qualquer:


calcula_carry_over_anual <- function(data) {
A <- dplyr::lag(RcppRoll::roll_meanr(data, n=12), n=12)
B <- dplyr::lag(data, n = 12)
C <- RcppRoll::roll_meanr(data, n=12)

carry_over <- (B-A)/A
cresc_real_do_periodo <- (C-B)/A

lista = data.frame(carry_over, cresc_real_do_periodo, carry_over+cresc_real_do_periodo)
return(lista)
}

 

Então, vamos fazer a decomposição da série de nível do IBC como exemplo:

library(BETS)
library(tidyverse)
library(ggplot2)
library(scales)

ibc = BETSget(24363, data.frame=TRUE)

tibble(ibc$date, calcula_carry_over_anual((ibc$value))*100) %>%
magrittr::set_colnames(c('date', 'carry_over', 'cresc_real', 'soma')) %>%
pivot_longer(-date, names_to = 'var', values_to = 'val') %>%
filter(date>as.Date('2018-01-01') & var != 'soma') %>%
mutate(idk = RcppRoll::roll_sumr(val, n=2),
idk = ifelse(rep(c(FALSE, TRUE), times = 39), idk, NA)) %>%
ggplot(aes(x=date, y = val, fill = var))+geom_bar(stat = 'identity')+
scale_x_date(breaks = date_breaks('3 months'),
labels = date_format("%b/%Y"))+
scale_fill_manual(labels = c('Carry over', 'Crescimento real'), values = c('#244747', '#9ae5de'))+
geom_line(aes(x=date,y=idk, color = 'Agregado'), size= 1.2, linetype='solid')+
scale_color_manual(values = c('Agregado' = '#e89835'))+
geom_hline(yintercept=0, colour='black', linetype='dashed')+
labs(title='Decomposição da variação do nível do IBC', y = '%',
caption='Fonte: IBGE')+
theme(panel.background = element_rect(colour = 'white', fill='white'),
legend.position = 'right',
strip.text = element_text(size=8, face='bold'),
axis.text.x = element_text(angle = 45, hjust=1),
plot.title = element_text(size=10, face='bold'),
legend.title = element_blank(),
plot.caption.position = 'plot',
axis.title.x = element_blank())

Relatório AM #10 - IPCA-15

By | Indicadores

No Relatório AM dessa semana, vamos comentar brevemente o comportamento do IPCA-15, divulgado na metade de todo mês. Esse índice é importante pois é um previsor do IPCA cheio, logo suas características tendem a ser semelhantes. Uma delas é a sazonalidade, como podemos ver no gráfico abaixo:

Como podemos ver, o componente sazonal nos meses da metade do ano é baixo, indicando que as variações recentes que estamos vendo podem ser atribuídas em quase sua totalidade a choques exógenos que estão afetando a economia. Existem diversas fontes para tais choques, mas a mais importante do momento é o setor de energia, com a iminente crise hídrica que estamos sofrendo. Com isso, o IPCA-15 do último mês fechou em forte alta:

 

Como estimar o índice de Gini no R

By | Data Science, Economia, Indicadores

O índice de Gini é uma medida de distribuição de renda muito interessante e conhecida, que tenta expressar em um valor único a desigualdade apresentada na curva de Lorenz. Neste exercício mostramos como podemos estimar essa medida facilmente no R.

O índice de Gini consiste em um número entre 0 e 1, onde 0 corresponde à completa igualdade e 1 corresponde à completa desigualdade e pode ser calculado com a fórmula de Brown abaixo:

Onde:

G = coeficiente de Gini
X = proporção acumulada da variável "população"
Y = proporção acumulada da variável "renda"

Para esse exercício usaremos os microdados da PNAD Contínua trimestral do IBGE, que possui a variável Rendimento mensal efetivo de todos os trabalhos (VD4020). E para tornar o exercício interessante faremos a estimação do índice de Gini por estado (UF) brasileiro.

Pacotes

Para a finalidade do exercício utilizaremos os seguintes pacotes do R, todos disponibilizados no CRAN:


# Instalar/carregar pacotes
if(!require("pacman")) install.packages("pacman")
pacman::p_load(
"PNADcIBGE",
"survey",
"convey",
"tidyverse"
)

Dados

O último trimestre da pesquisa disponível na data deste exercício é referente ao 1º trimestre de 2021. Apontamos esse período na função get_pnadc(), especificamos as variáveis desejadas para coleta e convertemos o objeto resultante para a classe convey para poder fazer a estimação:


pnadc_0121 <- PNADcIBGE::get_pnadc(year = 2021, quarter = 1, vars = c("UF", "VD4020")) %>%
convey::convey_prep()

Estimar índice de Gini

Para estimar o índice de Gini o pacote convey oferece a função svygini, bastando especificar a variável de renda desejada. Como queremos a estimação por estado, usaremos também a função svyby do pacote survey, que serve justamente para calcularmos estatísticas por grupos dos nossos dados, nesse caso a UF.


gini_uf <- survey::svyby(
~VD4020,
by = ~UF,
design = pnadc_0121,
FUN = convey::svygini,
na.rm = TRUE
)

Visualizar os resultados

Por fim, vamos fazer um gráfico simples para visualizar o resultado da estimação, usando o ggplot2:


gini_uf %>%
dplyr::as_tibble() %>%
dplyr::mutate(UF = forcats::fct_reorder(UF, VD4020)) %>%
ggplot2::ggplot(ggplot2::aes(x = VD4020, y = UF)) +
ggplot2::geom_col(fill = "darkblue") +
ggplot2::theme_classic() +
ggplot2::labs(
title = "Índice de Gini por Estado",
subtitle = "Dados do 1º trimestre de 2021",
x = NULL,
y = NULL,
caption = "Fonte: Microdados PNADC-T/IBGE"
)

O que você achou do resultado? Surpreendente? Deixe suas impressões!

Quer aprender mais sobre utilização de microdados? Inscreva-se no curso de R e Python para Economistas.

 

Pirâmide etária da vacinação

By | Notas

A vacinação tem avançado de forma acelerada no país. Em muitos regiões, já estamos vacinando as pessoas na casa dos 30 anos. Apesar disso, ainda deve levar algum tempo até uma parcela considerável da população estar plenamente imunizada.

Inspirado em um tweet que vi a algum tempo, resolvi replicar essa visualização bastante interessante, que mostra a quantidade de vacinados em relação à população por meio de uma pirâmide etária.

Como os microdados da vacinação para todo o país geram um arquivo muito grande, optei por selecionar apenas dos dados do Paraná. Os dados de vacinação estão disponíveis aqui. Seguindo o código abaixo, é possível construir o mesmo gráfico para qualquer um dos estados.


library(tidyverse)
library(sidrar)

dados = read.csv2("vacina_PR.csv")

Os dados de população serão obtidos por meio da PNAD contínua anual de 2019, com o pacote *sidrar*. Para o gráfico ficar no formato de pirâmide, iremos alterar o valor das mulheres para ser negativo.

piramide = get_sidra(6706, 
variable = 606, 
geo = "State",
geo.filter = list("State" = 41)) %>%
filter(Sexo != "Total") %>%
mutate(`Grupo de idade` = ifelse(`Grupo de idade` == "5 a 9 anos",
"05 a 9 anos",
`Grupo de idade`),
Valor = ifelse(Sexo == "Mulheres",-Valor, Valor))

 

Primeiramente, iremos selecionar apenas as primeiras doses. A mesma visualização pode ser feita para a segunda dose, alterando essa linha. Como os dados de população estão disponíveis por faixas de renda, iremos colocar os dados de vacinação nestas faixas. Assim, agrupamos por faixa de idade e gênero para contar as pessoas. Da mesma forma que fizemos para a população, também iremos tornar a quantidade de mulheres negativa.


dados_g = dados %>%
filter(str_replace_all(vacina_descricao_dose, "[^[:alnum:]]", " ") == "1ª Dose") %>%
mutate(`Grupo de idade` = case_when(paciente_idade < 5 ~ "0 a 4 anos",
paciente_idade > 4 & paciente_idade < 10 ~ "05 a 9 anos",
paciente_idade > 9 & paciente_idade < 15 ~ "10 a 14 anos",
paciente_idade > 14 & paciente_idade < 20 ~ "15 a 19 anos",
paciente_idade > 19 & paciente_idade < 25 ~ "20 a 24 anos",
paciente_idade > 24 & paciente_idade < 30 ~ "25 a 29 anos",
paciente_idade > 29 & paciente_idade < 35 ~ "30 a 34 anos",
paciente_idade > 34 & paciente_idade < 40 ~ "35 a 39 anos",
paciente_idade > 39 & paciente_idade < 45 ~ "40 a 44 anos",
paciente_idade > 44 & paciente_idade < 50 ~ "45 a 49 anos",
paciente_idade > 49 & paciente_idade < 55 ~ "50 a 54 anos",
paciente_idade > 54 & paciente_idade < 60 ~ "55 a 59 anos",
paciente_idade > 59 & paciente_idade < 65 ~ "60 a 64 anos",
paciente_idade > 64 & paciente_idade < 70 ~ "65 a 69 anos",
paciente_idade > 69 & paciente_idade < 75 ~ "70 a 74 anos",
paciente_idade > 74 & paciente_idade < 80 ~ "75 a 79 anos",
paciente_idade > 79 ~ "80 anos ou mais")) %>%
group_by(`Grupo de idade`, paciente_enumsexobiologico) %>%
tally() %>%
filter(paciente_enumsexobiologico %in% c("M", "F")) %>%
mutate(n = ifelse(paciente_enumsexobiologico == "F", -n, n),
paciente_enumsexobiologico = ifelse(paciente_enumsexobiologico == "F",
"Mulheres",
"Homens"))

Assim, basta fazer o gráfico de pirâmide etária com dois "geom_bar". O primeiro é o de população, que fica no fundo. Já o segundo é o de dados da vacinação.

ggplot() +
geom_bar(data=piramide, 
aes(y=Valor*1000, x = `Grupo de idade`, fill = Sexo), 
stat="identity", alpha = 0.4, width=1) +
geom_bar(data=dados_g, 
aes(y=n, x = `Grupo de idade`, fill = paciente_enumsexobiologico),
stat="identity", width=1) +
scale_fill_manual(values = c("#c23c3c", "#276f9c", "#276f9c", "#c23c3c")) +
coord_flip() +
theme_classic() +
ggtitle("Pirâmide etária de vacinados com a primeira dose no Paraná",
subtitle="Com dados do dia 14/07/2021") +
scale_y_continuous(labels = abs) +
ylab("Número de pessoas") 

 

Receba diretamente em seu e-mail gratuitamente nossas promoções especiais
e conteúdos exclusivos sobre Análise de Dados!

Assinar Gratuitamente
{"cart_token":"","hash":"","cart_data":""}