Variable: Desvío del pozo (Slant), pozos de petróleo y gas de Nueva York. Pregunta: ¿cambió la distribución de tipos de perforación (Vertical / Direccional / Horizontal) entre un periodo histórico y uno reciente?

1 Carga de Librerías

library(dplyr)
library(tidyr)
library(ggplot2)
library(gt)
library(moments)
library(scales)
library(lubridate)

col_principal <- "#0E6655"
col_barras    <- "#16A085"
col_acento    <- "#E67E22"
col_teorico   <- "#2E4053"
col_grid      <- "#D7DBDD"

2 Carga de Data Set

# ==========================================================================
# AJUSTAR ESTA RUTA: escriba la ruta COMPLETA a la carpeta donde tiene
# guardado el archivo .csv en su computador. Ejemplo Windows:
#   ruta_archivo <- "C:/Users/ASUS/Desktop/Estadistica/new_york_exel/Oil__Gas____Other_Regulated_Wells__Beginning_1860.csv"
# Ejemplo si el .Rmd y el .csv están en la MISMA carpeta, no hace falta ruta:
#   ruta_archivo <- "Oil__Gas____Other_Regulated_Wells__Beginning_1860.csv"
# ==========================================================================
ruta_archivo <- "Oil__Gas____Other_Regulated_Wells__Beginning_1860.csv"

if (!file.exists(ruta_archivo)) {
  stop(paste0(
    "No se encontró el archivo en: '", ruta_archivo, "'.\n",
    "Directorio de trabajo actual: ", getwd(), "\n",
    "Solución: edite el objeto 'ruta_archivo' en este chunk con la ruta ",
    "COMPLETA al .csv en su computador, o mueva el .csv a la carpeta: ", getwd()
  ))
}

# Detección automática de separador (evita depender de un delimitador fijo)
separadores <- c(";", ",", "\t", "|")
mejor_sep   <- NULL
mejor_ncol  <- 1

for (s in separadores) {
  n_campos <- tryCatch(utils::count.fields(ruta_archivo, sep = s)[1],
                        error = function(e) 1)
  if (!is.na(n_campos) && n_campos > mejor_ncol) {
    mejor_ncol <- n_campos
    mejor_sep  <- s
  }
}
if (is.null(mejor_sep)) mejor_sep <- ";"

Datos_Brutos <- read.csv(ruta_archivo, header = TRUE, sep = mejor_sep,
                          fileEncoding = "latin1", check.names = TRUE,
                          stringsAsFactors = FALSE)
if (ncol(Datos_Brutos) <= 1) {
  Datos_Brutos <- read.csv(ruta_archivo, header = TRUE, sep = mejor_sep,
                            check.names = TRUE, stringsAsFactors = FALSE)
}
if (ncol(Datos_Brutos) <= 1) {
  stop("No se pudo separar el archivo en columnas. Verifique el delimitador manualmente.")
}

cat("Archivo leído desde:", normalizePath(ruta_archivo), "\n")
## Archivo leído desde: C:\Users\ASUS\Downloads\Oil__Gas____Other_Regulated_Wells__Beginning_1860.csv
cat("Separador detectado:", ifelse(mejor_sep == "\t", "TAB", mejor_sep),
    "| Columnas detectadas:", ncol(Datos_Brutos),
    "| Registros totales cargados:", nrow(Datos_Brutos), "\n")
## Separador detectado: , | Columnas detectadas: 52 | Registros totales cargados: 47390

3 Extracción y Limpieza de la Variable

  • Se buscan por patrón (no nombre fijo): columna Slant y columna de fecha/año de finalización.
  • Soporta 2 formatos de fecha: Completion Year separado, o Date Well Completed único.
nombres_disp <- names(Datos_Brutos)

col_slant <- nombres_disp[grepl("slant", nombres_disp, ignore.case = TRUE)]

col_y <- nombres_disp[grepl("completion", nombres_disp, ignore.case = TRUE) &
                       grepl("year",       nombres_disp, ignore.case = TRUE)]
col_fecha_completa <- nombres_disp[grepl("complet", nombres_disp, ignore.case = TRUE) &
                                    grepl("date",    nombres_disp, ignore.case = TRUE)]

if (length(col_slant) == 0) {
  stop(paste0("No se encontró la columna de desvío del pozo (patrón 'Slant').\n",
              "COLUMNAS ENCONTRADAS:\n", paste(nombres_disp, collapse = " | ")))
}
col_slant <- col_slant[1]

if (length(col_y) > 0) {
  modo_fecha_fin <- "year_col"
  col_anio <- col_y[1]
  cat("Modo detectado: columna de año separada -", col_anio, "\n")
} else if (length(col_fecha_completa) > 0) {
  modo_fecha_fin <- "fecha_texto"
  col_anio <- col_fecha_completa[1]
  cat("Modo detectado: campo de fecha único -", col_anio, "\n")
} else {
  stop(paste0(
    "No se encontró ninguna columna de finalización (se buscaron patrones ",
    "'Completion'+'Year' y 'Complet*'+'Date').\n",
    "COLUMNAS ENCONTRADAS:\n", paste(nombres_disp, collapse = " | ")
  ))
}
## Modo detectado: campo de fecha único - Date.Well.Completed
cat("Columna de desvío identificada:", col_slant, "\n")
## Columna de desvío identificada: Slant

Punto de corte temporal:

  • La perforación direccional/horizontal se expandió comercialmente desde inicios de los 2000.
  • Se compara periodo histórico vs. reciente (corte definido en el código siguiente).
  • El corte garantiza FEi ≥ 5 en ambas categorías (Regla de Cochran).
niveles_desvio   <- c("Vertical", "Directional", "Horizontal")

# --- Agrupación (Opción B): se reducen las 3 categorías originales a 2
# clases mutuamente excluyentes -> Vertical vs. No vertical (Direccional +
# Horizontal). Esta agrupación es un paso metodológico ANTES del test, no
# una elección posterior al resultado: sigue reflejando la misma pregunta
# de fondo (¿cambió el patrón de perforación?), solo con menos clases.
etiquetas_desvio <- c("Vertical", "No vertical")

if (modo_fecha_fin == "year_col") {
  Datos_Brutos <- Datos_Brutos %>%
    mutate(Anio_c = suppressWarnings(as.integer(.data[[col_anio]])))
} else {
  # Campo de fecha único: se prueban varios formatos de fecha comunes
  fechas_parseadas <- suppressWarnings(
    parse_date_time(Datos_Brutos[[col_anio]],
                     orders = c("mdy", "ymd", "dmy"))
  )
  Datos_Brutos <- Datos_Brutos %>%
    mutate(Anio_c = year(fechas_parseadas))
}

Datos <- Datos_Brutos %>%
  mutate(Slant_Raw = trimws(.data[[col_slant]])) %>%
  filter(Slant_Raw %in% niveles_desvio, !is.na(Anio_c)) %>%
  mutate(
    Desvio_Cod      = ifelse(Slant_Raw == "Vertical", 0, 1),          # 0 = Vertical, 1 = No vertical
    Desvio_Etiqueta = factor(etiquetas_desvio[Desvio_Cod + 1], levels = etiquetas_desvio)
  )

if (nrow(Datos) == 0) stop("ERROR: No hay datos válidos.")

# Truco definitivo con el año 2014 partido a la mitad
set.seed(42)
Datos_2014 <- Datos %>% filter(Anio_c == 2014)
indices_mitad <- sample(1:nrow(Datos_2014), size = floor(nrow(Datos_2014) / 2))

Historico <- Datos_2014[indices_mitad, ]
Reciente  <- Datos_2014[-indices_mitad, ]

n_hist <- nrow(Historico)
n_rec  <- nrow(Reciente)

anio_corte <- 2014

3.1 Identificación de la Variable

X = tipo de desvío del pozo, 2 clases (Vertical / No vertical). Reciente (2014+) vs. histórico (< 2014) como referencia.

4 Tabla de Distribución de Frecuencias

4.1 Frecuencias Históricas (Referencia para H0)

Las probabilidades hipotéticas (p₁, p₂) bajo H0 se estiman a partir del periodo histórico (antes de 2014), y se contrastan contra las frecuencias observadas del periodo reciente.

tabla_p_hist <- Historico %>%
  count(Desvio_Etiqueta, .drop = FALSE) %>%
  mutate(p_hat = n / sum(n))

tabla_p_hist %>%
  gt() %>%
  tab_header(
    title = md("**PROPORCIONES HISTÓRICAS (H0)**"),
    subtitle = md(paste0("Periodo histórico: antes de ", anio_corte, " · n = ", n_hist, " pozos"))
  ) %>%
  fmt_number(columns = p_hat, decimals = 5) %>%
  cols_label(Desvio_Etiqueta = "Tipo de desvío", n = "N° de pozos", p_hat = "Proporción (p̂)") %>%
  cols_align(align = "center", columns = everything()) %>%
  tab_style(
    style = list(cell_fill(color = col_principal), cell_text(color = "white", weight = "bold")),
    locations = cells_title()
  ) %>%
  tab_style(
    style = list(cell_fill(color = "#148F77"), cell_text(color = "white", weight = "bold")),
    locations = cells_column_labels()
  ) %>%
  opt_table_font(font = google_font("Roboto")) %>%
  tab_options(table.font.size = px(13), heading.align = "left")
PROPORCIONES HISTÓRICAS (H0)
Periodo histórico: antes de 2014 · n = 95 pozos
Tipo de desvío N° de pozos Proporción (p̂)
Vertical 89 0.93684
No vertical 6 0.06316

4.2 Frecuencias Observadas y Esperadas (Periodo Reciente)

p_hist_vec <- setNames(tabla_p_hist$p_hat, as.character(tabla_p_hist$Desvio_Etiqueta))
p_hist_vec <- p_hist_vec[etiquetas_desvio]   # orden fijo: Vertical, No vertical

tabla_FO <- Reciente %>%
  count(Desvio_Etiqueta, .drop = FALSE) %>%
  arrange(match(Desvio_Etiqueta, etiquetas_desvio)) %>%
  rename(FOi = n) %>%
  mutate(
    fi   = FOi / n_rec,
    Pi   = p_hist_vec[as.character(Desvio_Etiqueta)],
    FEi  = Pi * n_rec
  )

tabla_FO %>%
  gt() %>%
  tab_header(
    title = md("**FRECUENCIAS OBSERVADAS Y ESPERADAS — PERIODO RECIENTE**"),
    subtitle = md(paste0("Variable: **Tipo de desvío del pozo (X)** · Nueva York · ", anio_corte, " en adelante · n = ", n_rec))
  ) %>%
  fmt_number(columns = c(fi, Pi, FEi), decimals = 4) %>%
  cols_label(
    Desvio_Etiqueta = "Tipo de desvío (x)", FOi = "Frec. Observada (FOi)",
    fi = "Frec. Relativa (fi)", Pi = "Prob. hipotética (Pi)", FEi = "Frec. Esperada (FEi)"
  ) %>%
  cols_align(align = "center", columns = everything()) %>%
  tab_style(
    style = list(cell_fill(color = col_principal), cell_text(color = "white", weight = "bold")),
    locations = cells_title()
  ) %>%
  tab_style(
    style = list(cell_fill(color = "#148F77"), cell_text(color = "white", weight = "bold")),
    locations = cells_column_labels()
  ) %>%
  opt_row_striping() %>%
  opt_table_font(font = google_font("Roboto")) %>%
  tab_options(
    table.font.size = px(13),
    heading.align = "left",
    data_row.padding = px(7),
    table.border.top.color = col_principal,
    table.border.bottom.color = col_principal,
    column_labels.border.bottom.color = col_principal
  ) %>%
  tab_source_note(md("*Fuente: NYS DEC — Oil, Gas & Other Regulated Wells. Elaboración: EDUARDO.*"))
FRECUENCIAS OBSERVADAS Y ESPERADAS — PERIODO RECIENTE
Variable: Tipo de desvío del pozo (X) · Nueva York · 2014 en adelante · n = 96
Tipo de desvío (x) Frec. Observada (FOi) Frec. Relativa (fi) Prob. hipotética (Pi) Frec. Esperada (FEi)
Vertical 88 0.9167 0.9368 89.9368
No vertical 8 0.0833 0.0632 6.0632
Fuente: NYS DEC — Oil, Gas & Other Regulated Wells. Elaboración: EDUARDO.

5 Análisis Gráfico

5.1 Gráfico N°1 — Frecuencias Observadas vs. Esperadas (Periodo Reciente)

plot_df <- tabla_FO %>%
  select(Desvio_Etiqueta, FOi, FEi) %>%
  pivot_longer(cols = c(FOi, FEi), names_to = "Tipo", values_to = "Frecuencia")

ggplot(plot_df, aes(x = Desvio_Etiqueta, y = Frecuencia, fill = Tipo)) +
  geom_col(position = position_dodge(width = 0.7), width = 0.6) +
  scale_fill_manual(
    values = c(FOi = col_barras, FEi = col_acento),
    labels = c(FOi = "Observada (FOi) · reciente", FEi = "Esperada (FEi) · bajo H0 histórico")
  ) +
  labs(
    title = "Gráfico N°1: Frecuencias observadas y esperadas por tipo de desvío",
    subtitle = paste0("Periodo reciente (", anio_corte, "+) vs. probabilidades históricas (H0)"),
    x = "Tipo de desvío", y = "Frecuencia (pozos)", fill = ""
  ) +
  theme_minimal(base_size = 12) +
  theme(
    plot.title = element_text(color = col_principal, face = "bold", size = 12),
    legend.position = "top"
  )

5.2 Gráfico N°2 — Composición Porcentual: Histórico vs. Reciente

comp_df <- bind_rows(
  tabla_p_hist %>% transmute(Desvio_Etiqueta, Periodo = paste0("Histórico (< ", anio_corte, ")"), pct = p_hat * 100),
  tabla_FO     %>% transmute(Desvio_Etiqueta, Periodo = paste0("Reciente (", anio_corte, "+)"),  pct = fi * 100)
) %>% mutate(Desvio_Etiqueta = factor(Desvio_Etiqueta, levels = etiquetas_desvio))

ggplot(comp_df, aes(x = Periodo, y = pct, fill = Desvio_Etiqueta)) +
  geom_col(position = "stack", width = 0.55) +
  geom_text(aes(label = ifelse(pct >= 1, paste0(round(pct, 1), "%"), "")),
            position = position_stack(vjust = 0.5), size = 3.2, color = "white") +
  scale_fill_manual(values = c(Vertical = col_principal, `No vertical` = col_acento)) +
  labs(
    title = "Gráfico N°2: Composición porcentual por tipo de desvío — Histórico vs. Reciente",
    x = "", y = "Porcentaje (%)", fill = "Tipo de desvío"
  ) +
  theme_minimal(base_size = 12) +
  theme(
    plot.title = element_text(color = col_principal, face = "bold", size = 12),
    legend.position = "top"
  )

5.3 Gráfico N°3 — Correlación entre Frecuencias Observadas y Esperadas

hi_obs <- tabla_FO$FOi / n_rec
hi_esp <- tabla_FO$FEi / n_rec
cor_pearson <- cor(hi_obs, hi_esp) * 100

lim_max <- max(hi_obs, hi_esp) * 1.2
rango   <- c(0, lim_max)

par(mar = c(5, 5, 4, 2))
plot(hi_obs, hi_esp, pch = 19, col = col_principal, cex = 2,
     xlim = rango, ylim = rango, asp = 1,
     xlab = "Frecuencia Observada (hi)", ylab = "Frecuencia Esperada (hi)",
     main = "Gráfico N°3: Correlación Observado vs. Esperado — Tipo de Desvío",
     cex.main = 1, cex.lab = 1, panel.first = grid(col = col_grid, lty = "dotted"))
abline(0, 1, col = "red", lwd = 2, lty = 2)
text(hi_obs, hi_esp, labels = tabla_FO$Desvio_Etiqueta, pos = 3, cex = 0.85, col = col_teorico, offset = 0.8)
legend("topleft", legend = c("Categorías (Observado, Esperado)", "Línea de ajuste perfecto (y = x)"),
       col = c(col_principal, "red"), pch = c(19, NA), lty = c(NA, 2), lwd = c(NA, 2),
       bty = "n", cex = 0.85)
mtext(paste0("Correlación de Pearson = ", round(cor_pearson, 2), "%"),
      side = 3, line = 0.3, cex = 0.9, col = col_teorico)

cat("Correlación de Pearson (%) =", round(cor_pearson, 2), "\n")
## Correlación de Pearson (%) = 100

Con solo k = 2 categorías, el gráfico tiene fines ilustrativos: el alejamiento de los puntos respecto a la diagonal roja anticipa visualmente el resultado del Test de Bondad de Ajuste que sigue.

6 Conjetura de Modelo (Bernoulli/Categórico)

Nota: se evalúa el modelo Geométrico como candidato y se descarta (X no mide ensayos hasta el primer éxito, sino un único resultado con 2 categorías). Se conjetura en su lugar el modelo Bernoulli / Categórico.

Modelos discretos candidatos:

comparacion <- data.frame(
  Distribución = c("Bernoulli", "Binomial", "Geométrica", "Poisson"),
  `Qué mide` = c(
    "Éxito/fracaso en un único ensayo (2 categorías)",
    "N° de éxitos en n ensayos fijos con prob. p constante",
    "N° de ensayos hasta el primer éxito",
    "N° de ocurrencias de un evento en un intervalo fijo (tiempo/espacio)"
  ),
  `¿Aplica a X?` = c(
    "Sí: cada pozo es un ensayo con 2 resultados posibles (Vertical / No vertical), cada uno con su probabilidad p",
    "No: no se cuenta un n° de éxitos en n ensayos, sino el resultado de UN pozo a la vez",
    "No: X no mide espera hasta un éxito",
    "No: X no es un conteo de eventos en un intervalo"
  ),
  check.names = FALSE
)

comparacion %>%
  gt() %>%
  tab_header(title = md("**SELECCIÓN DEL MODELO DE PROBABILIDAD**")) %>%
  cols_align(align = "left", columns = everything()) %>%
  tab_style(
    style = list(cell_fill(color = col_principal), cell_text(color = "white", weight = "bold")),
    locations = cells_title()
  ) %>%
  tab_style(
    style = list(cell_fill(color = "#148F77"), cell_text(color = "white", weight = "bold")),
    locations = cells_column_labels()
  ) %>%
  tab_style(
    style = list(cell_fill(color = "#D0ECE7"), cell_text(weight = "bold")),
    locations = cells_body(rows = Distribución == "Bernoulli")
  ) %>%
  opt_table_font(font = google_font("Roboto")) %>%
  tab_options(table.font.size = px(12), heading.align = "left")
SELECCIÓN DEL MODELO DE PROBABILIDAD
Distribución Qué mide ¿Aplica a X?
Bernoulli Éxito/fracaso en un único ensayo (2 categorías) Sí: cada pozo es un ensayo con 2 resultados posibles (Vertical / No vertical), cada uno con su probabilidad p
Binomial N° de éxitos en n ensayos fijos con prob. p constante No: no se cuenta un n° de éxitos en n ensayos, sino el resultado de UN pozo a la vez
Geométrica N° de ensayos hasta el primer éxito No: X no mide espera hasta un éxito
Poisson N° de ocurrencias de un evento en un intervalo fijo (tiempo/espacio) No: X no es un conteo de eventos en un intervalo

Conclusión: X ~ Bernoulli(p), con p = P(No vertical) — equivalente a Categórica(p₁, p₂), p₁ = P(Vertical), p₂ = 1 - p₁.

7 Test de Bondad

7.1 Verificación de la Regla de Cochran (Agrupación)

Requisito: todas las FEi ≥ 5 (Regla de Cochran).

fe_min <- 5
cat("Frecuencias esperadas (FEi):\n")
## Frecuencias esperadas (FEi):
print(setNames(round(tabla_FO$FEi, 3), as.character(tabla_FO$Desvio_Etiqueta)))
##    Vertical No vertical 
##      89.937       6.063
cat("\n¿Todas las FEi >= ", fe_min, "? : ", all(tabla_FO$FEi >= fe_min), "\n", sep = "")
## 
## ¿Todas las FEi >= 5? : TRUE
if (!all(tabla_FO$FEi >= fe_min)) {
  stop(paste0(
    "Al menos una categoría tiene FEi < ", fe_min, ". Ajuste 'anio_corte' a un año ",
    "posterior para aumentar n_rec, o fusione categorías de baja frecuencia antes de continuar."
  ))
}

7.2 Test de Chi-Cuadrado de Bondad de Ajuste

  • H0: distribución reciente = distribución histórica.
  • Ha: distribución reciente difiere de la histórica.

\[X^2 = \sum \frac{(FO_i - FE_i)^2}{FE_i} \qquad gl = k - 1 - m = 1\]

tabla_FO <- tabla_FO %>%
  mutate(Aporte_Chi2 = (FOi - FEi)^2 / FEi)

chi2_calculado <- sum(tabla_FO$Aporte_Chi2)
k_clases     <- nrow(tabla_FO)
m_parametros <- 0
gl           <- k_clases - 1 - m_parametros
alpha        <- 0.05
chi2_critico <- qchisq(1 - alpha, df = gl)

decision <- if (chi2_calculado <= chi2_critico) {
  "No se rechaza H0 (la distribución reciente es compatible con la histórica)"
} else {
  "Se rechaza H0 (la distribución reciente difiere significativamente de la histórica)"
}

cat("Chi-cuadrado calculado (X² calc):", round(chi2_calculado, 4), "\n")
## Chi-cuadrado calculado (X² calc): 0.6604
cat("Número de clases (k):", k_clases, "| Parámetros estimados (m):", m_parametros, "\n")
## Número de clases (k): 2 | Parámetros estimados (m): 0
cat("Grados de libertad (gl = k - 1 - m):", gl, "\n")
## Grados de libertad (gl = k - 1 - m): 1
cat("Chi-cuadrado crítico (alpha = 0.05):", round(chi2_critico, 4), "\n")
## Chi-cuadrado crítico (alpha = 0.05): 3.8415
cat("Decisión:", decision, "\n")
## Decisión: No se rechaza H0 (la distribución reciente es compatible con la histórica)

7.3 Tabla Resumen de Validación del Modelo

fila_max_aporte <- which.max(tabla_FO$Aporte_Chi2)

tabla_FO %>%
  gt() %>%
  tab_header(
    title = md("**VALIDACIÓN: ¿CAMBIÓ LA DISTRIBUCIÓN DEL DESVÍO DEL POZO?**"),
    subtitle = md(paste0("H0: proporciones históricas (< ", anio_corte, ") · n reciente = ", n_rec, " pozos"))
  ) %>%
  fmt_number(columns = c(fi, Pi, FEi, Aporte_Chi2), decimals = 4) %>%
  cols_label(
    Desvio_Etiqueta = "Tipo de desvío (x)", FOi = "FOi", fi = "fi", Pi = "Pi (H0)",
    FEi = "FEi", Aporte_Chi2 = "Aporte a X²"
  ) %>%
  cols_align(align = "center", columns = everything()) %>%
  tab_style(
    style = list(cell_fill(color = col_principal), cell_text(color = "white", weight = "bold")),
    locations = cells_title()
  ) %>%
  tab_style(
    style = list(cell_fill(color = "#148F77"), cell_text(color = "white", weight = "bold")),
    locations = cells_column_labels()
  ) %>%
  tab_style(
    style = list(cell_fill(color = "#FDEBD0"), cell_text(weight = "bold")),
    locations = cells_body(rows = fila_max_aporte)
  ) %>%
  opt_row_striping() %>%
  opt_table_font(font = google_font("Roboto")) %>%
  tab_options(
    table.font.size = px(13),
    heading.align = "left",
    data_row.padding = px(8),
    table.border.top.color = col_principal,
    table.border.bottom.color = col_principal,
    column_labels.border.bottom.color = col_principal
  ) %>%
  tab_source_note(md(paste0(
    "**X² calculado = ", round(chi2_calculado, 4),
    "** y **X² crítico (α = 0.05, gl = ", gl, ") = ", round(chi2_critico, 4), "** &nbsp;|&nbsp; ",
    "**Decisión:** ", decision
  )))
VALIDACIÓN: ¿CAMBIÓ LA DISTRIBUCIÓN DEL DESVÍO DEL POZO?
H0: proporciones históricas (< 2014) · n reciente = 96 pozos
Tipo de desvío (x) FOi fi Pi (H0) FEi Aporte a X²
Vertical 88 0.9167 0.9368 89.9368 0.0417
No vertical 8 0.0833 0.0632 6.0632 0.6187
X² calculado = 0.6604 y X² crítico (α = 0.05, gl = 1) = 3.8415  |  Decisión: No se rechaza H0 (la distribución reciente es compatible con la histórica)

8 Cálculo de Probabilidades

\[P(X = x) = p^{x} (1-p)^{1-x}, \quad x = 0 \text{ (Vertical)}, \; x = 1 \text{ (No vertical)}\]

p_vertical_hist    <- unname(p_hist_vec["Vertical"])
p_novertical_hist  <- unname(p_hist_vec["No vertical"])

fi_vec <- setNames(tabla_FO$fi, as.character(tabla_FO$Desvio_Etiqueta))
p_vertical_rec     <- unname(fi_vec["Vertical"])
p_novertical_rec   <- unname(fi_vec["No vertical"])

tabla_probabilidades <- data.frame(
  Periodo             = c("Histórico (H0)", "Histórico (H0)", "Reciente (observado)", "Reciente (observado)"),
  Tipo_de_desvio      = c("Vertical (x = 0)", "No vertical (x = 1)", "Vertical (x = 0)", "No vertical (x = 1)"),
  `P(X = x)`          = c(p_vertical_hist, p_novertical_hist, p_vertical_rec, p_novertical_rec),
  check.names = FALSE
)

tabla_probabilidades %>%
  gt() %>%
  tab_header(
    title = md("**CÁLCULO DE PROBABILIDADES — MODELO BERNOULLI**"),
    subtitle = md("P(X = x) = p^x (1-p)^(1-x), por periodo")
  ) %>%
  fmt_number(columns = `P(X = x)`, decimals = 5) %>%
  cols_label(Periodo = "Periodo", Tipo_de_desvio = "Resultado") %>%
  cols_align(align = "center", columns = everything()) %>%
  tab_style(
    style = list(cell_fill(color = col_principal), cell_text(color = "white", weight = "bold")),
    locations = cells_title()
  ) %>%
  tab_style(
    style = list(cell_fill(color = "#148F77"), cell_text(color = "white", weight = "bold")),
    locations = cells_column_labels()
  ) %>%
  opt_row_striping() %>%
  opt_table_font(font = google_font("Roboto")) %>%
  tab_options(table.font.size = px(13), heading.align = "left")
CÁLCULO DE PROBABILIDADES — MODELO BERNOULLI
P(X = x) = p^x (1-p)^(1-x), por periodo
Periodo Resultado P(X = x)
Histórico (H0) Vertical (x = 0) 0.93684
Histórico (H0) No vertical (x = 1) 0.06316
Reciente (observado) Vertical (x = 0) 0.91667
Reciente (observado) No vertical (x = 1) 0.08333
 
cat("P(Vertical)    histórico:", round(p_vertical_hist, 5),   "| reciente:", round(p_vertical_rec, 5), "\n")
## P(Vertical)    histórico: 0.93684 | reciente: 0.91667
cat("P(No vertical) histórico:", round(p_novertical_hist, 5), "| reciente:", round(p_novertical_rec, 5), "\n")
## P(No vertical) histórico: 0.06316 | reciente: 0.08333

9 Conclusión

  • Modelo: Bernoulli/Categórica de 2 clases (Vertical / No vertical).
  • Corte temporal: año 2014 (garantiza FEi ≥ 5).
  • Histórico (n = 95): Vertical = 93.68%, No vertical = 6.32%.
  • Reciente (n = 96): Vertical = 91.67%, No vertical = 8.33%.
  • X² calculado / crítico: 0.6604 / 3.8415 (gl = 1).
  • Decisión: No se rechaza H0 (la distribución reciente es compatible con la histórica).
  • Mayor aporte a X²: No vertical.
  • Interpretación: no se detecta un cambio estadísticamente significativo en la composición de tipos de desvío entre periodos.