1. Librerias

suppressMessages(suppressWarnings({
  library(dplyr)
  library(knitr)
  library(kableExtra)
  library(DT)
  library(plotly)
}))

2. Carga de datos

# usar Run All, no Knit
ruta <- file.choose()
datos <- read.csv(ruta, stringsAsFactors = FALSE)

datos <- datos %>%
  mutate(
    CUMULATIVE_PRODUCTION = as.numeric(CUMULATIVE_PRODUCTION),
    AVG_PRODUCTION        = as.numeric(AVG_PRODUCTION),
    YEARS_ACTIVE          = as.numeric(YEARS_ACTIVE)
  ) %>%
  select(CUMULATIVE_PRODUCTION, AVG_PRODUCTION, YEARS_ACTIVE) %>%
  na.omit() %>%
  filter(CUMULATIVE_PRODUCTION > 0, AVG_PRODUCTION > 0, YEARS_ACTIVE > 0)

3. Seleccion de variables

  • Variable Dependiente (Y): Producción Promedio. Se seleccionó como variable de respuesta porque representa la producción promedio estimada de un pozo, es decir, su desempeño operativo actual.

  • Variable Independiente 1 (X₁): Años Activo. Indica el tiempo que el pozo ha estado en operación. Esta variable refleja la madurez del pozo y su experiencia operativa . Es la variable dominante, ya que determina de manera significativa la capacidad de producción promedio actual.

  • Variable Independiente 2 (X₂): Producción Acumulada. Representa la producción acumulada del pozo desde su inicio. Esta variable cuantifica el efecto histórico del rendimiento del pozo sobre su producción promedio actual, aportando información complementaria sobre su desempeño operativo.

La Producción Promedio (Y) es, esencialmente, el resultado de la interacción entre la antigüedad del pozo (X₁) y su producción acumulada (X₂). Ambos factores interactúan para determinar la producción promedio estimada, proporcionando un panorama integral del desempeño del pozo.

4. Tabla de pares de valores

4.1 Tripleta de valores originales

condensado <- datos %>%
  group_by(YEARS_ACTIVE) %>%
  summarise(
    n_valores = n(),
    valores_avg = paste(AVG_PRODUCTION, collapse = ", "),
    valores_cum = paste(CUMULATIVE_PRODUCTION, collapse = ", "),
    .groups = "drop"
  ) %>%
  arrange(YEARS_ACTIVE)

condensado_html <- condensado %>%
  mutate(
    AVG_PRODUCTION = paste0("<details><summary>", n_valores, " valor(es)</summary>", valores_avg, "</details>"),
    CUMULATIVE_PRODUCTION = paste0("<details><summary>", n_valores, " valor(es)</summary>", valores_cum, "</details>")
  ) %>%
  select(YEARS_ACTIVE, AVG_PRODUCTION, CUMULATIVE_PRODUCTION)

nrow(datos)
## [1] 47757
datatable(condensado_html, escape = FALSE, rownames = FALSE,
          colnames = c("Años Activos (X2)", "Producción Promedio (valores, X1)", "Producción Acumulada (valores, Y)"),
          class = "stripe hover compact",
          options = list(pageLength = 10))
plot_ly(
  data = datos,
  x = ~AVG_PRODUCTION, y = ~YEARS_ACTIVE, z = ~CUMULATIVE_PRODUCTION,
  type = "scatter3d", mode = "markers",
  marker = list(size = 2, opacity = 0.4)
) %>%
  layout(scene = list(
    xaxis = list(title = "X1 AVG_PRODUCTION"),
    yaxis = list(title = "X2 YEARS_ACTIVE"),
    zaxis = list(title = "Y CUMULATIVE_PRODUCTION")
  ))
Estrategia
  1. Ordenar valores de X1, X2 e Y (un X para cada Y)
  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 tripletas (X1, X2, Y) usando la mediana

tabla_mediana <- datos %>%
  group_by(YEARS_ACTIVE) %>%
  summarise(
    AVG_PRODUCTION        = median(AVG_PRODUCTION),
    CUMULATIVE_PRODUCTION = median(CUMULATIVE_PRODUCTION),
    .groups = "drop"
  ) %>%
  arrange(YEARS_ACTIVE)

nrow(tabla_mediana)
## [1] 89
datatable(tabla_mediana, rownames = FALSE,
          colnames = c("Años Activos (X2)", "Mediana Producción Promedio (X1)", "Mediana Producción Acumulada (Y)"),
          class = "stripe hover compact",
          options = list(pageLength = 10)) %>%
  formatRound(columns = c("AVG_PRODUCTION", "CUMULATIVE_PRODUCTION"), digits = 2)
write.csv(condensado[, c("YEARS_ACTIVE", "valores_avg", "valores_cum")], "tabla_condensada_por_anios_oilgas.csv", row.names = FALSE)
write.csv(tabla_mediana, "tabla_mediana_por_anios_oilgas.csv", row.names = FALSE)

5. Grafica de dispersion

plot_ly(
  data = tabla_mediana,
  x = ~AVG_PRODUCTION, y = ~YEARS_ACTIVE, z = ~CUMULATIVE_PRODUCTION,
  type = "scatter3d", mode = "markers"
) %>%
  layout(scene = list(
    xaxis = list(title = "X1 AVG_PRODUCTION"),
    yaxis = list(title = "X2 YEARS_ACTIVE"),
    zaxis = list(title = "Y CUMULATIVE_PRODUCTION")
  ))

6. Conjetura

Se conjetura que la producción acumulada (Y) mantiene una relación lineal positiva con la producción promedio (X1) y con los años de actividad (X2). En el gráfico se observa que, a medida que aumenta la producción promedio y el tiempo de operación, la producción acumulada tiende a incrementarse. Aunque existen algunos valores atípicos, por lo que se ha deicidio usar un modelo de regresión lineal múltiple.

\[y = b_0 + b_1 x_1 + b_2 x_2\]

7. Calculo de parametros

modelo_lineal <- lm(CUMULATIVE_PRODUCTION ~ AVG_PRODUCTION + YEARS_ACTIVE, data = tabla_mediana)
invisible(summary(modelo_lineal))

b0 <- unname(coef(modelo_lineal)[1])
b1 <- unname(coef(modelo_lineal)[2])
b2 <- unname(coef(modelo_lineal)[3])

c(b0 = b0, b1 = b1, b2 = b2)
##            b0            b1            b2 
## -145857.37971      42.29425    3920.59028
kable(data.frame(Parámetro = c("b0 (intercepto)", "b1 (pendiente X1)", "b2 (pendiente X2)"),
                  Valor = round(c(b0, b1, b2), 4)),
      align = "lr") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = FALSE)
Parámetro Valor
b0 (intercepto) -145857.3797
b1 (pendiente X1) 42.2943
b2 (pendiente X2) 3920.5903

8. Realidad y modelo

x1_seq <- seq(min(tabla_mediana$AVG_PRODUCTION), max(tabla_mediana$AVG_PRODUCTION), length.out = 30)
x2_seq <- seq(min(tabla_mediana$YEARS_ACTIVE), max(tabla_mediana$YEARS_ACTIVE), length.out = 30)
z_plano <- outer(x2_seq, x1_seq, function(x2, x1) b0 + b1 * x1 + b2 * x2)

plot_ly() %>%
  add_markers(data = tabla_mediana, x = ~AVG_PRODUCTION, y = ~YEARS_ACTIVE, z = ~CUMULATIVE_PRODUCTION) %>%
  add_surface(x = x1_seq, y = x2_seq, z = z_plano, opacity = 0.6, showscale = FALSE) %>%
  layout(scene = list(
    xaxis = list(title = "X1 AVG_PRODUCTION"),
    yaxis = list(title = "X2 YEARS_ACTIVE"),
    zaxis = list(title = "Y CUMULATIVE_PRODUCTION")
  ))

Los puntos quedan cerca del plano, buen ajuste.

9. Test

r_x1 <- cor(tabla_mediana$AVG_PRODUCTION, tabla_mediana$CUMULATIVE_PRODUCTION)
r_x2 <- cor(tabla_mediana$YEARS_ACTIVE, tabla_mediana$CUMULATIVE_PRODUCTION)
r2 <- summary(modelo_lineal)$r.squared
p_modelo <- pf(summary(modelo_lineal)$fstatistic[1],
               summary(modelo_lineal)$fstatistic[2],
               summary(modelo_lineal)$fstatistic[3], lower.tail = FALSE)

resultado_x1 <- ifelse(abs(r_x1) > 0.7, "Aceptado", "Rechazado")
resultado_x2 <- ifelse(abs(r_x2) > 0.7, "Aceptado", "Rechazado")
resultado_modelo <- ifelse(p_modelo < 0.05, "Significativo", "No significativo")

tabla_test <- data.frame(
  Test = c("Pearson X1", "Pearson X2", "F-test modelo"),
  Valor = c(round(r_x1, 4), round(r_x2, 4), format(p_modelo, scientific = TRUE, digits = 4)),
  Resultado = c(resultado_x1, resultado_x2, resultado_modelo)
)

kable(tabla_test, align = "lrc") %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = FALSE)
Test Valor Resultado
Pearson X1 0.7735 Aceptado
Pearson X2 0.8576 Aceptado
F-test modelo 4.105e-48 Significativo 
cat("R² del modelo:", round(r2, 4), "\n")
## R² del modelo: 0.9209

Las dos variables pasan el test de Pearson (|r| > 0.7) y el modelo explica mas del 90% (R2 > 0.9).

10. Restricciones

El modelo presenta restricciones en valores para

  • X1 (AVG_PRODUCTION) desde -∞ hasta 3499

  • X2 (YEARS_ACTIVE) desde -∞ hasta 37.30

Advertencia : La inferencia de la muestra a la población en los reales positivos se respalda en el R2 > 90% dentro del rango observado, pero entre mas se aleje la predicción de ese rango, menor es la confiabilidad de la estimación.

11. Estimacion

nuevos_x1 <- c(3499, 5000, 10000, 25000, 50000, 100000)
nuevos_x2 <- c(38, 40, 45, 50, 60, 70)

estimaciones <- expand.grid(AVG_PRODUCTION = nuevos_x1, YEARS_ACTIVE = nuevos_x2)
estimaciones$Produccion_Acumulada_Estimada <- predict(modelo_lineal, newdata = estimaciones)

# solo estimaciones positivas
estimaciones <- estimaciones %>% filter(Produccion_Acumulada_Estimada > 0)
estimaciones <- estimaciones %>% arrange(AVG_PRODUCTION, YEARS_ACTIVE)

datatable(estimaciones, rownames = FALSE,
          colnames = c("Producción Promedio (X1)", "Años Activos (X2)", "Producción Acumulada Estimada"),
          class = "stripe hover compact",
          options = list(pageLength = 10)) %>%
  formatRound(columns = "Produccion_Acumulada_Estimada", digits = 2)

¿Cuál es la producción acumulada estimada para un pozo con una producción promedio de 5,000 y unos años activos de 40?

Explicación: predict(modelo, newdata = ...) evalúa la ecuación del plano ajustado en un punto específico (X1, X2) que no necesariamente está en los datos originales.

nuevo_dato <- data.frame(AVG_PRODUCTION = 5000, YEARS_ACTIVE = 40)
estimacion_puntual <- predict(modelo_lineal, newdata = nuevo_dato)

cat("## Para una Producción Promedio de 5,000 y unos Años Activos de 40, la Producción Acumulada estimada es:",
    format(round(estimacion_puntual, 2), big.mark = ","))
## ## Para una Producción Promedio de 5,000 y unos Años Activos de 40, la Producción Acumulada estimada es: 222,437.5

12. Conclusion

Entre la producción promedio (X₁), los años de actividad (X₂) y la producción acumulada (Y) existe una relación lineal múltiple cuya ecuación matemática es:

\[Y = -145857.4 + 42.29425\,X_1 + 3920.59\,X_2\]

Siendo X₁ la producción promedio, medida en las unidades correspondientes de producción; X₂ los años de actividad del pozo; y Y la producción acumulada. Con un coeficiente de correlación múltiple; El modelo explica 92.1% de la variabilidad de Y (R2), el resto (7.9%) es por otros factores, evidenciando una relación lineal positiva y de alta intensidad entre las variables.