No Hackeando o R de hoje, vamos falar sobre o pacote tidymodels, que tem como objetivo inserir dentro do tidyverse a modelagem estatística. É muito comum aprendermos a fazer a extração e a manipulação de dados utilizando o tidyverse, e então rodar modelos com funções de outras fontes, para voltar ao tidyverse para fazer a visualização dos resultados utilizando o ggplot2. Isso é de certa forma ineficiente, pois precisamos cuidar para que as funções sejam compatíveis e que a própria escrita do código seja legível para outras pessoas. Ao invés de se preocupar com isso, vamos mostrar como fazer os modelos diretamente no tidyverse, podendo assim maanter a consistência de código e facilitando a integração.
Primeiramente, vamos carregar alguns pacotes utilizados:
library(tidymodels) library(nycflights13)
Iremos trabalhar com os dados de voos do pacote nycflights13data. Antes de introduzirmos elementos específicos da modelagem, vamos manipular os dados:
set.seed(123) dados = flights %>% mutate( arr_delay = ifelse(arr_delay >= 30, "atraso", "sem_atraso"), arr_delay = factor(arr_delay), date = lubridate::as_date(time_hour) ) %>% inner_join(weather, by = c("origin", "time_hour")) %>% select(dep_time, flight, origin, dest, air_time, distance, carrier, date, arr_delay, time_hour) %>% na.omit() %>% mutate(across(where(is.character), as.factor))
A ideia de utilizar o tidymodels é implementar um workflow. Essencialmente, criamos componentes que são inteiramente integrados no tidyverse, e os juntamos dentro desse tipo de objeto, de modo que qualquer mudança nos dados, na metodologia de estimação, ou no modelo em si, é facilitada pelo fato de que conseguimos observar claramente cada módulo do modelo. No post de hoje, trabalharemos com um modelo simples, logo teremos duas partes: uma recipe e um model.
Para começar, vamos criar a recipe. Esse objeto engloba todo o pré-processamento dos dados, e algumas outras funções. Inicialmente, adicionamos a ela o modelo em si (em notação de modelo) e os dados que iremos utilizar. Após isso, nós podemos utilizar a função update_role() para dar um papel diferente a algumas variáveis. No caso, o modelo inicial considera todas as colunas como previsores de arr_delay, porém queremos que flight e time_hour não sejam previsores, e sim apenas indicadores de cada observação.
receita_flights = recipe(arr_delay ~ ., data = treino) %>% update_role(flight, time_hour, new_role = "ID")
A partir daí, podemos fazer transformações facilmente nos dados que melhoram nosso dataset. Vamos adicionar colunas de dia da semana e mês dos voos com a função step_date(), e adicionar uma coluna para feriados nos EUA.
receita_flights = receita_flights %>% step_date(date, features = c("dow", "month")) %>% step_holiday(date, holidays = timeDate::listHolidays("US"), keep_original_cols = FALSE) %>% step_dummy(all_nominal_predictors()) %>% step_zv(all_predictors())
Note que há uma outra função acima, chamada de step_dummy(). Essencialmente, ela transforma variáveis categóricas em dummies, para podermos fazer a regressão. Em outros métodos o R faz sozinho essa transformação, porém aqui devemos fazer ela manualmente pois há modelos que não dependem de valores numéricos, logo é melhor manter a estrutura original dos dados caso necessário. Para finalizar, adicionamos a step_zv, que elimina valores que ocorrem uma única vez (que podem causar problemas pois ou só estão no treino, ou só estão no teste).
Após criar a recipe, gerar um modelo é muito simples. Iniciamos com o tipo de regressão que queremos, e adicionamos a engine que queremos usar para calculá-lo. Abaixo, geramos uma regressão logística usando o glm:
modelo = logistic_reg() %>% set_engine("glm")
Após isso, basta gerar nosso workflow. Abaixo, criamos ele, adicionamos a recipe, o modelo, e então treinamos ele com dados de treino. Finalizamos extraindo o resultado com a extract_fit_parsnip().
workflow() %>% add_model(modelo) %>% add_recipe(receita_flights) %>% fit(data = treino) %>% extract_fit_parsnip() %>% tidy() term estimate std.error statistic p.value <chr> <dbl> <dbl> <dbl> <dbl> 1 (Intercept) 5.99 2.74 2.19 2.88e- 2 2 dep_time -0.00167 0.0000141 -118. 0 3 air_time -0.0441 0.000562 -78.5 0 4 distance 0.00565 0.00151 3.75 1.80e- 4 5 date_USChristmasDay 1.24 0.168 7.36 1.78e-13 6 date_USColumbusDay 0.801 0.179 4.48 7.40e- 6 7 date_USCPulaskisBirthday 0.836 0.140 5.97 2.31e- 9 8 date_USDecorationMemorialDay 0.436 0.113 3.87 1.08e- 4 9 date_USElectionDay 0.792 0.175 4.52 6.15e- 6 10 date_USGoodFriday 1.36 0.173 7.85 4.30e-15
________________________
(*) Para entender mais sobre a linguagem R e suas ferramentas, confira nosso Curso de Introdução ao R para análise de dados.