library(fitdistrplus)  # Cargar la librería fitdistrplus para ajustar distribuciones
fit_normal<-fitdist(data = modelo$residuals,distr = "norm")  # Ajustar una distribución normal a los residuos del modelo
plot(fit_normal)  # Graficar el ajuste de la distribución normal

summary(fit_normal)  # Obtener el resumen del ajuste de la distribución normal

library(tseries)  # Cargar la librería tseries para realizar la prueba de normalidad
salida_JB<-jarque.bera.test(modelo$residuals)  # Aplicar la prueba de Jarque-Bera a los residuos del modelo
salida_JB  # Mostrar los resultados de la prueba

library(normtest)  # Cargar la librería normtest para realizar pruebas de normalidad
jb.norm.test(modelo$residuals)  # Aplicar la prueba de Jarque-Bera a los residuos del modelo

# Evitar la notación científica en la salida
options(scipen = 999999)

# Cargar la librería fastGraph para graficar distribuciones estadísticas
library(fastGraph)

# Nivel de significancia para la prueba
alpha_sig <- 0.05

# Obtener el estadístico de Jarque-Bera de la prueba de normalidad
JB <- salida_JB$statistic

# Grados de libertad de la prueba de Jarque-Bera
gl <- salida_JB$parameter

# Calcular el valor crítico de la distribución chi-cuadrado
VC <- qchisq(1 - alpha_sig, gl, lower.tail = TRUE)

# Graficar la distribución chi-cuadrado con sombreado en la región de interés
shadeDist(JB, ddist = "dchisq",
          parm1 = gl,
          lower.tail = FALSE, xmin = 0,
          sub = paste("VC:", round(VC, 2), " ", "JB:", round(JB, 2)))

library(dplyr)  # Carga la librería dplyr para manipulación de datos
library(gt)  # Carga la librería gt para crear tablas de datos
library(gtExtras)  # Carga la librería gtExtras para agregar funcionalidades a las tablas creadas con gt

residuos<-modelo$residuals  # Crea un vector con los residuos del modelo estimado
residuos %>%  # Utiliza el operador %>% para encadenar las operaciones siguientes al vector residuos
  as_tibble() %>%  # Convierte el vector residuos en una tibble (tabla) de una columna
  mutate(posicion=row_number()) %>%  # Agrega una columna llamada "posicion" con el número de fila
  arrange(value) %>%  # Ordena la tabla por los valores de residuos en orden ascendente
  mutate(dist1=row_number()/n()) %>%  # Agrega una columna "dist1" con los percentiles según su posición en la tabla (usando la función row_number() y n() para obtener el número de filas)
  mutate(dist2=(row_number()-1)/n()) %>%  # Agrega una columna "dist2" con los percentiles según su posición en la tabla, pero ajustando en una unidad para evitar problemas con los extremos de la distribución
  mutate(zi=as.vector(scale(value,center=TRUE))) %>%  # Agrega una columna "zi" con los valores de residuos escalados para tener media cero y varianza uno
  mutate(pi=pnorm(zi,lower.tail = TRUE)) %>%  # Agrega una columna "pi" con los valores de la función de distribución acumulada (CDF) de una distribución normal estándar evaluada en los valores de zi
  mutate(dif1=abs(dist1-pi)) %>%  # Agrega una columna "dif1" con las diferencias absolutas entre los percentiles según la posición y los valores de pi
  mutate(dif2=abs(dist2-pi)) %>%  # Agrega una columna "dif2" con las diferencias absolutas entre los percentiles ajustados según la posición y los valores de pi
  rename(residuales=value) -> tabla_KS  # Renombra la columna "value" como "residuales" y asigna la tabla resultante a la variable tabla_KS


#Formato
 tabla_KS %>%  # Utiliza el operador %>% para encadenar las operaciones siguientes a la tabla tabla_KS
  gt() %>%  # Crea una tabla con la función gt()
  tab_header("Tabla para calcular el Estadistico KS") %>%  # Agrega un encabezado a la tabla
  tab_source_note(source_note = "Fuente: Elaboración propia.") %>%  # Agrega una nota de fuente a la tabla
  tab_style(  # Cambia el estilo de algunas celdas de la tabla
    style = list(
      cell_fill(color = "#A569BD"),  # Cambia el color de fondo de las celdas a un tono de morado
      cell_text(style = "italic")  # Cambia el estilo de texto de las celdas a itálico
      ),
    locations = cells_body(  # Aplica el estilo a las celdas del cuerpo de la tabla que cumplan las siguientes condiciones:
      columns = dif1,  # Que pertenezcan a la columna "dif1"
      rows = dif1==max(dif1)  # Que pertenezcan a la fila donde el valor de "dif1" es máximo
    )) %>%
   tab_style(  # Cambia el estilo de algunas celdas de la tabla
    style = list(
      cell_fill(color = "#3498DB"),  # Cambia el color de fondo de las celdas a un tono de azul
      cell_text(style = "italic")  # Cambia el estilo de texto de las celdas a itálico
      ),
    locations = cells_body(  # Aplica el estilo a las celdas del cuerpo de la tabla que cumplan las siguientes condiciones:
      columns = dif2,  # Que pertenezcan a la columna "dif2"
      rows = dif2==max(dif2)  # Que pertenezcan a la fila donde el valor de "dif2" es máximo
    ))

D<-max(max(tabla_KS$dif1),max(tabla_KS$dif2))  # Calcular el valor máximo entre las diferencias dif1 y dif2 en tabla_KS
print(D)  # Imprimir el valor calculado

library(nortest)  # Cargar la librería nortest para realizar pruebas de normalidad
prueba_KS<-lillie.test(modelo$residuals)  # Aplicar la prueba de Lilliefors a los residuos del modelo
prueba_KS  # Mostrar los resultados de la prueba

p.value<-prueba_KS$p.value  # Extraer el valor p de la prueba de Lilliefors

library(dplyr)   # Carga la librería dplyr, útil para la manipulación de datos en forma de "tibbles".
library(gt)      # Carga la librería gt, que permite crear tablas formateadas de manera atractiva.
residuos <- modelo$residuals # Asigna los residuos del objeto 'modelo' (presumiblemente un modelo estadístico) a la variable 'residuos'.
residuos %>%
  as_tibble() %>%   # Convierte el vector de residuos en un tibble (un tipo de data frame mejorado).
  rename(residuales = value) %>% # Renombra la columna 'value' del tibble a 'residuales'.
  arrange(residuales) %>%   # Ordena las filas del tibble en función de los valores de la columna 'residuales' de forma ascendente.
  mutate(pi = (row_number() - 0.375) / (n() + 0.25)) %>% # Calcula la probabilidad acumulada empírica (plotting position) 'pi' para cada residuo.
  mutate(mi = qnorm(pi, lower.tail = TRUE)) %>% # Calcula los cuantiles teóricos de una distribución normal estándar ('mi') correspondientes a las probabilidades 'pi'.
  mutate(ai = 0) -> tabla_SW  # Crea una nueva columna 'ai' en el tibble 'tabla_SW' y la inicializa con el valor 0 para todas las filas.

m <- sum(tabla_SW$mi^2) # Calcula la suma de los cuadrados de los valores en la columna 'mi' y la asigna a la variable 'm'.
n <- nrow(hprice1)      # Obtiene el número de filas del data frame 'hprice1' y lo asigna a la variable 'n'. Se asume que 'hprice1' contiene los datos originales del modelo.
theta <- 1 / sqrt(n)    # Calcula un valor 'theta' basado en el tamaño de la muestra '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 el valor de 'ai' para la última fila (n) utilizando una fórmula específica que involucra 'theta', 'mi' y 'm'.
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 la penúltima fila (n-1) utilizando una fórmula similar a la anterior.
tabla_SW$ai[1] <- -tabla_SW$ai[n]     # Asigna el valor negativo de 'ai' de la última fila a la primera fila.
tabla_SW$ai[2] <- -tabla_SW$ai[n - 1] # Asigna el valor negativo de 'ai' de la penúltima fila a la segunda 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 el valor de 'omega' utilizando 'm' y los valores de 'mi' y 'ai' de las dos últimas filas.
tabla_SW$ai[3:(n - 2)] <- tabla_SW$mi[3:(n - 2)] / sqrt(omega) # Calcula los valores de 'ai' para las filas intermedias (desde la tercera hasta la penúltima) dividiendo los valores correspondientes de 'mi' por la raíz cuadrada de 'omega'.

tabla_SW %>%
  mutate(ai_ui = ai * residuales, ui2 = residuales^2) -> tabla_SW # Crea dos nuevas columnas: 'ai_ui' (producto de 'ai' y 'residuales') y 'ui2' (cuadrado de 'residuales').

tabla_SW %>%
  gt() %>%   # Convierte el tibble 'tabla_SW' en un objeto 'gt_tbl' para crear una tabla formateada.
  tab_header("Tabla para calcular el Estadistico W") %>%       # Agrega un encabezado a la tabla con el texto especificado.
  tab_source_note(source_note = "Fuente: Elaboración propia.") # Agrega una nota al pie de la tabla indicando la fuente de los datos.

W<-(sum(tabla_SW$ai_ui)^2)/sum(tabla_SW$ui2)  # Calcula el estadístico de Shapiro-Wilk (W) elevando al cuadrado la suma de la columna 'ai_ui' y dividiéndola por la suma de la columna 'ui2'.
print(W) # Imprime el valor calculado del estadístico W en la consola.

mu<-0.0038915*log(n)^3-0.083751*log(n)^2-0.31082*log(n)-1.5861  # Calcula el valor de 'mu' utilizando una fórmula que depende del logaritmo del tamaño de la muestra 'n'. Esta fórmula es una aproximación para la media de la distribución del logaritmo de (1-W) bajo la hipótesis nula de normalidad.
sigma<-exp(0.0030302*log(n)^2-0.082676*log(n)-0.4803)  # Calcula el valor de 'sigma' utilizando una fórmula exponencial que depende del logaritmo del tamaño de la muestra 'n'. Esta fórmula es una aproximación para la desviación estándar de la distribución del logaritmo de (1-W) bajo la hipótesis nula de normalidad.
Wn<-(log(1-W)-mu)/sigma  # Calcula el estadístico ajustado 'Wn' restando la media aproximada ('mu') del logaritmo de (1 - W) y dividiendo el resultado por la desviación estándar aproximada ('sigma'). Este estadístico transformado sigue aproximadamente una distribución normal estándar bajo la hipótesis nula.
print(Wn)  # Imprime el valor calculado del estadístico ajustado 'Wn' en la consola.

p.value<-pnorm(Wn,lower.tail = FALSE)  # Calcula el valor p (p-value) asociado al estadístico ajustado 'Wn'. La función 'pnorm' calcula la probabilidad acumulada de una distribución normal estándar. Al establecer 'lower.tail = FALSE', se obtiene la probabilidad de observar un valor mayor que 'Wn', que es el valor p para una prueba de una cola derecha. En este contexto, se utiliza para evaluar la significancia estadística de la desviación de la normalidad.
print(p.value)  # Imprime el valor p calculado en la consola. Un valor p pequeño (típicamente menor que un nivel de significancia predefinido, como 0.05) sugiere evidencia en contra de la hipótesis nula de normalidad.

library(fastGraph)  # Carga la librería 'fastGraph', que proporciona funciones para visualizaciones estadísticas, incluyendo el sombreado de áreas bajo curvas de distribución.
shadeDist(Wn,ddist = "dnorm",lower.tail = FALSE)   # Genera un gráfico de la distribución normal estándar (especificado por 'ddist = "dnorm"'). La función sombrea el área bajo la curva que corresponde a la probabilidad de observar un valor mayor que 'Wn' (indicado por 'lower.tail = FALSE'). Esta área sombreada representa visualmente el valor p calculado anteriormente.

salida_SW<-shapiro.test(modelo$residuals)  # Realiza la prueba de Shapiro-Wilk directamente sobre los residuos del objeto 'modelo' utilizando la función 'shapiro.test()' que viene incorporada en R. Esta función devuelve un objeto con los resultados de la prueba.
print(salida_SW)  # Imprime los resultados de la prueba de Shapiro-Wilk en la consola. Estos resultados típicamente incluyen el estadístico W y el valor p asociado, así como el nombre de la prueba y los datos utilizados. Esta es una forma más directa de realizar la prueba de normalidad de Shapiro-Wilk en R.

Wn_salida<-qnorm(salida_SW$p.value,lower.tail = FALSE)  # Calcula el cuantil de la distribución normal estándar (valor Z) correspondiente al valor p obtenido de la función 'shapiro.test()'. Al establecer 'lower.tail = FALSE', se busca el valor Z tal que el área a su derecha sea igual al valor p. Este valor Z puede ser útil para comparar con un nivel de significancia o para otras interpretaciones relacionadas con la normalidad de los residuos.
print(Wn_salida)  # Imprime el valor Z calculado en la consola.