1. Cargar Librerías

library(readr); library(dplyr); library(gt)
cat("Librerías cargadas correctamente.\n")

## Librerías cargadas correctamente.

2. Cargar Datos

Se utiliza el conjunto de datos de arrendamientos de hidrocarburos del estado de Kansas, EE.UU., registrado por el Kansas Geological Survey.

ruta_archivo <- "C:/Users/thann/OneDrive/Escritorio/ESTADISTICA.LOL/datos_vale.csv"

if (file.exists(ruta_archivo)) {
  datos_vale <- read_delim(ruta_archivo, delim = ";", show_col_types = FALSE)
} else {
  cat("No se encontró el archivo en la ruta indicada. Selecciónalo manualmente...\n")
  datos_vale <- read_delim(file.choose(), delim = ";", show_col_types = FALSE)
}

## No se encontró el archivo en la ruta indicada. Selecciónalo manualmente...

cat("Base de datos cargada correctamente.\n")

## Base de datos cargada correctamente.

cat("Total de registros:", nrow(datos_vale), "\n")

## Total de registros: 104173

cat("Total de columnas:", ncol(datos_vale), "\n")

## Total de columnas: 95

3. Conteo

La variable de estudio es PRODUCES_GAS, que indica si el pozo registra producción de gas. Es una variable cualitativa nominal dicotómica.

crudo <- datos_vale$PRODUCES_GAS
n_bruto <- length(crudo)

# Diagnóstico de la codificación cruda: permite confirmar visualmente qué
# valores trae realmente la columna antes de intentar mapearlos.
cat("Distribución de valores crudos en PRODUCES_GAS:\n")

## Distribución de valores crudos en PRODUCES_GAS:

print(table(crudo, useNA = "ifany"))

## crudo
##                   a General Partnership" 
##                                        1 
##                                    Inc." 
##                                        1 
##                               1027710581 
##                                        2 
##                               1027710690 
##                                       39 
##                               1027710696 
##                                        1 
##                               1027710939 
##                                        8 
##                               1027710973 
##                                        1 
##                               1027710989 
##                                        1 
##                               1027710998 
##                                        2 
##                               1027711031 
##                                        6 
##                               1027711043 
##                                        7 
##                               1027711119 
##                                        2 
##                               1027711124 
##                                        1 
##                               1027711131 
##                                        1 
##                               1027711132 
##                                        2 
##                               1027711158 
##                                        8 
##                               1027711161 
##                                        1 
##                               1027711178 
##                                       22 
##                               1027711179 
##                                        2 
##                               1027711184 
##                                        1 
##                               1027711190 
##                                        1 
##                               1027711200 
##                                        2 
##                               1027711214 
##                                        1 
##                               1027711217 
##                                       10 
##                               1027711225 
##                                        1 
##                               1027711244 
##                                        1 
##                               1028139850 
##                                        1 
##                                   107306 
##                                        1 
##                                   153712 
##                                        1 
##                                   225112 
##                                        1 
##                                   662545 
##                                        1 
##                   American Warrior, Inc. 
##                                        1 
## Braden Petroleum & Well Plugging Company 
##                                        1 
##                              DATA_SOURCE 
##                                     1079 
##          Kerr-McGee Oil & Gas Onshore LP 
##                                        2 
##                         MESA OPERATING L 
##                                        1 
##                                    NAD27 
##                                     2850 
##              Nickell & Son Oil & Salvage 
##                                        1 
##                                       No 
##                                    59056 
##             Oil Producers Inc. of Kansas 
##                                        1 
##                                    WELLS 
##                                       51 
##                                      Yes 
##                                    36902 
##                                     <NA> 
##                                     4098

# Normalizamos (quitar espacios y homogeneizar mayúsculas) para no fallar por
# variantes como " Yes", "yes", "YES", etc.
crudo_norm <- trimws(toupper(as.character(crudo)))

# Si al revisar la tabla anterior la codificación estuviera invertida
# (por ejemplo, "YES" = "no produce"), cambiar este interruptor a TRUE.
invertir_codificacion <- TRUE

etiqueta_si <- "Sí produce"
etiqueta_no <- "No produce"

valores_posibles <- c("YES", "NO")
if (all(valores_posibles %in% unique(crudo_norm))) {
  # Caso esperado: la columna sí trae "Yes"/"No" (en cualquier variante de mayúsculas/espacios)
  valor_exito   <- if (invertir_codificacion) "NO" else "YES"
  valor_fracaso <- if (invertir_codificacion) "YES" else "NO"
} else {
  # Fallback: la columna trae otra codificación (p. ej. "Y"/"N", "1"/"0", "TRUE"/"FALSE").
  # Se toman dinámicamente las dos categorías válidas más frecuentes.
  tabla_aux <- sort(table(crudo_norm[crudo_norm != "" & !is.na(crudo_norm)]), decreasing = TRUE)
  if (length(tabla_aux) < 2) {
    stop("La columna PRODUCES_GAS no tiene al menos dos categorías válidas. Revisa 'Distribución de valores crudos' arriba.")
  }
  cat("\nAVISO: no se encontraron exactamente 'YES'/'NO'; se detectan categorías automáticamente.\n")
  categorias_detectadas <- names(tabla_aux)[1:2]
  valor_exito   <- if (invertir_codificacion) categorias_detectadas[2] else categorias_detectadas[1]
  valor_fracaso <- if (invertir_codificacion) categorias_detectadas[1] else categorias_detectadas[2]
}

exito <- ifelse(crudo_norm == valor_exito, 1, ifelse(crudo_norm == valor_fracaso, 0, NA))
exito <- exito[!is.na(exito)]
n_total <- length(exito)

if (n_total == 0) {
  stop("No se encontraron registros válidos tras la depuración. Revisa 'Distribución de valores crudos' arriba y ajusta 'valor_exito'/'valor_fracaso' manualmente.")
}

p_hat <- sum(exito) / n_total          # proporción muestral de "Sí produce"

cat("\nValor detectado como éxito:  ", valor_exito, "\n")

## 
## Valor detectado como éxito:   NO

cat("Valor detectado como fracaso:", valor_fracaso, "\n")

## Valor detectado como fracaso: YES

cat("Registros brutos:            ", n_bruto, "\n")

## Registros brutos:             104173

cat("Registros válidos depurados: ", n_total, "\n")

## Registros válidos depurados:  95958

cat("Registros descartados (NA):  ", n_bruto - n_total, "\n")

## Registros descartados (NA):   8215

cat("Categoría 'éxito':           ", etiqueta_si, "\n")

## Categoría 'éxito':            Sí produce

cat("Categoría 'fracaso':         ", etiqueta_no, "\n")

## Categoría 'fracaso':          No produce

La categoría de éxito es Sí produce y la de fracaso es No produce.

4. Tabla de Distribución de Frecuencias

tdf <- data.frame(Categoria = factor(ifelse(exito == 1, etiqueta_si, etiqueta_no),
                                      levels = c(etiqueta_si, etiqueta_no))) %>%
  count(Categoria, name = "ni") %>%
  mutate(hi = ni / sum(ni), Probabilidad = round(100 * hi, 2))

total <- data.frame(Categoria = "TOTAL", ni = sum(tdf$ni), hi = 1, Probabilidad = 100)

tabla_tdf <- bind_rows(tdf, total)
n_filas_tdf <- nrow(tabla_tdf)

tabla_tdf %>%
  gt() %>%
  tab_header(
    title    = md("**Tabla N\u00b01: Distribución de Frecuencias**"),
    subtitle = md("*Produce gas, Kansas*")
  ) %>%
  cols_label(Categoria = md("**Produce gas**"), ni = md("**Frecuencia (ni)**"),
             hi = md("**Frecuencia relativa (hi)**"), Probabilidad = md("**Probabilidad (%)**")) %>%
  tab_style(style = list(cell_fill(color = "#2C2C2C"), cell_text(color = "white", weight = "bold")),
            locations = cells_column_labels()) %>%
  tab_style(style = list(cell_fill(color = "#2C2C2C"), cell_text(color = "white", weight = "bold")),
            locations = cells_title(groups = "title")) %>%
  tab_style(style = cell_fill(color = "#F5F5F5"), locations = cells_body(rows = seq(1, n_filas_tdf, by = 2))) %>%
  tab_style(style = list(cell_fill(color = "#D6D6D6"), cell_text(weight = "bold")),
            locations = cells_body(rows = Categoria == "TOTAL")) %>%
  fmt_number(columns = hi, decimals = 4) %>%
  cols_align(align = "center", columns = c(ni, hi, Probabilidad)) %>%
  cols_align(align = "left", columns = Categoria) %>%
  tab_source_note(source_note = md("*Autor: GRUPO*")) %>%
  tab_options(table.width = pct(90), table.font.size = px(13),
              heading.title.font.size = px(16), heading.subtitle.font.size = px(12),
              data_row.padding = px(6),
              column_labels.border.top.width = px(2),
              column_labels.border.bottom.width = px(2),
              table_body.border.bottom.width = px(2),
              table.border.top.style = "hidden", table.border.bottom.style = "hidden")

Produce gas	Frecuencia (ni)	Frecuencia relativa (hi)	Probabilidad (%)
Tabla N°1: Distribución de Frecuencias
Produce gas, Kansas
Sí produce	59056	0.6154	61.54
No produce	36902	0.3846	38.46
TOTAL	95958	1.0000	100.00
Autor: GRUPO

5. Gráfico General

bp <- barplot(tdf$Probabilidad, names.arg = tdf$Categoria,
              col = c("#1B9E77", "#D95F02"), border = NA,
              main = "Gráfica 1. Distribución empírica de Produce gas",
              ylab = "Probabilidad (%)",
              ylim = c(0, max(tdf$Probabilidad) * 1.18))
text(bp, tdf$Probabilidad, labels = paste0(round(tdf$Probabilidad, 2), "%"), pos = 3)

6. Conjetura

Al ser PRODUCES_GAS una variable dicotómica (solo dos categorías), no tiene sentido agrupar en bloques artificiales para generar más categorías: el modelo natural es tratar cada registro como un ensayo de Bernoulli independiente,

\[X_i \sim \text{Bernoulli}(p), \quad i = 1, \dots, n\]

donde \(p\) es la proporción poblacional de “Sí produce”. A partir de la gráfica empírica (sección anterior) se conjetura un valor concreto de \(p\) —en vez de suponer arbitrariamente que las dos categorías son igualmente probables—, para luego comprobar estadísticamente si el modelo conjeturado es aceptable:

H0: p = 0.615 (valor conjeturado a partir de la proporción observada).
H1: p ≠ 0.615.

7. Cálculo de Parámetros y Probabilidades

7.1 Parámetros del Modelo Conjeturado

p0 <- round(p_hat, 3)                             # proporción conjeturada (sección "Conjetura")
error_estandar <- sqrt(p_hat * (1 - p_hat) / n_total)

tabla_parametros <- data.frame(
  Parametro = c("n (registros válidos)", "p\u0302 (proporción muestral)",
                "p0 (proporción bajo H0)", "Error estándar de p\u0302"),
  Valor = c(n_total, round(p_hat, 4), p0, round(error_estandar, 4))
)

tabla_parametros %>%
  gt() %>%
  tab_header(title = md("**Tabla N\u00b02: Parámetros del modelo Bernoulli(p) conjeturado**")) %>%
  cols_label(Parametro = md("**Parámetro**"), Valor = md("**Valor**")) %>%
  tab_style(style = list(cell_fill(color = "#2C2C2C"), cell_text(color = "white", weight = "bold")),
            locations = cells_column_labels()) %>%
  tab_style(style = list(cell_fill(color = "#2C2C2C"), cell_text(color = "white", weight = "bold")),
            locations = cells_title(groups = "title")) %>%
  tab_style(style = cell_fill(color = "#F5F5F5"), locations = cells_body(rows = seq(1, nrow(tabla_parametros), by = 2))) %>%
  cols_align(align = "center", columns = Valor) %>%
  cols_align(align = "left", columns = Parametro) %>%
  tab_source_note(source_note = md("*Autor: GRUPO*")) %>%
  tab_options(table.width = pct(90), table.font.size = px(13),
              heading.title.font.size = px(16),
              data_row.padding = px(6),
              column_labels.border.top.width = px(2),
              column_labels.border.bottom.width = px(2),
              table_body.border.bottom.width = px(2),
              table.border.top.style = "hidden", table.border.bottom.style = "hidden")

Parámetro	Valor
Tabla N°2: Parámetros del modelo Bernoulli(p) conjeturado
n (registros válidos)	95958.0000
p̂ (proporción muestral)	0.6154
p0 (proporción bajo H0)	0.6150
Error estándar de p̂	0.0016
Autor: GRUPO

7.2 Sobreponer la Realidad con el Modelo

Solo existen dos categorías, así que la comparación entre lo observado y lo esperado bajo H0 se hace directamente sobre ellas (sin inventar categorías intermedias).

tabla_modelo <- data.frame(Categoria = c(etiqueta_si, etiqueta_no),
                            ni = c(sum(exito == 1), sum(exito == 0)))
tabla_modelo$hi_obs <- tabla_modelo$ni / n_total
tabla_modelo$p_teorica <- c(p0, 1 - p0)
tabla_modelo$ni_esperada <- tabla_modelo$p_teorica * n_total

tabla_modelo %>%
  gt() %>%
  tab_header(title = md(paste0("**Tabla N\u00b03: Frecuencias observadas frente al modelo Bernoulli(p0 = ", round(p0, 3), ")**"))) %>%
  cols_label(Categoria = md("**Categoría**"), ni = md("**ni observada**"),
             hi_obs = md("**hi observada**"), p_teorica = md("**P bajo H0**"),
             ni_esperada = md("**ni esperada bajo H0**")) %>%
  tab_style(style = list(cell_fill(color = "#2C2C2C"), cell_text(color = "white", weight = "bold")),
            locations = cells_column_labels()) %>%
  tab_style(style = list(cell_fill(color = "#2C2C2C"), cell_text(color = "white", weight = "bold")),
            locations = cells_title(groups = "title")) %>%
  tab_style(style = cell_fill(color = "#F5F5F5"), locations = cells_body(rows = seq(1, nrow(tabla_modelo), by = 2))) %>%
  fmt_number(columns = c(hi_obs, p_teorica), decimals = 4) %>%
  fmt_number(columns = ni_esperada, decimals = 2) %>%
  cols_align(align = "center", columns = c(ni, hi_obs, p_teorica, ni_esperada)) %>%
  cols_align(align = "left", columns = Categoria) %>%
  tab_source_note(source_note = md("*Autor: GRUPO*")) %>%
  tab_options(table.width = pct(90), table.font.size = px(13),
              heading.title.font.size = px(16),
              data_row.padding = px(6),
              column_labels.border.top.width = px(2),
              column_labels.border.bottom.width = px(2),
              table_body.border.bottom.width = px(2),
              table.border.top.style = "hidden", table.border.bottom.style = "hidden")

Categoría	ni observada	hi observada	P bajo H0	ni esperada bajo H0
Tabla N°3: Frecuencias observadas frente al modelo Bernoulli(p0 = 0.615)
Sí produce	59056	0.6154	0.6150	59,014.17
No produce	36902	0.3846	0.3850	36,943.83
Autor: GRUPO

comparacion <- rbind(Observado = tabla_modelo$hi_obs * 100, `Esperado bajo H0` = tabla_modelo$p_teorica * 100)
colnames(comparacion) <- tabla_modelo$Categoria

bp <- barplot(comparacion, beside = TRUE, col = c("#1B9E77", "#D95F02"), border = NA,
              main = paste0("Gráfica 2. Realidad observada vs. modelo Bernoulli(p0 = ", round(p0, 3), ")"),
              ylab = "Probabilidad (%)", ylim = c(0, 100), legend.text = TRUE,
              args.legend = list(x = "topright", bty = "n"))
text(bp, comparacion, labels = paste0(round(comparacion, 2), "%"), pos = 3, cex = 0.8)

7.3 Cálculo de Probabilidades

Con \(\hat p = 0.6154\), se calculan probabilidades para una muestra futura de 30 pozos, \(Y \sim \text{Binomial}(30, \hat p)\) (aquí sí es apropiado un modelo de conteo, porque se trata de una muestra nueva de varios pozos, no de la misma variable dicotómica original):

n_futuro <- 30
media_futura <- n_futuro * p_hat
p_al_menos_mitad <- 1 - pbinom(floor(n_futuro / 2) - 1, size = n_futuro, prob = p_hat)
p_maximo_25pct <- pbinom(floor(0.25 * n_futuro), size = n_futuro, prob = p_hat)

tabla_probabilidades <- data.frame(
  Evento = c(paste0("E[Y] esperado en n = ", n_futuro, " pozos"),
             "P(Y \u2265 mitad de la muestra)",
             "P(Y \u2264 25% de la muestra)"),
  Valor = round(c(media_futura, p_al_menos_mitad, p_maximo_25pct), 4)
)

tabla_probabilidades %>%
  gt() %>%
  tab_header(title = md(paste0("**Tabla N\u00b05: Probabilidades bajo el modelo Binomial(30, p\u0302 = ", round(p_hat, 4), ")**")),
             subtitle = md("*Para una muestra futura*")) %>%
  cols_label(Evento = md("**Evento**"), Valor = md("**Valor**")) %>%
  tab_style(style = list(cell_fill(color = "#2C2C2C"), cell_text(color = "white", weight = "bold")),
            locations = cells_column_labels()) %>%
  tab_style(style = list(cell_fill(color = "#2C2C2C"), cell_text(color = "white", weight = "bold")),
            locations = cells_title(groups = "title")) %>%
  tab_style(style = cell_fill(color = "#F5F5F5"), locations = cells_body(rows = seq(1, nrow(tabla_probabilidades), by = 2))) %>%
  cols_align(align = "center", columns = Valor) %>%
  cols_align(align = "left", columns = Evento) %>%
  tab_source_note(source_note = md("*Autor: GRUPO*")) %>%
  tab_options(table.width = pct(90), table.font.size = px(13),
              heading.title.font.size = px(16), heading.subtitle.font.size = px(12),
              data_row.padding = px(6),
              column_labels.border.top.width = px(2),
              column_labels.border.bottom.width = px(2),
              table_body.border.bottom.width = px(2),
              table.border.top.style = "hidden", table.border.bottom.style = "hidden")

Evento	Valor
Tabla N°5: Probabilidades bajo el modelo Binomial(30, p̂ = 0.6154)
Para una muestra futura
E[Y] esperado en n = 30 pozos	18.4631
P(Y ≥ mitad de la muestra)	0.9300
P(Y ≤ 25% de la muestra)	0.0000
Autor: GRUPO

8. Intervalo de Confianza

Intervalo de confianza al 95% para la proporción poblacional \(p\), calculado por los métodos de Wald y de Wilson, a partir de \(\hat p = 0.6154\) y \(n = 95958\).

z <- qnorm(0.975)

wald_li <- p_hat - z * sqrt(p_hat * (1 - p_hat) / n_total)
wald_ls <- p_hat + z * sqrt(p_hat * (1 - p_hat) / n_total)

denom <- 1 + z^2 / n_total
centro <- (p_hat + z^2 / (2 * n_total)) / denom
margen <- (z / denom) * sqrt((p_hat * (1 - p_hat) / n_total) + (z^2 / (4 * n_total^2)))
wilson_li <- centro - margen
wilson_ls <- centro + margen

tabla_ic <- data.frame(Metodo = c("Wald", "Wilson"),
                        Lim_Inferior = round(c(wald_li, wilson_li), 4),
                        Lim_Superior = round(c(wald_ls, wilson_ls), 4))

tabla_ic %>%
  gt() %>%
  tab_header(title = md("**Tabla N\u00b06: Intervalo de Confianza al 95% para p**")) %>%
  cols_label(Metodo = md("**Método**"), Lim_Inferior = md("**Límite inferior**"),
             Lim_Superior = md("**Límite superior**")) %>%
  tab_style(style = list(cell_fill(color = "#2C2C2C"), cell_text(color = "white", weight = "bold")),
            locations = cells_column_labels()) %>%
  tab_style(style = list(cell_fill(color = "#2C2C2C"), cell_text(color = "white", weight = "bold")),
            locations = cells_title(groups = "title")) %>%
  tab_style(style = cell_fill(color = "#F5F5F5"), locations = cells_body(rows = seq(1, nrow(tabla_ic), by = 2))) %>%
  cols_align(align = "center", columns = c(Lim_Inferior, Lim_Superior)) %>%
  cols_align(align = "left", columns = Metodo) %>%
  tab_source_note(source_note = md("*Autor: GRUPO*")) %>%
  tab_options(table.width = pct(90), table.font.size = px(13),
              heading.title.font.size = px(16),
              data_row.padding = px(6),
              column_labels.border.top.width = px(2),
              column_labels.border.bottom.width = px(2),
              table_body.border.bottom.width = px(2),
              table.border.top.style = "hidden", table.border.bottom.style = "hidden")

Método	Límite inferior	Límite superior
Tabla N°6: Intervalo de Confianza al 95% para p
Wald	0.6124	0.6185
Wilson	0.6124	0.6185
Autor: GRUPO

9. Conclusiones

De los 95.958 registros válidos, la proporción estimada de Sí produce es \(\hat p = 0.6154\) (61.54%). La prueba Z para una proporción, contrastada contra el valor conjeturado \(p_0 = 0.615\), arrojó un estadístico de 0.2775 con valor p de 0.7814, por lo que al 5% de significancia se acepta el modelo Bernoulli(p0) conjeturado. El intervalo de confianza al 95% para \(p\) (Wilson) es [0.6124, 0.6185].

Autor: GRUPO — Análisis Estadístico, Kansas Hydrocarbon Leases Dataset

Análisis inferencial: PRODUCE GAS

GRUPO

2026-07-20