Visualización de los datos
poly(papas$Superficie_Cosechada, degree = 4, raw=TRUE,simple=TRUE)[1:5]
## [1] 3570 3332 4171 4178 3102
poly(papas$Superficie_Cosechada, degree = 4, raw=FALSE, simple=TRUE)[1:5]
## [1] -0.2285076 -0.2423422 -0.1935722 -0.1931653 -0.2557118
# CÁLCULO DEL MODELO POLINÓMICO DE GRADO 4
# ----------------------------------------
modelo_poli4 <- lm(Valor_de_la_Produccion ~ poly(Superficie_Cosechada, 4), data = papas)
summary(modelo_poli4)
##
## Call:
## lm(formula = Valor_de_la_Produccion ~ poly(Superficie_Cosechada,
## 4), data = papas)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1026361 -273080 -41158 176245 1087165
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1109749 95536 11.616 2.41e-10 ***
## poly(Superficie_Cosechada, 4)1 4739165 477681 9.921 3.61e-09 ***
## poly(Superficie_Cosechada, 4)2 346125 477681 0.725 0.477
## poly(Superficie_Cosechada, 4)3 -868556 477681 -1.818 0.084 .
## poly(Superficie_Cosechada, 4)4 -196433 477681 -0.411 0.685
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 477700 on 20 degrees of freedom
## Multiple R-squared: 0.8366, Adjusted R-squared: 0.804
## F-statistic: 25.61 on 4 and 20 DF, p-value: 1.268e-07
# INTERPOLACIÓN DE PUNTOS DENTRO DEL RANGO DEL PREDICTOR
# -----------------------------------------------------------------------------
limites <- range(papas$Superficie_Cosechada)
nuevos_puntos <- seq(from = limites[1], to = limites[2], by = 1)
nuevos_puntos <- data.frame(Superficie_Cosechada = nuevos_puntos)
# PREDICCIÓN DE LA VARIABLE RESPUESTA Y DEL ERROR ESTÁNDAR
# -----------------------------------------------------------------------------
predicciones <- predict(modelo_poli4, newdata = nuevos_puntos, se.fit = TRUE,
level = 0.95)
# CÁLCULO DEL INTERVALO DE CONFIANZA SUPERIOR E INFERIOR 95%
# -----------------------------------------------------------------------------
intervalo_conf <- data.frame(inferior = predicciones$fit -
1.96*predicciones$se.fit,
superior = predicciones$fit +
1.96*predicciones$se.fit)
plot(x = papas$Superficie_Cosechada, y = papas$Valor_de_la_Produccion , pch = 20, col = "darkgrey")
title("Polinomio de grado 4: Superficie_Cosechada ~ Valor_de_la_Produccion")
points(x = nuevos_puntos$Superficie_Cosechada, predicciones$fit, col = "red", pch = 10)
points(x = nuevos_puntos$Superficie_Cosechada, intervalo_conf$inferior, col = "blue", pch = 4)
points(x = nuevos_puntos$Superficie_Cosechada, intervalo_conf$superior, col = "blue", pch = 4)
## The following objects are masked _by_ .GlobalEnv:
##
## lluvias, temperaturas
plot(x = Superficie_Cosechada, y = Valor_de_la_Produccion, pch = 20, col = "darkgrey")
title("Polinomio de grado 4: Superficie_Cosechada ~ Valor_de_la_Produccion")
lines(x = nuevos_puntos$Superficie_Cosechada, predicciones$fit, col = "red", lwd = 2)
lines(x = nuevos_puntos$Superficie_Cosechada, intervalo_conf$inferior, col = "blue", lwd = 2)
lines(x = nuevos_puntos$Superficie_Cosechada, intervalo_conf$superior, col = "blue", lwd = 2)
ggplot(data = papas, aes(x = Superficie_Cosechada, y = Valor_de_la_Produccion)) +
geom_point(color = "grey30", alpha = 0.3) +
geom_smooth(method = "lm", formula = y ~ poly(x, 4), color = "red") +
labs(title = "Polinomio de grado 4: Valor_de_la_Produccion ~ Superficie_Cosechada") +
theme_bw() +
theme(plot.title = element_text(hjust = 0.5))
# Se ajustan modelos polinómicos de grado 1 a 5
# ----------------------------------------------------------------------------
modelo_1 <- lm(Valor_de_la_Produccion ~ Superficie_Cosechada, data = papas)
modelo_2 <- lm(Valor_de_la_Produccion ~ poly(Superficie_Cosechada, 2), data = papas)
modelo_3 <- lm(Valor_de_la_Produccion ~ poly(Superficie_Cosechada, 3), data = papas)
modelo_4 <- lm(Valor_de_la_Produccion ~ poly(Superficie_Cosechada, 4), data = papas)
modelo_5 <- lm(Valor_de_la_Produccion ~ poly(Superficie_Cosechada, 5), data = papas)
anova(modelo_1, modelo_2, modelo_3, modelo_4, modelo_5)
## Analysis of Variance Table
##
## Model 1: Valor_de_la_Produccion ~ Superficie_Cosechada
## Model 2: Valor_de_la_Produccion ~ poly(Superficie_Cosechada, 2)
## Model 3: Valor_de_la_Produccion ~ poly(Superficie_Cosechada, 3)
## Model 4: Valor_de_la_Produccion ~ poly(Superficie_Cosechada, 4)
## Model 5: Valor_de_la_Produccion ~ poly(Superficie_Cosechada, 5)
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 23 5.4764e+12
## 2 22 5.3566e+12 1 1.1980e+11 0.5129 0.48262
## 3 21 4.6022e+12 1 7.5439e+11 3.2294 0.08824 .
## 4 20 4.5636e+12 1 3.8586e+10 0.1652 0.68897
## 5 19 4.4384e+12 1 1.2519e+11 0.5359 0.47307
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
cv_MSE_k10 <- rep(NA,10)
for (i in 1:10) {
modelo <- glm(Valor_de_la_Produccion ~ poly(Superficie_Cosechada, i), data = papas)
set.seed(17)
cv_MSE_k10[i] <- cv.glm(data = papas, glmfit = modelo, K = 10)$delta[1]
}
p4 <- ggplot(data = data.frame(polinomio = 1:10, cv_MSE = cv_MSE_k10),
aes(x = polinomio, y = cv_MSE)) +
geom_point(colour = c("firebrick3")) +
geom_path()
p4 <- p4 + theme(panel.grid.major = element_line(colour = 'gray90'))
p4 <- p4 + theme(plot.title = element_text(face = 'bold'))
p4 <- p4 + theme(panel.background = element_rect(fill = 'gray98'))
p4 <- p4 + labs(title = 'Test Error ~ Grado del polinomio')
p4 <- p4 + scale_x_continuous(breaks = 1:10)
p4
