Análise temporal e de perfil da abertura e fechamento das empresas no Maranhão - 2014 a 2024 - Dados Receita Federal do Brasil

Carregando o banco de dados

dados_parquet_BRF = arrow::read_parquet("BRF_FINAL_1.parquet")
#dados_parquet_depara_CNAESEC = arrow::read_parquet("CNAESEC.parquet")
dados_parquet_SOCIOS = arrow::read_parquet("Socios.parquet")

dados_parquet_BRF_filtrado = dados_parquet_BRF %>%  select ( MATRIZ,
                            NOME_F,
                            SIT_CAD,
                            DATA_SIT_CAD,
                            MOT_SIT_CAD,
                            DATA_INICIO_ATV,
                            CNAE_PRINC,
                            CNAE_SEC,
                            MUNIC,
                            SIT_ESP,
                            DATA_SIT_ESP,
                            CNPJ,
                            RAZAO_SOCIAL,
                            NATUREZA_JURIDICA,
                            CAPITAL_SOCIAL,
                            PORTE,
                            OP_SIMPLES,
                            DT_OP_SIMPLES,
                            OP_MEI,
                            DT_OP_MEI,
                            DT_EX_SIMPLES,
                            DT_EX_MEI,
                            code_muni,
                            PUBLICO_SEBRAE
                            )

dados_parquet_SOCIOS_filtrado = dados_parquet_SOCIOS %>% select( CNPJ_BAS,
                                                              IDENTIFICADOR_DE_SÓCIO,
                                                           DATA_DE_ENTRADA_SOCIEDADE)

#transformando o campo CNPJ de 8 dígitos
dados_parquet_BRF_filtrado$CNPJ_BAS = str_sub(dados_parquet_BRF_filtrado$CNPJ,1,8)

Transformação das variáveis de data e rótulos

dados_parquet_BRF_filtrado$ANO_DATA_INICIO_ATV = Year(dados_parquet_BRF_filtrado$DATA_INICIO_ATV)

dados_parquet_BRF_filtrado$MES_DATA_INICIO_ATV = Month(dados_parquet_BRF_filtrado$DATA_INICIO_ATV)

dados_parquet_BRF_filtrado$DIA_DATA_INICIO_ATV = day(dados_parquet_BRF_filtrado$DATA_INICIO_ATV)

dados_parquet_BRF_filtrado = dados_parquet_BRF_filtrado %>% mutate(SIT_CAD_ROTULADA = case_when( 
  
     SIT_CAD == 1 ~ "Nula",
                                                                    SIT_CAD == 2 ~ "Ativa",
                                                                    SIT_CAD == 3 ~ "Suspensa",
                                                                    SIT_CAD == 4 ~ "Inapta",
                                                                    SIT_CAD == 8 ~ "Baixada"
))

dados_parquet_BRF_filtrado = dados_parquet_BRF_filtrado %>% mutate(MATRIZ_ROTULADA = case_when( 
  
     MATRIZ == 1 ~ "Matriz",
                                                                           MATRIZ == 2 ~ "Filial"
                                                                  
))

dados_parquet_BRF_filtrado = dados_parquet_BRF_filtrado %>% mutate(OP_SIMPLES_ROTULADA = case_when( 
  
     OP_SIMPLES == 0 ~ "Não",
                                                                         OP_SIMPLES == 1 ~ "Sim"
     
                                                                  
))

dados_parquet_BRF_filtrado = dados_parquet_BRF_filtrado %>% mutate(OP_MEI_ROTULADA = case_when( 
  
     OP_MEI == 0 ~ "Não",
                                                                                   OP_MEI == 1 ~ "Sim"
     
                                                                  
))

dados_parquet_BRF_filtrado = dados_parquet_BRF_filtrado %>% mutate(PUBLICO_SEBRAE_ROTULADA = case_when( 
  
     PUBLICO_SEBRAE == 0 ~ "Não",
                                                                                   PUBLICO_SEBRAE == 1 ~ "Sim"
     
                                                                  
))

dados_parquet_BRF_filtrado = dados_parquet_BRF_filtrado %>% mutate(PORTE_ROTULADA = case_when( 
  
     PORTE == 1 ~ "MEI",
                                                                           
     PORTE == 2 ~ "ME",
     
     PORTE == 3 ~ "EPP",
     
     PORTE == 4 ~ "Demais"
     
                                                                  
))


dados_parquet_BRF_filtrado$CAPITAL_SOCIAL_NUMERO = as.numeric(gsub(",",".",dados_parquet_BRF_filtrado$CAPITAL_SOCIAL))

#PEGAR AS CNAES PRINCIPAIS E CLASSIFICAR 

dados_parquet_BRF_filtrado$CNAE_PRINC_2DIGITOS = as.numeric(str_sub(dados_parquet_BRF_filtrado$CNAE_PRINC,1,2)) 

dados_parquet_BRF_filtrado = dados_parquet_BRF_filtrado %>% mutate(CNAE_PRINC_2DIGITOS_ROTULADA = case_when( 
  
     (CNAE_PRINC_2DIGITOS >= 1 & CNAE_PRINC_2DIGITOS <= 3) ~ "Agricultura, pecuária, produção florestal, pesca e aqüicultura",
                                                                           
     (CNAE_PRINC_2DIGITOS >= 5 & CNAE_PRINC_2DIGITOS <= 9) ~ "Indústrias extrativas",
     
     (CNAE_PRINC_2DIGITOS >= 10 & CNAE_PRINC_2DIGITOS <= 33) ~ "Indústrias de transformação",
     
     CNAE_PRINC_2DIGITOS == 35 ~ "Eletricidade e gás",
      
     (CNAE_PRINC_2DIGITOS >= 36 & CNAE_PRINC_2DIGITOS <= 39) ~ "Água, esgoto, atividades de gestão de resíduos e descontaminação",
      
     (CNAE_PRINC_2DIGITOS >= 41 & CNAE_PRINC_2DIGITOS <= 43) ~ "Construção",
     
     (CNAE_PRINC_2DIGITOS >= 45 & CNAE_PRINC_2DIGITOS <= 47) ~ "Comércio, reparação de veículos automotores e motocicletas",
     
     (CNAE_PRINC_2DIGITOS >= 49 & CNAE_PRINC_2DIGITOS <= 53) ~ "Transporte, armazenagem e correio",
     
     (CNAE_PRINC_2DIGITOS >= 55 & CNAE_PRINC_2DIGITOS <= 56) ~ "Alojamento e alimentação",
     
     (CNAE_PRINC_2DIGITOS >= 58 & CNAE_PRINC_2DIGITOS <= 63) ~ "Informação e comunicação",
     
     (CNAE_PRINC_2DIGITOS >= 64 & CNAE_PRINC_2DIGITOS <= 66) ~ "Atividades financeiras, de seguros e serviços relacionados",
     
     CNAE_PRINC_2DIGITOS == 68  ~ "Atividades imobiliárias",
     
     (CNAE_PRINC_2DIGITOS >= 69 & CNAE_PRINC_2DIGITOS <= 75) ~  "Atividades profissionais, científicas e técnicas",
     
     (CNAE_PRINC_2DIGITOS >= 77 & CNAE_PRINC_2DIGITOS <= 82) ~  "Atividades administrativas e serviços complementares",
       
      CNAE_PRINC_2DIGITOS == 84 ~ "Administração pública, defesa e seguridade social",
      
      CNAE_PRINC_2DIGITOS == 85 ~ "Educação",
      
     (CNAE_PRINC_2DIGITOS >= 86 & CNAE_PRINC_2DIGITOS <= 88) ~ "Saúde humana e serviços sociais",
      
     (CNAE_PRINC_2DIGITOS >= 90 & CNAE_PRINC_2DIGITOS <= 93) ~ "Artes, cultura, esporte e recreação",
      
     (CNAE_PRINC_2DIGITOS >= 94 & CNAE_PRINC_2DIGITOS <= 96) ~ "Outras atividades de serviços",
      
     CNAE_PRINC_2DIGITOS == 97 ~ "Serviços domésticos", 
     
     CNAE_PRINC_2DIGITOS == 99 ~ "Organismos internacionais e outras instituições extraterritoriais"  
))

1) Análise temporal das aberturas de empresas

Número de empresas matrizes que iniciaram as atividades no Maranhão por mês de 2014 a 2024

dados_parquet_BRF_filtrado_agregado = dados_parquet_BRF_filtrado %>% 
  filter(MATRIZ == 1) %>% 
  mutate(
    mes_atividade = floor_date(DATA_INICIO_ATV, "month"),
    ano_atividade = floor_date(DATA_INICIO_ATV, "year")
  ) %>% 
  filter(DATA_INICIO_ATV > "2014-01-01",
         DATA_INICIO_ATV < "2024-08-10") %>% 
  count(mes_atividade)

grafico_linha_inicio_atividades =  
  
  dados_parquet_BRF_filtrado_agregado %>% 
  #mutate(mes_atividade = ano_atividade) %>% 
  filter(year(mes_atividade) > 2014,
         year(mes_atividade) < 2024) %>% 
  ggplot(aes(x = mes_atividade, y = n)) +
  geom_line(size = 1) +
  theme_minimal(14) +
  scale_x_date(breaks = scales::date_breaks("3 month"),
               labels = scales::date_format("%Y-%m")) +
  scale_y_continuous(labels = scales::label_number(big.mark = "."),breaks = seq(0,35000,5000 ) )+
  ggtitle("Quantidade de empresas que iniciaram atividade por mês 
no Maranhão - 2014 a 2024"
        ) + 
  #scale_y_continuous(limits = c(0, 800)) +
  labs(x = "Ano de início das atividades",
       y = "Quantidade de empresas") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        axis.text = element_text(size = 9),
        panel.grid.minor.x = element_blank()) 

grafico_linha_inicio_atividades

Podemos já começar a vislumbrar como está a tendência e a sazonalidade ao longo dos anos por trimestre, mas ainda é importante decompor a série para melhor avaliarmos.

Plotando sazonalidades ao longo dos meses

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))

dados_ts1 %>% as_tsibble() %>% gg_season(value, labels="both")+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico sazonal da abertura de empresas - 2014 a 2024")

Podemos notar que os meses de agosto e setembro tem uma sazonalidade muito forte, destacando o pico de 2015, 2023 e abril de 2024 também se destacando frente ao mesmo mês nos outros anos.

Plotando médias mensais e observando sazonalidades ao longo dos meses

dados_ts1 %>% as_tsibble() %>% gg_subseries(value)+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico das médias mensais e sazonalidade da abertura de empresas 
- 2014 a 2024")

Por meio desse gráfico, podemos ver como a abertura das empresas variam ao longo de cada mês do ano. Podemos notar que os meses de agosto e setembro possuem maiores aberturas de empresas em média, ao passo que dezembro apresenta a menor média histórica.

Decomposição da série temporal em suas componentes para analisar tendência e sazonalidade

dados_ts2 = dts2(dados_ts1)
dados_ts2$date = as.Date(dados_ts2$date)

ggdecompose(dados_ts2) +
  xlab("Date")+
  ylab("Quantidade de empresas")+
  scale_x_date(date_breaks = "1 year", date_labels = "%Y-%m")+
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Tendência: Ao decompor a série mensal da abertura das empresas, podemos notar que no ano de 2014 a 2016 parece mais estacionário e logo após 2016, nos primeiros meses temos uma tendência de crescimento e depois volta a estacionaridade até 2018, onde começa a ter tendência crescente novamente, até chegar em 2021 e repetir um padrão parecido até 2022.
Sazonalidade: Há uma sazonalidade demarcada muito clara de forte alta no penúltimo trimestre de cada ano.

Analisando a decomposição por porte das empresas - MEI

dados_parquet_BRF_filtrado_agregado = dados_parquet_BRF_filtrado %>% 
  filter(MATRIZ == 1 & PORTE_ROTULADA == "MEI") %>% 
  mutate(
    mes_atividade = floor_date(DATA_INICIO_ATV, "month"),
    ano_atividade = floor_date(DATA_INICIO_ATV, "year")
  ) %>% 
  filter(DATA_INICIO_ATV > "2014-01-01",
         DATA_INICIO_ATV < "2024-08-10") %>% 
  count(mes_atividade)

grafico_linha_inicio_atividades =  
  
  dados_parquet_BRF_filtrado_agregado %>% 
  #mutate(mes_atividade = ano_atividade) %>% 
  filter(year(mes_atividade) > 2014,
         year(mes_atividade) < 2024) %>% 
  ggplot(aes(x = mes_atividade, y = n)) +
  geom_line(size = 1) +
  theme_minimal(14) +
  scale_x_date(breaks = scales::date_breaks("3 month"),
               labels = scales::date_format("%Y-%m")) +
  scale_y_continuous(labels = scales::label_number(big.mark = "."),breaks = seq(0,3000,600 ) )+
  ggtitle("Quantidade de empresas MEI que iniciaram atividade por mês 
no Maranhão - 2014 a 2024"
        ) + 
  #scale_y_continuous(limits = c(0, 800)) +
  labs(x = "Ano de início das atividades",
       y = "Quantidade de empresas") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        axis.text = element_text(size = 9),
        panel.grid.minor.x = element_blank()) 

grafico_linha_inicio_atividades

Plotando sazonalidades ao longo dos meses das MEI

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))

dados_ts1 %>% as_tsibble() %>% gg_season(value, labels="both")+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico sazonal da abertura de empresas - 2014 a 2024")

Podemos destacar que de fevereiro para março em todos os anos houve tendência crescente, menos os anos de 2019 e 2014 que foi descrescente ou estacionário, cabe destacar também o decréscimo acentuado de julho a agosto de 2024 frente aos demais anos.

Plotando médias mensais e observando sazonalidades ao longo dos meses das MEI

dados_ts1 %>% as_tsibble() %>% gg_subseries(value)+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico das médias mensais e sazonalidade da abertura de empresas 
- 2014 a 2024")

Por meio desse gráfico, podemos ver como a abertura das empresas variam ao longo de cada mês do ano. Podemos notar que os meses de janeiro e julho possuem maiores aberturas de empresas em média, ao passo que dezembro apresenta a menor média histórica, igual acontece de forma ampla considerando todas os portes.

Decomposição da série temporal em suas componentes para analisar tendência e sazonalidade das empresas MEI

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))
dados_ts2 = dts2(dados_ts1)
dados_ts2$date = as.Date(dados_ts2$date)

ggdecompose(dados_ts2) +
  xlab("Date")+
  ylab("Quantidade de empresas")+
  scale_x_date(date_breaks = "1 year", date_labels = "%Y-%m")+
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Tendência: Ao decompor a série mensal da abertura das empresas para as empresas de porte MEI, podemos notar que diferentemente do comportamento geral apresenta uma tendência crescente do ano de 2015 a 2016 e depois fica estacionária até 2017, onde retorna o padrão de crescimento até 2021, logo depois estaciona no período da pandemia e retoma o crescimento em 2023 a 2024.
Sazonalidade: Há uma sazonalidade demarcada muito clara de forte queda no último trimestre de cada ano.

Analisando a decomposição por porte das empresas - ME

dados_parquet_BRF_filtrado_agregado = dados_parquet_BRF_filtrado %>% 
  filter(MATRIZ == 1 & PORTE_ROTULADA == "ME") %>% 
  mutate(
    mes_atividade = floor_date(DATA_INICIO_ATV, "month"),
    ano_atividade = floor_date(DATA_INICIO_ATV, "year")
  ) %>% 
  filter(DATA_INICIO_ATV > "2014-01-01",
         DATA_INICIO_ATV < "2024-08-10") %>% 
  count(mes_atividade)

grafico_linha_inicio_atividades =  
  
  dados_parquet_BRF_filtrado_agregado %>% 
  #mutate(mes_atividade = ano_atividade) %>% 
  filter(year(mes_atividade) > 2014,
         year(mes_atividade) < 2024) %>% 
  ggplot(aes(x = mes_atividade, y = n)) +
  geom_line(size = 1) +
  theme_minimal(14) +
  scale_x_date(breaks = scales::date_breaks("3 month"),
               labels = scales::date_format("%Y-%m")) +
  scale_y_continuous(labels = scales::label_number(big.mark = "."),breaks = seq(0,3000,600 ) )+
  ggtitle("Quantidade de empresas ME que iniciaram atividade por mês 
no Maranhão - 2014 a 2024"
        ) + 
  #scale_y_continuous(limits = c(0, 800)) +
  labs(x = "Ano de início das atividades",
       y = "Quantidade de empresas") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        axis.text = element_text(size = 9),
        panel.grid.minor.x = element_blank()) 

grafico_linha_inicio_atividades

Plotando sazonalidades ao longo dos meses das ME

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))

dados_ts1 %>% as_tsibble() %>% gg_season(value, labels="both")+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico sazonal da abertura de empresas - 2014 a 2024")

Podemos destacar que de fevereiro para março em todos os anos houve tendência crescente, menos os anos de 2021 que foi descrescente, o crescimento acentuado e constante de maio, junho e julho de 2020 e cabe destacar também o decréscimo acentuado de julho a agosto de 2024 frente aos demais anos, assim como as empresas de porte MEI.

Plotando médias mensais e observando sazonalidades ao longo dos meses das ME

dados_ts1 %>% as_tsibble() %>% gg_subseries(value)+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico das médias mensais e sazonalidade da abertura de empresas 
- 2014 a 2024")

Por meio desse gráfico, podemos ver como a abertura das empresas variam ao longo de cada mês do ano. Podemos notar que os meses de janeiro , março, maio e julho possuem maiores aberturas de empresas em média, ao passo que novembro e dezembro apresenta a menor média histórica.

Decomposição da série temporal em suas componentes para analisar tendência e sazonalidade das empresas ME

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))
dados_ts2 = dts2(dados_ts1)
dados_ts2$date = as.Date(dados_ts2$date)

ggdecompose(dados_ts2) +
  xlab("Date")+
  ylab("Quantidade de empresas")+
  scale_x_date(date_breaks = "1 year", date_labels = "%Y-%m")+
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Tendência: Ao decompor a série mensal da abertura das empresas para as empresas de porte ME, podemos notar um comportamento de queda de 2014 a 2016 e logo depois alternâncias entre quedas e subidas, até 2020 onde tem uma tendência crescente até 2021 e depois descresce durante a pandemia em 2022 até 2023 e depois caminha a estacionariedade.
Sazonalidade: Há uma sazonalidade demarcada muito clara de forte queda no último trimestre de cada ano, assim como as MEI.

Analisando a decomposição por porte das empresas - EPP

dados_parquet_BRF_filtrado_agregado = dados_parquet_BRF_filtrado %>% 
  filter(MATRIZ == 1 & PORTE_ROTULADA == "EPP") %>% 
  mutate(
    mes_atividade = floor_date(DATA_INICIO_ATV, "month"),
    ano_atividade = floor_date(DATA_INICIO_ATV, "year")
  ) %>% 
  filter(DATA_INICIO_ATV > "2014-01-01",
         DATA_INICIO_ATV < "2024-08-10") %>% 
  count(mes_atividade)

grafico_linha_inicio_atividades =  
  
  dados_parquet_BRF_filtrado_agregado %>% 
  #mutate(mes_atividade = ano_atividade) %>% 
  filter(year(mes_atividade) > 2014,
         year(mes_atividade) < 2024) %>% 
  ggplot(aes(x = mes_atividade, y = n)) +
  geom_line(size = 1) +
  theme_minimal(14) +
  scale_x_date(breaks = scales::date_breaks("3 month"),
               labels = scales::date_format("%Y-%m")) +
  scale_y_continuous(labels = scales::label_number(big.mark = "."),breaks = seq(0,200,50 ) )+
  ggtitle("Quantidade de empresas EPP que iniciaram atividade por mês 
no Maranhão - 2014 a 2024"
        ) + 
  #scale_y_continuous(limits = c(0, 800)) +
  labs(x = "Ano de início das atividades",
       y = "Quantidade de empresas") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        axis.text = element_text(size = 9),
        panel.grid.minor.x = element_blank()) 

grafico_linha_inicio_atividades

Plotando sazonalidades ao longo dos meses das EPP

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))

dados_ts1 %>% as_tsibble() %>% gg_season(value, labels="both")+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico sazonal da abertura de empresas - 2014 a 2024")

Podemos destacar que de fevereiro para março em todos os anos houve tendência crescente, menos os anos de 2024 que foi descrescente, o crescimento acentuado e constante de junho, julho e agosto de 2020, o decréscimo acentuado e constante de agosto a outubro de 2016 e cabe destacar também o decréscimo acentuado de julho a agosto de 2024 frente aos demais anos, assim como as empresas de porte MEI.

Plotando médias mensais e observando sazonalidades ao longo dos meses das EPP

dados_ts1 %>% as_tsibble() %>% gg_subseries(value)+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico das médias mensais e sazonalidade da abertura de empresas 
- 2014 a 2024")

Por meio desse gráfico, podemos ver como a abertura das empresas variam ao longo de cada mês do ano. Podemos notar que os meses de maio e julho possuem maiores aberturas de empresas em média, ao passo que novembro e dezembro apresenta a menor média histórica, assim como quando se considera todas os portes de empresas.

Decomposição da série temporal em suas componentes para analisar tendência e sazonalidade das empresas EPP

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))
dados_ts2 = dts2(dados_ts1)
dados_ts2$date = as.Date(dados_ts2$date)

ggdecompose(dados_ts2) +
  xlab("Date")+
  ylab("Quantidade de empresas")+
  scale_x_date(date_breaks = "1 year", date_labels = "%Y-%m")+
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Tendência: Ao decompor a série trimestral da abertura das empresas para as empresas de porte EPP, podemos notar que há alternância entre subidas e descidas de 2014 a 2020, onde tem uma tendência crescente até 2021 e logo depois permanece estável durante a pandemia e em 2023 retoma o seu crescimento.
Sazonalidade: Há uma sazonalidade demarcada muito clara demarcada nos primeiros meses com leve alta e depois baixa nos últimos meses.

2) Análise temporal do fechamento das empresas

Número de empresas matrizes que encerraram as atividades no Maranhão por trimestre de 2019 a 2024

dados_parquet_BRF_filtrado_agregado = dados_parquet_BRF_filtrado %>% 
  filter(MATRIZ == 1 & SIT_CAD == 8) %>% 
  mutate(
    mes_atividade = floor_date(DATA_SIT_CAD, "month"),
    ano_atividade = floor_date(DATA_SIT_CAD, "year")
  ) %>% 
  filter(DATA_SIT_CAD > "2014-01-01",
         DATA_SIT_CAD < "2024-08-10") %>% 
  count(mes_atividade)

grafico_linha_inicio_atividades =  
  
  dados_parquet_BRF_filtrado_agregado %>% 
  #mutate(mes_atividade = ano_atividade) %>% 
  filter(year(mes_atividade) > 2014,
         year(mes_atividade) < 2024) %>% 
  ggplot(aes(x = mes_atividade, y = n)) +
  geom_line(size = 1) +
  theme_minimal(14) +
  scale_x_date(breaks = scales::date_breaks("3 month"),
               labels = scales::date_format("%Y-%m")) +
  scale_y_continuous(labels = scales::label_number(big.mark = "."),breaks = seq(0,35000,5000 ) )+
  ggtitle("Quantidade de empresas que encerraram atividade por mês 
no Maranhão - 2014 a 2024"
        ) + 
  #scale_y_continuous(limits = c(0, 800)) +
  labs(x = "Ano de início das atividades",
       y = "Quantidade de empresas") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        axis.text = element_text(size = 9),
        panel.grid.minor.x = element_blank()) 

grafico_linha_inicio_atividades

Podemos já começar a vislumbrar como está a tendência e a sazonalidade ao longo dos anos por trimestre, mas ainda é importante decompor a série para melhor avaliarmos.

Plotando sazonalidades ao longo dos meses

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))

dados_ts1 %>% as_tsibble() %>% gg_season(value, labels="both")+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico sazonal da abertura de empresas - 2014 a 2024")

Podemos destacar que o pico de alta de fechamento em fevereiro dos anos de 2016,2018,2022 e outro pico de alta em dezembro de 2016 e queda de fechamento de julho a agosto de 2024.

Plotando médias mensais e observando sazonalidades ao longo dos meses

dados_ts1 %>% as_tsibble() %>% gg_subseries(value)+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico das médias mensais e sazonalidade da abertura de empresas 
- 2014 a 2024")

Por meio desse gráfico, podemos ver como o fechamento das empresas variam ao longo de cada mês do ano. Podemos notar que os meses de fevereiro e dezembro possuem maiores fechamentos de empresas em média e os outros meses apresentam uma regularidade nas médias de fechamento ao longo dos anos.

Decomposição da série temporal em suas componentes para analisar tendência e sazonalidade

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))
dados_ts2 = dts2(dados_ts1)
dados_ts2$date = as.Date(dados_ts2$date)

ggdecompose(dados_ts2) +
  xlab("Date")+
  ylab("Quantidade de empresas")+
  scale_x_date(date_breaks = "1 year", date_labels = "%Y-%m")+
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Tendência: Ao decompor a série mensal do encerramento das empresas, podemos notar que há um padrão de alternância entre tendências crescentes, estacionariedades e tendências descrescentes e cabe destacar a estacionaridade duradoura de final de 2018 e início de 2021 e crescimento duradouro de 2022 a 2024.
Sazonalidade: Há uma sazonalidade demarcada muito clara de alta nos primeiros meses de cada ano.

Analisando a decomposição por porte das empresas - MEI (Não houve baixas de empresas matrizes MEI para o recorte temporal estudado)

dados_parquet_BRF_filtrado_agregado = dados_parquet_BRF_filtrado %>% 
  filter(MATRIZ == 1 & SIT_CAD == 8 & PORTE_ROTULADA=="MEI") %>% 
  mutate(
    mes_atividade = floor_date(DATA_SIT_CAD, "month"),
    ano_atividade = floor_date(DATA_SIT_CAD, "year")
  ) %>% 
  filter(DATA_SIT_CAD > "2014-01-01",
         DATA_SIT_CAD < "2024-08-10") %>% 
  count(mes_atividade)

grafico_linha_inicio_atividades =  
  
  dados_parquet_BRF_filtrado_agregado %>% 
  #mutate(mes_atividade = ano_atividade) %>% 
  filter(year(mes_atividade) > 2014,
         year(mes_atividade) < 2024) %>% 
  ggplot(aes(x = mes_atividade, y = n)) +
  geom_line(size = 1) +
  theme_minimal(14) +
  scale_x_date(breaks = scales::date_breaks("3 month"),
               labels = scales::date_format("%Y-%m")) +
  scale_y_continuous(labels = scales::label_number(big.mark = "."),breaks = seq(0,35000,5000 ) )+
  ggtitle("Quantidade de empresas MEI que encerraram atividade por mês 
no Maranhão - 2014 a 2024"
        ) + 
  #scale_y_continuous(limits = c(0, 800)) +
  labs(x = "Ano de encerramento das atividades",
       y = "Quantidade de empresas") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        axis.text = element_text(size = 9),
        panel.grid.minor.x = element_blank()) 

grafico_linha_inicio_atividades

Analisando a decomposição por porte das empresas - ME

dados_parquet_BRF_filtrado_agregado = dados_parquet_BRF_filtrado %>% 
  filter(MATRIZ == 1 & SIT_CAD == 8 & PORTE_ROTULADA=="ME") %>% 
  mutate(
    mes_atividade = floor_date(DATA_SIT_CAD, "month"),
    ano_atividade = floor_date(DATA_SIT_CAD, "year")
  ) %>% 
  filter(DATA_SIT_CAD > "2014-01-01",
         DATA_SIT_CAD < "2024-08-10") %>% 
  count(mes_atividade)

grafico_linha_inicio_atividades =  
  
  dados_parquet_BRF_filtrado_agregado %>% 
  #mutate(mes_atividade = ano_atividade) %>% 
  filter(year(mes_atividade) > 2014,
         year(mes_atividade) < 2024) %>% 
  ggplot(aes(x = mes_atividade, y = n)) +
  geom_line(size = 1) +
  theme_minimal(14) +
  scale_x_date(breaks = scales::date_breaks("3 month"),
               labels = scales::date_format("%Y-%m")) +
  scale_y_continuous(labels = scales::label_number(big.mark = "."),breaks = seq(0,35000,5000 ) )+
  ggtitle("Quantidade de empresas ME que encerraram atividade por mês 
no Maranhão - 2014 a 2024"
        ) + 
  #scale_y_continuous(limits = c(0, 800)) +
  labs(x = "Ano de encerramento das atividades",
       y = "Quantidade de empresas") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        axis.text = element_text(size = 9),
        panel.grid.minor.x = element_blank()) 

grafico_linha_inicio_atividades

Plotando sazonalidades ao longo dos meses das ME

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))

dados_ts1 %>% as_tsibble() %>% gg_season(value, labels="both")+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico sazonal da abertura de empresas - 2014 a 2024")

Podemos destacar que o pico de alta de fechamento em fevereiro dos anos de 2018, 2016, 2024 e queda de fechamento de julho a agosto de 2024.

Plotando médias mensais e observando sazonalidades ao longo dos meses das ME

dados_ts1 %>% as_tsibble() %>% gg_subseries(value)+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico das médias mensais e sazonalidade da abertura de empresas 
- 2014 a 2024")

Por meio desse gráfico, podemos ver como o fechamento das empresas variam ao longo de cada mês do ano. Podemos notar que os meseS de fevereiro possuem maiores fechamentos de empresas em média e os outros meses apresentam uma regularidade nas médias de fechamento ao longo dos anos, assim como quando considerado todos os portes.

Decomposição da série temporal em suas componentes para analisar tendência e sazonalidade das empresas ME

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))
dados_ts2 = dts2(dados_ts1)
dados_ts2$date = as.Date(dados_ts2$date)

ggdecompose(dados_ts2) +
  xlab("Date")+
  ylab("Quantidade de empresas")+
  scale_x_date(date_breaks = "1 year", date_labels = "%Y-%m")+
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Tendência: Ao decompor a série mensal do fechamento das empresas para as empresas de porte ME, podemos notar alternâncias entre altas, baixas e estacionariedades até final de 2018, onde fica estaionário até final de 2020 e depois tem uma tendência crescente até 2024.
Sazonalidade: Há uma sazonalidade demarcada muito clara de alta nos três primeiros meses de cada ano.

Analisando a decomposição por porte das empresas - EPP

dados_parquet_BRF_filtrado_agregado = dados_parquet_BRF_filtrado %>% 
  filter(MATRIZ == 1 & SIT_CAD == 8 & PORTE_ROTULADA=="EPP") %>% 
  mutate(
    mes_atividade = floor_date(DATA_SIT_CAD, "month"),
    ano_atividade = floor_date(DATA_SIT_CAD, "year")
  ) %>% 
  filter(DATA_SIT_CAD > "2014-01-01",
         DATA_SIT_CAD < "2024-08-10") %>% 
  count(mes_atividade)

grafico_linha_inicio_atividades =  
  
  dados_parquet_BRF_filtrado_agregado %>% 
  #mutate(mes_atividade = ano_atividade) %>% 
  filter(year(mes_atividade) > 2014,
         year(mes_atividade) < 2024) %>% 
  ggplot(aes(x = mes_atividade, y = n)) +
  geom_line(size = 1) +
  theme_minimal(14) +
  scale_x_date(breaks = scales::date_breaks("3 month"),
               labels = scales::date_format("%Y-%m")) +
  scale_y_continuous(labels = scales::label_number(big.mark = "."),breaks = seq(0,50,5 ) )+
  ggtitle("Quantidade de empresas EPP que encerraram atividade por mês 
no Maranhão - 2014 a 2024"
        ) + 
  #scale_y_continuous(limits = c(0, 800)) +
  labs(x = "Ano de encerramento das atividades",
       y = "Quantidade de empresas") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        axis.text = element_text(size = 9),
        panel.grid.minor.x = element_blank()) 

grafico_linha_inicio_atividades

Plotando sazonalidades ao longo dos meses das EPP

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))

dados_ts1 %>% as_tsibble() %>% gg_season(value, labels="both")+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico sazonal da abertura de empresas - 2014 a 2024")

Podemos destacar que o pico de alternancia entra alta e queda de fechamento consecutivos nos meses de junho a julho e julho a agosto, respectivamente no ano de 2024 e queda de agosto a setembro 2023, seguido de alta acentuada de setembro a dezembro do mesmo ano.

Plotando médias mensais e observando sazonalidades ao longo dos meses das EPP

dados_ts1 %>% as_tsibble() %>% gg_subseries(value)+
  labs(x="Meses",
       y="Quantidade de empresas",
       title = "Gráfico das médias mensais e sazonalidade da abertura de empresas 
- 2014 a 2024")

Por meio desse gráfico, podemos ver como o fechamento das empresas variam ao longo de cada mês do ano. Podemos notar que os meseS de junho e julho possuem maiores fechamentos de empresas em média e menor média de queda de fechamento em jan, fev e setembro.

Decomposição da série temporal em suas componentes para analisar tendência e sazonalidade das empresas EPP

dados_ts1 = ts(dados_parquet_BRF_filtrado_agregado$n,frequency=12, start=c(2014,1))
dados_ts2 = dts2(dados_ts1)
dados_ts2$date = as.Date(dados_ts2$date)

ggdecompose(dados_ts2) +
  xlab("Date")+
  ylab("Quantidade de empresas")+
  scale_x_date(date_breaks = "1 year", date_labels = "%Y-%m")+
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Tendência: Ao decompor a série mensal do fechamento das empresas para as empresas de porte EPP, o comportamento foi de tendência crescente ao longo de todos os anos.
Sazonalidade: Há uma sazonalidade demarcada muito clara demarcada de altas e quedas ao longo de todos os meses do ano.

3) Análise do perfil das empresas baixadas

Quantidade de empresas baixadas por porte - 2014 a 2024

dados_parquet_BRF_filtrado_porte = dados_parquet_BRF_filtrado %>% 
  filter(MATRIZ==1 & SIT_CAD == 8 & 
         DATA_SIT_CAD > "2014-01-01" &
         DATA_SIT_CAD < "2024-08-10")  %>% count(PORTE_ROTULADA) %>% arrange(desc(n)) 
dados_parquet_BRF_filtrado_porte$PORTE_ROTULADA_ORD = factor(dados_parquet_BRF_filtrado_porte$PORTE_ROTULADA, levels =c("ME","Demais","EPP")
)

ggplot(dados_parquet_BRF_filtrado_porte, aes(x = PORTE_ROTULADA_ORD, y = n)) +
  geom_bar(aes(fill = PORTE_ROTULADA_ORD),stat = "identity", show.legend = F) +
  xlab("Ano") + 
  ylab("Quantidade de empresas") +
  theme(axis.text.y = element_blank())+
  #scale_y_continuous(labels = scales::dollar_format(prefix="R$",big.mark = "."))+
  ggtitle("Quantidade de empresas baixadas por porte - 2014 a 2024") + 
  geom_text(aes(label = prettyNum(paste(round(n, 2),"(",round(n/sum(n)*100, 2),"% )"),big.mark = ".")), vjust=-0.2, color="black") +
  scale_fill_brewer(palette="BrBG")

A grande maior das parte das empresas baixadas são do porte ME (71,67%) e nenhuma delas eram MEI.

Top 5 Quantidade de empresas baixadas por CNAE - 2014 a 2024

dados_parquet_BRF_filtrado_CNAE = dados_parquet_BRF_filtrado %>% 
  filter(MATRIZ==1 & SIT_CAD == 8 & 
         DATA_SIT_CAD > "2014-01-01" &
         DATA_SIT_CAD < "2024-08-10") %>% count(CNAE_PRINC_2DIGITOS_ROTULADA) %>% arrange(desc(n)) %>% head(5)

dados_parquet_BRF_filtrado_CNAE$CNAE_PRINC_2DIGITOS_ROTULADA_ORD = factor(dados_parquet_BRF_filtrado_CNAE$CNAE_PRINC_2DIGITOS_ROTULADA, levels =c("Atividades profissionais, científicas e técnicas","Indústrias de transformação", "Alojamento e alimentação","Outras atividades de serviços","Comércio, reparação de veículos automotores e motocicletas"
  )
)

ggplot(dados_parquet_BRF_filtrado_CNAE, aes(x = CNAE_PRINC_2DIGITOS_ROTULADA_ORD, y = n)) +
  geom_bar(aes(fill = CNAE_PRINC_2DIGITOS_ROTULADA_ORD),stat = "identity", show.legend = F) +
  xlab("Ano") + 
  ylab("Quantidade de empresas") +
  #theme(axis.text.y = element_blank())+
  #scale_y_continuous(labels = scales::dollar_format(prefix="R$",big.mark = "."))+
  ggtitle("Top 5 Quantidade de empresas baixadas 
por CNAE - 2014 a 2024") + 
  geom_text(aes(label = prettyNum(paste(round(n, 2),"\n(",round(n/sum(n)*100, 2),"% )"),big.mark = ".")), vjust=-0.2, color="black",size=2.5) +
  scale_fill_brewer(palette="BrBG")+
  coord_flip()

A maior parte das empresas baixadas são de CNAE de comércio, reparação de veículos automotores e motocicletas (51,23%) e dentre as 5 maiores também se encontram as empresas cujas CNAES são: Outras atividades de serviços, alojamento e alimentação, indústrias de transformação e atividades profissionais, científicas e técnicas.

Quantidade de empresas baixadas por faixas de capital social - 2014 a 2024

minimo=min(dados_parquet_BRF_filtrado$CAPITAL_SOCIAL_NUMERO)

q2=quantile(dados_parquet_BRF_filtrado$CAPITAL_SOCIAL_NUMERO,0.5)

q3=quantile(dados_parquet_BRF_filtrado$CAPITAL_SOCIAL_NUMERO,0.75)

maximo=max(dados_parquet_BRF_filtrado$CAPITAL_SOCIAL_NUMERO)

dados_parquet_BRF_filtrado$capital_social_classes= 
  cut(dados_parquet_BRF_filtrado$CAPITAL_SOCIAL_NUMERO, breaks = c(minimo,q2,q3,maximo),
      labels = c("0 -|1500", "1500-|12000", "12000-|2.000.000.000"),
      include.lowest = T)

grafico_empilhado = ggplot(dados_parquet_BRF_filtrado %>% filter(MATRIZ==1 & SIT_CAD == 8 & 
         DATA_SIT_CAD > "2014-01-01" &
         DATA_SIT_CAD < "2024-08-10"),aes(PORTE_ROTULADA,fill=capital_social_classes)) + 
  geom_bar(position = 'fill') +
  ggtitle("Quantidade de empresas baixadas por faixas de capital social - 2014 a 2024")+
  scale_y_continuous(labels=percent) +
  labs(fill="Faixas de capital social")+
  scale_fill_brewer(palette="BrBG")+
  xlab("Porte das empresas") +
  ylab("Quantidade de empresas (%)")

ggplotly(grafico_empilhado)

Podemos ver que as empresas EPP possuem maior concentração do capital social na faixa de R$12.000,00 a R$2.000.000.000 (81,03%), o que não acontece com as ME, tendo capital social se concentrando mais na faixa de R$1.500,00 a R$12.000,00 (39,72%) e as demais se concentrando de R$0 a R$1.500,00 (96,00%).

4) Análise do perfil das empresas ativas

dados_parquet_BRF_filtrado_porte = dados_parquet_BRF_filtrado %>% 
  filter(MATRIZ==1 & SIT_CAD == 2 & 
         DATA_SIT_CAD > "2014-01-01" &
         DATA_SIT_CAD < "2024-08-10")  %>% count(PORTE_ROTULADA) %>% arrange(desc(n)) 
dados_parquet_BRF_filtrado_porte$PORTE_ROTULADA_ORD = factor(dados_parquet_BRF_filtrado_porte$PORTE_ROTULADA, levels =c("MEI","ME","Demais","EPP")
)

ggplot(dados_parquet_BRF_filtrado_porte, aes(x = PORTE_ROTULADA_ORD, y = n)) +
  geom_bar(aes(fill = PORTE_ROTULADA_ORD),stat = "identity", show.legend = F) +
  xlab("Ano") + 
  ylab("Quantidade de empresas") +
  theme(axis.text.y = element_blank())+
  #scale_y_continuous(labels = scales::dollar_format(prefix="R$",big.mark = "."))+
  ggtitle("Quantidade de empresas ativas por porte - 2014 a 2024") + 
  geom_text(aes(label = prettyNum(paste(round(n, 2),"(",round(n/sum(n)*100, 2),"% )"),big.mark = ".")), vjust=-0.2, color="black") +
  scale_fill_brewer(palette="BrBG")

A grande maior das parte das empresas ativas são do porte MEI (49,2%), em seguida as ME (34,99%).

Top 5 Quantidade de empresas ativas por CNAE - 2014 a 2024

dados_parquet_BRF_filtrado_CNAE = dados_parquet_BRF_filtrado %>% 
  filter(MATRIZ==1 & SIT_CAD == 2 & 
         DATA_SIT_CAD > "2014-01-01" &
         DATA_SIT_CAD < "2024-08-10") %>% count(CNAE_PRINC_2DIGITOS_ROTULADA) %>% arrange(desc(n)) %>% head(5)

dados_parquet_BRF_filtrado_CNAE$CNAE_PRINC_2DIGITOS_ROTULADA_ORD = factor(dados_parquet_BRF_filtrado_CNAE$CNAE_PRINC_2DIGITOS_ROTULADA, levels =c("Indústrias de transformação","Atividades profissionais, científicas e técnicas", "Alojamento e alimentação","Outras atividades de serviços","Comércio, reparação de veículos automotores e motocicletas"
  )
)

ggplot(dados_parquet_BRF_filtrado_CNAE, aes(x = CNAE_PRINC_2DIGITOS_ROTULADA_ORD, y = n)) +
  geom_bar(aes(fill = CNAE_PRINC_2DIGITOS_ROTULADA_ORD),stat = "identity", show.legend = F) +
  xlab("Ano") + 
  ylab("Quantidade de empresas") +
  #theme(axis.text.y = element_blank())+
  #scale_y_continuous(labels = scales::dollar_format(prefix="R$",big.mark = "."))+
  ggtitle("Top 5 Quantidade de empresas ativas 
por CNAE - 2014 a 2024") + 
  geom_text(aes(label = prettyNum(paste(round(n, 2),"\n(",round(n/sum(n)*100, 2),"% )"),big.mark = ".")), vjust=-0.2, color="black",size=2.5) +
  scale_fill_brewer(palette="BrBG")+
  coord_flip()

A maior parte das empresas ativas são de CNAE de comércio, reparação de veículos automotores e motocicletas (56,75%) e dentre as 5 maiores também se encontram as empresas cujas CNAES são: Outras atividades de serviços, alojamento e alimentação, atividades profissionais, científicas e técnicas e indústrias de transformação.

Quantidade de empresas ativas por faixas de capital social - 2014 a 2024

minimo=min(dados_parquet_BRF_filtrado$CAPITAL_SOCIAL_NUMERO)

q2=quantile(dados_parquet_BRF_filtrado$CAPITAL_SOCIAL_NUMERO,0.5)

q3=quantile(dados_parquet_BRF_filtrado$CAPITAL_SOCIAL_NUMERO,0.75)

maximo=max(dados_parquet_BRF_filtrado$CAPITAL_SOCIAL_NUMERO)

dados_parquet_BRF_filtrado$capital_social_classes= 
  cut(dados_parquet_BRF_filtrado$CAPITAL_SOCIAL_NUMERO, breaks = c(minimo,q2,q3,maximo),
      labels = c("0 -|1500", "1500-|12000", "12000-|2.000.000.000"),
      include.lowest = T)

grafico_empilhado = ggplot(dados_parquet_BRF_filtrado %>% filter(MATRIZ==1 & SIT_CAD == 2 & 
         DATA_SIT_CAD > "2014-01-01" &
         DATA_SIT_CAD < "2024-08-10"),aes(PORTE_ROTULADA,fill=capital_social_classes)) + 
  geom_bar(position = 'fill') +
  ggtitle("Quantidade de empresas ativas por faixas de capital social - 2014 a 2024")+
  scale_y_continuous(labels=percent) +
  labs(fill="Faixas de capital social")+
  scale_fill_brewer(palette="BrBG")+
  xlab("Porte das empresas") +
  ylab("Quantidade de empresas (%)")

ggplotly(grafico_empilhado)

Podemos ver que as empresas EPP possuem maior concentração do capital social na faixa de R$12.000,00 a R$2.000.000.000 (86,32%), o que também acontece com as ME (60,34%), jas as MEI tendo capital social se concentrando mais na faixa de R$1.500,00 a R$12.000,00 (60,39%) e as demais se concentrando de R$0 a R$1.500,00 (84,07%).

6) Modelo preditivo de classificação de empresas baixadas - Utilizando 1980 - 2024 (44 anos de base histórica)

Transformando para classe dicotômica

dados_parquet_BRF_filtrado_modelo = dados_parquet_BRF_filtrado %>% filter(MATRIZ==1 & 
         DATA_INICIO_ATV > "1980-01-01" &
         DATA_INICIO_ATV < "2024-08-10") %>% mutate(SIT_CAD_SIM_NAO = case_when( 
  
     SIT_CAD == 8  ~ "Sim",
                                                                                  TRUE ~ "Não"
                                                                  
))

Separando em treino e teste

#selectionando treino e teste 80%,20%
data_split <- initial_split(dados_parquet_BRF_filtrado_modelo, 
                            prop = 0.8
                            )

train_data <- training(data_split)
test_data <- testing(data_split)

Contando os exemplos da classe no treino e teste

#Quantos exemplos de cada classe na variável resposta no treino?
round(prop.table(table(train_data$SIT_CAD_SIM_NAO)),4)*100

## 
##   Não   Sim 
## 54.72 45.28

#Quantos exemplos de cada classe na variável no teste?
round(prop.table(table(test_data$SIT_CAD_SIM_NAO)),4)*100

## 
##   Não   Sim 
## 54.64 45.36

Podemos ver que a nossa variável target não sofre de um problema de desbalanceamento muito grande, o que pode facilitar a nossa abordagem de modelagem de classificação.

# treino
#x_train = subset(train_data,select = -c(SIT_CAD_SIM_NAO,SIT_CAD))
#y_train = train_data$SIT_CAD_SIM_NAO

# teste
#X_test = subset(test_data,select = -c(SIT_CAD_SIM_NAO,SIT_CAD))
#y_test = test_data$SIT_CAD_SIM_NAO

Criando o classificador da árvore de decisão

dt_model = decision_tree(tree_depth = 30) %>%
  set_engine("rpart") %>%
  set_mode("classification") #%>%
#translate()

Selecionando variáveis para compor o modelo

train_data_filter = train_data %>%  select(MATRIZ,
                                #MOT_SIT_CAD,
                                CNAE_PRINC_2DIGITOS,
                                NATUREZA_JURIDICA,
                                PORTE,
                                PUBLICO_SEBRAE,
                                CAPITAL_SOCIAL_NUMERO,
                                SIT_CAD_SIM_NAO
                                )

Criando o modelo na base de treino

dt_fit <- dt_model %>% 
  fit(factor(SIT_CAD_SIM_NAO) ~ .,
      data = train_data_filter)

Resumo do modelo

dt_fit

## parsnip model object
## 
## n= 700523 
## 
## node), split, n, loss, yval, (yprob)
##       * denotes terminal node
## 
##   1) root 700523 317174 Não (0.54723257 0.45276743)  
##     2) PORTE=1,3 169923   2097 Não (0.98765912 0.01234088) *
##     3) PORTE=2,4 530600 215523 Sim (0.40618734 0.59381266)  
##       6) CAPITAL_SOCIAL_NUMERO>=9999.995 162289  46821 Não (0.71149616 0.28850384) *
##       7) CAPITAL_SOCIAL_NUMERO< 9999.995 368311 100055 Sim (0.27165901 0.72834099)  
##        14) NATUREZA_JURIDICA< 4052 291502  97677 Sim (0.33508175 0.66491825)  
##          28) NATUREZA_JURIDICA>=2139 63518  25429 Não (0.59965679 0.40034321) *
##          29) NATUREZA_JURIDICA< 2139 227984  59588 Sim (0.26136922 0.73863078)  
##            58) PORTE=2 160141  53310 Sim (0.33289414 0.66710586)  
##             116) CAPITAL_SOCIAL_NUMERO< 0.86 48936  20515 Não (0.58077898 0.41922102) *
##             117) CAPITAL_SOCIAL_NUMERO>=0.86 111205  24889 Sim (0.22381188 0.77618812)  
##               234) NATUREZA_JURIDICA< 2098.5 6720   1253 Não (0.81354167 0.18645833) *
##               235) NATUREZA_JURIDICA>=2098.5 104485  19422 Sim (0.18588314 0.81411686) *
##            59) PORTE=4 67843   6278 Sim (0.09253718 0.90746282) *
##        15) NATUREZA_JURIDICA>=4052 76809   2378 Sim (0.03095991 0.96904009) *

#%>% tidy()

test_data_filter = test_data %>%  select(MATRIZ,
                                #MOT_SIT_CAD,
                                CNAE_PRINC_2DIGITOS,
                                NATUREZA_JURIDICA,
                                PORTE,
                                PUBLICO_SEBRAE,
                                CAPITAL_SOCIAL_NUMERO,
                                SIT_CAD_SIM_NAO
                                )

Predizer categorias da target

class_preds <- predict(dt_fit, new_data = test_data_filter,
                       type = 'class')

Predizer probabilidades da target

prob_preds <- predict(dt_fit, new_data = test_data_filter, 
                      type = 'prob')

Juntando teste e resultados

model_results <- test_data %>% 
  select(SIT_CAD_SIM_NAO) %>% 
  bind_cols(class_preds, prob_preds)

Visualizando resultados em uma tabela

model_results %>%
  head()

## # A tibble: 6 × 4
##   SIT_CAD_SIM_NAO .pred_class .pred_Não .pred_Sim
##   <chr>           <fct>           <dbl>     <dbl>
## 1 Sim             Não            0.600      0.400
## 2 Sim             Não            0.581      0.419
## 3 Sim             Não            0.581      0.419
## 4 Sim             Sim            0.0925     0.907
## 5 Não             Não            0.600      0.400
## 6 Não             Não            0.600      0.400

transformando em fator a classe

model_results$SIT_CAD_SIM_NAO = factor(model_results$SIT_CAD_SIM_NAO)

Criando a matriz de confusão

conf_mat(model_results, truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class)

##           Truth
## Prediction   Não   Sim
##        Não 88498 24289
##        Sim  7198 55146

calculando a acurácia

yardstick::accuracy(model_results, truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class)

## # A tibble: 1 × 3
##   .metric  .estimator .estimate
##   <chr>    <chr>          <dbl>
## 1 accuracy binary         0.820

calculando a sensibilidade

yardstick::sens(model_results, truth = SIT_CAD_SIM_NAO,
     estimate = .pred_class)

## # A tibble: 1 × 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 sens    binary         0.925

Calculando a especificidade

yardstick::spec(model_results, truth = SIT_CAD_SIM_NAO,
     estimate = .pred_class)

## # A tibble: 1 × 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 spec    binary         0.694

Customizando as métricas

resultados_metricas <- yardstick::metric_set(yardstick::accuracy, yardstick::sens, yardstick::spec)

Calculando metricas de forma customizada

resultados_metricas(model_results, 
                    truth = SIT_CAD_SIM_NAO,
                    estimate = .pred_class)

## # A tibble: 3 × 3
##   .metric  .estimator .estimate
##   <chr>    <chr>          <dbl>
## 1 accuracy binary         0.820
## 2 sens     binary         0.925
## 3 spec     binary         0.694

Criando a Matriz de confusão

conf_mat(model_results,
         truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class) %>% 
  summary()

## # A tibble: 13 × 3
##    .metric              .estimator .estimate
##    <chr>                <chr>          <dbl>
##  1 accuracy             binary         0.820
##  2 kap                  binary         0.631
##  3 sens                 binary         0.925
##  4 spec                 binary         0.694
##  5 ppv                  binary         0.785
##  6 npv                  binary         0.885
##  7 mcc                  binary         0.644
##  8 j_index              binary         0.619
##  9 bal_accuracy         binary         0.810
## 10 detection_prevalence binary         0.644
## 11 precision            binary         0.785
## 12 recall               binary         0.925
## 13 f_meas               binary         0.849

Plotando a matriz de confusão

conf_mat(model_results,
         truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class)  %>% 
  # criando o heat map
  autoplot(type = "heatmap")

Calculando métricas usando os limiares

threshold_df <- model_results %>% 
  roc_curve(truth = SIT_CAD_SIM_NAO, .pred_Não )

Visualizando os limiares

threshold_df %>%
  head()

## # A tibble: 6 × 3
##   .threshold specificity sensitivity
##        <dbl>       <dbl>       <dbl>
## 1  -Inf            0           1    
## 2     0.0310       0           1    
## 3     0.0925       0.234       0.994
## 4     0.186        0.429       0.976
## 5     0.581        0.694       0.925
## 6     0.600        0.760       0.851

Plotando a curva ROC

threshold_df %>% 
  autoplot()

Calculando a área embaixo da curva

roc_auc(model_results, truth = SIT_CAD_SIM_NAO, .pred_Não)

## # A tibble: 1 × 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 roc_auc binary         0.893

Criando o classificador da Floresta aleatória

rf_model = rand_forest(mode = "classification", trees = 100) %>%
  set_engine("ranger") %>% #,seed = 63233 
  set_mode("classification")

Selecionando variáveis para compor o modelo

train_data_filter = train_data %>%  select(MATRIZ,
                                #MOT_SIT_CAD,
                                CNAE_PRINC_2DIGITOS,
                                NATUREZA_JURIDICA,
                                PORTE,
                                PUBLICO_SEBRAE,
                                CAPITAL_SOCIAL_NUMERO,
                                SIT_CAD_SIM_NAO
                                )

Criando o modelo na base de treino

rf_fit <- rf_model %>% 
  fit(factor(SIT_CAD_SIM_NAO) ~ .,
      data = train_data_filter)

Resumo do modelo

#rf_fit %>% tidy()

test_data_filter = test_data %>%  select(MATRIZ,
                                #MOT_SIT_CAD,
                                CNAE_PRINC_2DIGITOS,
                                NATUREZA_JURIDICA,
                                PORTE,
                                PUBLICO_SEBRAE,
                                CAPITAL_SOCIAL_NUMERO,
                                SIT_CAD_SIM_NAO
                                )

Predizer categorias da target

class_preds <- predict(rf_fit, new_data = test_data_filter,
                       type = 'class')

Predizer probabilidades da target

prob_preds <- predict(rf_fit, new_data = test_data_filter, 
                      type = 'prob')

Juntando teste e resultados

model_results <- test_data %>% 
  select(SIT_CAD_SIM_NAO) %>% 
  bind_cols(class_preds, prob_preds)

Visualizando resultados em uma tabela

model_results %>%
  head()

## # A tibble: 6 × 4
##   SIT_CAD_SIM_NAO .pred_class .pred_Não .pred_Sim
##   <chr>           <fct>           <dbl>     <dbl>
## 1 Sim             Sim             0.285     0.715
## 2 Sim             Não             0.550     0.450
## 3 Sim             Não             0.553     0.447
## 4 Sim             Sim             0.162     0.838
## 5 Não             Sim             0.445     0.555
## 6 Não             Não             0.570     0.430

transformando em fator a classe

model_results$SIT_CAD_SIM_NAO = factor(model_results$SIT_CAD_SIM_NAO)

Criando a matriz de confusão

conf_mat(model_results, truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class)

##           Truth
## Prediction   Não   Sim
##        Não 82939 17269
##        Sim 12757 62166

calculando a acurácia

yardstick::accuracy(model_results, truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class)

## # A tibble: 1 × 3
##   .metric  .estimator .estimate
##   <chr>    <chr>          <dbl>
## 1 accuracy binary         0.829

calculando a sensibilidade

yardstick::sens(model_results, truth = SIT_CAD_SIM_NAO,
     estimate = .pred_class)

## # A tibble: 1 × 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 sens    binary         0.867

Calculando a especificidade

yardstick::spec(model_results, truth = SIT_CAD_SIM_NAO,
     estimate = .pred_class)

## # A tibble: 1 × 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 spec    binary         0.783

Customizando as métricas

resultados_metricas <- yardstick::metric_set(yardstick::accuracy, yardstick::sens, yardstick::spec)

Calculando metricas de forma customizada

resultados_metricas(model_results, 
                    truth = SIT_CAD_SIM_NAO,
                    estimate = .pred_class)

## # A tibble: 3 × 3
##   .metric  .estimator .estimate
##   <chr>    <chr>          <dbl>
## 1 accuracy binary         0.829
## 2 sens     binary         0.867
## 3 spec     binary         0.783

Criando a Matriz de confusão

conf_mat(model_results,
         truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class) %>% 
  summary()

## # A tibble: 13 × 3
##    .metric              .estimator .estimate
##    <chr>                <chr>          <dbl>
##  1 accuracy             binary         0.829
##  2 kap                  binary         0.652
##  3 sens                 binary         0.867
##  4 spec                 binary         0.783
##  5 ppv                  binary         0.828
##  6 npv                  binary         0.830
##  7 mcc                  binary         0.653
##  8 j_index              binary         0.649
##  9 bal_accuracy         binary         0.825
## 10 detection_prevalence binary         0.572
## 11 precision            binary         0.828
## 12 recall               binary         0.867
## 13 f_meas               binary         0.847

Plotando a matriz de confusão

conf_mat(model_results,
         truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class)  %>% 
  # criando o heat map
  autoplot(type = "heatmap")

Calculando métricas usando os limiares

threshold_df <- model_results %>% 
  roc_curve(truth = SIT_CAD_SIM_NAO, .pred_Não )

Visualizando os limiares

threshold_df %>%
  head()

## # A tibble: 6 × 3
##   .threshold specificity sensitivity
##        <dbl>       <dbl>       <dbl>
## 1  -Inf        0               1    
## 2     0.0618   0               1    
## 3     0.0653   0.0000252       1    
## 4     0.0753   0.0000378       1    
## 5     0.134    0.233           0.994
## 6     0.135    0.234           0.994

Plotando a curva ROC

threshold_df %>% 
  autoplot()

Calculando a área embaixo da curva

roc_auc(model_results, truth = SIT_CAD_SIM_NAO, .pred_Não)

## # A tibble: 1 × 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 roc_auc binary         0.919

Criando o classificador dos k vizinhos mais próximos (knn)

knn_model = nearest_neighbor(neighbors = 5) %>%
  set_engine("kknn") %>% #,seed = 63233 
  set_mode("classification") #%>%

Selecionando variáveis para compor o modelo

train_data_filter = train_data %>%  select(MATRIZ,
                                #MOT_SIT_CAD,
                                CNAE_PRINC_2DIGITOS,
                                NATUREZA_JURIDICA,
                                PORTE,
                                PUBLICO_SEBRAE,
                                CAPITAL_SOCIAL_NUMERO,
                                SIT_CAD_SIM_NAO
                                )

Criando o modelo na base de treino

knn_fit <- knn_model %>% 
  fit(factor(SIT_CAD_SIM_NAO) ~ .,
      data = train_data_filter)

Resumo do modelo

#rf_fit %>% tidy()

test_data_filter = test_data %>%  select(MATRIZ,
                                #MOT_SIT_CAD,
                                CNAE_PRINC_2DIGITOS,
                                NATUREZA_JURIDICA,
                                PORTE,
                                PUBLICO_SEBRAE,
                                CAPITAL_SOCIAL_NUMERO,
                                SIT_CAD_SIM_NAO
                                )

Predizer categorias da target

class_preds <- predict(knn_fit, new_data = test_data_filter,
                       type = 'class')

Predizer probabilidades da target

prob_preds <- predict(knn_fit, new_data = test_data_filter, 
                      type = 'prob')

Juntando teste e resultados

model_results <- test_data %>% 
  select(SIT_CAD_SIM_NAO) %>% 
  bind_cols(class_preds, prob_preds)

Visualizando resultados em uma tabela

model_results %>%
  head()

## # A tibble: 6 × 4
##   SIT_CAD_SIM_NAO .pred_class .pred_Não .pred_Sim
##   <chr>           <fct>           <dbl>     <dbl>
## 1 Sim             Sim             0        1     
## 2 Sim             Sim             0.295    0.705 
## 3 Sim             Sim             0.357    0.643 
## 4 Sim             Sim             0        1     
## 5 Não             Sim             0.118    0.882 
## 6 Não             Não             0.972    0.0280

transformando em fator a classe

model_results$SIT_CAD_SIM_NAO = factor(model_results$SIT_CAD_SIM_NAO)

Criando a matriz de confusão

conf_mat(model_results, truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class)

##           Truth
## Prediction   Não   Sim
##        Não 77948 14974
##        Sim 17748 64461

calculando a acurácia

yardstick::accuracy(model_results, truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class)

## # A tibble: 1 × 3
##   .metric  .estimator .estimate
##   <chr>    <chr>          <dbl>
## 1 accuracy binary         0.813

calculando a sensibilidade

yardstick::sens(model_results, truth = SIT_CAD_SIM_NAO,
     estimate = .pred_class)

## # A tibble: 1 × 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 sens    binary         0.815

Calculando a especificidade

yardstick::spec(model_results, truth = SIT_CAD_SIM_NAO,
     estimate = .pred_class)

## # A tibble: 1 × 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 spec    binary         0.811

Customizando as métricas

resultados_metricas <- yardstick::metric_set(yardstick::accuracy, yardstick::sens, yardstick::spec)

Calculando metricas de forma customizada

resultados_metricas(model_results, 
                    truth = SIT_CAD_SIM_NAO,
                    estimate = .pred_class)

## # A tibble: 3 × 3
##   .metric  .estimator .estimate
##   <chr>    <chr>          <dbl>
## 1 accuracy binary         0.813
## 2 sens     binary         0.815
## 3 spec     binary         0.811

Criando a Matriz de confusão

conf_mat(model_results,
         truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class) %>% 
  summary()

## # A tibble: 13 × 3
##    .metric              .estimator .estimate
##    <chr>                <chr>          <dbl>
##  1 accuracy             binary         0.813
##  2 kap                  binary         0.624
##  3 sens                 binary         0.815
##  4 spec                 binary         0.811
##  5 ppv                  binary         0.839
##  6 npv                  binary         0.784
##  7 mcc                  binary         0.624
##  8 j_index              binary         0.626
##  9 bal_accuracy         binary         0.813
## 10 detection_prevalence binary         0.531
## 11 precision            binary         0.839
## 12 recall               binary         0.815
## 13 f_meas               binary         0.827

Plotando a matriz de confusão

conf_mat(model_results,
         truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class)  %>% 
  # criando o heat map
  autoplot(type = "heatmap")

Calculando métricas usando os limiares

threshold_df <- model_results %>% 
  roc_curve(truth = SIT_CAD_SIM_NAO, .pred_Não )

Visualizando os limiares

threshold_df %>%
  head()

## # A tibble: 6 × 3
##   .threshold specificity sensitivity
##        <dbl>       <dbl>       <dbl>
## 1  -Inf            0           1    
## 2     0            0           1    
## 3     0.0280       0.237       0.974
## 4     0.0902       0.257       0.970
## 5     0.118        0.292       0.966
## 6     0.166        0.363       0.906

Plotando a curva ROC

threshold_df %>% 
  autoplot()

Calculando a área embaixo da curva

roc_auc(model_results, truth = SIT_CAD_SIM_NAO, .pred_Não)

## # A tibble: 1 × 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 roc_auc binary         0.861

Criando o classificador dos Gradient boosting

gb_model = boost_tree(mode = "classification", trees = 100) %>%
  set_engine("xgboost") %>% #,seed = 63233 
  set_mode("classification") #%>%

Selecionando variáveis para compor o modelo

train_data_filter = train_data %>%  select(MATRIZ,
                                #MOT_SIT_CAD,
                                CNAE_PRINC_2DIGITOS,
                                NATUREZA_JURIDICA,
                                PORTE,
                                PUBLICO_SEBRAE,
                                CAPITAL_SOCIAL_NUMERO,
                                SIT_CAD_SIM_NAO
                                )

Criando o modelo na base de treino

gb_fit <- gb_model %>% 
  fit(factor(SIT_CAD_SIM_NAO) ~ .,
      data = train_data_filter)

Resumo do modelo

#rf_fit %>% tidy()

test_data_filter = test_data %>%  select(MATRIZ,
                                #MOT_SIT_CAD,
                                CNAE_PRINC_2DIGITOS,
                                NATUREZA_JURIDICA,
                                PORTE,
                                PUBLICO_SEBRAE,
                                CAPITAL_SOCIAL_NUMERO,
                                SIT_CAD_SIM_NAO
                                )

Predizer categorias da target

class_preds <- predict(gb_fit, new_data = test_data_filter,
                       type = 'class')

Predizer probabilidades da target

prob_preds <- predict(gb_fit, new_data = test_data_filter, 
                      type = 'prob')

Juntando teste e resultados

model_results <- test_data %>% 
  select(SIT_CAD_SIM_NAO) %>% 
  bind_cols(class_preds, prob_preds)

Visualizando resultados em uma tabela

model_results %>%
  head()

## # A tibble: 6 × 4
##   SIT_CAD_SIM_NAO .pred_class .pred_Não .pred_Sim
##   <chr>           <fct>           <dbl>     <dbl>
## 1 Sim             Sim            0.150      0.850
## 2 Sim             Não            0.589      0.411
## 3 Sim             Não            0.618      0.382
## 4 Sim             Sim            0.0154     0.985
## 5 Não             Não            0.507      0.493
## 6 Não             Não            0.745      0.255

transformando em fator a classe

model_results$SIT_CAD_SIM_NAO = factor(model_results$SIT_CAD_SIM_NAO)

Criando a matriz de confusão

conf_mat(model_results, truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class)

##           Truth
## Prediction   Não   Sim
##        Não 87589 20697
##        Sim  8107 58738

calculando a acurácia

yardstick::accuracy(model_results, truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class)

## # A tibble: 1 × 3
##   .metric  .estimator .estimate
##   <chr>    <chr>          <dbl>
## 1 accuracy binary         0.836

calculando a sensibilidade

yardstick::sens(model_results, truth = SIT_CAD_SIM_NAO,
     estimate = .pred_class)

## # A tibble: 1 × 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 sens    binary         0.915

Calculando a especificidade

yardstick::spec(model_results, truth = SIT_CAD_SIM_NAO,
     estimate = .pred_class)

## # A tibble: 1 × 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 spec    binary         0.739

Customizando as métricas

resultados_metricas <- yardstick::metric_set(yardstick::accuracy, yardstick::sens, yardstick::spec)

Calculando metricas de forma customizada

resultados_metricas(model_results, 
                    truth = SIT_CAD_SIM_NAO,
                    estimate = .pred_class)

## # A tibble: 3 × 3
##   .metric  .estimator .estimate
##   <chr>    <chr>          <dbl>
## 1 accuracy binary         0.836
## 2 sens     binary         0.915
## 3 spec     binary         0.739

Criando a Matriz de confusão

conf_mat(model_results,
         truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class) %>% 
  summary()

## # A tibble: 13 × 3
##    .metric              .estimator .estimate
##    <chr>                <chr>          <dbl>
##  1 accuracy             binary         0.836
##  2 kap                  binary         0.664
##  3 sens                 binary         0.915
##  4 spec                 binary         0.739
##  5 ppv                  binary         0.809
##  6 npv                  binary         0.879
##  7 mcc                  binary         0.671
##  8 j_index              binary         0.655
##  9 bal_accuracy         binary         0.827
## 10 detection_prevalence binary         0.618
## 11 precision            binary         0.809
## 12 recall               binary         0.915
## 13 f_meas               binary         0.859

Plotando a matriz de confusão

conf_mat(model_results,
         truth = SIT_CAD_SIM_NAO,
         estimate = .pred_class)  %>% 
  # criando o heat map
  autoplot(type = "heatmap")

Calculando métricas usando os limiares

threshold_df <- model_results %>% 
  roc_curve(truth = SIT_CAD_SIM_NAO, .pred_Não )

Visualizando os limiares

threshold_df %>%
  head()

## # A tibble: 6 × 3
##    .threshold specificity sensitivity
##         <dbl>       <dbl>       <dbl>
## 1 -Inf          0                   1
## 2    0.000197   0                   1
## 3    0.000265   0.0000126           1
## 4    0.000329   0.0000252           1
## 5    0.000426   0.0000378           1
## 6    0.000486   0.0000881           1

Plotando a curva ROC

threshold_df %>% 
  autoplot()

Calculando a área embaixo da curva

roc_auc(model_results, truth = SIT_CAD_SIM_NAO, .pred_Não)

## # A tibble: 1 × 3
##   .metric .estimator .estimate
##   <chr>   <chr>          <dbl>
## 1 roc_auc binary         0.926

Tabela de comparação de métricas dos modelos de classificação

tabela_metricas = data.frame( Modelo= c("Árvore de decisão","Floresta Aleatória","K vizinhos mais próximos","Gradient boosting"),
Acuracia = c("0,823","0,831","0,820","0,838"), F1=c("0,851","0,849","0,843","0,861"),Precisão =c("0,788","0,830","0,809","0,812"),Revocação=c("0,926","0,869","0,880","0,916"),Espeficidade=c("0,697","0,784","0,784","0,743"),Roc_Auc=c("0,895","0,920","0,894","0,927"))

kableExtra::kable(tabela_metricas,align = "lccrr",caption = "Tabela de comparação de métricas dos modelos de classificação")

Tabela de comparação de métricas dos modelos de classificação
Modelo	Acuracia	F1	Precisão	Revocação	Espeficidade	Roc_Auc
Árvore de decisão	0,823	0,851	0,788	0,926	0,697	0,895
Floresta Aleatória	0,831	0,849	0,830	0,869	0,784	0,920
K vizinhos mais próximos	0,820	0,843	0,809	0,880	0,784	0,894
Gradient boosting	0,838	0,861	0,812	0,916	0,743	0,927

Analisando os modelos e suas respectivas métricas, cabe destacar que o modelo que obteve a melhor perfomance geral foi o Gradient Boosting (GB), sendo superior nas métricas:
F1 Score GB (86,1%);
Acurácia GB (83,8%);
Roc AUC GB (92,7%);
F1-Score : É a média harmônica da precisão e do recall. É útil quando se deseja um bom equilíbrio entre precisão e recall.
Interpretação: Um F1-score de 86,1% indica um bom equilíbrio entre precisão e recall, sugerindo que o modelo tem um bom desempenho geral.
Acurácia : Indica a proporção total de previsões corretas. Em outras palavras, o modelo acertou 82,19% das vezes.
Interpretação: É um bom indicador geral de desempenho, mas pode ser enganoso em casos de desbalanceamento de classes, que não é o caso da análise.
Roc: A curva ROC pode fornecer uma visão mais completa do desempenho do modelo, especialmente em cenários com diferentes limiares de classificação. Ela plota a taxa de verdadeiros positivos (sensibilidade) contra a taxa de falsos positivos (1 - especificidade) em vários pontos de corte.
Área Sob a Curva ROC: A área sob a curva ROC é um único valor que resume o desempenho global do modelo. Um AUC de 1 indica um classificador perfeito, enquanto um AUC de 0.5 indica um classificador aleatório, o nosso foi 0,927, indicando um desempenho bom do modelo.

Observações: O modelo perdeu nos seguintes quesitos:

Precisão (Perdeu para floresta aleatória)
Especificidade (Perdeu para floresta aleatória e k vizinhos mais próximos)
Revocação (Perdeu para Árvore de decisão)
Precisão: Indica a proporção de previsões positivas que são realmente positivas, ou seja quantas empresas que o modelo previu como baixadas e que realmente foram baixadas e o mesmo estava correto em 81,2% das vezes. É mais importante quando o custo de um falso positivo é alto.
Especificidade: : Indica a proporção de exemplos negativos que foram corretamente identificados, ou seja quantas empresas que de fato não foram baixadas e o modelo previu como não baixadas e estava correto em 74,3%.
Revocação: Indicam a proporção de exemplos positivos que foram corretamente identificados, ou seja quantas empresas que de fato foram baixadas e o modelo previu que foram baixadas e o mesmo estava correto em 91,6% das vezes, mais importante quando o custo de um falso negativo é alto.

Considerações finais

Com o objetivo de compreender os ciclos de vida das empresas maranhenses e aprimorar a previsão de encerramentos, analisamos uma década de dados da Receita Federal. Identificamos padrões sazonais nas aberturas e fechamentos, caracterizamos os diferentes perfis empresariais e desenvolvemos um modelo de classificação com alta precisão para prever quais empresas encerrariam suas atividades. Os resultados indicam que o modelo proposto é capaz de identificar corretamente 82,19% dos casos, com um F1-score de 86,1% e uma AUC de 92,7%.