#Configuración

{r setup, include=FALSE}

knitr::opts_chunk$set(echo = TRUE, eval = FALSE)

require(pacman)
pacman::p_load(dplyr, readr, knitr, ggplot2, corrplot, flextable, readxl,
GGally, lmtest, car, stargazer, rsm, PerformanceAnalytics, openxlsx)

Producción de madera aserrada en un período de 10 años (1990-1999) expresada en miles de m³.

df <- as.data.frame(read.xlsx("E:/[FAT32]/ACTIVIDAD DE TITULACION I (2023)/CLASE 10 REGRESION LINEAL/madera.xlsx"))

#Mostrar tabla con datos de las variables

flextable(df)

#Estructura del set de datos

str(df)

#Fijar nombre de variables en la memoria

attach(df)

#Resumen del set de datos (datos numéricos de las variables)

summary(df)

#Gráfico de dispersión del conjunto de variables

ggplot(data=df, mapping= aes(x= Años, y= Madera)) +
  geom_point(color = "blue", size= 2) +
  labs(title= "Diagrama de Dispersión", x= "pH", y= "Biomasa") +
  theme_bw() + 
  theme(plot.title = element_text(hjust = 0.5))

#Gráfico de correlaciones

chart.Correlation(df)
cor(df)  

#Test de correlación (“pearson”, “kendall”, “spearman”)

cor.test(x = Años, y = Madera, method = "pearson") #cuantitativas con comportamiento normal
cor.test(x = Años, y = Madera, method = "kendall") #variables cualitativas de tipo ordinal
cor.test(x = Años, y = Madera, method = "spearman") #cuantitativas no siguen un comportamiento normal

#Otra forma de aplicar la Correlación de Pearson entre 2 variables

cor(Años, Madera)
cor(Madera, Años)

#Modelo de regresión lineal simple

modelo.lineal <- lm(Madera ~ Años, data = df)
summary(modelo.lineal)

Nota:

Modelo 1: \[\hat{y}=38.681{x}-76560.279\]

#Contraste de hipótesis (normalidad de los residuos)

shapiro.test(modelo.lineal$residuals)

Comentarios:

#Test de Breush-Pagan (homocedasticidad de los residuos)

bptest(modelo.lineal)

Comentarios:

#Prueba de Durbin – Watson Usando libreria “lmtest”

• \(H_0\): No existe correlación entre los residuales

• \(H_1\): Los resisuales están autocorrelacionados

dwtest(modelo.lineal,alternative ="two.sided",iterations = 1000)

Comentarios:

#Residuos estudentizados

df$studentized_residual <- rstudent(modelo.lineal)
ggplot(data = df, aes(x = predict(modelo.lineal), y = abs(df$studentized_residual))) +
  geom_hline(yintercept = 3, color = "grey", linetype = "dashed") +
  geom_point(aes(color = ifelse(abs(df$studentized_residual) > 3, 'red', 'black'))) +
  scale_color_identity() +
  labs(title = "Distribución de los residuos studentized",
       x = "predicción modelo") + 
  theme_bw() + theme(plot.title = element_text(hjust = 0.5))
which(abs(df$studentized_residual) > 3)

No aparecen valores de color rojo, es decir, no existen valores atípicos o extremos

#Valores influyentes

summary(influence.measures(model = modelo.lineal))

Nota:

#Gráfico de valores influyentes

influencePlot(model = modelo.lineal)

#Recta de regresión con y sin valor de mayor influencia

ggplot(data = df, mapping = aes(x = Años, y = Madera)) +
  geom_point(color = "grey50", size = 2) +
  geom_smooth(method = "lm", se = FALSE, color = "black") +
  geom_point(data = df[10, ], color = "red", size = 2) +
  geom_smooth(data = df[-10, ], method="lm", se =FALSE,color = "blue") +
  labs(x = "Años", y = "Madera Aserrada") +   
  theme_bw() + 
  theme(plot.title = element_text(hjust = 0.5), legend.position = "none")

#Modelo eliminando el dato número 10

modelo.lineal2 <- lm(Madera ~ Años, data = df[-10,])
summary(modelo.lineal2)

Modelo 2: \[\hat{y}=36.40{x}-72008.94\]

• Volumen de Madera Aserrada PROMEDIO que esperaríamos para los años 2000, 2001 y 2002

#Predecir nuevos valores creando un data frame

nuevos.Valores <- data.frame(Años = c(2000, 2001, 2002))
nuevos.Valores

#Ejecución de la predicción

prediccion <- predict(modelo.lineal2, newdata = nuevos.Valores, interval="prediction")
prediccion

#Liberar nombres de variables de la memoria

detach(df)

#Limpiar Ambiente de Trabajo

rm(list = ls())