sociology rural

1 1. Cargar librerías necesarias

if (!require("readxl")) install.packages("readxl")

## Cargando paquete requerido: readxl

if (!require("lmtest")) install.packages("lmtest")

## Cargando paquete requerido: lmtest

## Cargando paquete requerido: zoo

## 
## Adjuntando el paquete: 'zoo'

## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric

if (!require("sandwich")) install.packages("sandwich")

## Cargando paquete requerido: sandwich

## Warning: package 'sandwich' was built under R version 4.6.1

if (!require("ggplot2")) install.packages("ggplot2")

## Cargando paquete requerido: ggplot2

library(ggplot2)
library(readxl)
library(lmtest)
library(sandwich)

2 2. Cargar datos y ajustar modelo MCO clásico

datos <- read_excel("C:/Users/Lenovo/Downloads/EJERCICIO.xlsx")
View(datos)

modelo <- lm(PR ~ AT, data = datos)

3 3. Crear el gráfico de dispersión con línea de tendencia

grafico_dispersion <- ggplot(datos, aes(x = AT, y = PR)) +
  # Puntos de dispersión
  geom_point(color = "navy", size = 3.5, alpha = 0.8) +
  # Línea de tendencia lineal (MCO) con intervalo de confianza
  geom_smooth(method = "lm", formula = y ~ x, color = "firebrick", fill = "firebrick", alpha = 0.15, se = TRUE) +
  # Etiqueta con el año sobre cada punto
  geom_text(aes(label = AÑO), vjust = -0.8, size = 3.5, color = "black") +
  # Etiqueta flotante con la ecuación del modelo y el R^2
  annotate("label", x = 45, y = 51, 
           label = "PR = 57.22 - 0.32 · AT\nR² = 0.9380", 
           fill = "#FFFDD0", color = "black", size = 4, fontface = "bold") +
  # Títulos y etiquetas de los ejes
  labs(
    title = "Diagrama de Dispersión y Línea de Tendencia",
    subtitle = "Influencia del Acceso a Tecnología (AT) sobre la Pobreza Rural (PR)",
    x = "Acceso a Tecnología (AT)",
    y = "Pobreza Rural (PR)"
  ) +
  # Estilo visual limpio
  theme_minimal(base_size = 12) +
  theme(
    plot.title = element_text(face = "bold", size = 14, hjust = 0.5),
    plot.subtitle = element_text(size = 11, hjust = 0.5, color = "gray30"),
    axis.title = element_text(face = "bold")
  )

4 4. Mostrar el gráfico en pantalla

print(grafico_dispersion)

5 5. MCO Clásico

summary(modelo)

## 
## Call:
## lm(formula = PR ~ AT, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1.1397 -0.8873 -0.2324  0.6234  2.0713 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 57.21870    1.19363   47.94 5.52e-09 ***
## AT          -0.32380    0.03398   -9.53 7.62e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.191 on 6 degrees of freedom
## Multiple R-squared:  0.938,  Adjusted R-squared:  0.9277 
## F-statistic: 90.82 on 1 and 6 DF,  p-value: 7.616e-05

6 6. Inferencia Robusta con Matriz HC3 (para muestras pequeñas)

coeftest(modelo, vcov = vcovHC(modelo, type = "HC3"))

## 
## t test of coefficients:
## 
##              Estimate Std. Error t value  Pr(>|t|)    
## (Intercept) 57.218705   1.707001 33.5200 4.693e-08 ***
## AT          -0.323805   0.047901 -6.7599 0.0005114 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

7 7. Intervalos de confianza robustos HC3 (al 95%)

coefci(modelo, vcov = vcovHC(modelo, type = "HC3"))

##                  2.5 %     97.5 %
## (Intercept) 53.0418246 61.3955847
## AT          -0.4410136 -0.2065956

8 8. Métricas del modelo y Pruebas diagnósticas

r_cuadrado <- summary(modelo)$r.squared
varianza_residual <- summary(modelo)$sigma^2

cat("\n--- MÉTRICAS DEL MODELO ---\n")

## 
## --- MÉTRICAS DEL MODELO ---

cat("R2:", r_cuadrado, "\n")

## R2: 0.9380323

cat("Varianza Residual (MSE):", varianza_residual, "\n\n")

## Varianza Residual (MSE): 1.419421

cat("--- PRUEBA DE NORMALIDAD (Shapiro-Wilk) ---\n")

## --- PRUEBA DE NORMALIDAD (Shapiro-Wilk) ---

print(shapiro.test(residuals(modelo)))

## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(modelo)
## W = 0.89391, p-value = 0.2543

qqnorm(residuals(modelo), main = "3. Gráfico Q-Q Normal de Residuos", pch = 19, col = "purple")
qqline(residuals(modelo), col = "firebrick", lwd = 2)

cat("--- PRUEBA DE HOMOCEDASTICIDAD (Breusch-Pagan) ---\n")

## --- PRUEBA DE HOMOCEDASTICIDAD (Breusch-Pagan) ---

print(bptest(modelo))

## 
##  studentized Breusch-Pagan test
## 
## data:  modelo
## BP = 0.97674, df = 1, p-value = 0.323

matplot(datos$AÑO, cbind(datos$PR, datos$AT), type = "b", pch = 19,
        col = c("firebrick", "navy"), lty = 1, lwd = 2,
        main = "4. Tendencia Temporal (2018–2025)",
        xlab = "Año", ylab = "Porcentaje (%)")
legend("topright", legend = c("Pobreza Rural (PR)", "Acceso Tecnología (AT)"),
       col = c("firebrick", "navy"), lty = 1, pch = 19)