1 Librerías

library(dplyr)
library(ggplot2)
library(gt)
library(stringr)
library(zoo)     # para la media móvil (rollmean)
library(DT)      # para la tabla interactiva del extracto de datos

2 Carga de Datos

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

# Detección automática del separador, a partir de la primera línea del archivo
primera_linea <- readLines(ruta_csv, n = 1, encoding = "Latin1")

candidatos <- c(";", ",", "\t", "|")
conteos <- sapply(candidatos, function(s) lengths(regmatches(primera_linea, gregexpr(s, primera_linea, fixed = TRUE))))
sep_detectado <- candidatos[which.max(conteos)]

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

cat("Separador detectado:", sep_detectado, "| Columnas detectadas:", ncol(Datos), "\n")
## Separador detectado: ; | Columnas detectadas: 55
cat("Número de registros:", nrow(Datos), "\n")
## Número de registros: 47407
cat("Número de variables:", ncol(Datos), "\n")
## Número de variables: 55

Extracto del dataset (55 variables, primeras 5 filas):

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 Profundidad (Depth Fee, USD) es la variable dependiente (y), pues se cobra en función de la profundidad alcanzada. A mayor TVD, mayor tarifa, pero con crecimiento decreciente, de ahí el modelo logarítmico.

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

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

cat("Registros luego de seleccionar y depurar las variables:", nrow(datos_raw), "\n")
## Registros luego de seleccionar y depurar las variables: 9931

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 Profundidad)\n")
## Tamaño muestral: N = 9931 pares de valores (TVD, Tarifa de Profundidad)
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 Profundidad (USD)") %>%
  cols_align(align = "center", columns = everything())
Tabla N°1: Primeras filas de los pares de valores (x, y)
TVD (ft) Tarifa de Profundidad (USD)
1800 760
1445 375
2070 625
6321 1625
1206 375
1590 375
9914 280
9674 5130
1529 950
4281 1125

5 Gráfica de Dispersión

plot(datos_raw$x_raw, datos_raw$y_raw,
     pch = 19,
     col = "darkorange",
     main = "Tarifa de Profundidad vs TVD",
     xlab = "TVD (ft)",
     ylab = "Tarifa de Profundidad (USD)")

La nube de puntos es muy densa y ruidosa: para el mismo TVD existen pozos con tarifas muy distintas, por lo que a simple vista no se aprecia con claridad ninguna tendencia.

6 Conjetura del Modelo de Regresión

6.1 Depuración: Reducción de la Cantidad de Puntos

Para que la tendencia sea visible, se depuran los datos reduciendo su cantidad: se ordenan por x y se calcula el promedio de y dentro de una ventana deslizante (media móvil), tras eliminar outliers globales en y (fuera del percentil 2.5%-97.5%). Esto conserva la posición real de cada x, pero sustituye miles de puntos dispersos por una nube mucho más compacta y representativa de la tendencia central.

# Orden de la ventana (impar, para tener un centro simétrico)
ventana <- 51

# 1) Eliminación de outliers globales en y
lim_y_raw <- quantile(datos_raw$y_raw, probs = c(0.025, 0.975))

datos_limpios <- datos_raw %>%
  filter(y_raw >= lim_y_raw[1] & y_raw <= lim_y_raw[2]) %>%
  arrange(x_raw)

# 2) Media móvil de y sobre los datos ordenados por x -> reduce la cantidad de puntos
pares_dep <- datos_limpios %>%
  mutate(
    x = x_raw,
    y = rollmean(y_raw, k = ventana, fill = NA, align = "center")
  ) %>%
  filter(!is.na(y), x > 0) %>%   # se excluye x = 0 (ln(0) no está definido)
  select(x, y) %>%
  distinct(x, .keep_all = TRUE)  # un único punto representativo por x

x <- pares_dep$x
y <- pares_dep$y

cat("N original:", nrow(datos_raw), " -> N depurado:", nrow(pares_dep), "pares (x, y)\n")
## N original: 9931  -> N depurado: 3467 pares (x, y)

6.2 Ajuste del Eje X para la Gráfica

El TVD real va de valores cercanos a cero hasta más de 10 000 ft. Aun así, si el logaritmo se evalúa demasiado lejos del origen la curvatura característica del logaritmo puede verse casi como una recta. Para que la curvatura sea claramente visible, se reasignan los valores de x únicamente para fines de graficación, buscando el menor desplazamiento (offset) que logre una correlación cor(log(x), y) > 0.75.

Importante: esta reasignación se usa solo en las gráficas de esta sección y de la Sección 8. El resto del análisis (parámetros del modelo, test, restricciones, predicción) se sigue calculando con los valores reales de TVD, ya que es la variable con significado físico real.

# Búsqueda del menor desplazamiento (offset) que logre correlación > 0.75
# en la versión graficada, manteniendo el eje lo más cercano posible a cero
buscar_offset <- function(x, y, objetivo = 0.75, max_offset = 500) {
  for (c in seq(1, max_offset, by = 1)) {
    xg <- x - min(x) + c
    r  <- cor(log(xg), y)
    if (r > objetivo) return(c)
  }
  return(NA)
}

offset_optimo <- buscar_offset(pares_dep$x, pares_dep$y)
cat("Menor desplazamiento que logra correlación > 0.75:", offset_optimo, "\n")
## Menor desplazamiento que logra correlación > 0.75: 1
pares_dep <- pares_dep %>%
  mutate(x_grafico = x - min(x) + offset_optimo)

cat("Primer TVD registrado:", min(pares_dep$x), "-> reasignado a:", min(pares_dep$x_grafico), "\n")
## Primer TVD registrado: 340 -> reasignado a: 1
cat("Último TVD registrado :", max(pares_dep$x), "-> reasignado a:", max(pares_dep$x_grafico), "\n")
## Último TVD registrado : 10260 -> reasignado a: 9921
r_grafico <- cor(log(pares_dep$x_grafico), pares_dep$y)
cat("Correlación de Pearson del eje reajustado    :", round(r_grafico, 4), "\n")
## Correlación de Pearson del eje reajustado    : 0.8899

6.3 Gráfica Simplificada

plot(pares_dep$x_grafico, pares_dep$y,
     pch = 20,
     col = rgb(0.1, 0.4, 0.5, 0.6),
     xlab = "TVD reasignado desde el primer registro (X)",
     ylab = "Tarifa de Profundidad (Y)",
     main = "Relación entre TVD Reasignado y Tarifa de Profundidad (Datos Depurados)")

Con la depuración y el reajuste del eje, el número de puntos se redujo considerablemente respecto a la gráfica original y ahora se aprecia con claridad que la tendencia no es una recta: el crecimiento de Y es más rápido en los primeros valores tras el reajuste y se va aplanando a medida que X aumenta, un comportamiento típico de una curva logarítmica.

6.4 Nueva Conjetura

Observando la gráfica depurada, se propone un Modelo de Regresión Logarítmico, ya que los datos muestran un crecimiento rápido al inicio que luego se desacelera (rendimientos decrecientes), sin cambios de dirección ni curvatura en “S” — justo la forma característica de una curva logarítmica. Este modelo tiene la forma:

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

7 Cálculo de Parámetros

El modelo se ajusta con los valores reales de TVD (no con el eje reasignado de la sección 6.2), ya que esa reasignación es solo un recurso visual.

modelo_log <- lm(y ~ log(x), data = pares_dep)

a_bin <- coef(modelo_log)[1]
b_bin <- coef(modelo_log)[2]

if (b_bin >= 0) {
  ecuacion <- paste0("y = ",
                      round(a_bin, 4),
                      " + ",
                      round(b_bin, 4),
                      " ln(x)")
} else {
  ecuacion <- paste0("y = ",
                      round(a_bin, 4),
                      " - ",
                      abs(round(b_bin, 4)),
                      " ln(x)")
}

cat("La ecuación estimada del modelo es:\n\n", ecuacion)
## La ecuación estimada del modelo es:
## 
##  y = -3688.3962 + 589.1185 ln(x)
r <- cor(log(x), y, use = "complete.obs")
r2 <- summary(modelo_log)$r.squared

tabla_resumen <- data.frame(
  Variable = c("TVD (ft)", "Tarifa de Profundidad (USD)"),
  Tipo = c("Independiente (x)", "Dependiente (y)"),
  R = c("", round(r, 2)),
  R2 = c("", round(r2, 2)),
  Intercepto_a = c("", round(a_bin, 4)),
  Pendiente_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 Logarítmica**")) %>%
  tab_source_note(source_note = "Autor: JENNY") %>%
  cols_align(align = "center", columns = everything())
Tabla N°2: Resumen del Modelo de Regresión Logarítmica
Variable Tipo R R2 Intercepto_a Pendiente_b Ecuación
TVD (ft) Independiente (x)
Tarifa de Profundidad (USD) Dependiente (y) 0.94 0.89 -3688.3962 589.1185 y = -3688.3962 + 589.1185 ln(x)
Autor: JENNY

8 Realidad y Modelo

Se presenta el ajuste del modelo logarítmico sobre la nube de puntos depurada, usando el mismo eje reasignado (x_grafico) de la Sección 6.2 para que la curvatura del modelo sea claramente visible junto a los puntos. La curva se evalúa con los coeficientes a y b calculados en la Sección 7 sobre el TVD real, revertiendo el desplazamiento para ubicarla correctamente sobre el eje graficado.

plot(pares_dep$x_grafico, pares_dep$y,
     pch = 19,
     col = "darkorange",
     main = "Tarifa de Profundidad vs TVD (Datos Depurados)",
     xlab = "TVD reasignado desde el primer registro (X)",
     ylab = "Tarifa de Profundidad (USD)")

min_x <- min(pares_dep$x)

curve(a_bin + b_bin * log(x - offset_optimo + min_x),
      add = TRUE,
      col = "darkslateblue",
      lwd = 2)

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.94

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.89

10 Restricciones

10.1 Dominio Teórico de las Variables

  • X — TVD (True Vertical Depth, ft): \(D = \{x \in \mathbb{R} \mid x \geq 0\}\), \(R = \{x \mid 0 \leq x \leq 27\,500\}\) ft.
  • Y — Depth Fee (Tarifa de Profundidad, USD): \(D = \{y \mid y \geq 0\}\), \(R = \{y \mid 0 \leq y \leq 10\,900\}\) USD.
# Verificación empírica contra el dataset real (antes y después del tratamiento)
cat("Dataset crudo -> TVD: [", round(min(datos_raw$x_raw), 1), ",", round(max(datos_raw$x_raw), 1), "] ft\n")
## Dataset crudo -> TVD: [ 8 , 11953 ] ft
cat("Dataset crudo -> Depth Fee: [", round(min(datos_raw$y_raw), 1), ",", round(max(datos_raw$y_raw), 1), "] USD\n\n")
## Dataset crudo -> Depth Fee: [ 90 , 10900 ] USD
cat("Datos depurados (Sección 6.1) -> TVD: [", round(min(pares_dep$x), 1), ",", round(max(pares_dep$x), 1), "] ft\n")
## Datos depurados (Sección 6.1) -> TVD: [ 340 , 10260 ] ft
cat("Datos depurados (Sección 6.1) -> Depth Fee: [", round(min(pares_dep$y), 1), ",", round(max(pares_dep$y), 1), "] USD\n")
## Datos depurados (Sección 6.1) -> Depth Fee: [ 272.5 , 1916.5 ] USD

El máximo real de Depth Fee en el dataset coincide exactamente con el límite superior de 10 900 USD documentado en la Tabla de Variables, confirmando que ese valor no es arbitrario sino el techo real del fenómeno.

10.2 Restricción Implícita del Modelo sobre X

El dominio matemático de X exige, además, \(x > 0\) (por \(\ln(x)\)). Dado que el modelo es una curva logarítmica monótona creciente (\(b > 0\)), existe la posibilidad de que, para valores de X fuera del rango de los datos observados, el valor estimado de Y se salga de su dominio real \([0,\ 10\,900]\) USD. Para comprobarlo, se resuelve la ecuación igualándola a los límites del dominio de Y, y se despeja X:

\[a + b\ln(x) = 10900 \qquad y \qquad a + b\ln(x) = 0\]

Despejando x en cada caso:

\[x = e^{\frac{10900 - a}{b}} \qquad y \qquad x = e^{\frac{0 - a}{b}}\]

x_lim_sup <- exp((10900 - a_bin) / b_bin)
x_lim_inf <- exp((0 - a_bin) / b_bin)

cat("Valor de X donde el modelo predice Y = 10900 USD:", round(x_lim_sup, 0), "ft\n")
## Valor de X donde el modelo predice Y = 10900 USD: 56816529793 ft
cat("Valor de X donde el modelo predice Y = 0 USD    :", round(x_lim_inf, 0), "ft\n")
## Valor de X donde el modelo predice Y = 0 USD    : 524 ft

Como el modelo es monótono (una sola curva sin cambios de dirección), el intervalo válido de X queda delimitado directamente por estos dos valores, y además debe mantenerse dentro del dominio teórico \([0,\ 27\,500]\) ft de la Tabla de Variables:

x_min_valido <- max(0, min(x_lim_inf, x_lim_sup))
x_max_valido <- min(27500, max(x_lim_inf, x_lim_sup))

cat("El modelo se mantiene dentro del dominio de Y [0, 10900] USD\n")
## El modelo se mantiene dentro del dominio de Y [0, 10900] USD
cat("y del dominio teórico de X [0, 27500] ft\n")
## y del dominio teórico de X [0, 27500] ft
cat("para valores de X (TVD) comprendidos entre:", round(x_min_valido, 0), "y", round(x_max_valido, 0), "ft\n")
## para valores de X (TVD) comprendidos entre: 524 y 27500 ft

Conclusión de la sección: el modelo sí presenta restricciones. Es válido únicamente para TVD comprendidos entre 524 y 2.75^{4} ft — un intervalo respaldado tanto por el dominio real y como por el rango efectivamente observado en el dataset. Fuera de este intervalo, el modelo produce valores de Tarifa de Profundidad que exceden el dominio documentado, por lo que no debe usarse para estimar fuera de ese rango.

11 Estimación del Modelo

Aprovechando la ecuación del modelo logarítmico, se realizan estimaciones dentro del rango válido determinado en la sección anterior (entre 524 y 2.75^{4} ft).

tvd_test <- 3000

if (tvd_test >= x_min_valido & tvd_test <= x_max_valido) {

  y_est <- predict(modelo_log, newdata = data.frame(x = tvd_test))

  cat("Estimación para un TVD de", tvd_test, "ft:\n")
  cat("Tarifa de Profundidad estimada:", round(y_est, 4), "USD\n")

} else {

  cat("El TVD de", tvd_test, "ft está fuera del rango válido del modelo",
      "[", round(x_min_valido, 0), ",", round(x_max_valido, 0), "] ft.\n")
  cat("No se realiza la estimación porque el resultado excedería el dominio real de Y.\n")

}
## Estimación para un TVD de 3000 ft:
## Tarifa de Profundidad estimada: 1028.303 USD

12 Conclusión

Entre el TVD (ft) y la Tarifa de Profundidad (USD) existe una relación de tipo logarítmica, con un coeficiente de determinación R² = 0.89, lo que indica un ajuste aceptable del modelo.

La ecuación estimada es: y = -3688.3962 + 589.1185 ln(x).

El modelo presenta restricciones: es válido únicamente para valores de TVD comprendidos entre 524 y 27500 ft. Fuera de ese intervalo, el modelo predice tarifas que exceden el dominio real de la variable Y [0, 10900] USD, por lo que no debe usarse para estimar fuera de ese rango.