Relatório de Modelo: Prophet

Forecasting Robusto com Detecção de Change Points

Authors

UIVS - Unidade de Inteligência em Vigilância em Saúde

Caio Sain Vallio

Published

11/01/2026

1. Introdução

1.1 Objetivo

Este relatório avalia o modelo Prophet (Meta/Facebook) para previsão de incidência de dengue. O Prophet é especialmente projetado para ser robusto a outliers, dados faltantes e mudanças de tendência.

1.2 Metodologia

Modelo: Prophet com sazonalidade multiplicativa
Estratégia: Rolling Replication (Janela deslizante)
Horizontes de Previsão: 4, 6 e 8 semanas
Comparativo: SNaive, SARIMA, SARIMAX, TBATS
Métricas: RMSE, MAE, MAPE, RMSLE, MASE

1.3 Sobre o Prophet

Prophet é um modelo de decomposição aditiva/multiplicativa:

\[y(t) = g(t) + s(t) + h(t) + \epsilon_t\]

Componentes do Prophet

\(g(t)\): Tendência (linear ou logística com changepoints)
\(s(t)\): Sazonalidade (Fourier series)
\(h(t)\): Efeitos de feriados/eventos
\(\epsilon_t\): Ruído

1.3.1 Vantagens para Séries Epidemiológicas

Robustez a Outliers: Surtos extremos não distorcem tanto as previsões
Change Points Automáticos: Detecta mudanças de regime na série
Sazonalidade Multiplicativa: Amplitude varia com o nível da série
Interpretabilidade: Componentes visualizáveis separadamente

2. Dados e Preparação

Code

df_raw <- load_raw_data()

#> [2026-01-11 15:19:56] INFO: Carregando dados de: /Users/caiosainvallio/ses-sp/forecast_dengue/data/raw/dengue.RData 
#> [2026-01-11 15:19:57] INFO: Dados carregados: 248865 linhas, 23 colunas

Code

df_state <- aggregate_state(df_raw)

#> [2026-01-11 15:19:57] INFO: Agregando dados por estado... 
#> [2026-01-11 15:19:57] INFO: Dados agregados: 626 semanas

Code

df <- preprocess_data(df_state)

#> [2026-01-11 15:19:57] INFO: ========== INICIO: Preprocessamento ========== 
#> [2026-01-11 15:19:57] INFO: Removidas 3 linhas da semana 53 
#> [2026-01-11 15:19:57] WARN: Semanas faltantes detectadas: 2 
#> [2026-01-11 15:19:57] WARN: ATENCAO: Alvo nao sera imputado (sera mantido como NA) 
#> [2026-01-11 15:19:57] INFO: Imputados 8 NAs em mean_temp (fill) 
#> [2026-01-11 15:19:57] INFO: Imputados 8 NAs em mean_max_temp (fill) 
#> [2026-01-11 15:19:57] INFO: Imputados 8 NAs em mean_min_temp (fill) 
#> [2026-01-11 15:19:57] INFO: Imputados 8 NAs em mean_precip (fill) 
#> [2026-01-11 15:19:57] INFO: Dados preprocessados: 625 linhas 
#> [2026-01-11 15:19:57] INFO: ========== FIM: Preprocessamento ==========

Code

df_base <- df |>
    mutate(data_iniSE = as.Date(data_iniSE)) |>
    arrange(data_iniSE)

# Carregar config de features
config_features <- tryCatch(load_config("features"), error = function(e) list())
if (length(config_features) == 0) {
    config_features <- list(temporal = list(enabled = TRUE), seasonal = list(enabled = TRUE))
}

cat("Dados carregados:", nrow(df_base), "observações\n")

#> Dados carregados: 625 observações

Code

cat("Período:", min(df_base$data_iniSE), "a", max(df_base$data_iniSE), "\n")

#> Período: 16075 a 20443

Code

glimpse(df_base |> select(data_iniSE, est_inc100k, semana, ano))

#> Rows: 625
#> Columns: 4
#> $ data_iniSE  <date> 2014-01-05, 2014-01-12, 2014-01-19, 2014-01-26, 2014-02-0…
#> $ est_inc100k <dbl> 4.688448, 5.657748, 6.589646, 8.467320, 10.189317, 14.3069…
#> $ semana      <dbl> 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18…
#> $ ano         <dbl> 2014, 2014, 2014, 2014, 2014, 2014, 2014, 2014, 2014, 2014…

3. Definição do Modelo

Code

model_name <- "prophet"
model <- get_model(model_name)

cat(sprintf(
    "Modelo: %s\nFamília: %s\nDescrição: %s\n",
    model$name, model$family, model$description
))

#> Modelo: prophet
#> Família: statistical
#> Descrição: Prophet - Modelo de forecasting robusto do Meta

Code

cat("\n### Configuração:\n")

#> 
#> ### Configuração:

Code

cat("- Sazonalidade anual: TRUE\n")

#> - Sazonalidade anual: TRUE

Code

cat("- Modo de sazonalidade: Multiplicativa\n")

#> - Modo de sazonalidade: Multiplicativa

Code

cat("- Changepoint prior scale: 0.05\n")

#> - Changepoint prior scale: 0.05

Code

cat("- Seasonality prior scale: 0.1 (prior forte → sazonalidade rígida)\n")

#> - Seasonality prior scale: 0.1 (prior forte → sazonalidade rígida)

4. Backtesting

Code

df_features <- make_features_for_prediction(df_base, config_features)

#> [2026-01-11 15:19:57] INFO: ========== INICIO: Feature Engineering ========== 
#> [2026-01-11 15:19:57] INFO: Adicionando features temporais... 
#> [2026-01-11 15:19:57] INFO: Adicionando features de sazonalidade... 
#> [2026-01-11 15:19:57] INFO: Adicionando features climaticas... 
#> [2026-01-11 15:19:57] INFO: Features criadas: 47 novas colunas 
#> [2026-01-11 15:19:57] INFO: ========== FIM: Feature Engineering ==========

Code

# Selecao de regressores (Climaticos)
selected_xregs <- c("temperature_lag4", "precipitation_lag4")

# Regressors disponiveis?
available_xregs <- intersect(selected_xregs, names(df_features))
if (length(available_xregs) < length(selected_xregs)) {
    warning("Alguns regressores nao encontrados: ", paste(setdiff(selected_xregs, available_xregs), collapse = ", "))
}
selected_xregs <- available_xregs

backtest_config <- list(
    backtest = list(
        horizons = c(4, 6, 8),
        initial_window = 52 * 5,
        step = 4
    ),
    data = list(
        date_col = "data_iniSE",
        target = "est_inc100k"
    )
)

model_config <- list(
    target = "est_inc100k",
    date_col = "data_iniSE",
    yearly.seasonality = TRUE,
    weekly.seasonality = FALSE,
    yearly.seasonality = TRUE,
    weekly.seasonality = FALSE,
    yearly.seasonality = TRUE,
    weekly.seasonality = FALSE,
    seasonality.mode = "multiplicative",
    changepoint.prior.scale = 0.5, # Retornando a 0.5 conforme pedido
    seasonality.prior.scale = 10.0, # Retornando a 10.0 conforme pedido
    xreg_cols = selected_xregs # Regressores externos
)

has_parallel <- requireNamespace("furrr", quietly = TRUE) &&
    requireNamespace("future", quietly = TRUE)

if (has_parallel) {
    cat("✅ Usando execução PARALELA\n")
    bt_result <- run_backtest_parallel(
        model_name = model_name,
        data_base = df_features, # Usar df com features
        config = backtest_config,
        model_config = model_config,
        verbose = TRUE
    )
} else {
    cat("⚠️ Usando execução SEQUENCIAL\n")
    bt_result <- run_backtest(
        model_name = model_name,
        data = df_features, # Usar df com features
        config = backtest_config,
        model_config = model_config,
        verbose = TRUE
    )
}

#> ✅ Usando execução PARALELA

#> [2026-01-11 15:19:57] INFO: ========== INICIO: Backtest PARALLEL: prophet ========== 
#> [2026-01-11 15:19:57] INFO: Origens: 90 | Horizontes: 4, 6, 8 | Workers: 15 
#> [2026-01-11 15:19:57] INFO: Modo: Features recalculadas por fold (sem leakage) 
#> [2026-01-11 15:20:19] INFO: Previsoes geradas: 270 
#> [2026-01-11 15:20:19] INFO: Erros: 0 
#> [2026-01-11 15:20:19] INFO: Tempo total: 21.6 segundos 
#> [2026-01-11 15:20:19] INFO: ========== FIM: Backtest PARALLEL: prophet (21.62s) ==========

Code

cat(sprintf("⏱️ Tempo total: %.1f segundos\n", bt_result$duration))

#> ⏱️ Tempo total: 21.6 segundos

Code

print(bt_result)

#> 
#> === Resultado de Backtest ===
#> 
#> Modelo: prophet 
#> Origens: 90 
#> Previsoes: 270 
#> Duracao: 21.6 s
#> 
#> Metricas por Horizonte:
#>  h      mae     rmse     mape    smape     mase    rmsle  n
#>  4 47.54691 117.4943 77.05757 84.85714 1.403894 1.177860 90
#>  6 50.35713 118.8359 83.75051 90.06271 1.486869 1.270754 90
#>  8 58.25221 136.7107 98.33974 96.17499 1.719984 1.404034 90

5. Avaliação de Desempenho

5.1 Métricas por Horizonte

Code

# Debug: verificar estrutura do resultado
cat("=== DEBUG ===\n")

#> === DEBUG ===

Code

cat("Predictions NULL?", is.null(bt_result$predictions), "\n")

#> Predictions NULL? FALSE

Code

cat("Metrics NULL?", is.null(bt_result$metrics), "\n")

#> Metrics NULL? FALSE

Code

cat("Número de erros:", length(bt_result$errors), "\n")

#> Número de erros: 0

Code

if (length(bt_result$errors) > 0) {
    cat("\nPrimeiros 5 erros:\n")
    for (i in seq_len(min(5, length(bt_result$errors)))) {
        cat(sprintf("  Origin %d: %s\n", bt_result$errors[[i]]$origin, bt_result$errors[[i]]$error))
    }
}

# Verificar se há metrics
if (is.null(bt_result$metrics) || nrow(bt_result$metrics) == 0) {
    stop("Métricas não foram geradas! Verifique os erros acima.")
}

bt_result$metrics |>
    mutate(
        across(where(is.numeric), \(x) round(x, 4)),
        h_desc = paste(h, "semanas")
    ) |>
    select(Horizonte = h_desc, RMSE = rmse, MAE = mae, MAPE = mape, RMSLE = rmsle, MASE = mase) |>
    kable(caption = "Métricas de Performance do Prophet") |>
    kable_styling(bootstrap_options = c("striped", "hover"))

Métricas de Performance do Prophet
Horizonte	RMSE	MAE	MAPE	RMSLE	MASE
4 semanas	117.4943	47.5469	77.0576	1.1779	1.4039
6 semanas	118.8359	50.3571	83.7505	1.2708	1.4869
8 semanas	136.7107	58.2522	98.3397	1.4040	1.7200

5.2 Comparação com Modelos Anteriores

Code

metrics_path <- file.path(PROJECT_ROOT, "data", "model_metrics.RData")

if (file.exists(metrics_path)) {
    load(metrics_path)

    comparison_data <- model_metrics |>
        group_by(model_name) |>
        slice(which.max(execution_date)) |>
        ungroup() |>
        select(model_name, horizon_4w_rmse, horizon_4w_mae, horizon_4w_mape, horizon_4w_rmsle)

    prophet_row <- data.frame(
        model_name = "prophet",
        horizon_4w_rmse = bt_result$metrics[bt_result$metrics$h == 4, "rmse"],
        horizon_4w_mae = bt_result$metrics[bt_result$metrics$h == 4, "mae"],
        horizon_4w_mape = bt_result$metrics[bt_result$metrics$h == 4, "mape"],
        horizon_4w_rmsle = bt_result$metrics[bt_result$metrics$h == 4, "rmsle"]
    )

    all_models <- comparison_data |>
        filter(model_name != "prophet") |>
        dplyr::bind_rows(prophet_row) |>
        arrange(horizon_4w_rmse)

    all_models |>
        mutate(across(where(is.numeric), \(x) round(x, 4))) |>
        rename(
            Modelo = model_name,
            RMSE = horizon_4w_rmse,
            MAE = horizon_4w_mae,
            MAPE = horizon_4w_mape,
            RMSLE = horizon_4w_rmsle
        ) |>
        kable(caption = "Comparativo de Modelos - Horizonte 4 semanas") |>
        kable_styling(bootstrap_options = c("striped", "hover")) |>
        row_spec(which(all_models$model_name == "prophet"), bold = TRUE, background = "#E8F4F8")
} else {
    cat("Arquivo de métricas anteriores não encontrado.")
}

Comparativo de Modelos - Horizonte 4 semanas
Modelo	RMSE	MAE	MAPE	RMSLE
tbats	20.4688	10.3990	23.2937	0.3153
sarima	36.3290	14.2967	25.1256	0.2930
sarimax	90.6956	34.7848	43.8269	0.7821
seasonal_mean	91.2846	39.0308	65.0013	1.0641
seasonal_naive	96.3101	42.7400	66.1433	0.9683
prophet	117.4943	47.5469	77.0576	1.1779

5.3 Análise Visual Multi-Horizonte

Code

preds_multi <- bind_rows(
    extract_predictions(bt_result, horizon = 4) |> mutate(horizon = "4 semanas"),
    extract_predictions(bt_result, horizon = 6) |> mutate(horizon = "6 semanas"),
    extract_predictions(bt_result, horizon = 8) |> mutate(horizon = "8 semanas")
)

p_backtest <- ggplot() +
    geom_line(data = preds_multi, aes(x = target_date, y = actual), color = "black", linewidth = 0.5) +
    geom_line(data = preds_multi, aes(x = target_date, y = predicted, color = horizon), linewidth = 0.5, alpha = 0.8) +
    scale_color_manual(values = c("4 semanas" = "#E74C3C", "6 semanas" = "#3498DB", "8 semanas" = "#2ECC71")) +
    labs(
        title = "Backtest Prophet: Previsões Multi-Horizonte",
        subtitle = "Sazonalidade multiplicativa com detecção de change points",
        x = "Data", y = "Incidência / 100k hab"
    ) +
    theme_minimal()

ggplotly(p_backtest) |>
    layout(legend = list(orientation = "h", x = 0.5, xanchor = "center", y = -0.2))

5.4 Erros por Regime Epidemiológico

Code

preds_all <- bt_result$predictions
preds_all$regime <- classify_regime(preds_all$actual)

metrics_regime <- compute_metrics_by_regime(preds_all)

metrics_regime |>
    mutate(across(where(is.numeric), \(x) round(x, 3))) |>
    datatable(options = list(dom = "t"), caption = "Métricas por Regime")

6. Previsão Futura

Code

# Definir atraso de notificação
delay_weeks <- 2
last_reliable_date <- max(df_base$data_iniSE) - (delay_weeks * 7)

# Filtrar dados para treino final (truncagem)
df_final <- df_base |>
    filter(data_iniSE <= last_reliable_date)

cat(sprintf("Ignorando últimas %d semanas (dados provisórios).\n", delay_weeks))

#> Ignorando últimas 2 semanas (dados provisórios).

Code

cat(sprintf("Treinando modelo com dados até: %s\n", last_reliable_date))

#> Treinando modelo com dados até: 2025-12-07

Code

cat(sprintf("Treinando modelo com dados até: %s\n", last_reliable_date))

#> Treinando modelo com dados até: 2025-12-07

Code

# Gerar features para o modelo final (importante para regressores)
df_final_features <- make_features_for_prediction(df_final, config_features)

#> [2026-01-11 15:20:19] INFO: ========== INICIO: Feature Engineering ========== 
#> [2026-01-11 15:20:19] INFO: Adicionando features temporais... 
#> [2026-01-11 15:20:19] INFO: Adicionando features de sazonalidade... 
#> [2026-01-11 15:20:19] INFO: Adicionando features climaticas... 
#> [2026-01-11 15:20:19] INFO: Features criadas: 47 novas colunas 
#> [2026-01-11 15:20:19] INFO: ========== FIM: Feature Engineering ==========

Code

model_fit <- model$fit(df_final_features, model_config)

params <- model$get_params(model_fit)

cat("### Especificação do Modelo Final\n")

#> ### Especificação do Modelo Final

Code

cat(sprintf("- **Modo de sazonalidade**: %s\n", params$seasonality_mode))

#> - **Modo de sazonalidade**: multiplicative

Code

cat(sprintf("- **Changepoint prior scale**: %.4f\n", params$changepoint_prior_scale))

#> - **Changepoint prior scale**: 0.5000

Code

cat(sprintf("- **Número de change points**: %d\n", params$n_changepoints))

#> - **Número de change points**: 25

Code

if (params$n_changepoints > 0) {
    cat("\n**Change points detectados** (primeiros 5):\n")
    print(head(params$changepoints, 5))
}

#> 
#> **Change points detectados** (primeiros 5):
#> [1] "2014-06-22 GMT" "2014-11-02 GMT" "2015-03-29 GMT" "2015-08-16 GMT"
#> [5] "2015-12-27 GMT"

Code

# Componentes do Prophet
decomp <- model$decompose(model_fit)

# Visualizar
decomp_df <- data.frame(
    t = seq_along(decomp$trend),
    trend = decomp$trend,
    yearly = if (!is.null(decomp$yearly)) decomp$yearly else 0
)

decomp_long <- decomp_df |>
    pivot_longer(cols = -t, names_to = "component", values_to = "value")

ggplot(decomp_long, aes(x = t, y = value)) +
    geom_line(color = "steelblue") +
    facet_wrap(~component, scales = "free_y", ncol = 1) +
    labs(title = "Decomposição Prophet", x = "Tempo", y = "Valor") +
    theme_minimal()

Code

last_date <- max(df_base$data_iniSE)
future_dates <- seq(from = last_date + 7, by = "week", length.out = 8)

# Importante: passar df_features como new_data para fornecer os regressores
forecast_h8 <- model$predict(model_fit, new_data = df_features, h = 8)
intervals <- model$predict_interval(model_fit, new_data = df_features, h = 8, level = 0.95)

forecast_df <- data.frame(
    Data = future_dates,
    Semana_Epidem = lubridate::epiweek(future_dates),
    Previsao = round(forecast_h8, 2),
    Inferior_95 = round(intervals$lower, 2),
    Superior_95 = round(intervals$upper, 2)
)

forecast_df |>
    kable(caption = "Previsão Prophet - Próximas 8 Semanas") |>
    kable_styling(full_width = FALSE, bootstrap_options = c("striped", "hover"))

Previsão Prophet - Próximas 8 Semanas
Data	Semana_Epidem	Previsao	Inferior_95	Superior_95
2025-12-28	53	2.27	0.00	33.87
2026-01-04	1	1.06	0.00	32.96
2026-01-11	2	3.85	0.00	36.01
2026-01-18	3	9.04	0.00	39.99
2026-01-25	4	15.73	0.00	45.28
2026-02-01	5	22.21	0.00	52.35
2026-02-08	6	27.98	0.48	58.64
2026-02-15	7	34.49	3.74	64.87

Code

# Historico usado no treino
history_train <- df_final |>
    tail(104) |>
    select(Date = data_iniSE, Value = est_inc100k) |>
    mutate(Type = "Histórico", Lower = NA, Upper = NA)

# Historico ignorado
history_ignored <- df_base |>
    filter(data_iniSE > last_reliable_date) |>
    select(Date = data_iniSE, Value = est_inc100k) |>
    mutate(Type = "Provisório (Ignorado)", Lower = NA, Upper = NA)

history <- bind_rows(history_train, history_ignored)

future_viz <- data.frame(
    Date = future_dates,
    Value = forecast_h8,
    Lower = intervals$lower,
    Upper = intervals$upper,
    Type = "Previsão (Prophet)"
)

viz_df <- bind_rows(history, future_viz)

p_future <- ggplot(viz_df, aes(x = Date, y = Value, color = Type)) +
    geom_line(linewidth = 0.8) +
    geom_ribbon(aes(ymin = Lower, ymax = Upper, fill = Type), alpha = 0.2, color = NA) +
    scale_color_manual(values = c("Histórico" = "black", "Previsão (Prophet)" = "darkorange", "Provisório (Ignorado)" = "gray")) +
    scale_fill_manual(values = c("Histórico" = NA, "Previsão (Prophet)" = "darkorange", "Provisório (Ignorado)" = NA)) +
    labs(
        title = "Previsão Prophet - Próximas 8 Semanas",
        subtitle = "Com intervalos de confiança (80%)",
        y = "Incidência / 100k hab"
    ) +
    theme_minimal()

ggplotly(p_future)

7. Salvamento e Registro

7.1 Salvar Modelo

Code

models_dir <- file.path(PROJECT_ROOT, "data", "models")
if (!dir.exists(models_dir)) dir.create(models_dir, recursive = TRUE)

saveRDS(model_fit, file.path(models_dir, "prophet.rds"))
cat("Modelo salvo em data/models/prophet.rds\n")

#> Modelo salvo em data/models/prophet.rds

7.2 Atualizar Leaderboard

Code

metrics_path <- file.path(PROJECT_ROOT, "data", "model_metrics.RData")

new_entry <- data.frame(
    model_name = "prophet",
    execution_date = Sys.time(),
    horizon_4w_rmse = bt_result$metrics[bt_result$metrics$h == 4, "rmse"],
    horizon_4w_mae = bt_result$metrics[bt_result$metrics$h == 4, "mae"],
    horizon_4w_mape = bt_result$metrics[bt_result$metrics$h == 4, "mape"],
    horizon_4w_rmsle = bt_result$metrics[bt_result$metrics$h == 4, "rmsle"],
    horizon_8w_rmse = bt_result$metrics[bt_result$metrics$h == 8, "rmse"],
    horizon_8w_mae = bt_result$metrics[bt_result$metrics$h == 8, "mae"],
    horizon_8w_mape = bt_result$metrics[bt_result$metrics$h == 8, "mape"],
    horizon_8w_rmsle = bt_result$metrics[bt_result$metrics$h == 8, "rmsle"],
    avg_rmse = mean(bt_result$metrics$rmse),
    avg_mape = mean(bt_result$metrics$mape)
)

if (file.exists(metrics_path)) {
    load(metrics_path)
    model_metrics <- dplyr::bind_rows(model_metrics, new_entry)
} else {
    model_metrics <- new_entry
}

save(model_metrics, file = metrics_path)

model_metrics |>
    group_by(model_name) |>
    slice(which.max(execution_date)) |>
    ungroup() |>
    select(Modelo = model_name, RMSE_4w = horizon_4w_rmse, MAPE_4w = horizon_4w_mape, RMSLE_4w = horizon_4w_rmsle) |>
    arrange(RMSE_4w) |>
    kable(caption = "Leaderboard de Modelos") |>
    kable_styling(bootstrap_options = c("striped", "hover"))

Leaderboard de Modelos
Modelo	RMSE_4w	MAPE_4w	RMSLE_4w
tbats	20.46878	23.29374	0.3153439
sarima	36.32899	25.12562	0.2929710
sarimax	90.69557	43.82694	0.7821497
seasonal_mean	91.28455	65.00128	1.0641488
seasonal_naive	96.31014	66.14327	0.9683374
prophet	117.49434	77.05757	1.1778603

8. Conclusão

O Prophet oferece uma alternativa robusta para séries com outliers e mudanças de regime.

Pontos de Destaque

Detecção de Change Points: Identifica automaticamente mudanças de tendência
Sazonalidade Multiplicativa: Amplitude varia proporcionalmente ao nível
Robustez: Menos sensível a picos extremos que ARIMA
Interpretabilidade: Componentes facilmente visualizáveis

Limitações

Pode subestimar picos em séries epidêmicas
Não incorpora variáveis exógenas nativamente (requer regressors adicionais)
Intervalos de confiança podem ser estreitos demais