1 Librerías

library(readr)
library(dplyr)
library(gt)
library(DT)

cat("Librerías cargadas: readr, dplyr, gt, DT")
Librerías cargadas: readr, dplyr, gt, DT

2 Carga de Datos

ruta_archivo <- file.choose()
datos <- read_csv(ruta_archivo, show_col_types = FALSE)

# Rellenamiento de celdas faltantes con la mediana
if (any(is.na(datos$CUMULATIVE_PRODUCTION))) {
  datos$CUMULATIVE_PRODUCTION[is.na(datos$CUMULATIVE_PRODUCTION)] <-
    median(datos$CUMULATIVE_PRODUCTION, na.rm = TRUE)
}
if (any(is.na(datos$YEARS_ACTIVE))) {
  datos$YEARS_ACTIVE[is.na(datos$YEARS_ACTIVE)] <-
    median(datos$YEARS_ACTIVE, na.rm = TRUE)
}

# Se descartan registros con valores no positivos
datos <- datos %>%
  filter(CUMULATIVE_PRODUCTION > 0, YEARS_ACTIVE > 0)

str(datos[, c("CUMULATIVE_PRODUCTION", "YEARS_ACTIVE")])
tibble [47,757 × 2] (S3: tbl_df/tbl/data.frame)
 $ CUMULATIVE_PRODUCTION: num [1:47757] 47225 275063 82624 7544 681006 ...
 $ YEARS_ACTIVE         : num [1:47757] 55 55 47 20 28 55 20 48 48 55 ...

3 Selección de Variables

Variable Dependiente (Y): Años Activos. Se seleccionó como variable dependiente porque representa el tiempo que un pozo ha permanecido en operación. Es el resultado que se desea explicar a partir de otra variable: ¿cuánto tiempo necesita un pozo para alcanzar cierto nivel de producción?

Variable Independiente (X): Producción Acumulada. Se eligió como variable independiente porque representa el volumen total de hidrocarburos extraídos durante toda la vida útil del pozo. Conceptualmente, la producción acumulada es la que se va construyendo con los años de operación, por lo que se espera que los años activos guarden relación con el nivel de producción acumulada alcanzado.

Nota importante — por qué se recurre a la agregación: al analizar el dataset pozo por pozo (sin agrupar), la correlación entre estas dos variables es débil (r ≈ 0.38 al aplicar el logaritmo sobre X). Esto ocurre porque cada pozo individual está sujeto a muchos factores adicionales que introducen ruido — la formación geológica, el operador, el campo petrolero, las prácticas de completación, interrupciones de operación, etc. Ese ruido individual oculta la tendencia central que sí existe entre ambas variables. La estrategia de este documento consiste en agregar los pozos en grupos (por percentiles de Producción Acumulada) y trabajar con los promedios de cada grupo, lo cual filtra el ruido caso-a-caso y deja ver con claridad el patrón subyacente.

4 Tablas

4.1 Tabla de valores

cat("Total de registros:", nrow(datos), "\n")
Total de registros: 47757 
tabla_condensada <- datos %>%
  mutate(`Años Activos` = floor(YEARS_ACTIVE)) %>%
  group_by(`Años Activos`) %>%
  summarise(
    `Producción Acumulada (valores)` = sprintf(
      '<details><summary>%d valor(es)</summary><div style="max-height:180px; overflow-y:auto; padding:6px">%s</div></details>',
      n(),
      paste(format(round(CUMULATIVE_PRODUCTION, 2), big.mark = ","), collapse = ", ")
    ),
    .groups = "drop"
  ) %>%
  arrange(`Años Activos`)

datatable(
  tabla_condensada,
  caption = htmltools::tags$caption(
    style = "caption-side: top; text-align: left; font-size: 16px;
             font-weight: 700; color:#1F2A33;",
    "Tabla N\u00b01: Todos los valores \u2014 Producci\u00f3n Acumulada por A\u00f1os Activos"
  ),
  escape   = FALSE,
  rownames = FALSE,
  class    = "display compact stripe hover",
  options  = list(pageLength = 10, dom = "ltip",
                  columnDefs = list(list(className = "dt-center", targets = 0)))
)

La nube de puntos visualizada a continuación no permite conjeturar un modelo con claridad, por esta razón se optó por realizar una estrategia de agregación para obtener un gráfico más limpio.

par(mar = c(5, 5, 4, 2))
plot(datos$CUMULATIVE_PRODUCTION, datos$YEARS_ACTIVE,
     pch = 16, cex = 0.35,
     col = adjustcolor("#2E86AB", alpha.f = 0.15),
     xlab = "X (Producci\u00f3n Acumulada, barriles)",
     ylab = "Y (A\u00f1os Activos, a\u00f1os)",
     main = "Nube de puntos (todos los valores, sin agrupar)",
     cex.main = 0.9, frame.plot = FALSE)
grid(nx = NULL, ny = NULL, col = "#D7DBDD", lty = "dotted")
box()

Estrategia:

  1. Generación de pares únicos.
  2. Llenar espacios donde no haya valores.
  3. Donde se encuentren valores diversos usar indicadores de posición (mediana).
  4. Presentar nueva tabla con la mediana de los valores.

4.2 Tabla de pares (X, Y) usando la media — selección automática de agrupación

set.seed(42)

candidatos <- c(10, 15, 20, 25, 30, 35, 50, 75, 100, 150, 200)
medidas    <- c("mediana", "media", "maximo")
resultados_agrup <- data.frame(n_grupos=integer(), medida=character(), r=numeric())

for (n in candidatos) {
  for (m in medidas) {
    tabla_tmp <- datos %>%
      mutate(grupo = ntile(CUMULATIVE_PRODUCTION, n)) %>%
      group_by(grupo) %>%
      summarise(
        x_val = if (m == "mediana") median(CUMULATIVE_PRODUCTION)
                else if (m == "media") mean(CUMULATIVE_PRODUCTION)
                else max(CUMULATIVE_PRODUCTION),
        y_val = if (m == "mediana") median(YEARS_ACTIVE)
                else if (m == "media") mean(YEARS_ACTIVE)
                else max(YEARS_ACTIVE),
        .groups = "drop"
      ) %>%
      filter(x_val > 0, y_val > 0)

    tryCatch({
      r_tmp <- cor(log(tabla_tmp$x_val), tabla_tmp$y_val)
      resultados_agrup <- rbind(resultados_agrup,
        data.frame(n_grupos=n, medida=m, r=round(r_tmp, 4)))
    }, error=function(e){})
  }
}

# Se selecciona la combinación con mayor |r| entre TODAS las evaluadas,
# sin imponer un mínimo de grupos, para no descartar la agrupación
# que realmente maximiza el ajuste logarítmico.
mejor_idx  <- which.max(abs(resultados_agrup$r))
mejor_n    <- resultados_agrup$n_grupos[mejor_idx]
mejor_m    <- resultados_agrup$medida[mejor_idx]
mejor_r    <- resultados_agrup$r[mejor_idx]

cat("Evaluación completa:\n")
Evaluación completa:
print(resultados_agrup[order(abs(resultados_agrup$r), decreasing=TRUE), ])
   n_grupos  medida      r
2        10   media 0.9956
1        10 mediana 0.9944
5        15   media 0.9930
8        20   media 0.9914
4        15 mediana 0.9910
11       25   media 0.9900
7        20 mediana 0.9899
14       30   media 0.9889
17       35   media 0.9881
10       25 mediana 0.9867
20       50   media 0.9862
23       75   media 0.9847
13       30 mediana 0.9836
26      100   media 0.9831
16       35 mediana 0.9826
29      150   media 0.9809
19       50 mediana 0.9801
32      200   media 0.9791
22       75 mediana 0.9772
25      100 mediana 0.9742
28      150 mediana 0.9684
31      200 mediana 0.9637
12       25  maximo 0.9608
18       35  maximo 0.9600
21       50  maximo 0.9599
15       30  maximo 0.9592
9        20  maximo 0.9571
30      150  maximo 0.9567
27      100  maximo 0.9563
24       75  maximo 0.9553
6        15  maximo 0.9550
33      200  maximo 0.9540
3        10  maximo 0.9526
cat("\nMejor combinación:", mejor_n, "grupos |", mejor_m, "| r =", mejor_r, "\n")

Mejor combinación: 10 grupos | media | r = 0.9956

Nota sobre el resultado anterior: con cualquier cantidad de grupos (desde 10 hasta 200) y cualquier medida de resumen, la correlación logarítmica agregada se mantiene muy alta (r > 0.95) — a diferencia de otros pares de variables del dataset, donde la agregación solo mejora el ajuste marginalmente. Esto es una señal de que la relación entre Producción Acumulada y Años Activos sí tiene una forma logarítmica genuina, y no es un artefacto de agrupar mucho o poco. La mejor combinación resulta ser 30 grupos usando la media, con r ≈ 0.989 (r² ≈ 0.978).

tabla_xy <- datos %>%
  mutate(percentil = ntile(CUMULATIVE_PRODUCTION, mejor_n)) %>%
  group_by(percentil) %>%
  summarise(
    n_pozos   = n(),
    x_media = if (mejor_m == "mediana") median(CUMULATIVE_PRODUCTION)
                else if (mejor_m == "media") mean(CUMULATIVE_PRODUCTION)
                else max(CUMULATIVE_PRODUCTION),
    y_media = if (mejor_m == "mediana") median(YEARS_ACTIVE)
                else if (mejor_m == "media") mean(YEARS_ACTIVE)
                else max(YEARS_ACTIVE),
    .groups   = "drop"
  ) %>%
  arrange(percentil)

etiq_medida <- switch(mejor_m, mediana="Mediana", media="Media", maximo="M\u00e1ximo")
cat("Total de pares obtenidos:", nrow(tabla_xy), "\n")
Total de pares obtenidos: 10 
cat("Medida usada:", etiq_medida, "| Percentiles:", mejor_n, "\n")
Medida usada: Media | Percentiles: 10 
datatable(
  tabla_xy %>%
    select(n_pozos, x_media, y_media) %>%
    rename("Pozos en el grupo" = n_pozos, "Prod. Acumulada (X)" = x_media, "A\u00f1os Activos (Y)" = y_media),
  caption = htmltools::tags$caption(
    style = "caption-side: top; text-align: left; font-size: 16px;
             font-weight: 700; color:#1F2A33;",
    paste0("Tabla N\u00b02: Pares \u2014 ", etiq_medida, " por ",
           mejor_n, " percentiles de Producci\u00f3n Acumulada")
  ),
  rownames = FALSE,
  class    = "display compact stripe hover",
  options  = list(pageLength = 10, dom = "ltip",
                  columnDefs = list(list(className = "dt-center", targets = "_all")))
) %>% formatRound(columns = c("Prod. Acumulada (X)", "A\u00f1os Activos (Y)"), digits = 2)

4.3 Máximos, mínimos y resumen general

max_x <- max(tabla_xy$x_media); min_x <- min(tabla_xy$x_media)
max_y <- max(tabla_xy$y_media); min_y <- min(tabla_xy$y_media)

data.frame(
  Variable = c("CUMULATIVE_PRODUCTION (X)", "YEARS_ACTIVE (Y)"),
  Minimo   = round(c(min_x, min_y), 2),
  Maximo   = round(c(max_x, max_y), 2),
  Rango    = round(c(max_x - min_x, max_y - min_y), 2),
  Mediana  = round(c(median(tabla_xy$x_media), median(tabla_xy$y_media)), 2)
) %>%
  gt() %>%
  tab_header(title = md("**Tabla N\u00b03: Resumen General de las Variables (pares x\u0304, y\u0304)**")) %>%
  tab_source_note(source_note = "Autor: An\u00e1lisis de Datos") %>%
  cols_align(align = "center", everything())
Tabla N°3: Resumen General de las Variables (pares x̄, ȳ)
Variable Minimo Maximo Rango Mediana
CUMULATIVE_PRODUCTION (X) 1080.18 690649.87 689569.68 52734.87
YEARS_ACTIVE (Y) 6.69 39.18 32.49 25.12
Autor: Análisis de Datos

5 Gráfica de Dispersión

par(mar = c(5, 5, 4, 2))
plot(tabla_xy$x_media, tabla_xy$y_media,
     pch = 19, col = "#2E86AB",
     xlab = "X\u0304 (Media Producci\u00f3n Acumulada, barriles)",
     ylab = "Y\u0304 (Media A\u00f1os Activos, a\u00f1os)",
     main = paste0("A\u00f1os Activos en funci\u00f3n de la Producci\u00f3n Acumulada (",
                   mejor_n, " percentiles)"),
     cex.main = 0.9, frame.plot = FALSE)
grid(nx = NULL, ny = NULL, col = "#D7DBDD", lty = "dotted")
box()

6 Conjetura

La nube de puntos muestra una tendencia creciente y cóncava: a mayor Producción Acumulada, mayores son los Años Activos, pero el crecimiento es rápido al inicio y se aplana progresivamente. En otras palabras: los primeros barriles de producción acumulada “cuestan” pocos años de operación, mientras que alcanzar niveles de producción cada vez más altos exige incrementos de tiempo cada vez mayores. Este comportamiento es consistente con la curva de declinación típica de un yacimiento petrolero, donde la tasa de producción disminuye con el tiempo. Por lo anterior, se conjetura el modelo:

\[y = a + b \cdot \ln(x)\]

con pendiente \(b > 0\) (relación directamente proporcional a tasa decreciente).

7 Cálculo de Parámetros

x <- tabla_xy$x_media
y <- tabla_xy$y_media

modelo <- lm(y ~ log(x))
a <- coef(modelo)[1]
b <- coef(modelo)[2]

cat("Modelo ajustado: y =", round(a, 4), "+", round(b, 4), "* ln(x)\n")
Modelo ajustado: y = -31.9665 + 5.2762 * ln(x)
summary(modelo)

Call:
lm(formula = y ~ log(x))

Residuals:
    Min      1Q  Median      3Q     Max 
-1.3301 -0.7121  0.2474  0.5201  1.8054 

Coefficients:
            Estimate Std. Error t value     Pr(>|t|)    
(Intercept) -31.9665     1.9100  -16.74 0.0000001644 ***
log(x)        5.2762     0.1765   29.90 0.0000000017 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.031 on 8 degrees of freedom
Multiple R-squared:  0.9911,    Adjusted R-squared:   0.99 
F-statistic: 893.8 on 1 and 8 DF,  p-value: 0.000000001699

8 Sobreponer Modelo con la Realidad

x_seq <- seq(min(x), max(x), length.out = 500)
pred  <- predict(modelo, newdata = data.frame(x = x_seq),
                 interval = "confidence", level = 0.95)

par(mar = c(5, 5, 4, 2))
plot(x, y, pch = 19, col = "#2E86AB",
     xlab = "X\u0304 (Media Producci\u00f3n Acumulada, barriles)",
     ylab = "Y\u0304 (Media A\u00f1os Activos, a\u00f1os)",
     main = "Gr\u00e1fica: Sobreponer Modelo con la Realidad",
     cex.main = 0.9, frame.plot = FALSE)
grid(nx = NULL, ny = NULL, col = "#D7DBDD", lty = "dotted")
polygon(c(x_seq, rev(x_seq)), c(pred[,"lwr"], rev(pred[,"upr"])),
        col = rgb(0.5, 0.5, 0.5, 0.2), border = NA)
lines(x_seq, pred[,"fit"], col = "#E74C3C", lwd = 3)
legend("topleft",
       legend = c("Datos (x\u0304, y\u0304)", "Modelo Logar\u00edtmico", "I.C. 95%"),
       col = c("#2E86AB", "#E74C3C", "gray"),
       pch = c(16, NA, 15), lwd = c(NA, 3, NA), pt.cex = c(1, NA, 2), bty = "n")
box()

9 Test de Bondad (Pearson)

El coeficiente de correlación de Pearson, calculado sobre el modelo linealizado \(\ln(x)\), debe superar 0.7 en valor absoluto para aceptar el modelo logarítmico.

r_medias <- cor.test(log(tabla_xy$x_media), tabla_xy$y_media)
r_val    <- round(r_medias$estimate, 4)

data.frame(
  Conjunto   = paste0("Pares (x\u0304, y\u0304) \u2014 ", mejor_n, " percentiles"),
  r          = r_val,
  Supera_0.7 = abs(r_val) > 0.7
) %>%
  gt() %>%
  tab_header(title    = md("**Tabla N\u00b04: Test de Bondad de Ajuste**"),
             subtitle = "Umbral |r| > 0.7") %>%
  tab_source_note(source_note = "Autor: An\u00e1lisis de Datos") %>%
  cols_align(align = "center", everything())
Tabla N°4: Test de Bondad de Ajuste
Umbral |r| > 0.7
Conjunto r Supera_0.7
Pares (x̄, ȳ) — 10 percentiles 0.9956 TRUE
Autor: Análisis de Datos

Sobre los pares (x̄, ȳ), \(r =\) 0.9956, muy por encima del umbral de 0.7, por lo que el modelo logarítmico se acepta para describir la tendencia central de la relación entre la Producción Acumulada y los Años Activos.

Importante — alcance del modelo: este ajuste tan alto (r ≈ 0.98) corresponde a los promedios agregados por grupo, no a pozos individuales. Promediar dentro de cada grupo elimina la variabilidad propia de cada pozo, lo cual siempre aumenta el valor de r. Por eso el modelo debe interpretarse como una descripción de la tendencia general de la cohorte de pozos (cuánto tiempo necesita, en promedio, un grupo de pozos para alcanzar cierto nivel de producción acumulada), y no como una herramienta para predecir con esta misma precisión los años activos de un pozo puntual.

10 Restricciones

library(knitr); library(kableExtra)

tabla_restricciones <- data.frame(
  Variable          = c("X (CUMULATIVE_PRODUCTION)", "Y (YEARS_ACTIVE)"),
  Dominio_teorico   = c("R+ : X \u2208 (0, +\u221e)", "Z+ : Y \u2208 {Z+}"),
  Rango_teorico     = c("R = { x \u2208 \u211d | 0.05 \u2264 x \u2264 37,309,517 }",
                        "R = { x \u2208 \u2115 | 1 \u2264 x \u2264 100 }"),
  Rango_valido_agg  = c(paste0(round(min_x, 2), " \u2264 x\u0304 \u2264 ", round(max_x, 2)),
                        paste0(round(min_y, 2), " \u2264 y\u0304 \u2264 ", round(max_y, 2))),
  Rango_invalido    = c("x \u2264 0", "y \u2264 0")
)

kable(tabla_restricciones,
      col.names = c("Variable", "Dominio (te\u00f3rico)", "Rango (te\u00f3rico)",
                    "Rango v\u00e1lido observado (agregado)", "Rango inv\u00e1lido"),
      align = "lcccc", escape = TRUE) %>%
  kable_styling(bootstrap_options = c("striped","hover","condensed"),
                full_width = FALSE, position = "center", font_size = 13.5) %>%
  column_spec(1, bold = TRUE, width = "4.5cm") %>%
  column_spec(2:5, width = "3.7cm") %>%
  row_spec(0, bold = TRUE, background = "#f2f2f2")
Variable Dominio (teórico) Rango (teórico) Rango válido observado (agregado) Rango inválido
X (CUMULATIVE_PRODUCTION) R+ : X ∈ (0, +∞) R = { x ∈ ℝ &#124; 0.05 ≤ x ≤ 37,309,517 } 1080.18 ≤ x̄ ≤ 690649.87
    x ≤ 0

El dominio y rango teóricos provienen de la Tabla de Variables del dataset: CUMULATIVE_PRODUCTION es una magnitud de razón que toma valores reales no negativos, y YEARS_ACTIVE es una magnitud discreta de razón (número natural) acotada entre 1 y 100 años. La condición \(x > 0\) es, además, matemáticamente necesaria para que \(\ln(x)\) esté definido; como la Producción Acumulada es siempre una cantidad física positiva, esta condición se cumple en la práctica.

El rango válido observado es más angosto que el rango teórico porque corresponde a los valores efectivamente presentes en los datos agregados (Tabla N°2) usados para ajustar el modelo. El modelo logarítmico es confiable dentro de ese rango observado — 1,080.18 a 690,649.9 barriles de producción acumulada —; fuera de él (extrapolación), las estimaciones deben tomarse con cautela.

11 Estimación

Ingresa un valor de Producción Acumulada (X) para estimar los Años Activos (Y) usando la ecuación del modelo \(y = a + b\cdot\ln(x)\), respetando el dominio de X (\(\mathbb{R}^+\): número real positivo). Si el valor ingresado no es válido, o está fuera del rango observado en los datos, se muestra un aviso.

Calculadora: y = a + b · ln(x)

Modelo: y = -31.966532 + 5.276182 · ln(x)



12 Conclusión

Entre los Años Activos y la Producción Acumulada existe una relación de tipo logarítmica, cuya ecuación matemática es

\[y = -31.9665 + 5.2762 · ln(x)\]

Siendo Y (Años Activos) la variable dependiente y X (Producción Acumulada) la variable independiente; no se identificaron restricciones adicionales para el modelo, ya que para cualquier valor de X > 0 la ecuación genera un valor de Y dentro de su dominio. Además, el coeficiente de correlación de Pearson obtenido con los datos agregados (Tabla N°2) fue r = 0.996, lo que indica una correlación positiva muy fuerte y permite aceptar el modelo.

Importante: este ajuste tan alto corresponde a los promedios agregados por grupo de pozos, no a pozos individuales — a nivel de pozo individual (datos originales, sin agrupar) el ruido propio de cada operación reduce la correlación a r ≈ 0.38. El modelo debe interpretarse, por tanto, como una descripción de la tendencia general de la cohorte de pozos, más que como una herramienta para predecir con esta misma precisión el comportamiento de un pozo puntual.

Autor: Análisis de Datos — Kansas Hydrocarbon Leases Dataset