library(readxl)
ventas_empresa <- read_excel("C:/Users/gisel/Downloads/ventas_empresa.xlsx")
View(ventas_empresa)
# 1. CARGA DE LIBRERÍAS Y DATOS
library(dplyr)      
library(gt)         
library(tseries)    # Jarque-Bera
library(nortest)    # Lilliefors
library(fastGraph)  # Gráficas
library(stargazer)  # Tablas de regresión
#ESTIMACIÓN DEL MODELO
modelo_ventas <- lm(V ~ C + P + M, data = ventas_empresa)

# Resultados del modelo
stargazer(modelo_ventas, title = "Modelo de Regresión de Ventas", type = "text")
## 
## Modelo de Regresión de Ventas
## ===============================================
##                         Dependent variable:    
##                     ---------------------------
##                                  V             
## -----------------------------------------------
## C                            0.923***          
##                               (0.223)          
##                                                
## P                            0.950***          
##                               (0.156)          
##                                                
## M                            1.298***          
##                               (0.431)          
##                                                
## Constant                    107.444***         
##                              (18.057)          
##                                                
## -----------------------------------------------
## Observations                    24             
## R2                             0.980           
## Adjusted R2                    0.977           
## Residual Std. Error       9.506 (df = 20)      
## F Statistic           323.641*** (df = 3; 20)  
## ===============================================
## Note:               *p<0.1; **p<0.05; ***p<0.01
# Extracción de residuos
residuos <- modelo_ventas$residuals
# 3. VERIFICACIÓN DEL SUPUESTO DE NORMALIDAD
# A) PRUEBA JARQUE-BERA (JB)
options(scipen = 999)
salida_JB <- jarque.bera.test(residuos)
print(salida_JB)
## 
##  Jarque Bera Test
## 
## data:  residuos
## X-squared = 1.4004, df = 2, p-value = 0.4965
# Gráfica JB (Distribución Chi-cuadrado)
alpha_sig <- 0.05
JB_stat <- salida_JB$statistic
gl_jb <- salida_JB$parameter
VC_jb <- qchisq(1 - alpha_sig, gl_jb)

shadeDist(JB_stat, ddist = "dchisq", parm1 = gl_jb, lower.tail = FALSE, xmin = 0,
          sub = paste("VC:", round(VC_jb, 2), " ", "JB:", round(JB_stat, 2)))

# B) PRUEBA LILLIEFORS (KS)
prueba_KS <- lillie.test(residuos)
print(prueba_KS)
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  residuos
## D = 0.13659, p-value = 0.2935
# Tabla de cálculo manual (Estructura solicitada)
tabla_KS <- residuos %>%
  as_tibble() %>%
  mutate(posicion = row_number()) %>%
  arrange(value) %>%
  mutate(dist1 = row_number() / n(),
         dist2 = (row_number() - 1) / n(),
         zi = as.vector(scale(value, center = TRUE)),
         pi = pnorm(zi),
         dif1 = abs(dist1 - pi),
         dif2 = abs(dist2 - pi))

# Mostrar tabla con destacados
tabla_KS %>%
  gt() %>%
  tab_header("Tabla para cálculo del Estadístico KS (Lilliefors)") %>%
  tab_style(style = cell_fill(color = "lightblue"),
            locations = cells_body(columns = dif1, rows = dif1 == max(dif1)))
Tabla para cálculo del Estadístico KS (Lilliefors)
value posicion dist1 dist2 zi pi dif1 dif2
-17.27950969 11 0.04166667 0.00000000 -1.949405386 0.02562352 0.0160431507 0.025623516
-15.50384667 10 0.08333333 0.04166667 -1.749082164 0.04013841 0.0431949244 0.001528258
-13.84401165 19 0.12500000 0.08333333 -1.561826195 0.05916447 0.0658355305 0.024168864
-9.91393148 5 0.16666667 0.12500000 -1.118450221 0.13168738 0.0349792889 0.006687378
-8.43534410 9 0.20833333 0.16666667 -0.951641887 0.17063932 0.0376940180 0.003972649
-7.90808316 24 0.25000000 0.20833333 -0.892158410 0.18615402 0.0638459847 0.022179318
-6.65267950 23 0.29166667 0.25000000 -0.750528776 0.22646815 0.0651985168 0.023531850
-3.33614422 8 0.33333333 0.29166667 -0.376370489 0.35332074 0.0199874078 0.061654074
-3.32226362 4 0.37500000 0.33333333 -0.374804535 0.35390292 0.0210970793 0.020569587
-2.43553165 3 0.41666667 0.37500000 -0.274766969 0.39174764 0.0249190286 0.016747638
-0.40251845 21 0.45833333 0.41666667 -0.045410527 0.48189005 0.0235567120 0.065223379
0.01280369 14 0.50000000 0.45833333 0.001444461 0.50057626 0.0005762565 0.042242923
3.09667909 20 0.54166667 0.50000000 0.349354988 0.63658859 0.0949219225 0.136588589
4.04562361 7 0.58333333 0.54166667 0.456411125 0.67595282 0.0926194830 0.134286150
4.98802235 12 0.62500000 0.58333333 0.562728793 0.71319021 0.0881902126 0.129856879
5.88324503 22 0.66666667 0.62500000 0.663724246 0.74656659 0.0798999242 0.121566591
6.28004614 18 0.70833333 0.66666667 0.708489765 0.76067942 0.0523460837 0.094012750
6.69809559 15 0.75000000 0.70833333 0.755652438 0.77507120 0.0250711962 0.066737863
6.78865323 16 0.79166667 0.75000000 0.765868789 0.77812281 0.0135438559 0.028122811
7.37251119 2 0.83333333 0.79166667 0.831737317 0.79722138 0.0361119495 0.005554717
8.70403920 6 0.87500000 0.83333333 0.981954997 0.83693899 0.0380610096 0.003605657
9.77138853 13 0.91666667 0.87500000 1.102369093 0.86484938 0.0518172880 0.010150621
10.67367778 1 0.95833333 0.91666667 1.204161768 0.88573647 0.0725968632 0.030930197
14.71907878 17 1.00000000 0.95833333 1.660547780 0.95159785 0.0484021522 0.006735485
 
# C) PRUEBA SHAPIRO-WILK (SW)
# ---------------------------------------------------------
# 1. Preparación de tabla manual y coeficientes ai
n_obs <- length(residuos)
residuos %>%    
  as_tibble() %>%  
  rename(residuales = value) %>%  
  arrange(residuales) %>%  
  mutate(pi = (row_number() - 0.375) / (n() + 0.25)) %>%  
  mutate(mi = qnorm(pi)) %>%   
  mutate(ai = 0) -> tabla_SW

# Cálculo de ai (Aproximación de Royston para n=24)
m_val <- sum(tabla_SW$mi^2)
theta <- 1 / sqrt(n_obs)

# Coeficientes extremos
tabla_SW$ai[n_obs] <- -2.706056*theta^5 + 4.434685*theta^4 - 2.071190*theta^3 - 0.147981*theta^2 + 0.2211570*theta + tabla_SW$mi[n_obs]/sqrt(m_val)
tabla_SW$ai[n_obs-1] <- -3.582633*theta^5 + 5.682633*theta^4 - 1.752461*theta^3 - 0.293762*theta^2 + 0.042981*theta + tabla_SW$mi[n_obs-1]/sqrt(m_val)
tabla_SW$ai[1] <- -tabla_SW$ai[n_obs]
tabla_SW$ai[2] <- -tabla_SW$ai[n_obs-1]

# Valores intermedios
omega <- (m_val - 2*tabla_SW$mi[n_obs]^2 - 2*tabla_SW$mi[n_obs-1]^2) / (1 - 2*tabla_SW$ai[n_obs]^2 - 2*tabla_SW$ai[n_obs-1]^2)
tabla_SW$ai[3:(n_obs-2)] <- tabla_SW$mi[3:(n_obs-2)] / sqrt(omega)

# Estadístico W manual
tabla_SW %>% mutate(ai_ui = ai * residuales, ui2 = residuales^2) -> tabla_SW
W_manual <- (sum(tabla_SW$ai_ui)^2) / sum(tabla_SW$ui2)

# Normalización Wn para la gráfica
mu_sw <- 0.0038915*log(n_obs)^3 - 0.083751*log(n_obs)^2 - 0.31082*log(n_obs) - 1.5861
sigma_sw <- exp(0.0030302*log(n_obs)^2 - 0.082676*log(n_obs) - 0.4803)
Wn_stat <- (log(1 - W_manual) - mu_sw) / sigma_sw

# Gráfica Shapiro-Wilk (Distribución Normal Estándar)
shadeDist(Wn_stat, ddist = "dnorm", lower.tail = FALSE, 
          sub = paste("Wn:", round(Wn_stat, 3), "W:", round(W_manual, 4)))

# Validación automática
shapiro.test(residuos)
## 
##  Shapiro-Wilk normality test
## 
## data:  residuos
## W = 0.95315, p-value = 0.3166

  1. Prueba de Jarque-Bera (JB)Este test analizó si la asimetría y la curtosis de los residuales coinciden con las de una distribución normal. El estadístico obtenido fue de 1.4004 con un p-valor de 0.4965. Al superar ampliamente el umbral crítico de 0.05, se mantiene la hipótesis nula . Visualmente, el estadístico se sitúa en la zona de confianza de una distribución chi2 con dos grados de libertad, confirmando que los errores siguen la morfología de la campana de Gauss.

  2. Prueba de Lilliefors (Corrección de Kolmogorov-Smirnov)Considerando que los parámetros poblacionales fueron estimados a partir de la muestra n=24, se empleó la corrección de Lilliefors. El análisis arrojó un estadístico D = 0.13659 y un p-valor de 0.2935. Dado que el p-valor es superior al nivel de significancia de 5%, se concluye que no existe una diferencia significativa entre la distribución observada y la teórica. Esto ratifica que las desviaciones máximas son despreciables y que los errores operan bajo condiciones de normalidad.

  3. Prueba de Shapiro-Wilk (SW)Dada la alta potencia de esta prueba en muestras pequeñas, se utilizó como validación final. El resultado mostró un valor W = 0.95315 y un p-valor de 0.3166. Al no encontrar evidencia suficiente para rechazar la normalidad (p > 0.05), se confirma que el estadístico normalizado W_n se ubica fuera de las regiones críticas de rechazo. Este hallazgo asegura la linealidad y el comportamiento estable de los términos de error en el modelo.

# Evaluación de Multicolinealidad
# ---

# 1. Carga de librerías necesarias
library(stargazer)
library(mctest)
library(psych)
library(car)
library(fastGraph)
# --- 1. Índice de Condición ---
X_mat <- model.matrix(modelo_ventas) # Matriz de diseño
XX_matrix <- t(X_mat) %*% X_mat

# Normalización
options(scipen = 999)
Sn <- solve(diag(sqrt(diag(XX_matrix))))
XX_norm <- (Sn %*% XX_matrix) %*% Sn

# Autovalores y Cálculo de K
lambdas <- eigen(XX_norm, symmetric = TRUE)
K <- sqrt(max(lambdas$values) / min(lambdas$values))
print(paste("Índice de Condición (K):", round(K, 4)))
## [1] "Índice de Condición (K): 71.1635"
# --- 2. Prueba de Farrar-Glauber (Cálculo Manual) ---
Zn <- scale(X_mat[,-1]) # Escalamiento de regresores
n <- nrow(Zn)
m <- ncol(Zn)
R <- (t(Zn) %*% Zn) * (1 / (n - 1)) # Matriz de Correlación R
det_R <- det(R)

# Estadístico Chi-cuadrado FG
chi_FG <- -(n - 1 - (2 * m + 5) / 6) * log(det_R)
gl_fg <- m * (m - 1) / 2
VC_fg <- qchisq(p = 0.95, df = gl_fg)
p_val_fg <- pchisq(chi_FG, df = gl_fg, lower.tail = FALSE)

# Gráfica FG
shadeDist(chi_FG, ddist = "dchisq", parm1 = gl_fg, lower.tail = FALSE, 
          sub = paste("VC:", round(VC_fg, 2), " | FG:", round(chi_FG, 2)))

# --- 3. Factores Inflacionarios de la Varianza (VIF) ---
VIFs_car <- vif(modelo_ventas)
print(VIFs_car)
##        C        P        M 
## 7.631451 3.838911 9.449210