library(stargazer)#carga la biblioteca 'stargazer', que se utiliza para crear tablas de resumen de modelos de regresión.
modelo <- lm(formula= price ~ lotsize + sqrft + bdrms, data = hprice1)#ajusta un modelo de regresión lineal utilizando 'lm'. El modelo predice 'price' a partir de 'lotsize', 'sqrft' y 'bdrms' usando el conjunto de datos 'hprice1'.
stargazer(modelo, title = "Modelo Estimacion", type = "text", digits = 6)#crea una tabla de resumen del modelo 'modelo' y la muestra en la consola. El título de la tabla es "Modelo Estimacion" y se ajustan seis dígitos de precisión.

library(fitdistrplus)#carga la biblioteca 'fitdistrplus', que se utiliza para ajustar distribuciones de probabilidad a datos y realizar análisis de bondad de ajuste.
fit_normal<-fitdist(data = modelo$residuals,distr = "norm")#ajusta una distribución normal a los residuos del modelo de regresión lineal 'modelo'.
plot(fit_normal)#Traza el resultado del ajuste de la distribución normal.

summary(fit_normal)#muestra un resumen de los parámetros estimados para la distribución normal ajustada.

library(tseries)#carga la biblioteca 'tseries', que proporciona funciones para el análisis de series temporales.
Salida_JB<-jarque.bera.test(modelo$residuals)#realiza el test de Jarque-Bera para evaluar si los residuos del modelo 'modelo' siguen una distribución normal.
print(Salida_JB)#imprime los resultados del test de Jarque-Bera en la consola.

library(fastGraph) #carga la biblioteca 'fastGraph', que proporciona funciones para gráficos rápidos.
sig<-0.05 # Nivel de significancia para el test de Jarque-Bera.
JB<-Salida_JB$statistic # Estadístico de prueba del test de Jarque-Bera.
gl<-Salida_JB$parameter # Grados de libertad del test de Jarque-Bera.
Vc<-qchisq(1-sig,gl,lower.tail = TRUE)# Valor crítico para el test de Jarque-Bera.
shadeDist(JB,ddist = "dchisq",
          parm1 = gl,
          lower.tail = FALSE,xmin = 0,
          sub=paste("Vc:",round(Vc,2)," ","JB:",round(JB,2)))# Grafica la distribución chi-cuadrado con el valor crítico y el estadístico de prueba del test de Jarque-Bera.

library(dplyr) # carga la biblioteca 'dplyr', que proporciona funciones para manipulación de datos.
library(gt) #carga la biblioteca 'gt', que se utiliza para crear tablas atractivas en R.
library(gtExtras) #carga la extensión 'gtExtras' para añadir estilos y notas a las tablas creadas con 'gt'.
residuos_modelo<-modelo$residuals # Almacena los residuos del modelo en 'residuos_modelo'.
residuos_modelo %>%
  as_tibble() %>% 
  mutate(Posicion = row_number()) %>% # Agrega una columna 'Posicion' que indica la posición de cada residuo.
  arrange(value) %>% # Ordena los residuos de menor a mayor.
  mutate(dist1 = row_number() / n()) %>% # Calcula la distribución acumulada empírica.
  mutate(dist2 = (row_number() - 1) / n()) %>% # Calcula la segunda distribución acumulada empírica.
  mutate(zi = as.vector(scale(value, center = TRUE))) %>% # Estandariza los residuos.
  mutate(pi = pnorm(zi, lower.tail = TRUE)) %>% # Calcula la probabilidad acumulada teórica.
  mutate(dif1 = abs(dist1 - pi)) %>% # Calcula la diferencia absoluta entre las dos distribuciones acumuladas.
  mutate(dif2 = abs(dist2 - pi)) %>% # Calcula la diferencia absoluta entre las dos distribuciones acumuladas.
  rename(residuales = value) -> tabla_KS # Renombra la columna de residuos y almacena el resultado en 'tabla_KS'.

# Formato de la Tabla

tabla_KS %>%  
  gt() %>%  
  tab_header("Tabla para calcular el Estadistico KS") %>% # Título de la tabla
  tab_source_note(source_note = "Fuente: Elaboración propia") %>%  # Nota de la fuente
  tab_style( style = list( cell_fill(color = "#4ED173"), # agrega un estilo a las celdas de la tabla, se colorea la celda de color verde claro y se escribe en cursiva el texto que contiene dicha celda. Además, se identifican las celdas que tienen el valor de maxima diferencia 1.
      cell_text(style = "italic")),locations = cells_body(columns = dif1,  
      rows = dif1==max(dif1))) %>%
   tab_style( style = list(cell_fill(color = "#3498DB"),  
      cell_text(style = "italic")),locations = cells_body(columns = dif2, 
      rows = dif2==max(dif2)))# agrega un estilo a las celdas de la tabla, se colorea la celda de color cian azul y se escribe en cursiva el texto que contiene dicha celda. Además, se identifican las celdas que tienen el valor de maxima diferencia 2.

D<- max(max(tabla_KS$dif1),max(tabla_KS$dif2))# Calcula el máximo de las diferencias máximas entre las dos distribuciones acumuladas.
print(D)# Imprime el valor de D, que representa el estadístico de Kolmogorov-Smirnov.

# n = 88 , α = 0.05
# las observaciones son mayores a 50
# 0.0875897/√88 = 0.093370
# Dado que 0.0754392 < 0.093371. No se rechaza la hipotesis nula por lo que los residuos siguen una distribución normal.

library(nortest)# carga la biblioteca 'nortest', que proporciona funciones para realizar pruebas de normalidad.
Salida_KS<-lillie.test(modelo$residuals) # Realiza la prueba de Lilliefors para evaluar si los residuos del modelo siguen una distribución normal.
print(Salida_KS) #Imprime los resultados de la prueba de Lilliefors en la consola.

library(dplyr) # Carga la biblioteca 'dplyr' para manipulación de datos.
library(gt) # Carga la biblioteca 'gt' para crear tablas atractivas en R.
library(forecast) # Carga la biblioteca 'forecast' para funciones relacionadas con pronósticos de series temporales.
residuos<-modelo$residuals #Crea un objeto llamado "residuos" que almacena los residuos del modelo lineal en la variable "modelo_lineal".
residuos %>%
  as_tibble() %>%
  rename(residuales=value) %>%
  arrange(residuales) %>%
  mutate(pi=(row_number()-0.375)/(n()+0.25)) %>%
  mutate(mi=qnorm(pi,lower.tail = TRUE)) %>% 
  mutate(ai=0)->tabla_SW #Toma los residuos del modelo, los convierte a un "tibble", los renombra como "residuales", los ordena de manera ascendente y luego se realizan una serie de cálculos para obtener las variables "pi", "mi", y "ai", que se almacenan en el objeto "tabla_SW".
m<-sum(tabla_SW$mi^2) # Calcula la suma de los cuadrados de los valores de mi en la tabla_SW y almacena el resultado en la variable 'm'.
n<-nrow(hprice1)  #Se calcula el número de filas de "hprice1" y se almacena en la variable "n".
theta<-1/sqrt(n) #Se calcula la variable "theta" como el inverso de la raíz cuadrada de "n".
tabla_SW$ai[n]<- -2.706056*theta^5+4.434685*theta^4-2.071190*theta^3-0.147981*theta^2+0.2211570*theta+tabla_SW$mi[n]/sqrt(m) # Calcula y asigna el valor de ai para la última fila de la tabla_SW utilizando la fórmula especificada.
tabla_SW$ai[n-1]<- -3.582633*theta^5+5.682633*theta^4-1.752461*theta^3-0.293762*theta^2+0.042981*theta+tabla_SW$mi[n-1]/sqrt(m) # Calcula el valor de ai para el penúltimo elemento de la variable "tabla_SW", según una fórmula dada.
tabla_SW$ai[1]<- -tabla_SW$ai[n] # Asigna el valor de ai para la primera fila de tabla_SW como el negativo del valor de ai de la última fila. 
tabla_SW$ai[2]<- -tabla_SW$ai[n-1] # Asigna el valor de ai para la segunda fila de tabla_SW como el negativo del valor de ai de la penúltima fila.
omega<-(m-2*tabla_SW$mi[n]^2-2*tabla_SW$mi[n-1]^2)/(1-2*tabla_SW$ai[n]^2-2*tabla_SW$ai[n-1]^2) # Calcula omega utilizando la fórmula dada, que implica los valores de m, mi y ai para las dos últimas filas de la tabla_SW.
tabla_SW$ai[3:(n-2)]<-tabla_SW$mi[3:(n-2)]/sqrt(omega) #  Calcula y asigna los valores de 'ai' para las filas 3 a (n-2) de tabla_SW.

tabla_SW %>% mutate(ai_ui=ai*residuales,ui2=residuales^2) ->tabla_SW # Calcula y agrega dos nuevas columnas a la tabla_SW:'ai_ui', que representa el producto de 'ai' y los residuales.'ui2', que representa los residuales al cuadrado.Estos cálculos son útiles para el cálculo del estadístico de Shapiro-Wilk.

tabla_SW %>% gt() %>% tab_header("Tabla para calcular el Estadistico W") %>%
  tab_source_note(source_note = "Fuente: Elaboración propia") # Crea una tabla utilizando la función 'gt' y luego agrega un encabezado con el título "Tabla para calcular el Estadistico W".También se añade una nota de fuente con la frase "Fuente: Elaboración propia".

W<-(sum(tabla_SW$ai_ui)^2)/sum(tabla_SW$ui2) #Calcula el estadístico de Shapiro-Wilk (W) utilizando la fórmula dada.
print(W) #imprime el valor de W.

mu<-0.0038915*log(n)^3-0.083751*log(n)^2-0.31082*log(n)-1.5861 # Calcula el valor esperado de W (mu) utilizando la fórmula dada, que depende del logaritmo del tamaño de muestra.
sigma<-exp(0.0030302*log(n)^2-0.082676*log(n)-0.4803) # Calcula la desviación estándar de W (sigma) utilizando la fórmula dada, que también depende del logaritmo del tamaño de muestra.
Wn<-(log(1-W)-mu)/sigma # Calcula el valor normalizado de W (Wn) utilizando la fórmula dada.
print(Wn) # Imprime el valor de Wn.

p.value<-pnorm(Wn,lower.tail = FALSE) # Calcula el valor p utilizando la función pnorm, que representa la probabilidad de observar un valor de Wn o mayor.
print(p.value) # Imprime el valor p.value.

library(fastGraph) # Carga la biblioteca 'fastGraph', que proporciona funciones para gráficos rápidos.
shadeDist(Wn,ddist = "dnorm",lower.tail = FALSE) # Grafica la distribución normal estándar sombreada desde Wn hacia la derecha. Esto muestra visualmente la probabilidad de observar un valor de Wn o mayor bajo la distribución normal estándar.

# En este caso dado que 0.0005937472 < 0.05 Se rechaza la Hipótesis Nula por lo que los residuos no siguen una distribución normal.

Salida_SW<-shapiro.test(modelo$residuals) # Realiza el test de Shapiro-Wilk para evaluar si los residuos del modelo siguen una distribución normal.
print(Salida_SW) # Imprime los resultados del test de Shapiro-Wilk en la consola.

Salida_Wn<-qnorm(Salida_SW$p.value,lower.tail = FALSE) # Calcula el valor normalizado de Wn correspondiente al valor p obtenido del test de Shapiro-Wilk.
print(Salida_Wn) # Imprime el valor de Wn obtenido del valor p del test de Shapiro-Wilk.