Relación Potencial: relación entre Liberación involuntaria (barriles) y Costos de pérdida de mercancía
Brandon
1 Librerías
El primer paso para el desarrollo del modelo es preparar nuestro entorno de trabajo cargando los paquetes necesarios para la manipulación de datos, generación de tablas formales y renderizado del documento. Evitaremos la notación científica para facilitar la lectura de los ejes financieros.
library(readr)
library(rmarkdown)
library(knitr)
library(kableExtra)2 Carga de datos
Realizamos la importación del conjunto de datos en bruto hacia nuestro entorno en R. Este paso nos permite acceder al historial documentado de incidentes en la infraestructura, estableciendo la base empírica de nuestro análisis estadístico.
library(readr)
datos <- read_csv("database-_1_.csv")## Warning: One or more parsing issues, call `problems()` on your data frame for details,
## e.g.:
## dat <- vroom(...)
## problems(dat)## Rows: 2795 Columns: 36
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (18): Accident Date/Time, Operator Name, Pipeline/Facility Name, Pipelin...
## dbl (18): Report Number, Supplemental Number, Accident Year, Operator ID, Ac...
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.3 Selección de variables
Para comprender el verdadero impacto de los incidentes en la infraestructura, este análisis se centra en relacionar la magnitud física de un accidente con su consecuencia financiera directa.
Justificación Causa-Efecto: Definimos los “Barriles Liberados” como nuestra variable independiente o causa (X), y al “Costo del Producto Perdido” como nuestra variable dependiente o efecto (Y). Asumimos lógicamente que el volumen físico derramado al medio ambiente es el principal detonante del impacto financiero tangible.
df_bruto <- data.frame(
X = datos$`Unintentional Release (Barrels)`,
Y = datos$`Lost Commodity Costs`
)4 Tabla de pares de valores
Verificamos el tamaño de la muestra sin procesar y presentamos una vista estructurada de los datos resultantes antes de la limpieza, indicando el tamaño muestral empírico original.
# Preparamos la tabla visual forzando los nombres de columnas y la numeración
df_visual <- df_bruto
colnames(df_visual) <- c("Barriles Liberados (X)", "Costos de pérdida de mercancía (Y)")
rownames(df_visual) <- 1:nrow(df_visual)
cat("El tamaño muestral inicial es de:", nrow(df_bruto), "observaciones.\n\n")## El tamaño muestral inicial es de: 2795 observaciones.cat("--- PRIMERAS FILAS DE LA TABLA DE PARES ---\n")## --- PRIMERAS FILAS DE LA TABLA DE PARES ---head(df_visual)## Barriles Liberados (X) Costos de pérdida de mercancía (Y)
## 1 21.00 1517
## 2 0.12 8
## 3 2.00 200
## 4 0.48 40
## 5 700.00 150
## 6 3784.00 167775paged_table(df_visual)5 Gráfica de dispersión
Construimos un diagrama de dispersión inicial (nube de puntos) para realizar un análisis exploratorio visual y observar de manera empírica cómo se distribuyen los datos crudos en el plano cartesiano.
plot(df_bruto$X, df_bruto$Y,
main="Gráfica N° 1: Dispersión Inicial (Datos Crudos)",
xlab="Barriles Liberados (X)",
ylab="Costo Producto Perdido ($)",
pch=16, col=rgb(0.1, 0.4, 0.8, 0.5))
6 Conjetura del modelo
La gráfica inicial presenta un agrupamiento denso y ruido provocado por valores nulos y ceros, lo que impide ver la tendencia real del costo. Por lo tanto, debemos proceder con un refinamiento crítico.
6.1 Tratamiento de los datos
Aplicamos un proceso riguroso de depuración que incluye la eliminación de valores nulos (NA), el descarte de registros en cero y la aplicación de filtros para excluir valores atípicos extremos.
# 1. OMITIMOA NAs de la base original
df_pot <- na.omit(df_bruto)
# 2. Rango operativo
df_pot <- df_pot[df_pot$X > 50 & df_pot$X < 5000, ]
df_pot <- df_pot[df_pot$Y > (df_pot$X * 10), ]
df_pot <- df_pot[!(df_pot$X < 1500 & df_pot$Y > 500000), ]
df_pot <- df_pot[!(df_pot$X < 1000 & df_pot$Y > 150000), ]
cat("Total de incidentes empíricos válidos tras la limpieza:", nrow(df_pot), "observaciones.\n")## Total de incidentes empíricos válidos tras la limpieza: 320 observaciones.6.2 Nueva gráfica de dispersión
Graficamos nuevamente los datos depurados.
plot(df_pot$X, df_pot$Y,
main="Gráfica N° 2: Nube de Puntos Tratada",
xlab="Barriles Liberados",
ylab="Costo Producto Perdido ($)",
pch=16, col=rgb(0.1, 0.4, 0.8, 0.5))
6.3 Nueva conjetura
A partir de la observación de la gráfica tratada, notamos que los costos financieros tienden a escalar de forma multiplicativa conforme aumenta el volumen. Justificamos la elección de un modelo potencial, asumiendo que la elasticidad del costo respecto al volumen derramado se mantiene proporcional.Se conjetura un modelo potencial, cuya ecuación matemática base es:\[Y = a \cdot X^b\]
7 Cálculo de parámetros
Procedemos a calcular los estimadores matemáticos. Para un modelo potencial, transformamos la ecuación a una forma lineal aplicando logaritmos a ambas variables (\(\ln(Y) = \ln(a) + b \cdot \ln(X)\)) y extraemos los coeficientes.
# Ojo: Aplicamos logaritmo a Y y a X
modelo_potencial <- lm(log(Y) ~ log(X), data = df_pot)
a <- exp(coef(modelo_potencial)[1]) # Intersección (Parámetro a)
b <- coef(modelo_potencial)[2] # Exponente (Parámetro b)
cat("Parámetro a (Multiplicador):", a, "\n")## Parámetro a (Multiplicador): 61.53506cat("Parámetro b (Exponente):", b, "\n")## Parámetro b (Exponente): 0.9781634cat("\nLa ecuación del modelo potencial es: Y =", round(a, 4), "* X ^", round(b, 6), "\n")##
## La ecuación del modelo potencial es: Y = 61.5351 * X ^ 0.9781638 Realidad y modelo
Para validar nuestro cálculo, superponemos la curva teórica generada por nuestra ecuación matemática sobre la nube de puntos reales.
plot(df_pot$X, df_pot$Y,
main="Gráfica N° 3: Relación entre la Realidad y Modelo Potencial",
xlab="Barriles Liberados",
ylab="Costo Producto Perdido ($)",
pch=16, col=rgb(0.1, 0.4, 0.8, 0.5))
# Creamos la línea que cruza la nube de puntos en formato de línea negra continua (lty=1)
x_seq <- seq(min(df_pot$X), max(df_pot$X), length.out=1000)
# Fórmula de la línea de predicción potencial: Y = a * X^b
y_pred <- a * (x_seq ^ b)
lines(x_seq, y_pred, col="red", lty=1, lwd=3)
9 Test de Pearson y Ajuste
Sometemos nuestro modelo a una evaluación estadística rigurosa cuantificando la fuerza y la dirección de la relación con el coeficiente de Pearson y el Coeficiente de Determinación (\(R^2\)).
# El coeficiente se calcula sobre ambas relaciones transformadas (log-log)
pearson_val <- cor(log(df_pot$X), log(df_pot$Y))
r2_val <- summary(modelo_potencial)$r.squared
cat("Coeficiente de correlación de Pearson (R):", round(pearson_val, 4), "\n")## Coeficiente de correlación de Pearson (R): 0.8634cat("Coeficiente de Determinación (R2):", round(r2_val, 4), "\n")## Coeficiente de Determinación (R2): 0.7454ecuacion_str <- paste0("y = ", round(a, 4), " * x^", round(b, 6))
tabla_resumen <- data.frame(
Variable = c("Barriles Liberados", "Costo del Producto"),
Tipo = c("Independiente (x)", "Dependiente (y)"),
R = c("", as.character(round(pearson_val, 4))),
R2 = c("", as.character(round(r2_val, 4))),
Parametro_a = c("", as.character(round(a, 4))),
Parametro_b = c("", as.character(round(b, 6))),
Ecuacion = c("", ecuacion_str)
)
colnames(tabla_resumen) <- c("Variable", "Tipo", "R", "R2", "Parámetro a", "Parámetro b", "Ecuación")
tabla_resumen %>%
kable(format = "html", caption = "<b>Tabla N°1: Resumen del Modelo Potencial</b>", align = "c", escape = FALSE) %>%
kable_styling(bootstrap_options = c("hover", "condensed"), full_width = FALSE, position = "center") %>%
row_spec(0, bold = TRUE, background = "#F2F2F2", color = "#333333") %>%
footnote(general = "Autor: Brandon", footnote_as_chunk = TRUE)| Variable | Tipo | R | R2 | Parámetro a | Parámetro b | Ecuación |
|---|---|---|---|---|---|---|
| Barriles Liberados | Independiente (x) | |||||
| Costo del Producto | Dependiente (y) | 0.8634 | 0.7454 | 61.5351 | 0.978163 | y = 61.5351 * x^0.978163 |
| Note: Autor: Brandon |
10 Restricciones
Restricciones El uso de este modelo matemático está sujeto a las siguientes restricciones:
Rango Operativo: El modelo pierde confiabilidad estadística si se utiliza para extrapolar incidentes menores a 50 barriles o superiores a 5000 barriles.
Dominio Natural: Debido al uso de logaritmos en su formulación, el modelo no admite el valor cero en ninguna de las dos variables. Un derrame de 0 barriles no puede calcularse directamente mediante esta formulación logarítmica.
11 Estimación
La verdadera utilidad de un modelo de regresión radica en su capacidad predictiva. En esta sección, pondremos a prueba nuestro modelo utilizando la ecuación matemática previamente definida (Y=a⋅Xb ) para realizar estimaciones. Calcularemos el valor esperado de la variable de respuesta (Y ) dados nuevos valores hipotéticos de la variable predictora (X ) que no necesariamente estaban incluidos en nuestra muestra original.
Ejemplo:
Supongamos un incidente donde se liberan 2,500 barriles, ¿cuál es el costo estimado del producto perdido?
val_test <- 2000
# Predicción potencial: Y = a * X^b
prediccion <- a * (val_test ^ b)
cat("\n--- ESTIMACIÓN ---\n")##
## --- ESTIMACIÓN ---cat("Para un derrame de 2,000 barriles, el valor estimado es: $", format(round(prediccion, 2), big.mark=","), "\n")## Para un derrame de 2,000 barriles, el valor estimado es: $ 104,248.412 Conclusión
Entre los Barriles Liberados y el Costo del Producto Perdido (USD) existe una relación directa de tipo potencial, explicada por un coeficiente de determinación R² \(\approx\) 74.54%. Esto indica que el 74.54% de la variabilidad del Costo del Producto Perdido (USD) puede ser explicada por la variación en los Barriles Liberados, lo que representa una relación estadísticamente significativa y de magnitud considerable. La ecuación matemática estimada del modelo es:\[y = 61.54 \cdot x^{0.98}\]Dado que el exponente b = 0.98 es positivo, la relación es creciente; y al ser un valor cercano a 1, el efecto de los Barriles Liberados sobre el Costo del Producto Perdido (USD) es fuerte y casi proporcional, indicando que por cada 1% que aumentan los Barriles Liberados, el costo se incrementa en aproximadamente un 0.98%.