##### UNIVERSIDAD CENTRAL DEL ECUADOR #####
#### AUTOR: ESTALIN ####
### CARRERA: INGENIERÍA EN PETRÓLEOS ###
#### MODELO DE REGRESIÓN EXPONENCIAL ####
## DATASET: Global Oil and Gas Extraction Tracker (GOGET) ##
library(readxl)
library(dplyr)
library(ggplot2)
library(gt)
# Selección del archivo Excel
ruta_archivo <- file.choose()
# guess_max se sube porque las primeras filas del archivo son solo metadatos
# (sin reservas ni unidades); con el valor por defecto (1000), readxl
# adivina mal el tipo de esas columnas y termina leyendo todo como NA.
Datos <- read_excel(ruta_archivo, guess_max = 50000)
## Estructura de los datos
str(Datos)## tibble [49,212 × 32] (S3: tbl_df/tbl/data.frame)
## $ Unit ID : chr [1:49212] "OG0000001" "OG0000002" "OG0000006" "OG0000007" ...
## $ Unit Name : chr [1:49212] "Matzen" "Abalone" "Aguilhada" "Agulha" ...
## $ Unit name local script : chr [1:49212] NA "Abalone" "Aguilhada" "Agulha" ...
## $ Fuel type : chr [1:49212] "oil and gas" "oil and gas" "oil and gas" "oil and gas" ...
## $ Unit type : chr [1:49212] "field" "field" "field" "field" ...
## $ Country : chr [1:49212] "Austria" "Brazil" "Brazil" "Brazil" ...
## $ Subnational unit (province, state): chr [1:49212] NA "Espírito Santo" "Sergipe" "Rio Grande do Norte" ...
## $ Latitude : num [1:49212] 48.4 -21.4 -10.7 -4.9 -22.1 ...
## $ Longitude : num [1:49212] 16.7 -39.6 -36.9 -36.3 -40 ...
## $ Location accuracy : chr [1:49212] "approximate" "exact" "exact" "exact" ...
## $ Status : chr [1:49212] "operating" "operating" "operating" "operating" ...
## $ Status year : num [1:49212] 2023 2022 2022 2022 2022 ...
## $ Discovery year : num [1:49212] 1949 2001 1966 1975 1984 ...
## $ FID Year : chr [1:49212] NA NA NA NA ...
## $ Production start year : chr [1:49212] "1951" "2009" "1969" "1979" ...
## $ Operator : chr [1:49212] "OMV" "Shell Brasil Petróleo Ltda." NA NA ...
## $ Owner : chr [1:49212] "OMV (100%)" "Shell Brasil (50%);ONGC Campos (27%);Qatarenergy (23%)" "Petrobras (100%)" "Petrobras (100%)" ...
## $ Parent : chr [1:49212] "OMV Aktiengesellschaft (100%)" "Shell plc (50%);Oil and Natural Gas Corporation (ONGC) (27%)" "Petróleo Brasileiro S.A. (100%)" "Petróleo Brasileiro S.A. (100%)" ...
## $ Basin : chr [1:49212] NA NA NA NA ...
## $ Concession / block : chr [1:49212] NA NA NA NA ...
## $ Project or complex : chr [1:49212] "Matzen" NA NA NA ...
## $ Government unit ID : chr [1:49212] NA NA NA NA ...
## $ Wiki URL : chr [1:49212] "https://www.gem.wiki/Matzen_Oil_and_Gas_Field_(Austria)" "https://www.gem.wiki/Abalone_Oil_and_Gas_Field_%28Esp%C3%ADrito_Santo%2C_Brazil%29" "https://www.gem.wiki/Aguilhada_Oil_and_Gas_Field_%28Sergipe%2C_Brazil%29" "https://www.gem.wiki/Agulha_Oil_and_Gas_Field_%28Rio_Grande_do_Norte%2C_Brazil%29" ...
## $ Unit name2 : chr [1:49212] NA NA NA NA ...
## $ Production/reserves : chr [1:49212] NA NA NA NA ...
## $ Fuel description : chr [1:49212] NA NA NA NA ...
## $ Reserves classification (original): chr [1:49212] NA NA NA NA ...
## $ Quantity (original) : num [1:49212] NA NA NA NA NA NA NA NA NA NA ...
## $ Units (original) : chr [1:49212] NA NA NA NA ...
## $ Data year : chr [1:49212] NA NA NA NA ...
## $ Quantity (converted) : num [1:49212] NA NA NA NA NA NA NA NA NA NA ...
## $ Units (converted) : chr [1:49212] NA NA NA NA ...
Se aplica un modelo de crecimiento/decrecimiento exponencial sobre los datos agrupados mediante técnica de binning, para reducir el ruido estadístico y capturar la tendencia general del comportamiento de las reservas petroleras a lo largo del tiempo.
Se estableció el Año de Descubrimiento (Discovery year) como variable independiente (x), por representar el momento histórico en que cada yacimiento fue identificado. Las Reservas (Quantity converted, en millones de barriles) actúan como variable dependiente (y), ya que expresan el volumen de hidrocarburo estimado para ese yacimiento.
Esta relación busca modelar cómo ha evolucionado el tamaño de los yacimientos descubiertos a lo largo de la historia de la exploración petrolera: la hipótesis de trabajo es que los campos gigantes fueron descubiertos preferentemente en las primeras décadas de exploración, mientras que los descubrimientos más recientes tienden a corresponder a yacimientos de menor tamaño remanente, lo que generaría una tendencia de tipo exponencial decreciente.
Nota sobre datos faltantes: por indicación del docente, los registros con valores faltantes no se eliminan. En su lugar, tanto el año de descubrimiento como la reserva se completan con la mediana de los valores disponibles, preservando así el tamaño de la muestra.
# Año de descubrimiento por yacimiento (dato único a nivel de campo)
disc_por_unidad <- Datos %>%
group_by(`Unit ID`) %>%
summarise(Descubrimiento = suppressWarnings(max(as.numeric(`Discovery year`), na.rm = TRUE))) %>%
mutate(Descubrimiento = ifelse(is.infinite(Descubrimiento), NA, Descubrimiento))
# Reservas máximas reportadas (millones de barriles) por yacimiento
reservas_por_unidad <- Datos %>%
filter(`Production/reserves` == "reserves", `Units (converted)` == "million bbl") %>%
group_by(`Unit ID`) %>%
summarise(Reservas = suppressWarnings(max(as.numeric(`Quantity (converted)`), na.rm = TRUE))) %>%
mutate(Reservas = ifelse(is.infinite(Reservas), NA, Reservas))
# Se parte de TODOS los yacimientos del dataset (no se elimina ninguno)
todos_yacimientos <- Datos %>% distinct(`Unit ID`)
datos_unidad <- todos_yacimientos %>%
left_join(disc_por_unidad, by = "Unit ID") %>%
left_join(reservas_por_unidad, by = "Unit ID")
# Imputación con la mediana (en lugar de eliminar registros con NA)
mediana_disc <- median(datos_unidad$Descubrimiento, na.rm = TRUE)
mediana_res <- median(datos_unidad$Reservas, na.rm = TRUE)
datos_raw <- datos_unidad %>%
mutate(
x_raw = ifelse(is.na(Descubrimiento), mediana_disc, Descubrimiento),
y_raw = ifelse(is.na(Reservas), mediana_res, Reservas)
) %>%
filter(y_raw > 0)
cat("Yacimientos totales:", nrow(todos_yacimientos), "\n")## Yacimientos totales: 8334
## Mediana usada para anio de descubrimiento: 1976
## Mediana usada para reservas (millones bbl): 74.51
## Registros finales tras imputacion: 8290
datos_plot <- datos_raw %>% filter(y_raw < quantile(y_raw, 0.99))
par(mar = c(5, 5, 4, 2))
color_trans <- rgb(0.2, 0.6, 0.86, 0.4)
plot(datos_plot$x_raw, datos_plot$y_raw,
main = "Gráfica N.1: Diagrama de Dispersión de las Reservas (millones bbl)\nen función del Año de Descubrimiento",
xlab = "Año de Descubrimiento",
ylab = "Reservas (millones bbl)",
col = color_trans, pch = 16, cex = 0.6, cex.main = 0.9, frame.plot = FALSE)
grid(nx = NULL, ny = NULL, col = "#D7DBDD", lty = "dotted")Debido a la alta variabilidad observada en la Gráfica N.1 y al efecto de la imputación por mediana, se implementa la técnica de binning con el objetivo de disminuir el ruido estadístico y visualizar con mayor claridad la tendencia general de los datos.
# Agrupación técnica (cada 10 años)
datos_agrupados <- datos_raw %>%
mutate(disc_bin = floor(x_raw / 10) * 10) %>%
group_by(disc_bin) %>%
summarise(
y = mean(y_raw, na.rm = TRUE),
conteo = n()
) %>%
rename(x = disc_bin) %>%
filter(conteo >= 5)
x_val <- datos_agrupados$x
y_val <- datos_agrupados$y
datos_agrupados## # A tibble: 15 × 3
## x y conteo
## <dbl> <dbl> <int>
## 1 1860 2513. 6
## 2 1880 74.5 6
## 3 1900 157. 53
## 4 1910 162. 87
## 5 1920 2514. 81
## 6 1930 795. 154
## 7 1940 1013. 292
## 8 1950 1191. 609
## 9 1960 1009. 666
## 10 1970 266. 4091
## 11 1980 408. 619
## 12 1990 360. 490
## 13 2000 333. 588
## 14 2010 562. 430
## 15 2020 271. 114
La ecuación es \(y = a \cdot e^{b \cdot x}\). Para linealizar, aplicamos logaritmo solo a la variable dependiente: \(\ln(y) = \ln(a) + b \cdot x\).
Se presenta el ajuste del modelo incluyendo la banda de incertidumbre estadística (Intervalo de Confianza del 95%).
par(mar = c(5, 5, 4, 2))
plot(x_val, y_val,
main = "Gráfica N.2: Modelo Exponencial de las Reservas (millones bbl)\nen función del Año de Descubrimiento",
xlab = "Año de Descubrimiento",
ylab = "Reservas (millones bbl)",
col = "#3498DB", pch = 16, cex = 1.0, cex.main = 0.9, frame.plot = FALSE)
grid(nx = NULL, ny = NULL, col = "#D7DBDD", lty = "dotted")
x_seq <- seq(min(x_val), max(x_val), length.out = 500)
pred_log <- predict(modelo_exponencial,
newdata = data.frame(x_val = x_seq),
interval = "confidence",
level = 0.95)
y_pred_fit <- exp(pred_log[, "fit"])
y_pred_lwr <- exp(pred_log[, "lwr"])
y_pred_upr <- exp(pred_log[, "upr"])
# Intervalo de Confianza
polygon(c(x_seq, rev(x_seq)),
c(y_pred_lwr, rev(y_pred_upr)),
col = rgb(0.5, 0.5, 0.5, 0.2), border = NA)
lines(x_seq, y_pred_fit, col = "#E74C3C", lwd = 3)
legend("topleft",
legend = c("Datos promediados (binning)",
"Modelo Exponencial",
"I.C. 95%"),
col = c("#3498DB", "#E74C3C", "gray"),
pch = c(16, NA, 15),
lwd = c(NA, 3, NA),
pt.cex = c(1, NA, 2),
bty = "n")## El coeficiente de correlacion es: -0.08
# Coeficientes del modelo
log_a <- coef(modelo_exponencial)[1]
b_param <- coef(modelo_exponencial)[2]
a_param <- exp(log_a)
# Construcción elegante de ecuación
ecuacion <- paste0("y = ", round(a_param, 4),
" * e^(", round(b_param, 6), "x)")
cat("La ecuacion estimada del modelo es:\n\n", "y =", round(a_param, 4), "* e^(", round(b_param, 6), "x)")## La ecuacion estimada del modelo es:
##
## y = 13996.99 * e^( -0.00172 x)
tabla_resumen <- data.frame(
Variable = c("Año de Descubrimiento", "Reservas (millones bbl)"),
Tipo = c("Independiente (X)", "Dependiente (Y)"),
R = c("", round(r, 2)),
R2 = c("", round(r2, 2)),
Parametro_a = c("", round(a_param, 4)),
Exponente_b = c("", round(b_param, 6)),
Ecuacion = c("", ecuacion)
)
tabla_resumen %>%
gt() %>%
tab_header(title = md("**Tabla N.1 del Resumen del Modelo de Regresión Exponencial**")) %>%
cols_label(Ecuacion = "Ecuación") %>%
tab_source_note(source_note = "Autor: Estalin — UCE") %>%
cols_align(align = "center", columns = everything())| Tabla N.1 del Resumen del Modelo de Regresión Exponencial | ||||||
| Variable | Tipo | R | R2 | Parametro_a | Exponente_b | Ecuación |
|---|---|---|---|---|---|---|
| Año de Descubrimiento | Independiente (X) | |||||
| Reservas (millones bbl) | Dependiente (Y) | -0.08 | 0.01 | 13996.9867 | -0.00172 | y = 13996.9867 * e^(-0.00172x) |
| Autor: Estalin — UCE | ||||||
¿Qué volumen de reservas se estima para un yacimiento descubierto en 1980?
anio_test <- 1980
res_est <- a_param * exp(b_param * anio_test)
cat("Para un anio de descubrimiento de", anio_test, ", las reservas estimadas son:", round(res_est, 4), "millones de bbl")## Para un anio de descubrimiento de 1980 , las reservas estimadas son: 464.329 millones de bbl
# Nota: se evitan tildes y simbolos especiales (n, 2 en vez de ny R2)
# dentro de este texto generado por codigo, porque cat()/paste0() puede
# corromper esos caracteres en sesiones de R que no usan locale UTF-8
# (el resto del documento, en Markdown normal, si usa tildes sin problema).
tendencia_txt <- if (b_param < 0) "decreciente" else "creciente"
interpretacion_txt <- if (b_param < 0) {
"lo cual es consistente con la teoria de agotamiento exploratorio: los yacimientos de mayor tamanio tienden a descubrirse en las primeras etapas de exploracion de una cuenca, mientras que los descubrimientos mas recientes corresponden, en promedio, a campos remanentes de menor volumen"
} else {
"lo cual sugiere que, en promedio, los descubrimientos mas recientes registran mayores volumenes de reservas, posiblemente asociado a mejoras tecnologicas de exploracion y estimacion"
}
cat(paste0(
"Entre el **Anio de Descubrimiento** y las **Reservas (millones bbl)** existe una relacion de tipo exponencial ", tendencia_txt,
", explicada por un coeficiente de determinacion R2 de ", round(r2 * 100, 2), "%, lo que indica que el modelo logra explicar una proporcion ",
if (r2 < 0.3) "modesta" else if (r2 < 0.6) "moderada" else "alta",
" de la variabilidad de las reservas en funcion del anio de descubrimiento.\n\n",
"El coeficiente de correlacion fue de ", round(r * 100, 2), "%, evidenciando una relacion ",
if (abs(r) < 0.3) "debil" else if (abs(r) < 0.6) "moderada" else "fuerte",
" entre las variables analizadas, ", interpretacion_txt, ".\n\n",
"La ecuacion matematica estimada del modelo es: y = ", round(a_param, 4), " * e^(", round(b_param, 6), "x)\n\n",
"Cabe senialar que, al haberse imputado con la mediana una parte considerable de los valores faltantes de anio de descubrimiento y de reservas (en lugar de eliminarlos, segun lo indicado por el docente), el ajuste del modelo puede verse atenuado respecto a un analisis realizado unicamente sobre los registros originales completos. Esto se evidencia en la concentracion de observaciones alrededor de los valores de la mediana utilizada para la imputacion."
))Entre el Anio de Descubrimiento y las Reservas (millones bbl) existe una relacion de tipo exponencial decreciente, explicada por un coeficiente de determinacion R2 de 0.66%, lo que indica que el modelo logra explicar una proporcion modesta de la variabilidad de las reservas en funcion del anio de descubrimiento.
El coeficiente de correlacion fue de -8.14%, evidenciando una relacion debil entre las variables analizadas, lo cual es consistente con la teoria de agotamiento exploratorio: los yacimientos de mayor tamanio tienden a descubrirse en las primeras etapas de exploracion de una cuenca, mientras que los descubrimientos mas recientes corresponden, en promedio, a campos remanentes de menor volumen.
La ecuacion matematica estimada del modelo es: y = 13996.9867 * e^(-0.00172x)
Cabe senialar que, al haberse imputado con la mediana una parte considerable de los valores faltantes de anio de descubrimiento y de reservas (en lugar de eliminarlos, segun lo indicado por el docente), el ajuste del modelo puede verse atenuado respecto a un analisis realizado unicamente sobre los registros originales completos. Esto se evidencia en la concentracion de observaciones alrededor de los valores de la mediana utilizada para la imputacion.