1 CARGAR LIBRERÍAS

library(tidyverse)
library(randomForest)
library(caret)
library(ranger)
library(Metrics)
library(ggplot2)
library(corrplot)
library(skimr)
library(knitr)
library(kableExtra)
library(rpart)
library(rpart.plot)
library(DT)
library(scales)

2 CARGAR DATASET

datos_originales <- read_csv2("Oil__Gas____Other_Regulated_Wells__Beginning_1860.csv")

dims <- dim(datos_originales)
cat("Dimensiones del dataset original:", dims[1], "filas y", dims[2], "columnas.\n")

## Dimensiones del dataset original: 47407 filas y 55 columnas.

datatable(
  head(datos_originales, 50),
  caption = htmltools::tags$caption(
    style = 'caption-side: top; text-align: center; font-weight: bold; color: #1b365d; font-size: 15px;',
    'Tabla 1: Visualización interactiva y estilizada de las primeras observaciones del dataset original'
  ),
  options = list(
    pageLength = 5,
    autoWidth = TRUE,
    scrollX = TRUE,
    dom = 'Bfrtip',
    language = list(url = '//cdn.datatables.net/plug-ins/1.10.25/i18n/Spanish.json')
  ),
  class = 'cell-border stripe hover compact'
)

3 PROCESAMIENTO DE DATOS

variables_requeridas <- c(
  "True Vertical Depth, ft",
  "Proposed Depth",
  "Elevation",
  "Surface Longitude",
  "Surface Latitude",
  "County",
  "Region",
  "Well Type",
  "Original Well Type",
  "Objective Formation",
  "Producing Formation",
  "Slant",
  "Spacing Acres"
)

datos_proc <- datos_originales %>%
  select(any_of(variables_requeridas))

names(datos_proc) <- c(
  "True_Vertical_Depth",
  "Proposed_Depth",
  "Elevation",
  "Surface_Longitude",
  "Surface_Latitude",
  "County",
  "Region",
  "Well_Type",
  "Original_Well_Type",
  "Objective_Formation",
  "Producing_Formation",
  "Slant",
  "Spacing_Acres"
)

datos_proc <- datos_proc %>% distinct()

numeric_vars <- c("True_Vertical_Depth", "Proposed_Depth", "Elevation", "Surface_Longitude", "Surface_Latitude", "Spacing_Acres")
for(v in numeric_vars) {
  if(v %in% names(datos_proc)) {
    datos_proc[[v]] <- as.numeric(as.character(datos_proc[[v]]))
  }
}

cat_vars <- c("County", "Region", "Well_Type", "Original_Well_Type", "Objective_Formation", "Producing_Formation", "Slant")
for(v in cat_vars) {
  if(v %in% names(datos_proc)) {
    datos_proc[[v]] <- addNA(as.factor(datos_proc[[v]]))
    levels(datos_proc[[v]])[is.na(levels(datos_proc[[v]]))] <- "Desconocido"
  }
}

datos_limpios <- datos_proc %>% 
  filter(!is.na(True_Vertical_Depth))

for(v in numeric_vars) {
  if(v %in% names(datos_limpios) && sum(is.na(datos_limpios[[v]])) > 0) {
    datos_limpios[[v]][is.na(datos_limpios[[v]])] <- median(datos_limpios[[v]], na.rm = TRUE)
  }
}

cat("Dimensiones finales del dataset limpio:", nrow(datos_limpios), "filas y", ncol(datos_limpios), "columnas.\n")

## Dimensiones finales del dataset limpio: 46365 filas y 11 columnas.

3.1 Distribución de la variable Objetivo

ggplot(datos_limpios, aes(x = True_Vertical_Depth)) +
  geom_histogram(bins = 40, fill = "steelblue", color = "black", alpha = 0.7) +
  theme_minimal() +
  labs(title = "Distribución de la Profundidad Vertical Verdadera", x = "Profundidad (ft)", y = "Frecuencia")

4 DIVISIÓN DE DATOS

set.seed(123)

if(nrow(datos_limpios) > 10) {
  train_index <- createDataPartition(datos_limpios$True_Vertical_Depth, p = 0.70, list = FALSE)
  train_set <- datos_limpios[train_index, ]
  test_set <- datos_limpios[-train_index, ]
} else {
  train_index <- sample(1:nrow(datos_limpios), size = 0.7 * nrow(datos_limpios))
  train_set <- datos_limpios[train_index, ]
  test_set <- datos_limpios[-train_index, ]
}

partition_summary <- data.frame(
  Conjunto = c("Entrenamiento (Train)", "Prueba (Test)", "Total"),
  Observaciones = c(nrow(train_set), nrow(test_set), nrow(datos_limpios)),
  Porcentaje = c("70%", "30%", "100%")
)

partition_summary %>%
  tabla_formato("Tabla 2: Resumen de la división del dataset")

Tabla 2: Resumen de la división del dataset
Conjunto	Observaciones	Porcentaje
Entrenamiento (Train)	32457	70%
Prueba (Test)	13908	30%
Total	46365	100%

5 CONSTRUCCIÓN DEL MODELO

set.seed(123)

p_preds <- ncol(train_set) - 1
mtry_opt <- floor(p_preds / 3)

rf_model_ranger <- ranger(
  formula = True_Vertical_Depth ~ .,
  data = train_set,
  num.trees = 500,
  mtry = mtry_opt,
  importance = "impurity",
  sample.fraction = 1.0,
  seed = 123
)

6 ENTRENAMIENTO DEL MODELO

6.1 Entrenamiento del Random Forest

print(rf_model_ranger)

## Ranger result
## 
## Call:
##  ranger(formula = True_Vertical_Depth ~ ., data = train_set, num.trees = 500,      mtry = mtry_opt, importance = "impurity", sample.fraction = 1,      seed = 123) 
## 
## Type:                             Regression 
## Number of trees:                  500 
## Sample size:                      32457 
## Number of independent variables:  10 
## Mtry:                             3 
## Target node size:                 5 
## Variable importance mode:         impurity 
## Splitrule:                        variance 
## OOB prediction error (MSE):       416766.3 
## R squared (OOB):                  0.8136499

6.2 Error Out-of-Bag (OOB)

cat("Error OOB (Prediction error / MSE aproximado):", rf_model_ranger$prediction.error, "\n")

## Error OOB (Prediction error / MSE aproximado): 416766.3

6.3 Importancia de las variables

importance_df <- data.frame(
  Variable = names(rf_model_ranger$variable.importance),
  IncNodePurity = as.numeric(rf_model_ranger$variable.importance)
) %>% arrange(desc(IncNodePurity))

head(importance_df, 10) %>%
  tabla_formato("Tabla 3: Top 10 Variables más influyentes")

Tabla 3: Top 10 Variables más influyentes
Variable	IncNodePurity
County	18339092467
Elevation	10325499181
Well_Type	10108900839
Original_Well_Type	8889413454
Proposed_Depth	8397100213
Region	2133786874
Objective_Formation	1902914164
Surface_Latitude	1589563748
Producing_Formation	1065163132
Surface_Longitude	0

ggplot(head(importance_df, 10), aes(x = reorder(Variable, IncNodePurity), y = IncNodePurity)) +
  geom_point(color = "steelblue", size = 3) +
  geom_segment(aes(x = Variable, xend = Variable, y = 0, yend = IncNodePurity), color = "steelblue", size = 1) +
  scale_x_discrete() +
  scale_y_continuous(labels = label_number(big.mark = ",")) +
  coord_flip() +
  theme_minimal() +
  labs(title = "Top 10 Variables de Importancia", x = "Variables", y = "Importancia")

7 PREDICCIONES

predictions <- predict(rf_model_ranger, data = test_set)$predictions

results_df <- data.frame(
  Valor_Real = test_set$True_Vertical_Depth,
  Valor_Predicho = predictions
) %>%
  mutate(
    Error_Absoluto = abs(Valor_Real - Valor_Predicho),
    Error_Relativo_Porc = ifelse(Valor_Real != 0, (Error_Absoluto / abs(Valor_Real)) * 100, 0)
  )

head(results_df, 10) %>%
  tabla_formato("Tabla 4: Comparación entre el Valor Real y el Predicho")

Tabla 4: Comparación entre el Valor Real y el Predicho
Valor_Real	Valor_Predicho	Error_Absoluto	Error_Relativo_Porc
2006	2423.9495	417.94954	20.834972
0	692.5764	692.57639	0.000000
6321	4751.1313	1569.86869	24.835765
928	863.7021	64.29795	6.928658
979	775.5399	203.46014	20.782445
1473	1665.7382	192.73823	13.084741
1608	1562.3738	45.62616	2.837448
5293	4917.5158	375.48424	7.093978
0	820.5711	820.57115	0.000000
1842	1682.2695	159.73053	8.671582

ggplot(results_df, aes(x = Valor_Real, y = Valor_Predicho)) +
  geom_point(alpha = 0.4, color = "darkblue") +
  geom_abline(intercept = 0, slope = 1, color = "red", linetype = "dashed", size = 1) +
  theme_minimal() +
  labs(title = "Valores Observados y Predichos", x = "Valor Real (ft)", y = "Valor Predicho (ft)")

results_df$Index <- 1:nrow(results_df)
n_muestras <- min(30, nrow(results_df))

ggplot(results_df[1:n_muestras, ], aes(x = Index)) +
  geom_line(aes(y = Valor_Real, color = "Real"), size = 1) +
  geom_line(aes(y = Valor_Predicho, color = "Predicho"), size = 1, linetype = "dashed") +
  theme_minimal() +
  labs(title = "Comparación Directa: Real y Predicho", x = "Índice de Pozo", y = "Profundidad (ft)", color = "Leyenda") +
  scale_color_manual(values = c("Real" = "darkblue", "Predicho" = "coral3"))

8 EVALUACIÓN DEL MODELO

y_true <- test_set$True_Vertical_Depth
y_pred <- predictions

mse_val <- mse(y_true, y_pred)
rmse_val <- rmse(y_true, y_pred)
mae_val <- mae(y_true, y_pred)
r2_val <- cor(y_true, y_pred)^2
pearson_corr <- cor(y_true, y_pred)
residuals_vec <- y_true - y_pred

8.1 Error cuadrático medio (MSE)

cat("MSE:", mse_val, "\n")

## MSE: 419711.2

8.2 Raíz del error cuadrático medio (RMSE)

cat("RMSE:", rmse_val, "ft\n")

## RMSE: 647.8512 ft

8.3 Error absoluto medio (MAE)

cat("MAE:", mae_val, "ft\n")

## MAE: 383.2403 ft

8.4 Coeficiente de determinación (R²)

cat("R²:", r2_val, "\n")

## R²: 0.8122774

8.5 Correlación de Pearson

cat("Correlación de Pearson:", pearson_corr, "\n")

## Correlación de Pearson: 0.9012644

8.6 Análisis de residuos

resumen_residuos <- data.frame(
  Métrica_Estadística = c("Mínimo", "1er Cuartil", "Mediana", "Media", "3er Cuartil", "Máximo", "Desviación Estándar"),
  Valor_ft = c(
    min(residuals_vec, na.rm = TRUE),
    quantile(residuals_vec, 0.25, na.rm = TRUE),
    median(residuals_vec, na.rm = TRUE),
    mean(residuals_vec, na.rm = TRUE),
    quantile(residuals_vec, 0.75, na.rm = TRUE),
    max(residuals_vec, na.rm = TRUE),
    sd(residuals_vec, na.rm = TRUE)
  )
)

resumen_residuos %>%
  tabla_formato("Tabla 5: Análisis estadístico de los residuos del modelo")

Tabla 5: Análisis estadístico de los residuos del modelo
Métrica_Estadística	Valor_ft
Mínimo	-7427.767290
1er Cuartil	-166.395098
Mediana	-4.793049
Media	5.329854
3er Cuartil	255.071845
Máximo	9746.930113
Desviación Estándar	647.852597

res_df <- data.frame(Residuo = residuals_vec)
ggplot(res_df, aes(x = Residuo)) +
  geom_histogram(bins = 40, fill = "skyblue", color = "black", alpha = 0.7) +
  theme_minimal() +
  labs(title = "Histograma de Residuos", x = "Residuo (ft)", y = "Densidad")

ggplot(res_df, aes(sample = Residuo)) +
  stat_qq(color = "darkblue", alpha = 0.5) +
  stat_qq_line(color = "red", size = 1) +
  theme_minimal() +
  labs(title = "Gráfico Q-Q Normal de los Residuos", x = "Cuantiles Teóricos", y = "Cuantiles Muestrales")

Nota analítica: Si bien el modelo presenta un excelente ajuste global (\(R^2\) de ~0.81), el gráfico Q-Q evidencia colas pesadas en los residuos, lo que indica la presencia de ciertos pozos con errores atípicos de predicción debido a la alta variabilidad geológica.

9 CONCLUSIÓN

El modelo de Random Forest funcionó muy bien para calcular la profundidad real de los pozos petroleros utilizando la información histórica disponible. Logró capturar con bastante precisión los factores clave del subsuelo y ofreció resultados muy cercanos a la realidad, aunque en ciertos pozos específicos se encontraron diferencias mayores debido a las complicaciones naturales del terreno y la geología de la zona.

Autor: Jennifer Cordones | Machine Learning de Regresión (Random Forest) — Oil, Gas & Other Regulated Wells - NY State

Modelo de Regresión con Random Forest: Predicción de la Profundidad Vertical Verdadera

Jennifer Cordones

2026-07-28