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)

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;',
    'Tabla 1: Visualización interactiva de las primeras observaciones del dataset original'
  ),
  options = list(
    pageLength = 5,
    autoWidth = TRUE,
    scrollX = TRUE,
    language = list(url = '//cdn.datatables.net/plug-ins/1.10.25/i18n/Spanish.json')
  ),
  class = 'display compact cell-border'
)

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.

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

7 PREDICCIONES

Se ejecutan las predicciones sobre el conjunto de prueba comparando los valores reales frente a los estimados.

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

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

9 VISUALIZACIÓN DE RESULTADOS

9.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")

9.2 Observados y Predichos

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 = "Real (ft)", y = "Predicho (ft)")

9.3 Importancia de variables

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) +
  coord_flip() +
  theme_minimal() +
  labs(title = "Top 10 Variables de Importancia", x = "Variables", y = "Importancia")

9.4 Residuos del modelo

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", y = "Densidad")

9.5 QQ Plot de los residuos

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")

9.6 Comparación de predicciones

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"))

10 CONCLUSIÓN

El modelo de regresión basado en Random Forest demostró ser una herramienta altamente confiable para predecir la profundidad vertical verdadera de los pozos petroleros. A partir del análisis de los datos históricos, se obtuvieron las siguientes conclusiones principales:

Factores determinantes: Se identificó que la profundidad propuesta y las características litoestratigráficas son las variables con mayor peso e influencia dentro del algoritmo, dictando en gran medida el comportamiento de la profundidad final.
Precisión y ajuste: Las evaluaciones de rendimiento mostraron un coeficiente de determinación (\(R^2\)) sobresaliente y una fuerte correlación de Pearson, lo que confirma que las estimaciones del modelo se alinean de manera muy cercana con los valores reales registrados en el campo.
Margen de error controlado: El análisis de los residuos y las métricas de error (como el RMSE y el MAE) evidenciaron desviaciones mínimas y una distribución normal centrada en cero, ratificando la estabilidad operativa del modelo frente a la variabilidad geológica de los pozos.

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-24