Modelo de Regresión Polinomial
CARGAR LIBRERÍAS Y DATOS
library(dplyr)## Warning: package 'dplyr' was built under R version 4.5.2##
## Adjuntando el paquete: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionsetwd("D:/Data")
datos <- read.csv("derrames_globales_.csv", header = TRUE, sep = ";", dec =".")1 Preparación de variables
x_raw <- as.numeric(as.character(datos[, 20])) # Volumen Derramado (Galones)## Warning: NAs introducidos por coercióny_raw <- as.numeric(as.character(datos[, 17])) # Máxima Liberación (Galones)
# Crear dataframe inicial y omitir nulos
datos_raw <- data.frame(x = x_raw, y = y_raw)
datos_clean <- na.omit(datos_raw)
# Filtramos valores mayores a cero
datos_finales <- subset(datos_clean, x > 0 & y > 0)2 Gráfica 1: Relación entre Volumen derramado y Máxima Liberación
plot(datos_finales$x, datos_finales$y,
col = rgb(0, 0, 1, 0.5),
pch = 16, cex = 0.8,
main = "Relación entre Volumen derramado y Máxima Liberación",
xlab = "Volumen derramado ",
ylab = "Máximo Liberación ")
abline(v=0, h=0, col="gray", lty=2) 
3 Conjetura de modelo Polinomial
\[ Y = \beta_0 + \beta_1 X + \beta_2 X^2 + \beta_3 X^3 \]
# Generar Modelo Polinómico de grado 3
# raw = TRUE usa los polinomios naturales, necesario para ver la ecuación clásica
modelo_poly <- lm(y ~ poly(x, 3, raw = TRUE), data = datos_finales)
# Resumen estadístico
res <- summary(modelo_poly)
res##
## Call:
## lm(formula = y ~ poly(x, 3, raw = TRUE), data = datos_finales)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2000359 -78807 -75032 -74062 29896272
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 7.431e+04 2.787e+04 2.666 0.007764 **
## poly(x, 3, raw = TRUE)1 1.900e+00 2.531e-01 7.504 1.12e-13 ***
## poly(x, 3, raw = TRUE)2 -1.118e-07 3.251e-08 -3.437 0.000606 ***
## poly(x, 3, raw = TRUE)3 2.158e-15 6.395e-16 3.375 0.000759 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1006000 on 1342 degrees of freedom
## Multiple R-squared: 0.6001, Adjusted R-squared: 0.5992
## F-statistic: 671.3 on 3 and 1342 DF, p-value: < 2.2e-16# Extracción de Coeficientes
coef_poli <- coef(modelo_poly)
a <- coef_poli[1] # Intercepto (beta 0)
b <- coef_poli[2] # Coeficiente x (beta 1)
c <- coef_poli[3] # Coeficiente x^2 (beta 2)
d <- coef_poli[4] # Coeficiente x^3 (beta 3)4 Gráfica 2: Relación entre Volumen derramado y Máxima Liberación (Modelo)
# Secuencia para la curva suave
x_seq <- seq(min(datos_finales$x), max(datos_finales$x), length.out = 500)
# Predicción: a + bx + cx^2 + dx^3
y_pred <- a + b*x_seq + c*x_seq^2 + d*x_seq^3
plot(datos_finales$x, datos_finales$y,
col = rgb(0, 0, 1, 0.5),
pch = 16, cex = 0.8,
main = "Relación entre Volumen derramado y Máxima Liberación",
xlab = "Volumen derramado ",
ylab = "Máximo Liberación ")
# Curva de regresión
lines(x_seq, y_pred, col = "red", lwd = 2)
5 Test de Pearson
pearson_r <- cor(datos_finales$x, datos_finales$y, method = "pearson")
cat("=== RESULTADO TEST DE PEARSON ===\n")## === RESULTADO TEST DE PEARSON ===cat("Coeficiente de Pearson (r):", round(pearson_r, 4), "\n")## Coeficiente de Pearson (r): 0.7722if(abs(pearson_r) > 0.7) {
cat("INTERPRETACIÓN: Correlación Lineal Fuerte.\n")
} else if(abs(pearson_r) > 0.4) {
cat("INTERPRETACIÓN: Correlación Lineal Moderada.\n")
} else {
cat("INTERPRETACIÓN: Correlación Lineal Débil.\n")
}## INTERPRETACIÓN: Correlación Lineal Fuerte.6 Relación en el intervalo [30100 - 35100]
# Filtramos estrictamente el intervalo óptimo
datos_intervalo <- subset(datos_finales, x >= 30100 & x <= 35100)
# Modelo Polinomial Grado 3 para el intervalo
modelo_int <- lm(y ~ poly(x, 3, raw = TRUE), data = datos_intervalo)
summary(modelo_int)##
## Call:
## lm(formula = y ~ poly(x, 3, raw = TRUE), data = datos_intervalo)
##
## Residuals:
## Min 1Q Median 3Q Max
## -8131.0 -2316.0 -791.4 -233.6 14627.9
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4.661e+08 4.200e+07 -11.10 1.07e-05 ***
## poly(x, 3, raw = TRUE)1 4.231e+04 3.824e+03 11.07 1.09e-05 ***
## poly(x, 3, raw = TRUE)2 -1.279e+00 1.159e-01 -11.04 1.11e-05 ***
## poly(x, 3, raw = TRUE)3 1.286e-05 1.169e-06 11.00 1.14e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 6681 on 7 degrees of freedom
## Multiple R-squared: 0.9522, Adjusted R-squared: 0.9318
## F-statistic: 46.53 on 3 and 7 DF, p-value: 5.436e-05# Coeficientes del intervalo
coef_int <- coef(modelo_int)
a_int <- coef_int[1] # beta 0
b_int <- coef_int[2] # beta 1
c_int <- coef_int[3] # beta 2
d_int <- coef_int[4] # beta 37 Gráfica 3: Relación entre Volumen derramado y Máxima Liberación (intervalo)
x_seq_int <- seq(min(datos_intervalo$x), max(datos_intervalo$x), length.out = 500)
y_pred_int <- a_int + b_int*x_seq_int + c_int*x_seq_int^2 + d_int*x_seq_int^3
plot(datos_intervalo$x, datos_intervalo$y,
col = rgb(0.8, 0.4, 0, 0.6), pch = 16, cex = 1.2,
main = "Modelo Polinomial Grado 3 [30100 - 35100] gal",
xlab = "Volumen Derramado (Galones)",
ylab = "Máxima Liberación (Galones)")
lines(x_seq_int, y_pred_int, col = "red", lwd = 3)
legend("bottomright", legend = "Modelo Intervalo", col = "red", lwd = 3, bty = "n")
8 Test de Pearson (Intervalo)
# Para polinomios, extraemos el R de la raíz del R cuadrado del modelo
r_poli_int <- sqrt(summary(modelo_int)$r.squared)
cat("=== COMPARACIÓN DE AJUSTE (R) ===\n")## === COMPARACIÓN DE AJUSTE (R) ===cat("Pearson Global:", round(pearson_r, 4), "\n")## Pearson Global: 0.7722cat("Ajuste Intervalo [30.1k - 35.1k]:", round(r_poli_int, 4), "\n\n")## Ajuste Intervalo [30.1k - 35.1k]: 0.9758if(r_poli_int > abs(pearson_r)) {
cat("CONCLUSIÓN: El ajuste MEJORA drásticamente en este intervalo (R ≈ 0.9758).\n")
} else {
cat("CONCLUSIÓN: El ajuste es débil en este intervalo.\n")
}## CONCLUSIÓN: El ajuste MEJORA drásticamente en este intervalo (R ≈ 0.9758).9 Ecuación del modelo polinomial de tercer orden
\[ Y = -466000000 + 42300 X - 1.28 X^2 + 0.0000129 X^3 \]
cat("\n=== ECUACIÓN FINAL (INTERVALO) ===\n")##
## === ECUACIÓN FINAL (INTERVALO) === cat("y =", format(a_int, scientific = TRUE, digits = 3),
"+", format(b_int, scientific = TRUE, digits = 3), "* x",
"+", format(c_int, scientific = TRUE, digits = 3), "* x^2",
"+", format(d_int, scientific = TRUE, digits = 3), "* x^3\n")## y = -4.66e+08 + 4.23e+04 * x + -1.28e+00 * x^2 + 1.29e-05 * x^310 Estimación de un punto
x_nueva <- 33000
y_pred_estimada <- a_int + b_int*x_nueva + c_int*x_nueva^2 + d_int*x_nueva^3
cat("\n=== PREDICCIÓN ===\n")##
## === PREDICCIÓN ===cat("Para un volumen derramado de", x_nueva, "galones:\n")## Para un volumen derramado de 33000 galones:cat("La máxima liberación estimada es =", round(y_pred_estimada, 2), "galones\n")## La máxima liberación estimada es = 69822.69 galonesConclusiones
Entre la variable independiente máxima liberación (X) y la variable dependiente respuesta actual (Y) existe una relación matemática de tipo regresión polinomial de tercer grado Esta relación se expresa mediante la fórmula del modelo:\(Y = -466000000 + 42300 X - 1.28 X^2 + 0.0000129 X^3\). El modelo permite realizar una estimación para un derrame de 33000, se estima una máxima liberación de 69822.69 galones.