library(wooldridge)
## Warning: package 'wooldridge' was built under R version 4.5.3
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.5.2
## 
## Adjuntando el paquete: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(nortest) # Útil para KS/Lilliefors
## Warning: package 'nortest' was built under R version 4.5.2
library(fastGraph)
## Warning: package 'fastGraph' was built under R version 4.5.3
data(hprice1)
head(force(hprice1), n=5)
##   price assess bdrms lotsize sqrft colonial   lprice  lassess llotsize   lsqrft
## 1   300  349.1     4    6126  2438        1 5.703783 5.855359 8.720297 7.798934
## 2   370  351.5     3    9903  2076        1 5.913503 5.862210 9.200593 7.638198
## 3   191  217.7     3    5200  1374        0 5.252274 5.383118 8.556414 7.225482
## 4   195  231.8     3    4600  1448        1 5.273000 5.445875 8.433811 7.277938
## 5   373  319.1     4    6095  2514        1 5.921578 5.765504 8.715224 7.829630

#Estimación del modelo

library(stargazer)
modelo_estimado <- lm(price ~ lotsize + sqrft + bdrms, data = hprice1)
stargazer(modelo_estimado,title = "MODELO_ESTIMADO_PRUEBAS",type = "text")
## 
## MODELO
## ===============================================
##                         Dependent variable:    
##                     ---------------------------
##                                price           
## -----------------------------------------------
## lotsize                      0.002***          
##                               (0.001)          
##                                                
## sqrft                        0.123***          
##                               (0.013)          
##                                                
## bdrms                         13.853           
##                               (9.010)          
##                                                
## Constant                      -21.770          
##                              (29.475)          
##                                                
## -----------------------------------------------
## Observations                    88             
## R2                             0.672           
## Adjusted R2                    0.661           
## Residual Std. Error      59.833 (df = 84)      
## F Statistic           57.460*** (df = 3; 84)   
## ===============================================
## Note:               *p<0.1; **p<0.05; ***p<0.01

#Ajuste de los residuos a la distribución normal

library(fitdistrplus)
## Warning: package 'fitdistrplus' was built under R version 4.5.3
## Cargando paquete requerido: MASS
## 
## Adjuntando el paquete: 'MASS'
## The following object is masked from 'package:dplyr':
## 
##     select
## The following object is masked from 'package:wooldridge':
## 
##     cement
## Cargando paquete requerido: survival
fit_normal<-fitdist(data = modelo_estimado$residuals,distr = "norm")
plot(fit_normal)

summary(fit_normal)
## Fitting of the distribution ' norm ' by maximum likelihood 
## Parameters : 
##           estimate Std. Error
## mean -2.321494e-15   6.231625
## sd    5.845781e+01   4.406424
## Loglikelihood:  -482.8775   AIC:  969.7549   BIC:  974.7096 
## Correlation matrix:
##      mean sd
## mean    1  0
## sd      0  1

#Prueba de normalidad de Jarque Bera ##Usando tseries

options(scipen = 9999)
library(tseries)
salida_JB<-jarque.bera.test(modelo_estimado$residuals)
salida_JB
## 
##  Jarque Bera Test
## 
## data:  modelo_estimado$residuals
## X-squared = 32.278, df = 2, p-value = 0.00000009794
#Resultado: Dado que el JB=32.278 y el nivel de signficancia es de 0.05 entonces 32.278>0.05 por lo tanto se rechaza la hipotesis nula.

#Calculo manual de la prueba Jarque Bera

library(dplyr)
library(gt)
## Warning: package 'gt' was built under R version 4.5.3
# 1. Extraer residuales y definir n
residuos<- modelo_estimado$residuals
n <- length(residuos)

# 2. Calcular los componentes de la prueba
# Asimetría (Skewness)
asimetria <- sum((residuos - mean(residuos))^3) / (n * (sd(residuos) * sqrt((n-1)/n))^3)

# Curtosis
curtosis <- sum((residuos - mean(residuos))^4) / (n * (sd(residuos) * sqrt((n-1)/n))^4)

# Estadístico Jarque-Bera (JB)
jb_est <- (n / 6) * (asimetria^2 + ((curtosis - 3)^2 / 4))

# Valor p (usa una distribución Chi-cuadrado con 2 grados de libertad)
p_val_jb <- pchisq(jb_est, df = 2, lower.tail = FALSE)

# 3. Crear la tabla de resultados
tabla_JB <- tibble(
  Estadistico = c("Asimetría", "Curtosis", "Jarque-Bera (JB)", "Valor p"),
  Valor = c(asimetria, curtosis, jb_est, p_val_jb),
  Interpretacion = c(
    ifelse(asimetria > 0, "Positiva", "Negativa"),
    ifelse(curtosis > 3, "Leptocúrtica", "Platicúrtica"),
    "Estadístico de prueba",
    ifelse(p_val_jb > 0.05, "No Rechaza Normalidad", "Rechaza Normalidad")
  )
)

# 4. Formatear con gt
tabla_JB %>%
  gt() %>%
  tab_header(title = "Resultados de la Prueba Jarque-Bera",
             subtitle = "Evaluación de Asimetría y Curtosis de los Residuales") %>%
  fmt_number(columns = Valor, decimals = 4) %>%
  cols_label(Estadistico = "Métrica", Valor = "Valor Calculado")
Resultados de la Prueba Jarque-Bera
Evaluación de Asimetría y Curtosis de los Residuales
Métrica Valor Calculado Interpretacion
Asimetría 0.9607 Positiva
Curtosis 5.2608 Leptocúrtica
Jarque-Bera (JB) 32.2779 Estadístico de prueba
Valor p 0.0000 Rechaza Normalidad

###Graficando la prueba JB

alpha_sig<-0.05
JB<-salida_JB$statistic
gl<-salida_JB$parameter
VC<-qchisq(1-alpha_sig,gl,lower.tail=TRUE)
shadeDist(JB,ddist = "dchisq",
          parm1=gl,
          lower.tail= FALSE, xmin=0,
          sub=paste("VC:",round(VC,2)," ","JB:",round(JB,2)))

#Calculo de la prueba de Kolmogorov Smimov-Liliefors

library(nortest)
prueba_KS<-lillie.test(modelo_estimado$residuals)
prueba_KS
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  modelo_estimado$residuals
## D = 0.075439, p-value = 0.2496
#Resultado: Dado que p_value>nivel de significancia, corresponde a 0.2496>0.05 entonces se afirma que no se rechaza la hipotesis nula.

#Calculo manual de la prueba Kolmogorov Smimov-Liliefors

library(dplyr)  
library(gt) 
library(gtExtras)
## Warning: package 'gtExtras' was built under R version 4.5.3
## 
## Adjuntando el paquete: 'gtExtras'
## The following object is masked from 'package:MASS':
## 
##     select
residuos<-modelo_estimado$residuals  
residuos |>  
  as_tibble() |>   
  mutate(posicion=row_number()) |>  
  arrange(value) |>   
  mutate(dist1=row_number()/n()) |> 
  mutate(dist2=(row_number()-1)/n()) |> 
  mutate(zi=as.vector(scale(value,center=TRUE))) |>  
  mutate(pi=pnorm(zi,lower.tail = TRUE)) |>   
  mutate(dif1=abs(dist1-pi)) |>  
  mutate(dif2=abs(dist2-pi)) |> 
  rename(residuales=value) ->tabla_KS
  
  #Formato
 tabla_KS |>  
  gt() |>  
  tab_header("Tabla para calcular el Estadistico KS") |> 
  tab_source_note(source_note = "Fuente: Elaboración propia") |> 
  tab_style(  
    style = list(
      cell_fill(color = "#A569BD"),  
      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)
    ))
Tabla para calcular el Estadistico KS
residuales posicion dist1 dist2 zi pi dif1 dif2
-120.026447 81 0.01136364 0.00000000 -2.041515459 0.02059981 0.0092361731 0.0205998094
-115.508697 77 0.02272727 0.01136364 -1.964673586 0.02472601 0.0019987418 0.0133623781
-107.080889 24 0.03409091 0.02272727 -1.821326006 0.03427866 0.0001877487 0.0115513850
-91.243980 48 0.04545455 0.03409091 -1.551957925 0.06033615 0.0148816002 0.0262452366
-85.461169 12 0.05681818 0.04545455 -1.453598781 0.07302879 0.0162106057 0.0275742421
-77.172687 32 0.06818182 0.05681818 -1.312620980 0.09465535 0.0264735301 0.0378371665
-74.702719 54 0.07954545 0.06818182 -1.270609602 0.10193378 0.0223883300 0.0337519664
-65.502849 39 0.09090909 0.07954545 -1.114130117 0.13261169 0.0417025941 0.0530662305
-63.699108 69 0.10227273 0.09090909 -1.083450505 0.13930425 0.0370315271 0.0483951634
-62.566594 83 0.11363636 0.10227273 -1.064187703 0.14362184 0.0299854747 0.0413491110
-59.845223 36 0.12500000 0.11363636 -1.017900230 0.15436269 0.0293626861 0.0407263225
-54.466158 13 0.13636364 0.12500000 -0.926408352 0.17711690 0.0407532663 0.0521169027
-54.300415 14 0.14772727 0.13636364 -0.923589260 0.17785010 0.0301228311 0.0414864675
-52.129801 15 0.15909091 0.14772727 -0.886669532 0.18762842 0.0285375141 0.0399011505
-51.441108 17 0.17045455 0.15909091 -0.874955638 0.19079902 0.0203444766 0.0317081129
-48.704980 47 0.18181818 0.17045455 -0.828417174 0.20371714 0.0218989601 0.0332625965
-48.350295 29 0.19318182 0.18181818 -0.822384375 0.20542908 0.0122472664 0.0236109028
-47.855859 11 0.20454545 0.19318182 -0.813974573 0.20782976 0.0032843043 0.0146479407
-45.639765 1 0.21590909 0.20454545 -0.776281294 0.21879146 0.0028823668 0.0142460032
-43.142550 9 0.22727273 0.21590909 -0.733806463 0.23153335 0.0042606233 0.0156242596
-41.749618 57 0.23863636 0.22727273 -0.710114247 0.23881665 0.0001802823 0.0115439187
-40.869022 27 0.25000000 0.23863636 -0.695136302 0.24348494 0.0065150566 0.0048485798
-37.749811 34 0.26136364 0.25000000 -0.642082009 0.26040997 0.0009536682 0.0104099682
-36.663785 71 0.27272727 0.26136364 -0.623609925 0.26644190 0.0062853771 0.0050782592
-36.646568 79 0.28409091 0.27272727 -0.623317083 0.26653809 0.0175528221 0.0061891857
-33.801248 37 0.29545455 0.28409091 -0.574921384 0.28267223 0.0127823120 0.0014186757
-29.766931 16 0.30681818 0.29545455 -0.506302171 0.30632227 0.0004959124 0.0108677240
-26.696234 22 0.31818182 0.30681818 -0.454073044 0.32488813 0.0067063089 0.0180699452
-24.271531 23 0.32954545 0.31818182 -0.412831567 0.33986501 0.0103195566 0.0216831929
-23.651448 86 0.34090909 0.32954545 -0.402284648 0.34373728 0.0028281851 0.0141918214
-19.683427 88 0.35227273 0.34090909 -0.334793052 0.36889060 0.0166178738 0.0279815102
-17.817835 10 0.36363636 0.35227273 -0.303061413 0.38092153 0.0172851663 0.0286488027
-16.762094 60 0.37500000 0.36363636 -0.285104441 0.38778206 0.0127820638 0.0241457002
-16.596960 21 0.38636364 0.37500000 -0.282295711 0.38885839 0.0024947507 0.0138583870
-16.271207 58 0.39772727 0.38636364 -0.276755010 0.39098411 0.0067431583 0.0046204781
-13.815798 56 0.40909091 0.39772727 -0.234991254 0.40710776 0.0019831485 0.0093804879
-13.462160 75 0.42045455 0.40909091 -0.228976273 0.40944368 0.0110108666 0.0003527698
-12.081520 4 0.43181818 0.42045455 -0.205493119 0.41859344 0.0132247451 0.0018611087
-11.629207 51 0.44318182 0.43181818 -0.197799788 0.42160086 0.0215809622 0.0102173258
-11.312669 74 0.45454545 0.44318182 -0.192415834 0.42370825 0.0308372092 0.0194735728
-8.236558 3 0.46590909 0.45454545 -0.140094626 0.44429261 0.0216164775 0.0102528411
-7.662789 70 0.47727273 0.46590909 -0.130335452 0.44815052 0.0291222111 0.0177585748
-6.752801 67 0.48863636 0.47727273 -0.114857588 0.45427900 0.0343573625 0.0229937262
-6.707262 31 0.50000000 0.48863636 -0.114083016 0.45458599 0.0454140074 0.0340503710
-6.402439 85 0.51136364 0.50000000 -0.108898313 0.45664157 0.0547220642 0.0433584278
-5.446904 82 0.52272727 0.51136364 -0.092645733 0.46309251 0.0596347676 0.0482711313
-3.537785 43 0.53409091 0.52272727 -0.060173762 0.47600862 0.0580822876 0.0467186512
-2.824941 61 0.54545455 0.53409091 -0.048049090 0.48083856 0.0646159857 0.0532523493
-2.745208 68 0.55681818 0.54545455 -0.046692922 0.48137899 0.0754391961 0.0640755598
-0.195089 65 0.56818182 0.55681818 -0.003318245 0.49867621 0.0695056040 0.0581419676
1.399296 55 0.57954545 0.56818182 0.023800450 0.50949411 0.0700513452 0.0586877088
5.363331 26 0.59090909 0.57954545 0.091224254 0.53634280 0.0545662924 0.0432026561
6.700640 53 0.60227273 0.59090909 0.113970383 0.54536936 0.0569033628 0.0455397265
7.386314 80 0.61363636 0.60227273 0.125632935 0.54998875 0.0636476093 0.0522839730
9.099900 41 0.62500000 0.61363636 0.154779103 0.56150227 0.0634977329 0.0521340965
12.433611 46 0.63636364 0.62500000 0.211481796 0.58374433 0.0526193043 0.0412556680
16.718018 62 0.64772727 0.63636364 0.284354766 0.61193074 0.0357965328 0.0244328965
18.093192 5 0.65909091 0.64772727 0.307744934 0.62086179 0.0382291219 0.0268654856
18.801816 38 0.67045455 0.65909091 0.319797835 0.62543921 0.0450153400 0.0336517036
19.168108 33 0.68181818 0.67045455 0.326028052 0.62779843 0.0540197476 0.0426561112
19.219211 72 0.69318182 0.68181818 0.326897255 0.62812720 0.0650546167 0.0536909803
20.334434 59 0.70454545 0.69318182 0.345865960 0.63527827 0.0692671805 0.0579035442
24.909926 78 0.71590909 0.70454545 0.423689939 0.66410402 0.0518050676 0.0404414312
26.236229 40 0.72727273 0.71590909 0.446248874 0.67229126 0.0549814685 0.0436178321
30.924022 25 0.73863636 0.72727273 0.525982978 0.70054998 0.0380863808 0.0267227444
32.253952 45 0.75000000 0.73863636 0.548603608 0.70836125 0.0416387548 0.0302751184
32.529367 49 0.76136364 0.75000000 0.553288104 0.70996693 0.0513967091 0.0400330727
32.675968 18 0.77272727 0.76136364 0.555781630 0.71081993 0.0619073452 0.0505437088
33.275839 20 0.78409091 0.77272727 0.565984762 0.71429793 0.0697929786 0.0584293423
36.031430 52 0.79545455 0.78409091 0.612854281 0.73001365 0.0654408934 0.0540772571
37.147186 84 0.80681818 0.79545455 0.631832029 0.73625168 0.0705665028 0.0592028664
40.320875 7 0.81818182 0.80681818 0.685812928 0.75358446 0.0645973596 0.0532337232
44.334467 30 0.82954545 0.81818182 0.754079634 0.77459930 0.0549461574 0.0435825211
46.907165 28 0.84090909 0.82954545 0.797838357 0.78751785 0.0533912405 0.0420276041
54.418366 87 0.85227273 0.84090909 0.925595465 0.82267187 0.0296008528 0.0182372164
55.091131 35 0.86363636 0.85227273 0.937038450 0.82563061 0.0380057535 0.0266421172
55.470305 44 0.87500000 0.86363636 0.943487765 0.82728426 0.0477157353 0.0363520989
62.939597 6 0.88636364 0.87500000 1.070532059 0.85781006 0.0285535797 0.0171899433
66.478628 50 0.89772727 0.88636364 1.130727018 0.87091500 0.0268122757 0.0154486394
67.426518 63 0.90909091 0.89772727 1.146849569 0.87427810 0.0348128083 0.0234491719
67.603959 19 0.92045455 0.90909091 1.149867648 0.87490081 0.0455537393 0.0341901029
69.707122 64 0.93181818 0.92045455 1.185640095 0.88211777 0.0497004123 0.0383367759
69.843246 8 0.94318182 0.93181818 1.187955411 0.88257451 0.0606073068 0.0492436705
74.848732 2 0.95454545 0.94318182 1.273093116 0.89850750 0.0560379553 0.0446743189
112.729191 66 0.96590909 0.95454545 1.917397313 0.97240626 0.0064971714 0.0178608078
163.795081 73 0.97727273 0.96590909 2.785970904 0.99733162 0.0200588896 0.0314225260
198.660139 42 0.98863636 0.97727273 3.378986513 0.99963623 0.0109998685 0.0223635048
209.375830 76 1.00000000 0.98863636 3.561248407 0.99981545 0.0001845478 0.0111790885
Fuente: Elaboración propia

##Cálculo del estadistico SK

D<-max(max(tabla_KS$dif1),max(tabla_KS$dif2))
print(D)
## [1] 0.0754392
##Resultado: dado que 0.0754392<0.094447917. Es decir que D<V.C entonces la hipótesis nula no se rechaza.

#Calculo manual de la prueba Shapiro-Wilk

library(dplyr) # Carga la librería para manipulación de datos (pipes, mutate, etc.)
library(gt)    # Carga la librería para la creación de tablas con formato profesional
residuos <-modelo_estimado$residuals # Extrae los residuales del objeto del modelo previamente ajustado
residuos |>   # Toma el vector de residuales y lo pasa al siguiente paso
  as_tibble() |>  # Convierte el vector en una estructura de datos tipo tabla (tibble)
  rename(residuales = value) |>   # Cambia el nombre de la columna por defecto 'value' a 'residuales'
  arrange(residuales) |>   # Ordena los residuales de menor a mayor (requisito de Shapiro-Wilk)
  mutate(pi = (row_number() - 0.375) / (n() + 0.25)) |>   # Calcula la posición de graficación (Blom's plotting position)
  mutate(mi = qnorm(pi, lower.tail = TRUE)) |>   # Calcula los cuantiles teóricos (esperados) de una normal estándar
  mutate(ai = 0) -> tabla_SW # Inicializa la columna 'ai' en cero y guarda el resultado en 'tabla_SW'

m <-sum(tabla_SW$mi^2) # Calcula la suma de los cuadrados de los cuantiles teóricos (m)
n <-nrow(tabla_SW) # Obtiene el tamaño de la muestra original
theta<-1/sqrt(n) # Define una variable auxiliar basada en la raíz inversa de n para las aproximaciones

# Calcula el coeficiente 'ai' para el último dato usando el polinomio de aproximación de Royston
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 coeficiente 'ai' para el penúltimo dato usando la fórmula polinómica correspondiente
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)

tabla_SW$ai[1] <- -tabla_SW$ai[n] # El primer coeficiente es el negativo del último (por simetría)
tabla_SW$ai[2] <- -tabla_SW$ai[n-1] # El segundo coeficiente es el negativo del penúltimo

# Calcula un factor de escala (omega) para los coeficientes restantes basándose en m y los extremos calculados
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)

# Asigna los coeficientes 'ai' intermedios escalándolos por la raíz de omega
tabla_SW$ai[3:(n-2)] <- tabla_SW$mi[3:(n-2)] / sqrt(omega)

# Calcula los productos (ai * residual) y los residuales al cuadrado para el estadístico W final
tabla_SW |>   
  mutate(ai_ui = ai * residuales, ui2 = residuales^2) -> tabla_SW # Actualiza la tabla con los cálculos de los productos y cuadrados
print(tabla_SW)
## # A tibble: 88 × 6
##    residuales      pi    mi     ai ai_ui    ui2
##         <dbl>   <dbl> <dbl>  <dbl> <dbl>  <dbl>
##  1     -120.  0.00708 -2.45 -0.286 34.3  14406.
##  2     -116.  0.0184  -2.09 -0.226 26.1  13342.
##  3     -107.  0.0297  -1.88 -0.202 21.6  11466.
##  4      -91.2 0.0411  -1.74 -0.186 17.0   8325.
##  5      -85.5 0.0524  -1.62 -0.173 14.8   7304.
##  6      -77.2 0.0637  -1.52 -0.163 12.6   5956.
##  7      -74.7 0.0751  -1.44 -0.154 11.5   5580.
##  8      -65.5 0.0864  -1.36 -0.146  9.55  4291.
##  9      -63.7 0.0977  -1.29 -0.138  8.82  4058.
## 10      -62.6 0.109   -1.23 -0.132  8.24  3915.
## # ℹ 78 more rows

##Calculo del estadístico W

W<-(sum(tabla_SW$ai_ui)^2)/sum(tabla_SW$ui2)
print(W)
## [1] 0.9413208

##Cálculo del Wn y su p value

mu<-0.0038915*log(n)^3-0.083751*log(n)^2-0.31082*log(n)-1.5861
sigma<-exp(0.0030302*log(n)^2-0.082676*log(n)-0.4803)
Wn<-(log(1-W)-mu)/sigma
print(Wn)
## [1] 3.241867
p.value<-pnorm(Wn,lower.tail = FALSE)
print(p.value)
## [1] 0.0005937472
library(fastGraph)
shadeDist(Wn,ddist = "dnorm",lower.tail = FALSE)

#Prueba Shapiro Wilk usando libreria stats

salida_SW<-shapiro.test(modelo_estimado$residuals)
print(salida_SW)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo_estimado$residuals
## W = 0.94132, p-value = 0.0005937
Wn_salida<-qnorm(salida_SW$p.value,lower.tail = FALSE)
print(Wn_salida)
## [1] 3.241867
#Conclusión: Dado que 0.0005937(p-value)<0.05 entonces se rechaza la hipótesis nula.

##Grafico de la prueba Shapiro Wilk

library(ggplot2)
library(patchwork) # Para juntar las gráficas
## 
## Adjuntando el paquete: 'patchwork'
## The following object is masked from 'package:MASS':
## 
##     area
# 1. Preparar los datos (usando tu tabla_SW ya creada)
df_grafico <- tabla_SW 

# 2. Gráfico Q-Q (El más importante para Shapiro-Wilk)
# Este gráfico compara tus residuales contra los cuantiles teóricos (mi)
p1 <- ggplot(df_grafico, aes(x = mi, y = residuales)) +
  geom_point(color = "steelblue", alpha = 0.7) +
  geom_smooth(method = "lm", color = "firebrick", se = FALSE) + # Línea de referencia
  labs(title = "Gráfico Q-Q de Residuales",
       subtitle = "Comparación con Distribución Normal",
       x = "Cuantiles Teóricos (mi)",
       y = "Residuales Ordenados") +
  theme_minimal()

# 3. Histograma con Curva Normal
p2 <- ggplot(df_grafico, aes(x = residuales)) +
  geom_histogram(aes(y = ..density..), bins = 15, fill = "gray80", color = "white") +
  stat_function(fun = dnorm, 
                args = list(mean = mean(df_grafico$residuales), sd = sd(df_grafico$residuales)),
                color = "darkblue", size = 1) +
  labs(title = "Distribución de Residuales",
       subtitle = "Histograma vs. Curva Normal",
       x = "Residuales",
       y = "Densidad") +
  theme_minimal()

# Juntar las dos gráficas en una sola imagen
(p1 | p2) + plot_annotation(title = "Diagnóstico de Normalidad (Shapiro-Wilk)")
## `geom_smooth()` using formula = 'y ~ x'