1 Propósito del análisis

En esta etapa se aplica Probabilidad a integridad de tuberías, ajustando distribuciones a variables numéricas para estimar probabilidades operativas y apoyar decisiones.

Se reporta:

  1. Estadística descriptiva (tabla + gráficos)
  2. Conjetura “a mano” (fórmulas + sustitución)
  3. Ajuste en R (solo los modelos definidos)
  4. Bondad de ajuste (Chi-cuadrado de Pearson) + AIC (desempate)
  5. Predicción probabilística

2 Preparación y validación de datos

# Nota: en R existe la función df() (densidad F). Por eso validamos tipo y no solo exists().
if (!exists("df", inherits = FALSE) || is.function(df) || !(is.data.frame(df) || tibble::is_tibble(df))) {
  df <- readxl::read_excel(path = params$data_path, sheet = params$sheet_name)
}
df <- as.data.frame(df)

kbl2(
  tibble(n = nrow(df), p = ncol(df), NA_totales = sum(is.na(df))),
  caption = "Resumen general del dataset"
)
Resumen general del dataset
n p NA_totales
1000 12 0 
model_map <- tibble::tribble(
  ~var_solicitada,       ~modelo,
  "Material_Loss_Norm",  "exp",
  "Material_Loss_Norm",  "log_normal",
  "Temperature",         "norm"
)

resolve_col <- function(v, cols) {
  hit <- cols[tolower(cols) == tolower(v)]
  if (length(hit) >= 1) return(hit[1])
  NA_character_
}

cols <- names(df)
model_map <- model_map |>
  dplyr::mutate(var_real = purrr::map_chr(var_solicitada, resolve_col, cols = cols))

if (any(is.na(model_map$var_real))) {
  faltan <- model_map |> dplyr::filter(is.na(var_real)) |> dplyr::pull(var_solicitada) |> unique()
  stop(
    "No se encontraron en df estas variables (revisar nombres): ",
    paste(faltan, collapse = ", "),
    "\nColumnas disponibles: ",
    paste(cols, collapse = ", ")
  )
}

vars_unicas <- unique(model_map$var_real)

# Validación: deben ser numéricas
no_num <- vars_unicas[!sapply(df[vars_unicas], is.numeric)]
if (length(no_num) > 0) stop("Estas variables no son numéricas: ", paste(no_num, collapse = ", "))

kbl2(model_map, caption = "Diccionario: variable y modelo(s) asignado(s)")
Diccionario: variable y modelo(s) asignado(s)
var_solicitada modelo var_real
Material_Loss_Norm exp Material_Loss_Norm
Material_Loss_Norm log_normal Material_Loss_Norm
Temperature norm Temperature 
thresholds <- purrr::map_dfr(vars_unicas, function(v) {
  x <- df[[v]]
  x <- x[is.finite(x)]
  tibble(
    variable = v,
    a   = as.numeric(stats::quantile(x, 0.25, na.rm = TRUE)),
    b   = as.numeric(stats::quantile(x, 0.75, na.rm = TRUE)),
    lim = as.numeric(stats::quantile(x, 0.10, na.rm = TRUE))
  )
}) |>
  dplyr::mutate(dplyr::across(c(a, b, lim), ~fmt_num(.x, 6)))

kbl2(thresholds, caption = "Ventana [a,b] y límite lim (cuantiles por defecto).")
Ventana [a,b] y límite lim (cuantiles por defecto).
variable a b lim
Material_Loss_Norm 0.048906 0.191279 0.016343
Temperature 13.400000 69.150000 -9.230000

3 Estadística descriptiva

desc_tbl <- purrr::map_dfr(vars_unicas, function(v) {
  x <- df[[v]]
  x <- x[is.finite(x)]
  tibble(
    variable   = v,
    etiqueta   = get_label(v),
    n          = length(x),
    media      = mean(x),
    mediana    = stats::median(x),
    varianza   = stats::var(x),
    sd         = stats::sd(x),
    min        = min(x),
    max        = max(x),
    asimetria  = psych::skew(x),
    curtosis   = psych::kurtosi(x)
  )
}) |>
  dplyr::mutate(dplyr::across(where(is.numeric), ~fmt_num(.x, 6)))

kbl2(desc_tbl, caption = "Tabla descriptiva")
Tabla descriptiva
variable etiqueta n media mediana varianza sd min max asimetria curtosis
Material_Loss_Norm Pérdida de material normalizada (Material_Loss_Norm) 1000 0.146445 0.099099 0.021386 0.146241 0 1.0 1.938118 4.604346
Temperature Temperatura (Temperature) 1000 42.595600 41.200000 1691.435977 41.127071 -50 149.7 0.268212 -0.268778 
# Funciones auxiliares para gráficas
eps_global <- 1e-6

make_nonneg <- function(x) {
  shift <- 0
  x2 <- x
  if (any(x2 < 0, na.rm = TRUE)) {
    shift <- abs(min(x2, na.rm = TRUE)) + eps_global
    x2 <- x2 + shift
  }
  x2[x2 == 0] <- x2[x2 == 0] + eps_global
  list(x = x2, shift = shift)
}

# Funciones PDF/CDF/QF (en escala original; shift interno si aplica)
make_norm_funs <- function(mu, sd) list(
  pdf = function(x) dnorm(x, mu, sd),
  cdf = function(q) pnorm(q, mu, sd),
  qf  = function(p) qnorm(p, mu, sd)
)

make_exp_funs <- function(rate, shift = 0) list(
  pdf = function(x) dexp(x + shift, rate = rate),
  cdf = function(q) pexp(q + shift, rate = rate),
  qf  = function(p) qexp(p, rate = rate) - shift
)

make_lnorm_funs <- function(meanlog, sdlog, shift = 0) list(
  pdf = function(x) dlnorm(x + shift, meanlog = meanlog, sdlog = sdlog),
  cdf = function(q) plnorm(q + shift, meanlog = meanlog, sdlog = sdlog),
  qf  = function(p) qlnorm(p, meanlog = meanlog, sdlog = sdlog) - shift
)

plot_hist_pdf <- function(x, v, th, pdf_list_named) {
  lab <- get_label(v)
  dfp <- tibble(x = x)
  xs <- seq(min(x), max(x), length.out = 500)

  df_lines <- purrr::imap_dfr(pdf_list_named, function(fn, nm) {
    tibble(xs = xs, y = fn(xs), modelo = nm)
  })

  ggplot(dfp, aes(x = x)) +
    geom_histogram(aes(y = after_stat(density)), bins = 30, color = "white") +
    annotate("rect", xmin = th$a, xmax = th$b, ymin = -Inf, ymax = Inf, alpha = 0.15) +
    geom_vline(xintercept = th$lim, linetype = "dashed", linewidth = 0.9) +
    geom_line(
      data = df_lines,
      aes(x = xs, y = y, color = modelo, linetype = modelo, group = modelo),
      linewidth = 0.9
    ) +
    scale_color_discrete(name = "Modelo") +
    scale_linetype_discrete(name = "Modelo") +
    labs(
      title = "Histograma + densidad teórica",
      subtitle = lab,
      x = lab,
      y = "Densidad",
      caption = paste("Sombreado: [a,b]. Línea punteada: lim.", CAP_ELAB, sep = "\n")
    )
}

plot_ecdf_vs_cdf <- function(x, v, th, cdf_list_named) {
  lab <- get_label(v)
  xs <- sort(unique(x))


  df_lines <- purrr::imap_dfr(cdf_list_named, function(fn, nm) {
    tibble(xs = xs, y = fn(xs), modelo = nm)
  })

  ggplot() +
    stat_ecdf(data = tibble(x = x), aes(x = x), linewidth = 0.9) +
    annotate("rect", xmin = th$a, xmax = th$b, ymin = -Inf, ymax = Inf, alpha = 0.10) +
    geom_vline(xintercept = th$lim, linetype = "dashed", linewidth = 0.9) +
    geom_line(
      data = df_lines,
      aes(x = xs, y = y, color = modelo, linetype = modelo, group = modelo),
      linewidth = 0.9
    ) +
    scale_color_discrete(name = "Modelo") +
    scale_linetype_discrete(name = "Modelo") +
    labs(
      title = "Comparación CDF: empírica vs teórica",
      subtitle = lab,
      x = lab,
      y = "Probabilidad acumulada",
      caption = paste("Escalón: ECDF (datos). Líneas: CDF teórica.", CAP_ELAB, sep = "\n")
    )
}

plot_qq <- function(x, v, qf, modelo_nombre) {
  lab <- get_label(v)
  x <- sort(x)
  n <- length(x)
  pp <- ppoints(n)
  theo <- qf(pp)
  ggplot(tibble(theo = theo, sample = x), aes(theo, sample)) +
    geom_point(alpha = 0.8) +
    geom_abline(intercept = 0, slope = 1, linewidth = 0.9) +
    labs(
      title = paste0("Q–Q plot (", modelo_nombre, ")"),
      subtitle = lab,
      x = "Cuantiles teóricos",
      y = "Cuantiles muestrales",
      caption = CAP_ELAB
    )
}

plot_box <- function(x, v) {
  lab <- get_label(v)

  x <- x[is.finite(x)]
  n <- length(x)
  mu <- mean(x); med <- median(x); s <- sd(x)

  qs <- quantile(x, probs = c(0, 0.25, 0.5, 0.75, 1), na.rm = TRUE)
  rng <- diff(range(x, na.rm = TRUE))
  pad <- ifelse(is.finite(rng) && rng > 0, 0.10 * rng, 1)

  stats_text <- paste0(
    "Min: ", round(qs[1], 4),
    "\nQ1: ", round(qs[2], 4),
    "\nMediana: ", round(qs[3], 4),
    "\nQ3: ", round(qs[4], 4),
    "\nMáx: ", round(qs[5], 4),
    "\nMedia (♦): ", round(mu, 4),
    "\nSD: ", round(s, 4)
  )

  dfp <- tibble(variable = lab, valor = x)

  ggplot(dfp, aes(y = variable, x = valor)) +
    geom_boxplot(width = 0.35, outlier.alpha = 0.35) +
    stat_summary(fun = mean, geom = "point", shape = 18, size = 3) +
    annotate(
      "label",
      y = lab,
      x = as.numeric(qs[5]) + pad,
      label = stats_text,
      hjust = 0,
      vjust = 0.5,
      size = 3,
      label.size = 0.25
    ) +
    coord_cartesian(clip = "off") +
    scale_x_continuous(expand = expansion(mult = c(0.05, 0.35))) +
    labs(
      title = paste0("Boxplot (horizontal) — ", lab),
      subtitle = paste0("n=", n, " | media=", round(mu, 4), " | mediana=", round(med, 4), " | sd=", round(s, 4), " | ♦ = media"),
      x = "Valor",
      y = lab,
      caption = CAP_ELAB
    ) +
    theme(
      # Evita redundancia: como solo hay 1 categoría, el texto del eje Y repite el nombre de la variable
      axis.text.y = element_blank(),
      axis.ticks.y = element_blank(),
      plot.margin = margin(5.5, 70, 5.5, 5.5)
    )
}


for (v in vars_unicas) {
  x <- df[[v]]
  x <- x[is.finite(x)]
  th <- thresholds |> dplyr::filter(variable == v) |> dplyr::slice(1)

  cat("\n\n## ", get_label(v), "\n\n", sep="")

  # Comentario compacto
  cat('<div class="kpi"><b>Resumen rápido:</b> ',
      'n=', length(x),
      ' | rango=[', round(min(x),4), ', ', round(max(x),4), ']',
      ' | asimetría=', round(psych::skew(x),4),
      '</div>', sep = "")

  if (tolower(v) == tolower("Material_Loss_Norm")) {
    nn <- make_nonneg(x)
    x2 <- nn$x; sh <- nn$shift
    rate0 <- 1/mean(x2)
    logx <- log(x2)
    meanlog0 <- mean(logx)
    sdlog0 <- sqrt(mean((logx - meanlog0)^2))

    p1 <- plot_hist_pdf(
      x, v, th,
      pdf_list_named = list(
        "Exponencial" = make_exp_funs(rate0, shift = sh)$pdf,
        "Lognormal"   = make_lnorm_funs(meanlog0, sdlog0, shift = sh)$pdf
      )
    )
    print(p1)

    p2 <- plot_ecdf_vs_cdf(
      x, v, th,
      cdf_list_named = list(
        "Exponencial" = make_exp_funs(rate0, shift = sh)$cdf,
        "Lognormal"   = make_lnorm_funs(meanlog0, sdlog0, shift = sh)$cdf
      )
    )
    print(p2)

  } else {
    mu0 <- mean(x); sd0 <- sd(x)
    p1 <- plot_hist_pdf(
      x, v, th,
      pdf_list_named = list("Normal" = make_norm_funs(mu0, sd0)$pdf)
    )
    print(p1)

    p2 <- plot_ecdf_vs_cdf(
      x, v, th,
      cdf_list_named = list("Normal" = make_norm_funs(mu0, sd0)$cdf)
    )
    print(p2)
  }

  print(plot_box(x, v))
}

3.1 Pérdida de material normalizada (Material_Loss_Norm)

Resumen rápido: n=1000 | rango=[0, 1] | asimetría=1.9381

3.2 Temperatura (Temperature)

Resumen rápido: n=1000 | rango=[-50, 149.7] | asimetría=0.2682

4 Conjetura del modelo y ajuste en R

fits <- list()

# -------- Material_Loss_Norm: Exponencial y Lognormal --------
v <- model_map$var_real[model_map$var_solicitada == "Material_Loss_Norm"][1]
x <- df[[v]]
x <- x[is.finite(x)]

cat("\n\n## ", get_label(v), " — Modelos: **Exponencial** y **Lognormal**\n\n", sep="")

4.1 Pérdida de material normalizada (Material_Loss_Norm) — Modelos: Exponencial y Lognormal

nn <- make_nonneg(x)
x2 <- nn$x
shift <- nn$shift

if (shift > 0) {
  cat("**Nota de soporte:** se aplicó desplazamiento \\(x' = x + \\delta\\) con \\(\\delta = ",
      round(shift,6),
      "\\) para asegurar \\(x'>0\\). Las probabilidades se reportan en la escala original.\n\n", sep="")
}

# Exponencial + MLE
rate_hat <- 1/mean(x2)
loglik_exp <- sum(dexp(x2, rate = rate_hat, log = TRUE))

cat("### Conjetura y parámetros — Exponencial\n\n")

4.1.1 Conjetura y parámetros — Exponencial

cat("- \\(\\hat\\lambda = 1/\\bar x' =\\) ", round(rate_hat,6), "\n\n", sep="")
  • \(\hat\lambda = 1/\bar x' =\) 6.82851
fits[[paste0(v,"__exp")]] <- list(
  variable = v, modelo = "exp", k = 1,
  pars = c(rate = rate_hat),
  loglik = loglik_exp,
  funs = make_exp_funs(rate_hat, shift = shift)
)

kbl2(
  tibble(parametro = "lambda", estimacion = fmt_num(rate_hat, 6)),
  caption = "Ajuste final: Exponencial"
)
Ajuste final: Exponencial
parametro estimacion
lambda 6.82851 
# Lognormal + MLE
logx <- log(x2)
meanlog_hat <- mean(logx)
sdlog_hat_mle <- sqrt(mean((logx - meanlog_hat)^2))  # MLE (denominador n)
sdlog_hat_s <- sd(logx)                              # s muestral

cat("### Conjetura y parámetros — Lognormal\n\n")

4.1.2 Conjetura y parámetros — Lognormal

cat("Sea \\(Y=\\ln(X')\\). Entonces:\n\n")

Sea \(Y=\ln(X')\). Entonces:

cat("- \\(\\hat\\mu_L = \\overline{\\ln x'} =\\) ", round(meanlog_hat,6), "\n", sep="")
  • \(\hat\mu_L = \overline{\ln x'} =\) -2.4821
cat("- \\(\\hat\\sigma_L \\approx sd(\\ln x') =\\) ", round(sdlog_hat_s,6), " (muestral)\n\n", sep="")
  • \(\hat\sigma_L \approx sd(\ln x') =\) 1.29581 (muestral)
loglik_lnorm <- sum(dlnorm(x2, meanlog = meanlog_hat, sdlog = sdlog_hat_mle, log = TRUE))

fits[[paste0(v,"__log_normal")]] <- list(
  variable = v, modelo = "log_normal", k = 2,
  pars = c(meanlog = meanlog_hat, sdlog = sdlog_hat_mle),
  loglik = loglik_lnorm,
  funs = make_lnorm_funs(meanlog_hat, sdlog_hat_mle, shift = shift)
)

kbl2(
  tibble(parametro = c("meanlog","sdlog (MLE)"),
         estimacion = fmt_num(c(meanlog_hat, sdlog_hat_mle), 6)),
  caption = "Ajuste final: Lognormal"
)
Ajuste final: Lognormal
parametro estimacion
meanlog -2.48210
sdlog (MLE) 1.29516 
# -------- Temperature: Normal --------
vT <- model_map$var_real[model_map$var_solicitada == "Temperature"][1]
xT <- df[[vT]]
xT <- xT[is.finite(xT)]

cat("\n\n## ", get_label(vT), " — Modelo: **Normal**\n\n", sep="")

4.2 Temperatura (Temperature) — Modelo: Normal

mu_hat <- mean(xT)
sd_hat <- sd(xT)
loglik_norm <- sum(dnorm(xT, mean = mu_hat, sd = sd_hat, log = TRUE))

cat("### Conjetura y parámetros — Normal\n\n")

4.2.1 Conjetura y parámetros — Normal

cat("- \\(\\hat\\mu = \\bar x =\\) ", round(mu_hat,6),
    ",  \\(\\hat\\sigma = s =\\) ", round(sd_hat,6), "\n\n", sep="")
  • \(\hat\mu = \bar x =\) 42.5956, \(\hat\sigma = s =\) 41.1271
fits[[paste0(vT,"__norm")]] <- list(
  variable = vT, modelo = "norm", k = 2,
  pars = c(mean = mu_hat, sd = sd_hat),
  loglik = loglik_norm,
  funs = make_norm_funs(mu_hat, sd_hat)
)

kbl2(
  tibble(parametro = c("mu","sigma"), estimacion = fmt_num(c(mu_hat, sd_hat), 6)),
  caption = "Ajuste final: Normal"
)
Ajuste final: Normal
parametro estimacion
mu 42.5956
sigma 41.1271

5 Bondad de ajuste (Chi-cuadrado de Pearson) + AIC

Interpretación del Chi-cuadrado:
- Si \(p\text{-valor} > 0.05\), no se rechaza que la distribución teórica se ajusta a los datos (ajuste aceptable).
- Si \(p\text{-valor} \le 0.05\), hay evidencia de desajuste.
Se usa además AIC como criterio de desempate (menor AIC es mejor).

make_bins <- function(x, rule = c("FD","Sturges")) {
  rule <- match.arg(rule)
  h <- hist(x, plot = FALSE, breaks = rule)
  br <- unique(h$breaks)
  br <- sort(br)
  if (length(br) < 5) br <- pretty(range(x), n = 6)
  br
}

combine_bins_expected_ge5 <- function(obs, exp, breaks) {
  i <- 1
  while (any(exp < 5) && length(exp) > 2) {
    idx <- which.min(exp)
    if (idx == 1) {
      obs[2] <- obs[1] + obs[2]
      exp[2] <- exp[1] + exp[2]
      obs <- obs[-1]; exp <- exp[-1]
      breaks <- breaks[-2]
    } else {
      obs[idx-1] <- obs[idx-1] + obs[idx]
      exp[idx-1] <- exp[idx-1] + exp[idx]
      obs <- obs[-idx]; exp <- exp[-idx]
      breaks <- breaks[-(idx+1)]
    }
    i <- i + 1
    if (i > 500) break
  }
  list(obs = obs, exp = exp, breaks = breaks)
}

chisq_gof <- function(x, funs, k_params, breaks_rule = "FD") {
  x <- x[is.finite(x)]
  n <- length(x)

  br <- make_bins(x, rule = breaks_rule)
  obs <- as.numeric(table(cut(x, breaks = br, include.lowest = TRUE, right = TRUE)))

  p_bin <- funs$cdf(br[-1]) - funs$cdf(br[-length(br)])
  exp <- n * p_bin

  comb <- combine_bins_expected_ge5(obs, exp, br)
  obs2 <- comb$obs; exp2 <- comb$exp

  chi <- sum((obs2 - exp2)^2 / exp2)
  k_bins <- length(obs2)
  gl <- k_bins - 1 - k_params

  if (gl <= 0) return(tibble(chi2=chi, gl=gl, p_value=NA_real_, bins=k_bins,
                             nota="gl<=0: demasiados parámetros o pocos bins."))

  pval <- stats::pchisq(chi, df = gl, lower.tail = FALSE)
  tibble(chi2=chi, gl=gl, p_value=pval, bins=k_bins, nota="")
}
gof_tbl <- purrr::map_dfr(names(fits), function(key) {
  fit <- fits[[key]]
  v <- fit$variable
  x <- df[[v]]
  x <- x[is.finite(x)]

  out <- chisq_gof(x, funs = fit$funs, k_params = fit$k, breaks_rule = "FD")
  aic <- 2*fit$k - 2*fit$loglik

  tibble(
    variable = v,
    etiqueta = get_label(v),
    modelo = fit$modelo,
    bins = out$bins,
    chi2 = out$chi2,
    gl = out$gl,
    p_value = out$p_value,
    AIC = aic,
    ajuste = case_when(
      is.na(out$p_value) ~ "NA",
      out$p_value > 0.05 ~ "Aceptable",
      TRUE ~ "Débil"
    ),
    nota = out$nota
  )
}) |>
  dplyr::mutate(
    chi2 = fmt_num(chi2, 6),
    p_value = fmt_num(p_value, 6),
    AIC = fmt_num(AIC, 3)
  ) |>
  dplyr::arrange(variable, AIC)

kbl2(gof_tbl, caption = "Bondad de ajuste: Chi-cuadrado (bins FD combinados hasta Exp≥5) + AIC.")
Bondad de ajuste: Chi-cuadrado (bins FD combinados hasta Exp≥5) + AIC.
variable etiqueta modelo bins chi2 gl p_value AIC ajuste nota
Material_Loss_Norm Pérdida de material normalizada (Material_Loss_Norm) exp 28 47.0266 26 0.006987 -1840.21 Débil
Material_Loss_Norm Pérdida de material normalizada (Material_Loss_Norm) log_normal 32 117.7290 29 0.000000 -1605.04 Débil
Temperature Temperatura (Temperature) norm 19 40.1191 16 0.000748 10274.21 Débil

6 Selección del modelo

# Regla:
# 1) Preferir mayor p_value (Chi-cuadrado)
# 2) Desempate por AIC mínimo
select_model <- function(tbl) {
  if (nrow(tbl) == 1) return(tbl[1, ])
  tbl2 <- tbl |>
    dplyr::mutate(p_ord = ifelse(is.na(p_value), -Inf, p_value)) |>
    dplyr::arrange(dplyr::desc(p_ord), AIC)
  tbl2[1, ] |> dplyr::select(-p_ord)
}

sel_tbl <- gof_tbl |>
  dplyr::group_by(variable) |>
  dplyr::group_modify(~select_model(.x)) |>
  dplyr::ungroup()

kbl2(sel_tbl, caption = "Modelo elegido por variable (mayor p-valor; desempate por AIC mínimo).")
Modelo elegido por variable (mayor p-valor; desempate por AIC mínimo).
variable etiqueta modelo bins chi2 gl p_value AIC ajuste nota
Material_Loss_Norm Pérdida de material normalizada (Material_Loss_Norm) exp 28 47.0266 26 0.006987 -1840.21 Débil
Temperature Temperatura (Temperature) norm 19 40.1191 16 0.000748 10274.21 Débil

7 6) Predicción probabilística (modelo elegido)

# N óptimo / tamaño operativo:
N_opt <- if (is.null(params$N_opt)) nrow(df) else as.numeric(params$N_opt)

# Mapa key elegido por variable
key_by_var <- purrr::map_chr(sel_tbl$variable, function(v) {
  m <- sel_tbl |> dplyr::filter(variable == v) |> dplyr::pull(modelo) |> dplyr::first()
  key <- names(fits)[purrr::map_lgl(names(fits), ~ {
    f <- fits[[.x]]
    tolower(f$variable) == tolower(v) && tolower(f$modelo) == tolower(m)
  })]
  key[1]
})
names(key_by_var) <- sel_tbl$variable

pred_tbl <- purrr::map_dfr(sel_tbl$variable, function(v) {
  key <- key_by_var[[v]]
  fit <- fits[[key]]
  funs <- fit$funs
  th <- thresholds |> dplyr::filter(variable == v) |> dplyr::slice(1)

  p_window <- funs$cdf(th$b) - funs$cdf(th$a)
  p_low <- funs$cdf(th$lim)

  tibble(
    variable = v,
    etiqueta = get_label(v),
    modelo_elegido = fit$modelo,
    parametros = paste(names(fit$pars), "=", fmt_num(fit$pars, 6), collapse = "; "),
    p_valor_chi2 = sel_tbl |> dplyr::filter(variable == v) |> dplyr::pull(p_value) |> dplyr::first(),
    P_a_b = p_window,
    P_menor_lim = p_low,
    Esperado_menor_lim_en_Nopt = N_opt * p_low
  )
}) |>
  dplyr::mutate(
    P_a_b = fmt_num(P_a_b, 6),
    P_menor_lim = fmt_num(P_menor_lim, 6),
    Esperado_menor_lim_en_Nopt = fmt_num(Esperado_menor_lim_en_Nopt, 6)
  )

kbl2(pred_tbl, caption = "Predicción probabilística con el modelo elegido")
Predicción probabilística con el modelo elegido
variable etiqueta modelo_elegido parametros p_valor_chi2 P_a_b P_menor_lim Esperado_menor_lim_en_Nopt
Material_Loss_Norm Pérdida de material normalizada (Material_Loss_Norm) exp rate = 6.828512 0.006987 0.445225 0.105597 105.597
Temperature Temperatura (Temperature) norm mean = 42.5956; sd = 41.127071 0.000748 0.501866 0.103811 103.811 
long_probs <- pred_tbl |>
  tidyr::pivot_longer(cols = c(P_a_b, P_menor_lim), names_to = "evento", values_to = "prob") |>
  mutate(
    evento = recode(evento, P_a_b = "P(a ≤ X ≤ b)", P_menor_lim = "P(X < lim)"),
    etiqueta = factor(etiqueta, levels = unique(pred_tbl$etiqueta))
  )

ggplot(long_probs, aes(x = etiqueta, y = prob)) +
  geom_col() +
  coord_flip() +
  facet_wrap(~evento, ncol = 1, scales = "free_x") +
  scale_y_continuous(labels = percent_format(accuracy = 1)) +
  labs(
    title = "Resumen visual de probabilidades (modelo elegido)",
    x = NULL,
    y = "Probabilidad",
    caption = CAP_ELAB
  )

8 Cálculo de Probabilidades y Toma de Decisiones

8.0.1 Interrogantes técnicas — Pérdida de material normalizada (Material_Loss_Norm)

Pregunta 1: ¿Cuál es la probabilidad de que, al seleccionar aleatoriamente un registro, el valor de Pérdida de material normalizada (Material_Loss_Norm) se encuentre entre 0.0489 y 0.1913?

Pregunta 2: Si se planifica una campaña sobre la muestra optimizada (N=1000), ¿cuántos registros se estima que presentarán Pérdida de material normalizada (Material_Loss_Norm) menor que 0.0163?

Pregunta 3: Aplicando el teorema del límite central, ¿en qué intervalo se estima la media poblacional de Pérdida de material normalizada (Material_Loss_Norm) con una confianza del 95%?

CONCLUSIONES: La variable Pérdida de material normalizada (Material_Loss_Norm) sigue un modelo de probabilidad exponencial, aprobando el test Chi-cuadrado de Pearson (p-valor = 0.007). De esta manera, logramos calcular probabilidades como, por ejemplo, que al seleccionar aleatoriamente un registro la probabilidad de que Pérdida de material normalizada (Material_Loss_Norm) se encuentre entre 0.0489 y 0.1913 es de 44.52%; y la probabilidad de que Pérdida de material normalizada (Material_Loss_Norm) sea menor que 0.0163 es de 10.56%, lo que implica un estimado de 106 registros en la muestra optimizada (N=1000). Finalmente, aplicando el teorema del límite central, la media aritmética poblacional se estima entre 0.14 y 0.16 con una confianza del 95%.

8.0.2 Interrogantes técnicas — Temperatura (Temperature)

Pregunta 1: ¿Cuál es la probabilidad de que, al seleccionar aleatoriamente un registro, el valor de Temperatura (Temperature) se encuentre entre 13.4 y 69.15?

Pregunta 2: Si se planifica una campaña sobre la muestra optimizada (N=1000), ¿cuántos registros se estima que presentarán Temperatura (Temperature) menor que -9.23?

Pregunta 3: Aplicando el teorema del límite central, ¿en qué intervalo se estima la media poblacional de Temperatura (Temperature) con una confianza del 95%?

CONCLUSIONES: La variable Temperatura (Temperature) sigue un modelo de probabilidad normal, aprobando el test Chi-cuadrado de Pearson (p-valor = 0.0007). De esta manera, logramos calcular probabilidades como, por ejemplo, que al seleccionar aleatoriamente un registro la probabilidad de que Temperatura (Temperature) se encuentre entre 13.4 y 69.15 es de 50.19%; y la probabilidad de que Temperatura (Temperature) sea menor que -9.23 es de 10.38%, lo que implica un estimado de 104 registros en la muestra optimizada (N=1000). Finalmente, aplicando el teorema del límite central, la media aritmética poblacional se estima entre 40.05 y 45.14 con una confianza del 95%.

9 Interpretación por variable

for (v in pred_tbl$variable) {
  th <- thresholds |> dplyr::filter(variable == v) |> dplyr::slice(1)
  fila <- pred_tbl |> dplyr::filter(variable == v) |> dplyr::slice(1)
  modelo <- fila$modelo_elegido

  cat("\n\n## ", fila$etiqueta, "\n\n", sep="")

  # KPI tipo tarjeta
  cat('<div class="kpi">',
      '<b>Modelo elegido:</b> ', modelo,
      '<br><span class="small-note">Parámetros: ', fila$parametros, '</span>',
      '<br><span class="small-note">p-valor (Chi²): ', fmt_num(fila$p_valor_chi2, 4), '</span>',
      '</div>', sep="")

  cat("**Interpretación operativa:**\n\n")
  cat("- **P(a ≤ X ≤ b)** con a=", round(th$a,4), ", b=", round(th$b,4),
      " es **", fmt_pct(fila$P_a_b, 2), "**.\n", sep="")
  cat("- **P(X < lim)** con lim=", round(th$lim,4),
      " es **", fmt_pct(fila$P_menor_lim, 2), "**.\n", sep="")
  cat("- Para **N_opt = ", N_opt, "**, el esperado con **X < lim** es **",
      round(fila$Esperado_menor_lim_en_Nopt,2), "**.\n\n", sep="")

  # Diagnóstico visual del modelo elegido (Q–Q + CDF vs ECDF)
  key <- key_by_var[[v]]
  fit <- fits[[key]]
  x <- df[[v]]
  x <- x[is.finite(x)]

  print(plot_ecdf_vs_cdf(x, v, th, cdf_list_named = list("Modelo elegido" = fit$funs$cdf)))
  print(plot_qq(x, v, fit$funs$qf, modelo_nombre = modelo))
}

9.1 Pérdida de material normalizada (Material_Loss_Norm)

Modelo elegido: exp
Parámetros: rate = 6.828512
p-valor (Chi²): 0.007

Interpretación operativa:

  • P(a ≤ X ≤ b) con a=0.0489, b=0.1913 es 44.52%.
  • P(X < lim) con lim=0.0163 es 10.56%.
  • Para N_opt = 1000, el esperado con X < lim es 105.6.

9.2 Temperatura (Temperature)

Modelo elegido: norm
Parámetros: mean = 42.5956; sd = 41.127071
p-valor (Chi²): 0.0007

Interpretación operativa:

  • P(a ≤ X ≤ b) con a=13.4, b=69.15 es 50.19%.
  • P(X < lim) con lim=-9.23 es 10.38%.
  • Para N_opt = 1000, el esperado con X < lim es 103.81.