1 Librerías

library(dplyr)
library(ggplot2)
library(gt)
library(stringr)

2 Carga de Datos

ruta_csv <- "C:/Users/PATRICIA/Desktop/pr-estadistica/Oil__Gas____Other_Regulated_Wells__Beginning_1860 (3).csv"

Datos <- read.csv(ruta_csv,
                   header = TRUE,
                   sep = ";",
                   dec = ".",
                   fileEncoding = "Latin1",
                   stringsAsFactors = FALSE)

cat("Dimensiones del dataset:", nrow(Datos), "filas x", ncol(Datos), "columnas\n")
## Dimensiones del dataset: 47407 filas x 55 columnas

3 Selección de Variables

El TVD (True Vertical Depth, ft) es la variable independiente (x): una magnitud física definida previamente, sin depender del costo. La Tarifa de Permiso (Permit Fee, USD) es la variable dependiente (y), pues se cobra en función de la profundidad alcanzada. A mayor TVD, mayor tarifa, con un crecimiento que se ajusta a una tasa proporcional constante (elasticidad), de ahí el modelo potencial.

  • Variable Independiente (X): TVD (ft).
  • Variable Dependiente (Y): Tarifa de Permiso (USD).
# Localizamos las columnas de forma robusta (los nombres pueden llegar
# alterados por R como "True.Vertical.Depth..ft" y "Permit.Fee")
col_tvd <- grep("True.?Vertical.?Depth", names(Datos), value = TRUE)[1]
col_fee <- grep("^Permit.?Fee$", names(Datos), value = TRUE)[1]

if (is.na(col_tvd)) stop("ERROR: No se encontró la columna de Profundidad Vertical Real (TVD).")
if (is.na(col_fee)) stop("ERROR: No se encontró la columna de Tarifa de Permiso (Permit Fee).")

# Selección de variables
datos_raw <- Datos %>%
  select(all_of(c(col_tvd, col_fee))) %>%
  setNames(c("tvd", "permit_fee")) %>%
  mutate(
    x_raw = abs(as.numeric(str_replace(as.character(tvd), ",", "."))),
    y_raw = abs(as.numeric(str_replace(as.character(permit_fee), ",", ".")))
  ) %>%
  filter(!is.na(x_raw) & !is.na(y_raw) & x_raw > 0 & y_raw > 0) %>%
  filter(x_raw <= 27500 & y_raw <= 11000)   # rango de dominio válido de cada variable

4 Tabla de Pares de Valores

Se indica el tamaño muestral obtenido y se muestran únicamente las primeras filas de los pares de valores (x, y).

cat("Tamaño muestral: N =", nrow(datos_raw), "pares de valores (TVD, Tarifa de Permiso)\n")
## Tamaño muestral: N = 11045 pares de valores (TVD, Tarifa de Permiso)
datos_raw %>%
  select(x_raw, y_raw) %>%
  head(10) %>%
  gt() %>%
  tab_header(title = md("**Tabla N°1: Primeras filas de los pares de valores (x, y)**")) %>%
  cols_label(x_raw = "TVD (ft)", y_raw = "Tarifa de Permiso (USD)") %>%
  cols_align(align = "center", columns = everything())
Tabla N°1: Primeras filas de los pares de valores (x, y)
TVD (ft) Tarifa de Permiso (USD)
1800 860
1445 475
2070 725
6321 1725
1206 475
1590 475
9914 380
9674 5230
1529 1050
4281 1225

5 Gráfica de Dispersión

par(mar = c(5, 5, 4, 2))

color_trans <- rgb(0.2, 0.6, 0.86, 0.4)

plot(datos_raw$x_raw, datos_raw$y_raw,
     main = "Gráfica N°1: Diagrama de Dispersión de la Tarifa de Permiso\nen función de la Profundidad Vertical Real (TVD)",
     xlab = "TVD (ft)",
     ylab = "Tarifa de Permiso (USD)",
     col = color_trans,
     pch = 16,
     cex = 0.6,
     cex.main = 0.9,
     frame.plot = FALSE)

grid(nx = NULL, ny = NULL, col = "#D7DBDD", lty = "dotted")

6 Conjetura

La nube de puntos observada en la Gráfica N°1 presenta una alta variabilidad, con un comportamiento disperso y ruidoso que dificulta identificar visualmente una tendencia clara. Al tratarse de una nube compleja/caótica, se procede con el tratamiento de los datos descrito a continuación.

6.1 Tratamiento de los Datos

Para reducir el ruido y resaltar la tendencia, se segmenta x en bins de 500 ft (único x), se promedia y por bin (único y), y se omiten outliers: bins con menos de 3 datos y valores de y fuera del percentil 5%-95%.

# Agrupamiento cada 500 unidades de TVD
datos_model <- datos_raw %>%
  mutate(x_bin = round(x_raw / 500) * 500) %>%
  group_by(x_bin) %>%
  summarise(
    y = mean(y_raw, na.rm = TRUE),
    conteo = n(),
    .groups = "drop"
  ) %>%
  rename(x = x_bin) %>%
  filter(conteo >= 3, x > 0)   # se excluye x = 0 (ln(0) no está definido)

# Limpieza de outliers
lim_y <- quantile(datos_model$y, probs = c(0.05, 0.95))

datos_model <- datos_model %>%
  filter(y >= lim_y[1] & y <= lim_y[2])

x <- datos_model$x
y <- datos_model$y

6.2 Nueva Gráfica de Dispersión

par(mar = c(5, 5, 4, 2))

plot(x, y,
     main = "Gráfica N°2: Dispersión de la Tarifa de Permiso promedio\npor intervalos de TVD (ft)",
     xlab = "TVD (ft)",
     ylab = "Tarifa de Permiso promedio (USD)",
     col = "#3498DB",
     pch = 16,
     cex = 0.9,
     cex.main = 0.9,
     frame.plot = FALSE)

grid(nx = NULL, ny = NULL, col = "#D7DBDD", lty = "dotted")

6.3 Nueva Conjetura

Tras el tratamiento de los datos, la nube de puntos muestra una tendencia creciente compatible con una tasa de crecimiento proporcional (elasticidad constante) entre las variables, patrón compatible con un modelo potencial:

\[y = a \cdot x^{b}\]

7 Cálculo de Parámetros

Para linealizar el modelo se aplica logaritmo natural a ambos lados de la ecuación:

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

modelo_pot <- lm(log(y) ~ log(x))

ln_a_bin <- coef(modelo_pot)[1]
b_bin <- coef(modelo_pot)[2]
a_bin <- exp(ln_a_bin)

ecuacion <- paste0("y = ",
                    round(a_bin, 4),
                    " * x^",
                    round(b_bin, 4))

cat("La ecuación estimada del modelo es:\n\n", ecuacion)
## La ecuación estimada del modelo es:
## 
##  y = 1.4542 * x^0.8175
r <- cor(log(x), log(y), use = "complete.obs")
r2 <- summary(modelo_pot)$r.squared

tabla_resumen <- data.frame(
  Variable = c("TVD (ft)", "Tarifa de Permiso (USD)"),
  Tipo = c("Independiente (x)", "Dependiente (y)"),
  R = c("", round(r, 2)),
  R2 = c("", round(r2, 2)),
  Coeficiente_a = c("", round(a_bin, 4)),
  Exponente_b = c("", round(b_bin, 4)),
  Ecuación = c("", ecuacion)
)

tabla_resumen %>%
  gt() %>%
  tab_header(title = md("**Tabla N°2: Resumen del Modelo de Regresión Potencial**")) %>%
  tab_source_note(source_note = "Autor: JENNY") %>%
  cols_align(align = "center", columns = everything())
Tabla N°2: Resumen del Modelo de Regresión Potencial
Variable Tipo R R2 Coeficiente_a Exponente_b Ecuación
TVD (ft) Independiente (x)
Tarifa de Permiso (USD) Dependiente (y) 0.96 0.92 1.4542 0.8175 y = 1.4542 * x^0.8175
Autor: JENNY

8 Realidad y Modelo

Se presenta el ajuste del modelo sobre los datos reales (agrupados), incluyendo la banda de incertidumbre estadística (Intervalo de Confianza del 95%). El intervalo se calcula sobre el modelo linealizado y se transforma de vuelta a la escala original con la función exponencial.

par(mar = c(5, 5, 4, 2))

plot(x, y,
     main = "Gráfica N°3: Modelo Potencial de la Tarifa de Permiso\nen función del TVD (ft)",
     xlab = "TVD (ft)",
     ylab = "Tarifa de Permiso (USD)",
     col = "#3498DB",
     pch = 16,
     cex = 1.0,
     cex.main = 0.9,
     frame.plot = FALSE)

grid(nx = NULL, ny = NULL, col = "#D7DBDD", lty = "dotted")

# Secuencia suave
x_seq <- seq(min(x), max(x), length.out = 500)

pred_pot <- predict(modelo_pot,
                     newdata = data.frame(x = x_seq),
                     interval = "confidence",
                     level = 0.95)

# Se transforma de la escala logarítmica a la escala original
pred_pot_original <- exp(pred_pot)

# Intervalo de confianza
polygon(c(x_seq, rev(x_seq)),
        c(pred_pot_original[, "lwr"], rev(pred_pot_original[, "upr"])),
        col = rgb(0.5, 0.5, 0.5, 0.2),
        border = NA)

# Línea ajustada
lines(x_seq, pred_pot_original[, "fit"], col = "#E74C3C", lwd = 3)

legend("topleft",
       legend = c("Datos promediados (binning)",
                  "Modelo Potencial",
                  "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")

9 Test

9.1 Coeficiente de Correlación del Modelo Linealizado

cat("El coeficiente de correlación es: ", round(r, 2))
## El coeficiente de correlación es:  0.96

9.2 Coeficiente de Determinación

cat(paste0("El coeficiente de determinación (R²) es: ", round(r2, 2)))
## El coeficiente de determinación (R²) es: 0.92

10 Restricciones

El modelo presenta como única condición matemática que x > 0 y y > 0, lo cual se cumple naturalmente en el contexto físico del TVD (ft) y de la Tarifa de Permiso (USD), por lo que no existen restricciones prácticas adicionales dentro del rango analizado.

11 Estimación

¿Cuál es la Tarifa de Permiso estimada para un pozo con un TVD de 3000 ft?

tvd_test <- 3000
ln_y_est <- predict(modelo_pot,
                     newdata = data.frame(x = tvd_test))
y_est <- exp(ln_y_est)

cat("Para un TVD de", tvd_test,
    "ft, la Tarifa de Permiso estimada es:",
    round(y_est, 4), "USD")
## Para un TVD de 3000 ft, la Tarifa de Permiso estimada es: 1012.127 USD

12 Conclusión

Entre el TVD (ft) y la Tarifa de Permiso (USD) existe una relación de tipo potencial, con un coeficiente de determinación R² = 0.92, lo que indica un ajuste excelente del modelo.

La ecuación estimada es: y = 1.4542 * x^0.8175.

El modelo presenta como única condición matemática que x > 0 y y > 0, lo cual se cumple naturalmente en el contexto físico del TVD (ft) y de la Tarifa de Permiso (USD), por lo que no existen restricciones prácticas adicionales dentro del rango analizado.