Red Flags

Carregando bibliotecas

library(tidyverse)
library(hrbrthemes)
library(ggbeeswarm)
library(readr)
library(scales)
library(solitude)
library(viridis)
library(kableExtra)
library(lubridate)
library(GGally)
library(ggfortify)
theme_set(theme_minimal())

Carregando dados

Backup: 22/10/2020

filtro <- "bd-merenda" #bd-covid

#Dados da Receita Federal
cnaes_info <- read_csv(here::here(paste0("data/", filtro, "/info_cnaes.csv")))
cnaes_secundarios_rs <- read_csv(here::here(paste0("data/",filtro,"/cnaes_secundarios.csv")))
dados_cnpjs_rs <- read_csv(here::here(paste0("data/", filtro, "/dados_cadastrais.csv")))
natureza_juridica <- read_csv(here::here(paste0("data/", filtro, "/natureza_juridica.csv")))
socios_rs <- read_csv(here::here(paste0("data/", filtro, "/socios.csv")))

#Dados do TCE RS
compras <- read_csv(here::here(paste0("data/", filtro, "/info_contrato.csv")))
empenhos <- read_csv(here::here(paste0("data/", filtro, "/info_empenhos.csv")), 
                          col_types = list(cnpj_cpf = col_character()))
fornecedores <- read_csv(here::here(paste0("data/", filtro, "/info_fornecedores_contrato.csv")))
orgaos <- read_csv(here::here(paste0("data/", filtro, "/info_orgaos.csv")))

Calculando features

# Remove os orgaos estaduais
ids_orgaos_municipais <- orgaos %>% 
  filter(esfera == "MUNICIPAL") %>% 
  select(id_orgao)

# Compras realizadas pelos municipios
compras_municipios <- compras %>% 
  inner_join(ids_orgaos_municipais)

# ----------- Calculo de features de contratos --------------------

# Contratos:  Totais
recebido_total <- compras_municipios %>% 
  filter(!is.na(vl_contrato)) %>% 
  group_by(cnpj_cpf = nr_documento_contratado) %>% 
  summarise(recebido_total = sum(vl_contrato), 
            n_compras_total = n(), 
            mediana_por_compra = median(vl_contrato),
            maior_compra = max(vl_contrato),
            desvio_maior_compra_mediana = maior_compra - mediana_por_compra) %>% 
  ungroup()

# Contratos: Anuais
recebido_anual_por_empresa <- compras_municipios %>% 
  filter(!is.na(vl_contrato)) %>% 
  group_by(cnpj_cpf = nr_documento_contratado, ano = ano_contrato) %>% 
  summarise(recebido = sum(vl_contrato), 
            n_compras = n()) %>% 
  ungroup()

# Orgaos: Totais
orgaos_por_empresa_total <- compras_municipios %>% 
  group_by(cnpj_cpf= nr_documento_contratado) %>% 
  summarise(total_orgaos_fornecidos = n_distinct(id_orgao)) 

# Orgaos: Anuais
orgaos_anuais_por_empresa <- compras_municipios %>% 
  group_by(cnpj_cpf= nr_documento_contratado, ano = ano_contrato) %>% 
  summarise(orgaos_fornecidos = n_distinct(id_orgao)) 

# Municipios: Totais
municipios_por_empresa_total <- compras_municipios %>% 
  left_join(orgaos %>%
              select(id_orgao, cd_municipio_ibge)) %>% 
  group_by(cnpj_cpf = nr_documento_contratado) %>% 
  summarise(total_municipios_fornecidos = n_distinct(cd_municipio_ibge)) 

# Municipios: Anuais
municipios_anuais_por_empresa <- compras_municipios %>% 
  left_join(orgaos %>%
              select(id_orgao, cd_municipio_ibge)) %>% 
  group_by(cnpj_cpf = nr_documento_contratado, ano = ano_contrato) %>% 
  summarise(municipios_fornecidos = n_distinct(cd_municipio_ibge)) 


# Ultimo contrato
ultimo_contrato <- compras_municipios %>% 
  select(cnpj_cpf = nr_documento_contratado, dt_inicio_vigencia) %>% 
  group_by(cnpj_cpf) %>% 
  summarise(data_ultimo_contrato = max(dt_inicio_vigencia)) %>% 
  ungroup() %>% 
  na.omit()

# Recebido nos últimos 12 meses (a partir da data do ultimo contrato da empresa)
recebido_ultimo_ano <- compras_municipios %>% 
  select(cnpj_cpf = nr_documento_contratado, vl_contrato, dt_inicio_vigencia) %>% 
  left_join(ultimo_contrato) %>% 
  mutate(data_limite_inferior = ymd(data_ultimo_contrato) - years(1)) %>% 
  filter(ymd(dt_inicio_vigencia) >= ymd(data_limite_inferior)) %>% 
  group_by(cnpj_cpf) %>% 
  summarise(recebido_ultimo_ano = sum(vl_contrato))

# ----------- Calculo de features de fornecedores --------------------

# Experiencias e eventos
experiencia <- fornecedores %>% 
  select(cnpj = nr_documento, data_primeiro_contrato) %>% 
  left_join(dados_cnpjs_rs) %>% 
  select(cnpj, data_primeiro_contrato, data_inicio_atividade) %>% 
  left_join(ultimo_contrato %>% 
              select(cnpj = cnpj_cpf, data_ultimo_contrato)) %>% 
  mutate(experiencia = interval(ymd(data_primeiro_contrato),
                                ymd(data_ultimo_contrato)) %/% days(1),
         tempo_primeiro_contrato = interval(ymd(data_inicio_atividade),
                                            ymd(data_primeiro_contrato)) %/% days(1),
         cnpj_cpf = cnpj)

# Porte
porte_empresas <- dados_cnpjs_rs %>%
  select(cnpj_cpf = cnpj, porte_empresa) %>% 
  na.omit() %>% 
  mutate(porte_empresa = if_else(porte_empresa == "Empresa de pequeno porte", 
                                 "Pequeno porte", 
                                 porte_empresa))

# ----------- Medianas e desvios -------------------- 

features_agrupamentos_anuais <- recebido_anual_por_empresa %>% 
  left_join(orgaos_anuais_por_empresa) %>% 
  left_join(municipios_anuais_por_empresa) %>% 
  group_by(cnpj_cpf) %>% 
  summarise(mediana_municipios_fornecidos = median(municipios_fornecidos),
            mediana_recebido = median(recebido),
            maior_desvio_recebido_anual = max(recebido) - mediana_recebido) %>% 
  ungroup()

features_totais <- recebido_total %>% 
  left_join(orgaos_por_empresa_total) %>% 
  left_join(municipios_por_empresa_total) %>% 
  left_join(ultimo_contrato) %>% 
  left_join(recebido_ultimo_ano) %>% 
  mutate(n_recompras = n_compras_total - total_municipios_fornecidos) %>% 
  na.omit()

features_merge <- features_totais %>% 
  left_join(features_agrupamentos_anuais) %>% 
  left_join(experiencia) %>% 
  left_join(porte_empresas)

Filtrando features

features_numericas <- features_merge %>% 
  select(-cnpj_cpf, -data_ultimo_contrato, 
         -cnpj, -data_primeiro_contrato, 
         -data_inicio_atividade, -porte_empresa) 

ggcorr(features_numericas,
    title = "Correlograma das features",
    palette = "RdYlBu",
    label = TRUE,
    label_round = 2,
    label_size = 3,
    hjust = 0.8,
    size = 3,
    color = "black") +
  hrbrthemes::theme_ipsum_rc()

Features selecionadas:

• experiencia:
  • Cálculo: Número de dias percorridos da data do primeiro contrato para a data do último;
  • Significado: Experiência do fornecedor com compras públicas.
• desvio_maior_compra_mediana:
  • Cálculo: (Valor da maior compra) - (Mediana das compras);
  • Significado: O quão destoante foi a compra de maior valor daquele fornecedor para suas demais compras.
• n_recompras:
  • Cálculo: (Nº total de compras) - (N° total de municípios fornecidos);
  • Significado:
• porte_empresa

Detectando os outliers

detect_outliers <- features_merge  %>% select(cnpj_cpf, 
                                              porte_empresa,
                                              desvio_maior_compra_mediana, 
                                              experiencia,
                                              n_recompras) %>% na.omit()

iforest <- isolationForest$new(nrow(detect_outliers))

iforest$fit(detect_outliers %>% select(-cnpj_cpf))

## INFO  [13:12:06.495] dataset has duplicated rows 
## INFO  [13:12:06.518] Building Isolation Forest ...  
## INFO  [13:12:07.379] done 
## INFO  [13:12:07.380] Computing depth of terminal nodes ...  
## INFO  [13:12:07.744] done 
## INFO  [13:12:07.786] Completed growing isolation forest

detect_outliers$pred <- iforest$predict(detect_outliers %>% select(-cnpj_cpf))

p <- detect_outliers %>% ggplot(aes("Fornecedores", pred$anomaly_score, label = "")) +
             theme_ipsum_rc() +
             theme(legend.position = "top") +
             geom_quasirandom(width = 0.4, color = "#CBA135") +
             labs(
              x = NULL, y = "Score", fill = NULL,
              title = "Score de anormalidade"
             ) 
p

Limiar escolhido: 0.7

detect_outliers$outlier <- as.factor(ifelse(detect_outliers$pred$anomaly_score >=0.7, "outlier", "normal"))

Comportamento

library(kableExtra)
difs.pca <- prcomp(detect_outliers %>% select(-cnpj_cpf, -porte_empresa, -outlier, -pred), scale. = T,  center = T)
difs.pca.df <- data.frame(difs.pca$rotation)
kable(difs.pca.df, format = 'html') %>%
  kable_styling(bootstrap_options = c('hover', 'striped'))

	PC1	PC2	PC3
desvio_maior_compra_mediana	0.2640893	-0.9618020	-0.0720680
experiencia	0.6757396	0.2378234	-0.6977220
n_recompras	0.6882099	0.1355618	0.7127343

detect_outliers$porte_empresa_f = factor(detect_outliers$porte_empresa, 
                                         levels=c('Microempresa','Pequeno porte','Demais'))

pca <- detect_outliers %>% mutate(score = pred$anomaly_score)

autoplot(difs.pca,
         data = pca,
         colour='score',
         loadings = TRUE,
         loadings.colour = 'grey',
         loadings.label.colour = 'grey50',
  loadings.label = TRUE,
  loadings.label.size = 3) +
  theme_ipsum_rc(grid=FALSE) +
  labs(title = "Distribuição nos componentes principais") +
  facet_wrap(~porte_empresa_f) +
  theme(legend.position = "right") +
  scale_color_fermenter(n.breaks = 5, palette = "PuOr")

Detalhes dos outliers

outliers <- detect_outliers %>% 
  group_by(porte_empresa) %>% 
  mutate(rank_desvio = n() + 1 - rank(desvio_maior_compra_mediana),
         rank_experiencia =  n() + 1 -  rank(experiencia),
         rank_recompras =  n() + 1 -  rank(n_recompras),
         mediana_desvio = quantile(desvio_maior_compra_mediana, .90),
         mediana_experiencia = quantile(experiencia, .90),
         mediana_recompras = quantile(n_recompras, .90),
         ) %>% 
  ungroup() %>% 
  left_join(dados_cnpjs_rs %>% select(cnpj_cpf = cnpj, razao_social)) %>% 
  distinct() %>% 
  mutate(link = paste0("https://tadepemerenda.transparencia.org.br/fornecedores/", cnpj_cpf)) %>% 
  filter(outlier == "outlier")

tabela <- outliers %>% select(`Fornecedor` = razao_social, Porte = porte_empresa, everything()) %>% 
  mutate(`Ranking de desvio da maior compra` = paste0(rank_desvio, '°'),
         `Desvio da maior compra (em reais)` = paste0(round(desvio_maior_compra_mediana/1e6, 2), "mi"),
         `90% do mesmo porte teve desvio menor que (em reais)` = paste0(round(mediana_desvio/1e6,2 ), "mi"),
         `Ranking de experiencia` = paste0(round(rank_experiencia), '°'),
         `Experiência (em dias)` = paste0(round(experiencia)),
         `90% do mesmo porte tem experiencia menor que (em dias)` = paste0(round(mediana_experiencia)),
         `Ranking de recompras` = paste0(round(rank_recompras), '°'),
         `Nº de recompras` = paste0(round(n_recompras)),
         `90% do mesmo porte realizou menos que (recompras)` = paste0(round(mediana_recompras)),
         Link = link) %>% 
  select(-rank_desvio, -rank_experiencia, -rank_recompras, 
         -mediana_desvio, -mediana_experiencia, -mediana_recompras, 
         -link, -pred, -n_recompras, -desvio_maior_compra_mediana, 
         -experiencia, -n_recompras, -cnpj_cpf, -outlier, -porte_empresa_f)

kable(tabela, format = 'html') %>%
  kable_styling(bootstrap_options = c('hover', 'striped'))

Fornecedor	Porte	Ranking de desvio da maior compra	Desvio da maior compra (em reais)	90% do mesmo porte teve desvio menor que (em reais)	Ranking de experiencia	Experiência (em dias)	90% do mesmo porte tem experiencia menor que (em dias)	Ranking de recompras	Nº de recompras	90% do mesmo porte realizou menos que (recompras)	Link
COOPERATIVA DE AGRICULTORES E AGROINDUSTRIAS FAMILIARES E CAXIAS DO SUL LTDA	Demais	2°	1.96mi	0.17mi	2°	1128	896	11°	28	11	https://tadepemerenda.transparencia.org.br/fornecedores/14169702000108
COOPERATIVA LANGUIRU LTDA.	Demais	5°	0.69mi	0.17mi	1°	5583	896	5°	45	11	https://tadepemerenda.transparencia.org.br/fornecedores/89774160000100
COOPERATIVA DOS SUINOCULTORES DO CAI SUPERIOR LTDA	Demais	1°	28.88mi	0.17mi	3°	1002	896	12°	26	11	https://tadepemerenda.transparencia.org.br/fornecedores/91360420000134

Um exemplo

Onde está o outlier?

• CNPJ: 89774160000100
  
• Razão social: COOPERATIVA LANGUIRU LTDA.
  
• Porte: Demais

features_merge$porte_empresa_f = factor(features_merge$porte_empresa, 
                                         levels=c('Microempresa','Pequeno porte','Demais'))

p <- features_merge %>% 
  na.omit() %>% 
  ggplot(aes("", log2(desvio_maior_compra_mediana), 
             label = "", color = if_else(cnpj_cpf == '89774160000100', TRUE , FALSE),
             size = if_else(cnpj_cpf == '89774160000100', 3.5, 3))) +
             theme_ipsum_rc() +
             theme(legend.position = "") +
             scale_y_continuous(breaks = c(0, 9.9, 14.2, 19.1, 24), 
                                labels = c(0, "1mil", "20mil", "600mil", "28mi")) +
             geom_quasirandom(width = 0.3) +
             labs(
              x = NULL, y = "Desvio (em reais)", fill = NULL,
              title = "Distribuição do desvio da maior compra"
             ) +
  scale_color_brewer(palette = "Dark2") +
  facet_wrap(~porte_empresa_f) 

p

p <- features_merge %>% 
  na.omit() %>% 
  filter(experiencia > 0) %>% 
  ggplot(aes("", experiencia, 
             label = "", color = if_else(cnpj_cpf == '89774160000100', TRUE , FALSE))) +
             theme_ipsum_rc() +
             theme(legend.position = "") +
             geom_quasirandom(width = 0.5) +
             scale_y_continuous(breaks = c(0, 1e3, 3e3, 5e3), 
                                labels = c(0, "1mil", "3mil", "5mil")) +
             labs(
              x = NULL, y = "Experiencia (em dias)", fill = NULL,
              title = "Distribuição das experiencias"
             ) +
  scale_color_brewer(palette = "Dark2") +
  facet_wrap(~porte_empresa_f) 

p

p <- features_merge %>% 
  na.omit() %>% 
  ggplot(aes("", n_recompras, 
             label = "", color = if_else(cnpj_cpf == '89774160000100', TRUE , FALSE))) +
             theme_ipsum_rc() +
             theme(legend.position = "") +
             geom_quasirandom(width = 0.5) +
             labs(
              x = NULL, y = "Nº recompras", fill = NULL,
              title = "Distribuição do número de recompras"
             ) +
  scale_color_brewer(palette = "Dark2") +
  facet_wrap(~porte_empresa_f) 

p