Nota metodológica: Este relatório analisa 5802 ataques mecânicos pesados registados na rede IP. A eficácia é avaliada por comparação directa entre inspecções de nivelamento antes e depois de cada ataque, usando os dados de inspecção do veículo de medição. Secções com dados insuficientes são assinaladas explicitamente.


Contexto e Pergunta Central

Os ataques mecânicos pesados são eficazes e justificam o seu uso?

Esta pergunta desdobra-se em quatro dimensões que este relatório analisa sequencialmente:

  1. Eficácia imediata — o ataque melhora o nivelamento?
  2. Duração do efeito — durante quanto tempo se mantém a melhoria?
  3. Relação eficácia-duração — ataques que melhoram mais duram mais?
  4. Qualidade dos dados — com que confiança podemos concluir?

1 Qualidade e Cobertura dos Dados

Quando é que NÃO devemos confiar nos resultados?

Antes de qualquer análise de eficácia, é fundamental perceber quantos ataques têm dados suficientes para ser avaliados — e quais as limitações dos restantes.

# Classificar cada ataque segundo a disponibilidade de dados
cobertura <- metricas %>%
  mutate(
    tem_pre  = !is.na(NIV_PRE),
    tem_pos  = !is.na(NIV_POS),
    tem_ambos = tem_pre & tem_pos,
    conf_class = case_when(
      # Dados completos + data directa = máxima confiança
      tem_ambos & !data_indireta                          ~ "Alta (dados completos, data directa)",
      # Dados completos mas data estimada = confiança média
      tem_ambos & data_indireta                           ~ "Média (dados completos, data indirecta)",
      # Só pós-ataque disponível
      !tem_pre & tem_pos                                  ~ "Baixa (sem dados pré-ataque)",
      # Só pré-ataque disponível
      tem_pre & !tem_pos                                  ~ "Baixa (sem dados pós-ataque)",
      # Sem dados de inspecção
      TRUE                                                ~ "Indisponível"
    )
  )

# Tabela de cobertura
tab_cobertura <- cobertura %>%
  count(conf_class) %>%
  mutate(
    pct   = n / sum(n),
    label = paste0(n, " (", percent(pct, accuracy = 0.1), ")")
  ) %>%
  arrange(desc(n))

tab_cobertura %>%
  gt() %>%
  cols_label(
    conf_class = "Classe de Confiança",
    n          = "Nº Ataques",
    pct        = "% Total",
    label      = "Resumo"
  ) %>%
  fmt_percent(columns = pct, decimals = 1) %>%
  tab_header(
    title    = "Cobertura de Dados por Ataque",
    subtitle = "Classificação segundo disponibilidade de inspecções pré e pós-ataque"
  ) %>%
  tab_style(
    style = cell_fill(color = "#d4edda"),
    locations = cells_body(rows = str_starts(conf_class, "Alta"))
  ) %>%
  tab_style(
    style = cell_fill(color = "#fff3cd"),
    locations = cells_body(rows = str_starts(conf_class, "Média"))
  ) %>%
  tab_style(
    style = cell_fill(color = "#f8d7da"),
    locations = cells_body(rows = str_starts(conf_class, "Baixa") | conf_class == "Indisponível")
  )
Cobertura de Dados por Ataque
Classificação segundo disponibilidade de inspecções pré e pós-ataque
Classe de Confiança Nº Ataques % Total Resumo
Alta (dados completos, data directa) 4362 75.2% 4362 (75.2%)
Média (dados completos, data indirecta) 1173 20.2% 1173 (20.2%)
Baixa (sem dados pré-ataque) 130 2.2% 130 (2.2%)
Baixa (sem dados pós-ataque) 75 1.3% 75 (1.3%)
Indisponível 62 1.1% 62 (1.1%)
# Distribuição da distância temporal entre ataque e inspecções mais próximas
dist_temp <- metricas %>%
  filter(!is.na(days_pre) | !is.na(days_pos)) %>%
  select(ataque_id, days_pre, days_pos, data_indireta) %>%
  pivot_longer(cols = c(days_pre, days_pos),
               names_to = "tipo", values_to = "dias") %>%
  filter(!is.na(dias)) %>%
  mutate(
    tipo  = if_else(tipo == "days_pre", "Antes do ataque", "Após o ataque"),
    dias  = abs(dias)
  )

p_dist <- dist_temp %>%
  ggplot(aes(x = dias, fill = tipo)) +
  geom_histogram(bins = 40, alpha = 0.8, colour = "white") +
  facet_wrap(~tipo, scales = "free_y") +
  scale_fill_manual(values = c("Antes do ataque" = COR_WARN,
                                "Após o ataque"  = COR_POS)) +
  labs(
    title    = "Distância temporal entre ataque e inspecção mais próxima",
    subtitle = "Valores elevados indicam maior incerteza na medição do efeito imediato",
    x        = "Dias",
    y        = "Nº ataques",
    fill     = NULL
  ) +
  theme(legend.position = "none")

ggplotly(p_dist, tooltip = c("x", "y")) %>%
  layout(showlegend = FALSE)
# Identificar ataques com distância temporal elevada (potencial contaminação)
limiar_dias_pre <- 365   # inspecção pré mais velha que 1 ano — baseline menos fiável
limiar_dias_pos <- 180   # inspecção pós mais de 6 meses depois — pode já ter degradado

n_dist_pre_alta <- sum(abs(metricas$days_pre) > limiar_dias_pre, na.rm = TRUE)
n_dist_pos_alta <- sum(metricas$days_pos > limiar_dias_pos, na.rm = TRUE)

⚠️ 41 ataques têm a inspecção pré-ataque a mais de 1 ano de distância (baseline possivelmente desactualizado). 163 ataques têm a primeira inspecção pós-ataque a mais de 6 meses (efeito imediato pode já incluir degradação posterior).


2 Eficácia Imediata

Os ataques funcionam? Melhoram quanto?

# Calcular estatísticas de eficácia imediata
n_total     <- nrow(m_validos)
n_melhora   <- sum(m_validos$delta_niv > 0, na.rm = TRUE)
n_piora     <- sum(m_validos$delta_niv < 0, na.rm = TRUE)
n_neutro    <- sum(m_validos$delta_niv == 0, na.rm = TRUE)
pct_melhora <- n_melhora / n_total
pct_piora   <- n_piora / n_total

De um total de 5535 ataques com dados suficientes para avaliação imediata:

  • 4739 (85.6%) melhoraram o nivelamento
  • 764 (13.8%) pioraram ou não produziram efeito
  • 32 sem variação detectável

2.1 Funcionam em média?

# Estatísticas descritivas do delta_niv
stats_delta <- m_validos %>%
  summarise(
    Media    = mean(delta_niv, na.rm = TRUE),
    Mediana  = median(delta_niv, na.rm = TRUE),
    P25      = quantile(delta_niv, 0.25, na.rm = TRUE),
    P75      = quantile(delta_niv, 0.75, na.rm = TRUE),
    DP       = sd(delta_niv, na.rm = TRUE),
    Min      = min(delta_niv, na.rm = TRUE),
    Max      = max(delta_niv, na.rm = TRUE)
  ) %>%
  pivot_longer(everything(), names_to = "Estatística", values_to = "Valor (mm)")

stats_delta %>%
  gt() %>%
  fmt_number(columns = `Valor (mm)`, decimals = 3) %>%
  tab_header(
    title    = "Eficácia Imediata — Δ Nivelamento (mm)",
    subtitle = "δ > 0 indica melhoria (nivelamento reduziu após ataque)"
  )
Eficácia Imediata — Δ Nivelamento (mm)
δ > 0 indica melhoria (nivelamento reduziu após ataque)
Estatística Valor (mm)
Media 0.455
Mediana 0.450
P25 0.170
P75 0.730
DP 0.518
Min −5.793
Max 3.821
# Histograma da distribuição de delta_niv
p_hist <- m_validos %>%
  ggplot(aes(x = delta_niv,
             fill  = delta_niv > 0,
             text  = paste0("Δ niv: ", round(delta_niv, 3), " mm"))) +
  geom_histogram(bins = 60, colour = "white", alpha = 0.85) +
  geom_vline(xintercept = 0,
             linetype = "dashed", colour = "grey30", linewidth = 0.8) +
  geom_vline(xintercept = median(m_validos$delta_niv, na.rm = TRUE),
             linetype = "solid", colour = "black", linewidth = 1) +
  scale_fill_manual(values = c("TRUE" = COR_POS, "FALSE" = COR_NEG),
                    labels = c("TRUE" = "Melhoria", "FALSE" = "Piora / sem efeito")) +
  annotate("text",
           x = median(m_validos$delta_niv, na.rm = TRUE) + 0.05,
           y = Inf, vjust = 2, hjust = 0, size = 3.5,
           label = paste0("Mediana: ",
                          round(median(m_validos$delta_niv, na.rm = TRUE), 3), " mm")) +
  labs(
    title    = "Distribuição da Eficácia Imediata (Δ Nivelamento)",
    subtitle = "Cada barra representa um ataque; linha vertical = mediana",
    x        = "Δ Nivelamento (mm) — positivo = melhoria",
    y        = "Nº ataques",
    fill     = NULL
  )

ggplotly(p_hist, tooltip = "text")

2.2 Melhoria em percentagem do baseline

# delta_niv_pct = (NIV_PRE - NIV_POS) / NIV_PRE * 100
p_pct <- m_validos %>%
  filter(!is.na(delta_niv_pct), is.finite(delta_niv_pct)) %>%
  ggplot(aes(x = delta_niv_pct,
             fill = delta_niv_pct > 0,
             text = paste0("Melhoria: ", round(delta_niv_pct, 1), "%"))) +
  geom_histogram(bins = 60, colour = "white", alpha = 0.85) +
  geom_vline(xintercept = 0, linetype = "dashed", colour = "grey30") +
  scale_fill_manual(values = c("TRUE" = COR_POS, "FALSE" = COR_NEG)) +
  scale_x_continuous(labels = label_percent(scale = 1)) +
  labs(
    title    = "Melhoria Relativa do Nivelamento (% do baseline pré-ataque)",
    subtitle = "Normaliza pelo estado inicial — comparável entre troços com níveis distintos",
    x        = "Melhoria relativa (%)",
    y        = "Nº ataques",
    fill     = NULL
  ) +
  theme(legend.position = "none")

ggplotly(p_pct, tooltip = "text")

2.3 Eficácia vs distância temporal à inspecção

# Verificar se a eficácia medida depende do tempo decorrido entre ataque
# e inspecção mais próxima — se sim, estamos a misturar efeito real com degradação
p_dist_ef <- m_validos %>%
  filter(!is.na(days_pos), days_pos >= 0) %>%
  ggplot(aes(x = days_pos, y = delta_niv,
             colour = data_indireta,
             text = paste0(
               "Ataque: ", ataque_id, "<br>",
               "Dias até inspecção: ", round(days_pos), "<br>",
               "Δ niv: ", round(delta_niv, 3), " mm"
             ))) +
  geom_point(alpha = 0.4, size = 1.5) +
  geom_smooth(method = "loess", se = TRUE, colour = "black", linewidth = 1) +
  geom_hline(yintercept = 0, linetype = "dashed", colour = COR_NEG) +
  scale_colour_manual(values = c("FALSE" = COR_POS, "TRUE" = COR_WARN),
                      labels = c("FALSE" = "Data directa", "TRUE" = "Data indirecta")) +
  labs(
    title    = "Eficácia imediata vs. tempo até inspecção pós-ataque",
    subtitle = "Se a linha de tendência descer com o tempo, estamos a medir degradação posterior",
    x        = "Dias entre ataque e primeira inspecção pós",
    y        = "Δ Nivelamento (mm)",
    colour   = "Confiança da data"
  )

ggplotly(p_dist_ef, tooltip = "text")

3 Duração do Efeito

Durante quanto tempo dura o efeito do ataque?

# Classificar ataques segundo o estado de observação da duração
dur_class <- metricas %>%
  mutate(
    estado_duracao = case_when(
      reverteu                        ~ "Reverteu (observado)",
      censored_dtr                    ~ "Censurado (novo ataque)",
      !reverteu & !censored_dtr &
        !is.na(days_to_revert)        ~ "Ainda eficaz (fim dos dados)",
      TRUE                            ~ "Sem dados suficientes"
    )
  )

tab_dur <- dur_class %>%
  count(estado_duracao) %>%
  mutate(pct = n / sum(n))

tab_dur %>%
  gt() %>%
  cols_label(estado_duracao = "Estado", n = "Nº Ataques", pct = "% Total") %>%
  fmt_percent(columns = pct, decimals = 1) %>%
  tab_header(
    title    = "Observabilidade da Duração do Efeito",
    subtitle = "Quantos ataques permitem medir quanto tempo o efeito durou"
  )
Observabilidade da Duração do Efeito
Quantos ataques permitem medir quanto tempo o efeito durou
Estado Nº Ataques % Total
Censurado (novo ataque) 1472 25.4%
Reverteu (observado) 2943 50.7%
Sem dados suficientes 1387 23.9%

3.1 Curva de sobrevivência — Kaplan-Meier

# Análise de sobrevivência: "sobreviver" = manter-se eficaz (não reverter ao baseline)
# Evento = reversão ao nível pré-ataque (reverteu = TRUE)
# Censura = novo ataque antes de reversão (censored_dtr = TRUE) ou fim dos dados

surv_data <- metricas %>%
  filter(!is.na(days_to_revert) | reverteu | censored_dtr) %>%
  mutate(
    # Para Surv(): time = dias observados; event = 1 se reverteu, 0 se censurado
    surv_time  = if_else(reverteu, days_to_revert,
                         if_else(!is.na(days_to_revert), days_to_revert, NA_real_)),
    surv_event = as.integer(reverteu)
  ) %>%
  filter(!is.na(surv_time), surv_time > 0)

if (nrow(surv_data) >= 10) {
  km_fit <- survfit(Surv(surv_time, surv_event) ~ 1, data = surv_data)

  # Usar ggsurvplot e converter para plotly
  p_km <- ggsurvplot(
    km_fit,
    data           = surv_data,
    conf.int       = TRUE,
    risk.table     = TRUE,
    palette        = COR_POS,
    xlab           = "Dias após ataque",
    ylab           = "Probabilidade de manter eficácia",
    title          = "Duração do Efeito — Curva de Kaplan-Meier",
    subtitle       = "Probabilidade de o nivelamento se manter abaixo do baseline pré-ataque",
    ggtheme        = theme_ip,
    risk.table.col = "strata"
  )
  print(p_km)

  # Mediana da sobrevivência
  km_median <- summary(km_fit)$table["median"]
  cat(sprintf("\nDuração mediana do efeito: %.0f dias (≈ %.1f meses)\n",
              km_median, km_median / 30.5))
} else {
  cat("Dados insuficientes para curva de Kaplan-Meier (< 10 ataques com duração observável).")
}

## 
## Duração mediana do efeito: 279 dias (≈ 9.1 meses)

3.2 Distribuição de days_to_revert

# Histograma dos ataques que reverteram (duração observada)
p_dtr <- metricas %>%
  filter(reverteu, !is.na(days_to_revert)) %>%
  ggplot(aes(x = days_to_revert,
             text = paste0("Dias: ", round(days_to_revert), " (≈",
                           round(days_to_revert/30.5, 1), " meses)"))) +
  geom_histogram(bins = 40, fill = COR_POS, colour = "white", alpha = 0.85) +
  geom_vline(xintercept = median(metricas$days_to_revert, na.rm = TRUE),
             linetype = "dashed", colour = "black") +
  scale_x_continuous(
    sec.axis = sec_axis(~ . / 30.5, name = "Meses")
  ) +
  labs(
    title    = "Distribuição da Duração do Efeito (ataques que reverteram)",
    subtitle = "Linha a tracejado = mediana",
    x        = "Dias até reversão ao baseline",
    y        = "Nº ataques"
  )

ggplotly(p_dtr, tooltip = "text")

3.3 Taxa de degradação pré vs pós-ataque

# Comparar declives pré e pós — se slope_pos < slope_pre, o ataque abrandou a degradação
p_slopes <- m_slopes %>%
  pivot_longer(cols = c(slope_pre, slope_pos),
               names_to = "periodo",
               values_to = "slope") %>%
  mutate(periodo = if_else(periodo == "slope_pre", "Pré-ataque", "Pós-ataque")) %>%
  ggplot(aes(x = slope * 30.5,   # converter para mm/mês
             fill = periodo,
             text = paste0(periodo, ": ", round(slope * 30.5, 4), " mm/mês"))) +
  geom_histogram(bins = 50, alpha = 0.7, colour = "white",
                 position = "identity") +
  scale_fill_manual(values = c("Pré-ataque" = COR_WARN, "Pós-ataque" = COR_POS)) +
  geom_vline(xintercept = 0, linetype = "dashed", colour = "grey30") +
  labs(
    title    = "Taxa de Degradação do Nivelamento — Pré vs Pós-ataque",
    subtitle = "Se a distribuição pós-ataque está mais à esquerda, o ataque abrandou a degradação",
    x        = "Taxa de degradação (mm/mês)",
    y        = "Nº ataques",
    fill     = NULL
  )

ggplotly(p_slopes, tooltip = "text")
# delta_slope = slope_pre - slope_pos
# positivo → ataque abrandou a degradação (bom)
# negativo → degradação acelerou após ataque (mau)

p_dslope <- m_slopes %>%
  ggplot(aes(x = delta_slope * 30.5,
             fill = delta_slope > 0,
             text = paste0("Δ slope: ", round(delta_slope * 30.5, 4), " mm/mês"))) +
  geom_histogram(bins = 50, colour = "white", alpha = 0.85) +
  geom_vline(xintercept = 0, linetype = "dashed", colour = "grey30") +
  scale_fill_manual(values = c("TRUE" = COR_POS, "FALSE" = COR_NEG),
                    labels = c("TRUE" = "Abrandou degradação",
                               "FALSE" = "Acelerou degradação")) +
  labs(
    title    = "Variação da Taxa de Degradação (Δ slope)",
    subtitle = "Positivo = ataque abrandou a degradação futura",
    x        = "Δ slope (mm/mês) — positivo = melhoria",
    y        = "Nº ataques",
    fill     = NULL
  )

ggplotly(p_dslope, tooltip = "text")

4 Relação Eficácia vs Duração

Ataques que melhoram mais duram mais?

# Scatter: delta_niv vs days_to_revert
p_scatter <- metricas %>%
  filter(!is.na(delta_niv), reverteu, !is.na(days_to_revert)) %>%
  ggplot(aes(
    x    = delta_niv,
    y    = days_to_revert,
    colour = MSI_categoria,
    text = paste0(
      "Ataque: ", ataque_id, "<br>",
      "Linha: ", LINHA_NOME, "<br>",
      "Δ niv: ", round(delta_niv, 3), " mm<br>",
      "Duração: ", round(days_to_revert), " dias<br>",
      "MSI: ", MSI_categoria
    )
  )) +
  geom_point(alpha = 0.5, size = 2) +
  geom_smooth(method = "lm", se = TRUE, colour = "black",
              linewidth = 0.8, linetype = "dashed") +
  scale_colour_manual(values = MSI_CORES, na.value = "grey70") +
  labs(
    title    = "Eficácia Imediata vs. Duração do Efeito",
    subtitle = "Cada ponto é um ataque; cor = classificação MSI",
    x        = "Eficácia imediata — Δ Nivelamento (mm)",
    y        = "Duração do efeito (dias)",
    colour   = "MSI"
  )

ggplotly(p_scatter, tooltip = "text")
# Correlação de Spearman (mais robusta a outliers que Pearson)
cor_data <- metricas %>%
  filter(!is.na(delta_niv), reverteu, !is.na(days_to_revert))

if (nrow(cor_data) >= 10) {
  cor_test <- cor.test(cor_data$delta_niv, cor_data$days_to_revert,
                       method = "spearman")
  cat(sprintf(
    "Correlação Spearman (eficácia vs duração): ρ = %.3f (p = %.4f, n = %d)\n",
    cor_test$estimate, cor_test$p.value, nrow(cor_data)
  ))
} else {
  cat("Dados insuficientes para correlação.")
}
## Correlação Spearman (eficácia vs duração): ρ = 0.673 (p = 0.0000, n = 2943)

4.1 Quadrantes de eficácia

# Dividir ataques em 4 quadrantes com base na mediana de delta_niv e days_to_revert
med_delta <- median(metricas$delta_niv, na.rm = TRUE)
med_dur   <- median(metricas$days_to_revert, na.rm = TRUE)

p_quad <- metricas %>%
  filter(!is.na(delta_niv), reverteu, !is.na(days_to_revert)) %>%
  mutate(
    quadrante = case_when(
      delta_niv >= med_delta & days_to_revert >= med_dur ~ "Alta eficácia / Longa duração",
      delta_niv >= med_delta & days_to_revert <  med_dur ~ "Alta eficácia / Curta duração",
      delta_niv <  med_delta & days_to_revert >= med_dur ~ "Baixa eficácia / Longa duração",
      TRUE                                               ~ "Baixa eficácia / Curta duração"
    )
  ) %>%
  ggplot(aes(
    x      = delta_niv,
    y      = days_to_revert,
    colour = quadrante,
    text   = paste0("Δ niv: ", round(delta_niv, 3),
                    " mm | Duração: ", round(days_to_revert), " dias")
  )) +
  geom_point(alpha = 0.5, size = 2) +
  geom_vline(xintercept = med_delta, linetype = "dashed", colour = "grey40") +
  geom_hline(yintercept = med_dur,   linetype = "dashed", colour = "grey40") +
  scale_colour_brewer(palette = "Set1") +
  annotate("text", x = max(metricas$delta_niv, na.rm=TRUE) * 0.9,
           y = max(metricas$days_to_revert, na.rm=TRUE) * 0.95,
           label = "Ideal", fontface = "bold", colour = "darkgreen") +
  annotate("text", x = min(metricas$delta_niv, na.rm=TRUE) * 0.9,
           y = min(metricas$days_to_revert, na.rm=TRUE) * 1.05,
           label = "Problemático", fontface = "bold", colour = "darkred") +
  labs(
    title    = "Quadrantes de Eficácia — Impacto vs. Duração",
    subtitle = "Linhas a tracejado = mediana; quadrante ideal: canto superior direito",
    x        = "Eficácia imediata — Δ Nivelamento (mm)",
    y        = "Duração do efeito (dias)",
    colour   = "Quadrante"
  )

ggplotly(p_quad, tooltip = "text")

5 MSI — Maintenance Success Indicator

Classificação integrada da eficácia de cada ataque

O MSI classifica cada ataque em quatro categorias, combinando a melhoria imediata (σt) com a variação da taxa de degradação (Δσt):

  • Cat. 1 — Efeitos Benéficos: abrandou a degradação E melhorou o nível
  • Cat. 2 — Efeitos de Retardo: melhorou o nível mas degradação continua alta
  • Cat. 3 — Sem Efeitos Significativos: sem impacto detectável
  • Cat. 4 — Efeitos Negativos: nivelamento piorou após o ataque
msi_counts <- metricas %>%
  filter(!is.na(MSI_categoria)) %>%
  count(MSI_categoria) %>%
  mutate(pct = n / sum(n))

# Gráfico de barras MSI
p_msi <- msi_counts %>%
  ggplot(aes(
    x    = reorder(MSI_categoria, -n),
    y    = n,
    fill = MSI_categoria,
    text = paste0(MSI_categoria, ": ", n, " ataques (", percent(pct, 0.1), ")")
  )) +
  geom_col(alpha = 0.9, width = 0.7) +
  geom_text(aes(label = paste0(n, "\n(", percent(pct, 0.1), ")")),
            vjust = -0.3, size = 3.5) +
  scale_fill_manual(values = MSI_CORES) +
  scale_y_continuous(expand = expansion(mult = c(0, 0.15))) +
  labs(
    title    = "Distribuição MSI — Classificação de Eficácia dos Ataques",
    subtitle = "Baseado em limiares empíricos (percentis da distribuição observada)",
    x        = NULL,
    y        = "Nº ataques",
    fill     = NULL
  ) +
  theme(legend.position = "none",
        axis.text.x = element_text(angle = 15, hjust = 1))

ggplotly(p_msi, tooltip = "text")
# Tabela resumo MSI com estatísticas por categoria
metricas %>%
  filter(!is.na(MSI_categoria)) %>%
  group_by(MSI_categoria) %>%
  summarise(
    N             = n(),
    `Δ niv médio (mm)` = round(mean(delta_niv,    na.rm = TRUE), 3),
    `Δ niv mediana`    = round(median(delta_niv,  na.rm = TRUE), 3),
    `Δ slope médio`    = round(mean(delta_slope,  na.rm = TRUE) * 30.5, 4),
    `Duração mediana (dias)` = round(median(days_to_revert, na.rm = TRUE), 0),
    `% data indirecta` = round(mean(data_indireta, na.rm = TRUE) * 100, 1)
  ) %>%
  gt() %>%
  tab_header(
    title    = "Resumo por Categoria MSI",
    subtitle = "Δ slope em mm/mês (positivo = abrandou degradação)"
  ) %>%
  tab_style(
    style = cell_fill(color = "#d4edda"),
    locations = cells_body(rows = str_starts(MSI_categoria, "1"))
  ) %>%
  tab_style(
    style = cell_fill(color = "#fff3cd"),
    locations = cells_body(rows = str_starts(MSI_categoria, "2"))
  ) %>%
  tab_style(
    style = cell_fill(color = "#f8d7da"),
    locations = cells_body(rows = str_starts(MSI_categoria, "3") |
                                  str_starts(MSI_categoria, "4"))
  )
Resumo por Categoria MSI
Δ slope em mm/mês (positivo = abrandou degradação)
MSI_categoria N Δ niv médio (mm) Δ niv mediana Δ slope médio Duração mediana (dias) % data indirecta
1 - Efeitos Beneficos 1055 0.759 0.693 0.0143 304 20.6
2 - Efeitos de Retardo 402 0.751 0.653 -0.0228 697 24.4
3 - Sem Efeitos Significativos 1077 0.245 0.250 -0.0028 330 17.0
4 - Efeitos Negativos 2500 0.392 0.444 -0.0099 191 19.0
# MSI por linha ferroviária — onde funcionam melhor?
metricas %>%
  filter(!is.na(MSI_categoria), !is.na(LINHA_NOME)) %>%
  count(LINHA_NOME, MSI_categoria) %>%
  group_by(LINHA_NOME) %>%
  mutate(pct = n / sum(n)) %>%
  ungroup() %>%
  ggplot(aes(
    x    = reorder(LINHA_NOME, -n),
    y    = pct,
    fill = MSI_categoria
  )) +
  geom_col(position = "stack", width = 0.7) +
  scale_fill_manual(values = MSI_CORES) +
  scale_y_continuous(labels = label_percent()) +
  coord_flip() +
  labs(
    title    = "Distribuição MSI por Linha Ferroviária",
    subtitle = "Proporção de cada categoria por linha — permite identificar linhas problemáticas",
    x        = NULL,
    y        = "% ataques",
    fill     = "MSI"
  )


6 Variabilidade Espacial

Onde funcionam melhor? Onde há padrões persistentes?

# Mapa leaflet com eficácia imediata por ataque
# Requer colunas Latitude e Latitude em metricas (do Bloco 8)
mapa_data <- metricas %>%
  filter(!is.na(Latitude), !is.na(Longitude), !is.na(delta_niv))

if (nrow(mapa_data) > 0) {

  # Paleta de cores: vermelho (piora) → verde (melhoria)
  pal_delta <- colorNumeric(
    palette = c(COR_NEG, "white", COR_POS),
    domain  = mapa_data$delta_niv,
    na.color = "grey"
  )

  leaflet(mapa_data) %>%
    addProviderTiles(providers$CartoDB.Positron) %>%
    addCircleMarkers(
      lng    = ~Longitude,
      lat    = ~Latitude,
      radius = ~pmax(3, pmin(12, abs(delta_niv) * 3)),
      color  = ~pal_delta(delta_niv),
      fillOpacity = 0.7,
      stroke = FALSE,
      popup  = ~paste0(
        "<b>Ataque ", ataque_id, "</b><br>",
        "Linha: ", LINHA_NOME, "<br>",
        "Data: ", DATA_ATAQUE, "<br>",
        "Δ Nivelamento: <b>", round(delta_niv, 3), " mm</b><br>",
        "MSI: ", MSI_categoria, "<br>",
        "Data indirecta: ", data_indireta
      )
    ) %>%
    addLegend(
      pal    = pal_delta,
      values = ~delta_niv,
      title  = "Δ Nivelamento (mm)",
      position = "bottomright"
    )

} else {
  cat("Sem coordenadas disponíveis para mapa. Verificar Bloco 8 do pipeline.")
}
# Eficácia média por segmento — identificar segmentos com padrão sistemático
metricas %>%
  filter(!is.na(delta_niv), !is.na(LINHA_INT)) %>%
  group_by(LINHA_INT, LINHA_NOME) %>%
  summarise(
    n_ataques    = n(),
    delta_medio  = mean(delta_niv, na.rm = TRUE),
    delta_dp     = sd(delta_niv, na.rm = TRUE),
    pct_melhora  = mean(delta_niv > 0, na.rm = TRUE),
    .groups = "drop"
  ) %>%
  filter(n_ataques >= 3) %>%   # só segmentos com pelo menos 3 ataques
  arrange(desc(delta_medio)) %>%
  gt() %>%
  cols_label(
    LINHA_INT   = "Segmento",
    LINHA_NOME  = "Linha",
    n_ataques   = "N ataques",
    delta_medio = "Δ médio (mm)",
    delta_dp    = "Desvio padrão",
    pct_melhora = "% melhora"
  ) %>%
  fmt_number(columns = c(delta_medio, delta_dp), decimals = 3) %>%
  fmt_percent(columns = pct_melhora, decimals = 1) %>%
  tab_header(
    title    = "Eficácia Média por Segmento (mínimo 3 ataques)",
    subtitle = "Ordenado por melhoria média decrescente"
  ) %>%
  data_color(
    columns = delta_medio,
    method  = "numeric",
    palette = c(COR_NEG, "white", COR_POS)
  )
Eficácia Média por Segmento (mínimo 3 ataques)
Ordenado por melhoria média decrescente
Segmento Linha N ataques Δ médio (mm) Desvio padrão % melhora
552 C. Bombel 6 0.899 0.263 100.0%
791 R. Neves Corvo 30 0.842 0.332 100.0%
521 C. Verride 4 0.776 0.507 100.0%
271 L. Leste 168 0.688 0.502 91.1%
222 R. Alfarelos 20 0.593 0.444 100.0%
281 L. Sintra 124 0.549 0.381 93.5%
242 R. Tomar 24 0.548 0.352 87.5%
292 C. Sete Rios 4 0.538 0.425 75.0%
451 L. Algarve 103 0.527 0.540 82.5%
61 L. Douro 134 0.519 0.486 87.3%
321 L. Cascais 183 0.517 0.468 87.4%
331 L. Vendas Novas 272 0.517 0.498 88.2%
12 L. Minho 333 0.514 0.556 85.9%
452 L. Algarve 179 0.484 0.452 89.4%
201 L. Beira Alta 257 0.470 0.515 84.0%
381 L. Sines 105 0.470 0.414 81.9%
231 L. Oeste 524 0.448 0.589 82.4%
81 L. Norte 1561 0.434 0.478 86.4%
252 L. Beira Baixa 445 0.434 0.382 87.6%
375 L. Sul 191 0.434 0.354 91.1%
372 L. Sul 57 0.425 0.695 86.0%
551 C. Bombel 3 0.414 0.113 100.0%
251 L. Beira Baixa 72 0.411 0.535 88.9%
374 L. Sul 177 0.399 0.545 80.2%
292 L. Cintura 75 0.395 0.613 78.7%
373 L. Sul 8 0.360 0.195 100.0%
371 L. Sul 44 0.313 0.646 77.3%
681 V. Alcácer 13 0.291 0.152 100.0%
341 L. Alentejo 391 0.281 0.715 76.7%
461 C. Poceirão 5 0.250 0.317 80.0%
901 R. Acesso ao Porto de Aveiro 3 0.017 0.090 33.3%
91 L. Guimarães 4 −0.055 0.255 50.0%

7 Factores Explicativos

Em que condições funcionam melhor? — preparação para fase seguinte (cadastro)

Esta secção analisa os factores disponíveis actualmente: linha, extensão do ataque, tipo de data. Os factores de cadastro (idade da via, tipo de travessas) serão integrados numa fase posterior.

# Eficácia vs extensão do ataque — ataques maiores são mais eficazes?
p_ext <- m_validos %>%
  filter(!is.na(EXT_M)) %>%
  ggplot(aes(
    x    = EXT_M / 1000,   # converter para km
    y    = delta_niv,
    colour = data_indireta,
    text = paste0(
      "Extensão: ", round(EXT_M/1000, 2), " km<br>",
      "Δ niv: ", round(delta_niv, 3), " mm"
    )
  )) +
  geom_point(alpha = 0.3, size = 1.5) +
  geom_smooth(method = "loess", se = TRUE, colour = "black") +
  geom_hline(yintercept = 0, linetype = "dashed", colour = COR_NEG) +
  scale_colour_manual(values = c("FALSE" = COR_POS, "TRUE" = COR_WARN),
                      labels = c("FALSE" = "Data directa", "TRUE" = "Data indirecta")) +
  labs(
    title    = "Eficácia vs. Extensão do Ataque",
    x        = "Extensão do ataque (km)",
    y        = "Δ Nivelamento (mm)",
    colour   = "Confiança da data"
  )

ggplotly(p_ext, tooltip = "text")
# Comparar eficácia entre ataques de data directa vs indirecta
# (controlo de qualidade — não deveria haver diferença sistemática)
m_validos %>%
  group_by(data_indireta) %>%
  summarise(
    N           = n(),
    media       = mean(delta_niv, na.rm = TRUE),
    mediana     = median(delta_niv, na.rm = TRUE),
    dp          = sd(delta_niv, na.rm = TRUE),
    pct_melhora = mean(delta_niv > 0, na.rm = TRUE)
  ) %>%
  mutate(tipo_data = if_else(data_indireta, "Indirecta", "Directa")) %>%
  select(tipo_data, N, media, mediana, dp, pct_melhora) %>%
  gt() %>%
  cols_label(
    tipo_data   = "Tipo de Data",
    N           = "N",
    media       = "Δ médio (mm)",
    mediana     = "Δ mediana (mm)",
    dp          = "Desvio padrão",
    pct_melhora = "% melhora"
  ) %>%
  fmt_number(columns = c(media, mediana, dp), decimals = 3) %>%
  fmt_percent(columns = pct_melhora, decimals = 1) %>%
  tab_header(
    title    = "Eficácia por Tipo de Data",
    subtitle = "Se houver diferença grande, a data indirecta está a enviesar os resultados"
  )
Eficácia por Tipo de Data
Se houver diferença grande, a data indirecta está a enviesar os resultados
Tipo de Data N Δ médio (mm) Δ mediana (mm) Desvio padrão % melhora
Directa 4362 0.459 0.450 0.506 87.8%
Indirecta 1173 0.438 0.444 0.559 77.5%

8 Decisão Operacional

O que fazer com estes resultados?

# Calcular indicadores chave para resumo executivo
n_analisados    <- nrow(m_validos)
pct_ef          <- mean(m_validos$delta_niv > 0, na.rm = TRUE)
delta_mediana   <- median(m_validos$delta_niv, na.rm = TRUE)
dur_mediana     <- median(metricas$days_to_revert, na.rm = TRUE)
pct_cat1        <- mean(metricas$MSI_num == 1, na.rm = TRUE)
pct_cat4        <- mean(metricas$MSI_num == 4, na.rm = TRUE)

8.1 Resumo executivo

# KPIs em formato de texto estruturado
cat(sprintf("
**Em %d ataques analisados:**

| Indicador | Valor |
|-----------|-------|
| %% ataques com melhoria de nivelamento | **%s** |
| Melhoria mediana de nivelamento | **%.3f mm** |
| Duração mediana do efeito | **%.0f dias (≈ %.1f meses)** |
| %% ataques MSI Categoria 1 (Efeitos Benéficos) | **%s** |
| %% ataques MSI Categoria 4 (Efeitos Negativos) | **%s** |
",
  n_analisados,
  percent(pct_ef, 0.1),
  delta_mediana,
  dur_mediana, dur_mediana / 30.5,
  percent(pct_cat1, 0.1),
  percent(pct_cat4, 0.1)
))

Em 5535 ataques analisados:

Indicador Valor
% ataques com melhoria de nivelamento 85.6%
Melhoria mediana de nivelamento 0.450 mm
Duração mediana do efeito 276 dias (≈ 9.0 meses)
% ataques MSI Categoria 1 (Efeitos Benéficos) 21.0%
% ataques MSI Categoria 4 (Efeitos Negativos) 49.7%

8.2 Frequência óptima de ataque

# Estimativa da frequência óptima: baseada na duração mediana do efeito
# Se a duração mediana é X dias, um ataque a cada X dias manteria a via
# continuamente abaixo do baseline — mas na prática há custos e logística

if (!is.na(dur_mediana)) {
  cat(sprintf(
    "Com base na duração mediana de %.0f dias:\n\n",
    dur_mediana
  ))
  cat(sprintf(
    "  · Frequência de ataque sugerida: a cada %.0f–%.0f dias (%.1f–%.1f meses)\n",
    dur_mediana * 0.8, dur_mediana * 1.0,
    dur_mediana * 0.8 / 30.5, dur_mediana / 30.5
  ))
  cat(sprintf(
    "  · Antecipar o ataque em %.0f dias face à reversão esperada\n",
    dur_mediana * 0.2
  ))
} else {
  cat("Dados insuficientes para estimar frequência óptima.")
}
## Com base na duração mediana de 276 dias:
## 
##   · Frequência de ataque sugerida: a cada 221–276 dias (7.2–9.0 meses)
##   · Antecipar o ataque em 55 dias face à reversão esperada
# Curvas de degradação médias antes e depois do ataque
# Normalizar days_since_attack para centrar no ataque (t=0)
# e agrupar por janelas de tempo para suavizar

curvas_data <- nivel_A %>%
  filter(!is.na(NIV_MEDIO_ATQ), !is.na(days_since_attack)) %>%
  # Criar janelas de 30 dias para suavizar
  mutate(janela_dias = round(days_since_attack / 30) * 30) %>%
  group_by(janela_dias) %>%
  summarise(
    niv_medio = mean(NIV_MEDIO_ATQ, na.rm = TRUE),
    niv_p25   = quantile(NIV_MEDIO_ATQ, 0.25, na.rm = TRUE),
    niv_p75   = quantile(NIV_MEDIO_ATQ, 0.75, na.rm = TRUE),
    n         = n(),
    .groups = "drop"
  ) %>%
  filter(n >= 5)    # só janelas com observações suficientes

if (nrow(curvas_data) >= 5) {
  p_curvas <- curvas_data %>%
    ggplot(aes(x = janela_dias / 30)) +   # converter para meses
    geom_ribbon(aes(ymin = niv_p25, ymax = niv_p75),
                fill = COR_POS, alpha = 0.2) +
    geom_line(aes(y = niv_medio), colour = COR_POS, linewidth = 1.2) +
    geom_vline(xintercept = 0, linetype = "solid",
               colour = "black", linewidth = 1) +
    annotate("text", x = 0.3, y = max(curvas_data$niv_p75, na.rm=TRUE),
             label = "Ataque", hjust = 0, fontface = "bold") +
    labs(
      title    = "Curva de Degradação Média — Antes e Depois do Ataque",
      subtitle = "Banda = percentis 25–75; linha = média; t=0 = data do ataque",
      x        = "Meses relativos ao ataque",
      y        = "Nivelamento médio (mm)",
      caption  = "Cada ponto agrega todas as inspecções em janelas de 30 dias"
    )

  ggplotly(p_curvas)
} else {
  cat("Dados insuficientes para curvas de degradação agregadas.")
}

Notas Metodológicas

Correspondência espacial: cada ataque foi cruzado com os segmentos de inspecção de 200m por sobreposição de intervalos de PK, com pesos proporcionais ao comprimento de sobreposição. O valor de nivelamento do ataque é uma média ponderada desses segmentos.

Baseline pré-ataque: definido como o valor da inspecção imediatamente anterior ao ataque no mesmo troço. A taxa de degradação pré usa todas as inspecções disponíveis antes do ataque (mínimo 3 observações).

Censura por re-ataque: quando um novo ataque ocorre antes de observarmos a reversão ao baseline, a duração do efeito é marcada como censurada (censored_by_reattack = TRUE) e tratada como censura à direita na análise de sobrevivência.

Limiares MSI: calculados empiricamente por percentis da distribuição observada. Quando disponíveis dados de velocidade de linha, os limiares serão substituídos pelos valores normativos da Tabela 21 (Classes de Qualidade TQC — EN 13848).

Dados de cadastro: a análise de factores explicativos (idade da via, tipo de travessas, balastro) será integrada numa fase posterior, após cruzamento com os dados de cadastro da rede.


Relatório gerado automaticamente a 20/03/2026 18:32 — Pipeline AMP v1.0 | DAM/AM-PR/AM-INF