#### UNIVERSIDAD CENTRAL DEL ECUADOR ######
#### AUTOR: JHONNY JOSÉ CATUCUAMBA CALCAN ####
#### CARRERA: INGENIERÍA AMBIENTAL ####
### VARIABLES: ####
### MODELO DE REGRESIÓN POLINÓMICO #####
setwd("/cloud/project")
datos <- read.csv("Conaf_Data_Chile_2017.csv")
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
###################### ESTADISTICA Multivariable ###############################
plot(1, type = "n", axes = FALSE, xlab = "", ylab = "") # Crear un gráfico vacío
text(x = 1, y = 1,
labels = "ESTADÍSTICA MULTIVARIABLE",
cex = 2, # Tamaño del texto (ajustable)
col = "blue", # Color del texto
font =6) #tipo
total_veg<- datos$total_veg
total_veg<- as.numeric(gsub(",", ".", total_veg))
sup_a <- datos$sup_t_a
sup_a <- as.numeric(gsub(",", ".", sup_a))
length(total_veg); length(sup_a)
## [1] 5234
## [1] 5234
tabla <- data.frame(total_veg,sup_a)
# Depuramos la tabla
tabla_limpia <- tabla[!apply(tabla, 1, function(fila) {
any(is.na(fila) | fila == 0 | fila == "")
}), ]
length(tabla_limpia$total_veg); length(tabla_limpia$sup_a)
## [1] 4549
## [1] 4549
######Regresión polinómica
# Filtrar según los mismos límites usados en el plot
tabla_filtrada <- tabla_limpia %>%
filter(total_veg >= 1000, total_veg <= 40000,
sup_a >= 0, sup_a <= 160000)
x <- tabla_filtrada$total_veg
y <- tabla_filtrada$sup_a
plot(total_veg,sup_a, ylim = c(0,160000),xlim = c(1000,40000),
col = 4, main ="Gráfica Nº 10.1: Diagrama de dispersión de Superficie Afectada
en función de la Vegetación Total",
xlab ="Vegetación Total (ha)",
ylab = "Sperficie Afectada Total (ha)")
####parametros
xcuad <- x^2
xcub <- x^3
xcuar <- x^4
xquin <- x^5
regrespoli <- lm(y~x+xcuad+xcub+xcuar+xquin)
regrespoli
##
## Call:
## lm(formula = y ~ x + xcuad + xcub + xcuar + xquin)
##
## Coefficients:
## (Intercept) x xcuad xcub xcuar xquin
## -7.128e+03 1.012e+01 -2.861e-03 3.631e-07 -1.668e-11 2.415e-16
summary(regrespoli)
##
## Call:
## lm(formula = y ~ x + xcuad + xcub + xcuar + xquin)
##
## Residuals:
## Min 1Q Median 3Q Max
## -7970.5 -2137.5 2.6 1794.0 9132.5
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -7.128e+03 3.375e+03 -2.112 0.040838 *
## x 1.012e+01 3.346e+00 3.025 0.004282 **
## xcuad -2.861e-03 9.731e-04 -2.940 0.005373 **
## xcub 3.631e-07 1.081e-07 3.357 0.001708 **
## xcuar -1.668e-11 4.651e-12 -3.586 0.000885 ***
## xquin 2.415e-16 6.473e-17 3.732 0.000578 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3412 on 41 degrees of freedom
## Multiple R-squared: 0.9821, Adjusted R-squared: 0.9799
## F-statistic: 449.5 on 5 and 41 DF, p-value: < 2.2e-16
beta0 <- regrespoli$coefficients[1]
beta1 <- regrespoli$coefficients[2]
beta2 <- regrespoli$coefficients[3]
beta3 <- regrespoli$coefficients[4]
beta4 <- regrespoli$coefficients[5]
beta5 <- regrespoli$coefficients[6]
a <- beta0
b <- beta1
c <-beta2
d <-beta3
e <-beta4
f <-beta5
a
## (Intercept)
## -7127.852
b
## x
## 10.12053
c
## xcuad
## -0.002861023
d
## xcub
## 3.630871e-07
e
## xcuar
## -1.667924e-11
f
## xquin
## 2.415449e-16
plot(x,y,
col = "deepskyblue",
main = "Gráfica Nº 10.2: Diagrama de dispersión con la línea de regresión",
xlab = "Vegetación total (ha)",
ylab = "Superficie total afectada (ha)",
pch = 16)
curve(a+(b*x)+(c*x^2)+(d*x^3)+(e*x^4)+(f*x^5),add = TRUE,col="blue")
# Crear el vector de coeficientes (orden: x^5 a x^0)
coef_pol <- c(f, e, d, c, b, a)
# Calcular raíces del polinomio
raices <- polyroot(coef_pol)
# Extraer raíces reales (Imaginaria ≈ 0)
raices_reales <- Re(raices[abs(Im(raices)) < 1e-6])
# Verlas
raices_reales
## [1] 2.923804e-05 4.193899e-05 1.093671e-03
# Test
r <- cor(x + xcuad + xcub + xcuar + xquin,y)
r
## [1] 0.9415804
### Restricciones
x_vals <- seq(4000, 40000, length.out = 1000)
y_vals <- a + b*x_vals + c*x_vals^2 + d*x_vals^3 + e*x_vals^4 + f*x_vals^5
# Filtrar donde y >= 0
x_validos <- x_vals[y_vals >= 0]
y_validos <- y_vals[y_vals >= 0]
# Graficar solo la parte válida del modelo
plot(x, y, col = "gray", pch = 16, main = "Gráfica Nº 10.3: Modelo polinómico restringido",
xlab = "Vegetación total (ha)",
ylab = "Superficie total afectada (ha)")
lines(x_validos, y_validos, col = "blue", lwd = 2)
abline(h = 0, col = "red", lty = 2)
# Estimaciones
sup_a10000 <- -7127.852 + 10.12053*10000 - 0.002861023*10000^2 +
3.630871e-07*10000^3 - 1.667924e-11*10000^4 + 2.415449e-16*10000^5
sup_a10000
## [1] 28424.34
plot(1, type = "n", axes = FALSE, xlab = "", ylab = "") # Crear un gráfico vacío
text(x = 1, y = 1,
labels = "¿Que cantidad de Superficie total
afectada se espera cuando se
tenga 10000(ha) de Vegetación total?
\n R= 28424.34 (ha)",
cex = 2, # Tamaño del texto (ajustable)
col = "blue", # Color del texto
font = 6)
# CONCLUSIONES
plot(1, type = "n", axes = FALSE, xlab = "", ylab = "") # Crear un gráfico vacío
text(x = 1, y = 1,
labels = "CONCLUSIONES",
cex = 2, # Tamaño del texto (ajustable)
col = "blue", # Color del texto
font = 6)
library(knitr)
Tabla_resumen <- data.frame(Variables = c("x", "y"),
Modelo = c("Vegetación total (ha)", "Superficie total afectada (ha)"),
Variab = c("y=-7127.85 + 10.12*x - 0.0028*x^2 + 3.6e-07*x^3 - 1.66e-11*x^4 + 2.41e-16*x^5", ""),
Restricciones = c(" x >= 4000 (ha)","y > 0 (ha)"),
Coeficiente_pearson = c("0.94 ", ""))
colnames(Tabla_resumen) <- c("Variables ", "Nombres ", " Ecuación ", "Restricciones",
"Coef. Pearson" )
kable(Tabla_resumen, align = 'c', caption = "Resumen del Modelo de Regresión Simple Polinómico")
Resumen del Modelo de Regresión Simple Polinómico
| x |
Vegetación total (ha) |
y=-7127.85 + 10.12x - 0.0028x^2 +
3.6e-07x^3 - 1.66e-11x^4 + 2.41e-16*x^5 |
x >= 4000 (ha) |
0.94 |
| y |
Superficie total afectada (ha) |
|
y > 0 (ha) |
|
###### Entre el total veg y la superficie total de los incendios de Chile en 2017 existe una relación polinónica
###### representada por y=-7127.85 + 10.12*x - 0.0028*x^2 + 3.6e-07*x^3 - 1.66e-11*x^4 + 2.41e-16*x^5
##### siendo X el total veg en (ha) Y la superficie total (ha) y aprueban el test de pearson con un 94% y con restricciones " x >= 4000 (ha)","y > 0 (ha)")