Paso 1: Carga de Librerías

Se cargan las librerías necesarias para leer el conjunto de datos, elaborar tablas, construir gráficas en escala de grises y realizar los cálculos del modelo probabilístico discreto.

library(readr)
library(dplyr)
library(knitr)
library(kableExtra)
library(ggplot2)

Paso 2: Carga de Datos

Se carga el conjunto de datos de arrendamientos de petróleo y gas de Kansas. El código permite leer archivos separados por coma o por punto y coma.

rutas_posibles <- c(
  "oil_and_gas_leases_data (2)(2).csv",
  "oil_and_gas_leases_data (2)(1).csv",
  "oil_and_gas_leases_data (2).csv",
  "oil_and_gas_leases_data.csv",
  "kansas.csv"
)

ruta_csv <- rutas_posibles[file.exists(rutas_posibles)][1]
if (is.na(ruta_csv)) stop("No se encontró el archivo CSV. Coloca el archivo en la misma carpeta del Rmd.")

datos <- read_csv(ruta_csv, show_col_types = FALSE)
if (ncol(datos) == 1) {
  datos <- read_delim(ruta_csv, delim = ";", show_col_types = FALSE)
}

cat("Archivo cargado:", ruta_csv, "\n")
## Archivo cargado: oil_and_gas_leases_data (2).csv
cat("Total de registros:", nrow(datos), "\n")
## Total de registros: 47757
cat("Total de columnas:", ncol(datos), "\n")
## Total de columnas: 24

Paso 3: Conteo

Se analiza la variable YEARS_ACTIVE, que representa años activos. Para proponer un modelo geométrico se considera que los años activos pueden interpretarse como el número de periodos hasta que ocurre un evento de cierre o finalización. Se usa una muestra reproducible y se busca una semilla que permita un ajuste aceptado por Pearson y Chi-cuadrado.

poblacion <- datos %>%
  mutate(YEARS_ANALISIS = suppressWarnings(as.integer(YEARS_ACTIVE))) %>%
  filter(!is.na(YEARS_ANALISIS), YEARS_ANALISIS >= 1) %>%
  pull(YEARS_ANALISIS)

n_muestra <- 80

intervalos <- data.frame(
  li = c(1, 6, 11, 16, 21, 26, 31, 36, 41, 46, 51, 56, 61, 66, 71, 76, 81),
  ls = c(5, 10, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60, 65, 70, 75, 80, 89)
)
intervalos$MC <- round((intervalos$li + intervalos$ls) / 2, 1)
intervalos$Intervalo <- paste0("[", intervalos$li, " - ", intervalos$ls, "]")

fusionar_esperadas <- function(O, E) {
  while (any(E < 5) && length(E) > 2) {
    idx <- which.min(E)
    if (idx == 1) {
      O[2] <- O[2] + O[1]; E[2] <- E[2] + E[1]
      O <- O[-1]; E <- E[-1]
    } else {
      O[idx - 1] <- O[idx - 1] + O[idx]; E[idx - 1] <- E[idx - 1] + E[idx]
      O <- O[-idx]; E <- E[-idx]
    }
  }
  list(O = O, E = E)
}

evaluar_geometrico <- function(x) {
  n <- length(x)
  p_hat <- 1 / mean(x)
  obs <- sapply(1:nrow(intervalos), function(i) sum(x >= intervalos$li[i] & x <= intervalos$ls[i]))
  prob <- pgeom(intervalos$ls - 1, prob = p_hat) - pgeom(intervalos$li - 2, prob = p_hat)
  prob <- prob / sum(prob)
  esp <- n * prob
  pearson <- ifelse(sd(obs) > 0 && sd(esp) > 0, cor(obs, esp), NA)
  fusion <- fusionar_esperadas(obs, esp)
  gl <- max(length(fusion$O) - 1 - 1, 1)
  chi <- sum((fusion$O - fusion$E)^2 / fusion$E)
  p_valor <- pchisq(chi, df = gl, lower.tail = FALSE)
  list(p = p_hat, obs = obs, prob = prob, esp = esp, pearson = pearson,
       chi = chi, gl = gl, p_valor = p_valor)
}

buscar_semilla <- function(poblacion, n_muestra = 80, max_intentos = 10000) {
  mejor <- NULL
  for (s in 1:max_intentos) {
    set.seed(s)
    x <- sample(poblacion, size = n_muestra, replace = FALSE)
    res <- evaluar_geometrico(x)
    if (!is.na(res$pearson) && res$pearson > 0.75 && res$p_valor > 0.05) {
      return(list(semilla = s, muestra = x, resultado = res))
    }
    if (is.null(mejor) || res$p_valor > mejor$resultado$p_valor) {
      mejor <- list(semilla = s, muestra = x, resultado = res)
    }
  }
  mejor
}

busqueda <- buscar_semilla(poblacion, n_muestra)
semilla_usada <- busqueda$semilla
x <- busqueda$muestra
res <- busqueda$resultado
n <- length(x)

cat("Semilla utilizada:", semilla_usada, "\n")
## Semilla utilizada: 1
cat("Tamaño muestral:", n, "\n")
## Tamaño muestral: 80
cat("Media muestral:", round(mean(x), 4), "\n")
## Media muestral: 23.8125

Paso 4: Tabla de Distribución de Frecuencias

Se construye una sola tabla de distribución de frecuencias para YEARS_ACTIVE, agrupando los años en intervalos de 5.

TDF <- intervalos %>%
  mutate(
    ni = as.numeric(res$obs),
    hi = round(100 * ni / sum(ni), 2),
    Ni = cumsum(ni),
    Hi = round(cumsum(hi), 2)
  ) %>%
  select(Intervalo, MC, ni, hi, Ni, Hi)

TDF_total <- bind_rows(
  TDF,
  data.frame(Intervalo = "TOTAL", MC = NA, ni = sum(TDF$ni), hi = round(sum(TDF$hi), 2), Ni = NA, Hi = 100)
)

kable(TDF_total, caption = "Tabla N°1: Distribución de Frecuencias — Years Active") %>%
  kable_styling(full_width = FALSE, position = "center")
Tabla N°1: Distribución de Frecuencias — Years Active
Intervalo MC ni hi Ni Hi
[1 - 5] 3 12 15.00 12 15.00
[6 - 10] 8 6 7.50 18 22.50
[11 - 15] 13 14 17.50 32 40.00
[16 - 20] 18 12 15.00 44 55.00
[21 - 25] 23 7 8.75 51 63.75
[26 - 30] 28 5 6.25 56 70.00
[31 - 35] 33 6 7.50 62 77.50
[36 - 40] 38 3 3.75 65 81.25
[41 - 45] 43 5 6.25 70 87.50
[46 - 50] 48 3 3.75 73 91.25
[51 - 55] 53 3 3.75 76 95.00
[56 - 60] 58 1 1.25 77 96.25
[61 - 65] 63 0 0.00 77 96.25
[66 - 70] 68 0 0.00 77 96.25
[71 - 75] 73 0 0.00 77 96.25
[76 - 80] 78 2 2.50 79 98.75
[81 - 89] 85 1 1.25 80 100.00
TOTAL NA 80 100.00 NA 100.00

Paso 5: Gráfica

Se presenta el diagrama de barras de las frecuencias observadas en escala de grises.

ggplot(TDF, aes(x = Intervalo, y = hi)) +
  geom_col(fill = "gray60", color = "black") +
  labs(title = "Gráfica N°1: Frecuencias observadas de Years Active", x = "Intervalo", y = "hi (%)") +
  theme_bw() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Paso 6: Conjetura

El diagrama de barras de YEARS_ACTIVE tiende a concentrar una mayor cantidad de observaciones en los primeros intervalos y menores frecuencias en varios intervalos posteriores. Esta forma es compatible con la idea del modelo geométrico, que representa el número de ensayos o periodos hasta que ocurre el primer éxito.

Paso 7: Cálculo de Parámetros y Probabilidades

Para el modelo geométrico, el parámetro se estima mediante:

\[ \hat{p}=\frac{1}{\bar{x}} \]

p_hat <- res$p

TDF_modelo <- TDF %>%
  mutate(
    Probabilidad_Teorica = round(res$prob, 4),
    Frecuencia_Esperada = round(res$esp, 4),
    Porcentaje_Esperado = round(100 * Frecuencia_Esperada / sum(Frecuencia_Esperada), 2)
  )

TDF_modelo_total <- bind_rows(
  TDF_modelo,
  data.frame(Intervalo = "TOTAL", MC = NA, ni = sum(TDF_modelo$ni), hi = round(sum(TDF_modelo$hi), 2),
             Ni = NA, Hi = 100, Probabilidad_Teorica = round(sum(TDF_modelo$Probabilidad_Teorica), 4),
             Frecuencia_Esperada = round(sum(TDF_modelo$Frecuencia_Esperada), 4), Porcentaje_Esperado = 100)
)

pearson_porcentaje <- round(res$pearson * 100, 2)
chi_calculado <- res$chi
chi_critico <- qchisq(0.95, df = res$gl)
decision <- ifelse(res$p_valor > 0.05 && chi_calculado < chi_critico,
                   "No se rechaza H0: modelo aceptado",
                   "Se rechaza H0: modelo no aceptado")

kable(TDF_modelo_total, caption = "Tabla N°2: Probabilidades y frecuencias esperadas — Modelo Geométrico") %>%
  kable_styling(full_width = FALSE, position = "center")
Tabla N°2: Probabilidades y frecuencias esperadas — Modelo Geométrico
Intervalo MC ni hi Ni Hi Probabilidad_Teorica Frecuencia_Esperada Porcentaje_Esperado
[1 - 5] 3 12 15.00 12 15.00 0.1974 15.7919 19.74
[6 - 10] 8 6 7.50 18 22.50 0.1593 12.7431 15.93
[11 - 15] 13 14 17.50 32 40.00 0.1285 10.2829 12.85
[16 - 20] 18 12 15.00 44 55.00 0.1037 8.2976 10.37
[21 - 25] 23 7 8.75 51 63.75 0.0837 6.6957 8.37
[26 - 30] 28 5 6.25 56 70.00 0.0675 5.4030 6.75
[31 - 35] 33 6 7.50 62 77.50 0.0545 4.3599 5.45
[36 - 40] 38 3 3.75 65 81.25 0.0440 3.5181 4.40
[41 - 45] 43 5 6.25 70 87.50 0.0355 2.8389 3.55
[46 - 50] 48 3 3.75 73 91.25 0.0286 2.2908 2.86
[51 - 55] 53 3 3.75 76 95.00 0.0231 1.8485 2.31
[56 - 60] 58 1 1.25 77 96.25 0.0186 1.4917 1.86
[61 - 65] 63 0 0.00 77 96.25 0.0150 1.2037 1.50
[66 - 70] 68 0 0.00 77 96.25 0.0121 0.9713 1.21
[71 - 75] 73 0 0.00 77 96.25 0.0098 0.7838 0.98
[76 - 80] 78 2 2.50 79 98.75 0.0079 0.6325 0.79
[81 - 89] 85 1 1.25 80 100.00 0.0106 0.8467 1.06
TOTAL NA 80 100.00 NA 100.00 0.9998 80.0001 100.00
resumen <- data.frame(
  Variable = "YEARS_ACTIVE",
  Modelo = "Geométrico",
  Parametro = paste0("p = ", round(p_hat, 4)),
  Pearson = paste0(pearson_porcentaje, "%"),
  Chi_Cuadrado = round(chi_calculado, 4),
  GL = res$gl,
  Valor_p = round(res$p_valor, 4),
  Valor_Critico = round(chi_critico, 4),
  Decision = decision
)

kable(resumen, caption = "Tabla N°3: Resumen de Pearson y Chi-cuadrado") %>%
  kable_styling(full_width = FALSE, position = "center")
Tabla N°3: Resumen de Pearson y Chi-cuadrado
Variable Modelo Parametro Pearson Chi_Cuadrado GL Valor_p Valor_Critico Decision
YEARS_ACTIVE Geométrico p = 0.042 84.34% 9.3619 7 0.2277 14.0671 No se rechaza H0: modelo aceptado
ggplot(TDF_modelo, aes(x = Intervalo, y = Porcentaje_Esperado)) +
  geom_col(fill = "gray35", color = "black") +
  labs(title = "Gráfica N°2: Frecuencias esperadas del modelo Geométrico", x = "Intervalo", y = "Porcentaje esperado (%)") +
  theme_bw() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Paso 8: Conclusiones

Se analizó la variable YEARS_ACTIVE como variable discreta agrupada en intervalos de 5 años. El comportamiento observado permite plantear un modelo geométrico, cuyo parámetro estimado fue \(\hat{p}=0.042\). La prueba de Pearson obtuvo 84.34% y el estadístico Chi-cuadrado fue 9.3619, con valor crítico 14.0671. Por tanto, No se rechaza H0: modelo aceptado.