TP3_Statistiques_spatiales
2026-03-04
# Chargement des packages
library(foreign)## Warning: package 'foreign' was built under R version 4.4.1library(sp)
library(sf)## Warning: package 'sf' was built under R version 4.4.1## Linking to GEOS 3.13.0, GDAL 3.8.5, PROJ 9.5.1; sf_use_s2() is TRUElibrary(spdep)## Warning: package 'spdep' was built under R version 4.4.3## Loading required package: spData## Warning: package 'spData' was built under R version 4.4.1## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`library(cartography)## Warning: package 'cartography' was built under R version 4.4.1## This project is in maintenance mode.
## Core functionalities of `cartography` can be found in `mapsf`.
## https://riatelab.github.io/mapsf/library(car)## Warning: package 'car' was built under R version 4.4.1## Loading required package: carData# Lecture du shapefile (les variables du .dbf sont déjà dedans)
south <- st_read("/Users/giacinto/Desktop/Master 2/Stat spatiale/south00-20260304/south00.shp",
quiet = TRUE)
# Option robuste : convertir en 'sp' pour spdep
south_sp <- as(south, "Spatial")
# Voisinage queen
nb_queen <- poly2nb(south_sp, queen = TRUE)
# Poids standardisés en ligne
lw_queen <- nb2listw(nb_queen, style = "W", zero.policy = TRUE)
# petit check
summary(nb_queen)## Neighbour list object:
## Number of regions: 1387
## Number of nonzero links: 7996
## Percentage nonzero weights: 0.4156424
## Average number of links: 5.76496
## Link number distribution:
##
## 1 2 3 4 5 6 7 8 9 10 11
## 3 20 46 150 313 460 294 83 13 4 1
## 3 least connected regions:
## 1050 1296 1347 with 1 link
## 1 most connected region:
## 266 with 11 linksQuestion 2 : Cartographie de la pauvreté infantile
library(classInt)
nclass <- 6
brks <- classIntervals(south$PPOV, n = nclass, style = "jenks")$brks
cols <- hcl.colors(nclass, "Reds")
cl <- cut(south$PPOV, breaks = brks, include.lowest = TRUE)
col_vec <- cols[as.integer(cl)]
plot(st_geometry(south), col = col_vec, border = NA,
main = "Pauvreté infantile (PPOV) — classes Jenks")
legend("topright",
legend = paste0(round(brks[-length(brks)], 3), "–", round(brks[-1], 3)),
fill = cols, bty = "n", cex = 0.8)
#La carte représentant la proportion d’enfants vivant dans la pauvreté (PPOV) met en évidence une structuration spatiale marquée de la pauvreté infantile dans le sud et le sud-est des États-Unis. Les classes ont été construites selon la méthode des « natural breaks » (Jenks), qui permet de regrouper les comtés présentant des niveaux de pauvreté similaires et de mieux faire ressortir les contrastes spatiaux.
#On observe des zones contiguës de forte pauvreté, notamment dans certaines régions du Sud-Est et du delta du Mississippi, ainsi que dans certaines zones des Appalaches. À l’inverse, quelques zones présentent des niveaux de pauvreté plus faibles. Cette organisation spatiale en clusters suggère que les comtés voisins ont tendance à présenter des niveaux de pauvreté similaires.
#Cette première exploration visuelle laisse donc supposer l’existence d’une autocorrélation spatiale positive de la pauvreté infantile, hypothèse qui sera vérifiée par les tests de Moran dans la suite de l’analyse.Question 3 : Tests pour étudier l’auto corrélation
moran_mc <- moran.mc(south$PPOV, lw_queen, nsim = 999)
moran_mc##
## Monte-Carlo simulation of Moran I
##
## data: south$PPOV
## weights: lw_queen
## number of simulations + 1: 1000
##
## statistic = 0.5893, observed rank = 1000, p-value = 0.001
## alternative hypothesis: greatermoran_test <- moran.test(south$PPOV, lw_queen)
moran_test##
## Moran I test under randomisation
##
## data: south$PPOV
## weights: lw_queen
##
## Moran I statistic standard deviate = 36.544, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.5892974827 -0.0007215007 0.0002606743#Le test de Moran permet d’étudier l’existence d’une autocorrélation spatiale globale de la proportion d’enfants vivant dans la pauvreté (PPOV). L’hypothèse nulle suppose l’absence d’autocorrélation spatiale, tandis que l’hypothèse alternative suppose l’existence d’une autocorrélation spatiale positive. Le test par permutation (moran.mc()) donne une statistique de Moran I égale à 0.589 avec une p-value de 0.001, ce qui conduit à rejeter l’hypothèse nulle au seuil de 5 %. Le test de Moran classique (moran.test()), basé sur une approximation asymptotique, confirme ce résultat avec une p-value très faible (< 2.2e−16). On conclut donc à l’existence d’une autocorrélation spatiale positive significative : les comtés proches géographiquement ont tendance à présenter des niveaux similaires de pauvreté infantile, ce qui indique la présence de clusters spatiaux de pauvreté.
# Tests d'association spatiale locale
lisa <- localmoran(south$PPOV, lw_queen)
# Ajouter les résultats au dataframe
south$Ii <- lisa[,1]
south$E.Ii <- lisa[,2]
south$Var.Ii <- lisa[,3]
south$Z.Ii <- lisa[,4]
south$pvalue <- lisa[,5]
south$pvalue_holm <- p.adjust(south$pvalue, method = "holm")
# variable centrée-réduite
pov_std <- as.numeric(scale(south$PPOV))
# Moran plot
moran.plot(pov_std, lw_queen)
# retard spatial
lag_pov <- lag.listw(lw_queen, south$PPOV)
# variable centrée-réduite
pov_std <- scale(south$PPOV)
lag_std <- scale(lag_pov)
quad <- rep(NA, length(pov_std))
quad[pov_std > 0 & lag_std > 0] <- "HH"
quad[pov_std < 0 & lag_std < 0] <- "LL"
quad[pov_std > 0 & lag_std < 0] <- "HL"
quad[pov_std < 0 & lag_std > 0] <- "LH"
south$quad <- quad
plot(south["quad"],
main="Clusters spatiaux locaux (LISA)")
alpha <- 0.05
south$quad_sig <- south$quad
south$quad_sig[south$pvalue_holm > alpha] <- "Non significatif"# palette de couleurs
cols <- c(
"HH" = "red",
"LL" = "blue",
"HL" = "orange",
"LH" = "lightblue",
"Non significatif" = "grey40"
)
plot(st_geometry(south),
col = cols[south$quad_sig],
border = NA,
main = "Clusters spatiaux locaux significatifs (LISA - correction Holm)")
legend("bottomleft",
legend = names(cols),
fill = cols,
bty = "n",
cex = 0.8)
table(south$quad_sig)##
## HH HL LH LL
## 48 1 1 3
## Non significatif
## 1334##Mettre en œuvre des tests d’association spatiale locale (fonction localmoran()). Expliquer pourquoi les p.values des tests doivent être corrigées et les corriger avec la méthode de Holm, que l’on expliquera, en utilisant la fonction p.adjust(). Faire un diagramme de Moran (fonctions scale() et moran.plot()), puis une carte permettant de situer les zones en fonction de leurs caractéristiques (HH, LL, LH, HL). Refaire cette même carte, mais en différenciant les slots non-significatifs (d’après les LISA) des autres. Commentez les sorties ; on étudiera en particulier les zones de non-stationnarité spatiale (ie qui ne suivent pas le processus global) et on statuera sur l’homogénéité de l’autocorrélation spatiale sur le territoire.
vars <- c("PPOV","PHSP","PFHH","PUNEM","PEXTR",
"PBLK","P65UP","METRO","PHSPLUS")
cor_matrix <- cor(st_drop_geometry(south)[,vars], use = "complete.obs")
cor_matrix## PPOV PHSP PFHH PUNEM PEXTR PBLK
## PPOV 1.0000000 0.172589119 0.55797094 0.72300591 0.3267298 0.40717204
## PHSP 0.1725891 1.000000000 -0.17651720 0.11568573 0.4020603 -0.23969955
## PFHH 0.5579709 -0.176517200 1.00000000 0.51592006 -0.2063338 0.87410559
## PUNEM 0.7230059 0.115685733 0.51592006 1.00000000 0.0642663 0.39821180
## PEXTR 0.3267298 0.402060306 -0.20633377 0.06426630 1.0000000 -0.15310572
## PBLK 0.4071720 -0.239699546 0.87410559 0.39821180 -0.1531057 1.00000000
## P65UP 0.1084852 -0.054428212 -0.05519853 -0.08620285 0.2312783 -0.18511442
## METRO -0.4243095 -0.061533645 -0.07996237 -0.25751787 -0.4112264 -0.01392120
## PHSPLUS -0.5301424 0.003830355 -0.08563649 -0.30602581 -0.2889653 -0.05416422
## P65UP METRO PHSPLUS
## PPOV 0.10848516 -0.42430952 -0.530142390
## PHSP -0.05442821 -0.06153365 0.003830355
## PFHH -0.05519853 -0.07996237 -0.085636490
## PUNEM -0.08620285 -0.25751787 -0.306025807
## PEXTR 0.23127827 -0.41122643 -0.288965318
## PBLK -0.18511442 -0.01392120 -0.054164223
## P65UP 1.00000000 -0.37354270 -0.180072534
## METRO -0.37354270 1.00000000 0.444050164
## PHSPLUS -0.18007253 0.44405016 1.000000000round(cor_matrix["PPOV",],3)## PPOV PHSP PFHH PUNEM PEXTR PBLK P65UP METRO PHSPLUS
## 1.000 0.173 0.558 0.723 0.327 0.407 0.108 -0.424 -0.530#L’étude des corrélations entre la variable expliquée PPOV et les variables explicatives montre que la proportion d’enfants vivant dans la pauvreté est fortement corrélée avec la proportion de chômeurs (PUNEM, 0.723) et la proportion de femmes chefs de famille (PFHH, 0.558). Elle est également positivement corrélée avec la proportion d’afro-américains (PBLK, 0.407) et l’emploi dans l’industrie extractive (PEXTR, 0.327). À l’inverse, la variable PHSPLUS (niveau d’éducation) est négativement corrélée avec PPOV (-0.530), ce qui suggère qu’un niveau d’éducation plus élevé est associé à une plus faible pauvreté infantile. Enfin, la variable METRO présente également une corrélation négative (-0.424), indiquant que les comtés métropolitains tendent à avoir des niveaux de pauvreté infantile plus faibles.
model_lm <- lm(PPOV ~ PHSP + PFHH + PUNEM + PEXTR +
PBLK + P65UP + METRO + PHSPLUS,
data = st_drop_geometry(south))
round(summary(model_lm)$coefficients,3)## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.053 0.010 5.527 0.000
## PHSP 0.048 0.010 5.005 0.000
## PFHH 0.748 0.039 19.387 0.000
## PUNEM 1.408 0.060 23.495 0.000
## PEXTR 0.391 0.025 15.904 0.000
## PBLK -0.111 0.015 -7.562 0.000
## P65UP 0.030 0.036 0.846 0.398
## METRO -0.008 0.003 -2.590 0.010
## PHSPLUS -0.243 0.013 -18.946 0.000summary(model_lm)##
## Call:
## lm(formula = PPOV ~ PHSP + PFHH + PUNEM + PEXTR + PBLK + P65UP +
## METRO + PHSPLUS, data = st_drop_geometry(south))
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.194206 -0.025234 -0.001538 0.023895 0.246325
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.053335 0.009651 5.527 3.90e-08 ***
## PHSP 0.048178 0.009625 5.005 6.30e-07 ***
## PFHH 0.748316 0.038598 19.387 < 2e-16 ***
## PUNEM 1.408049 0.059929 23.495 < 2e-16 ***
## PEXTR 0.391433 0.024612 15.904 < 2e-16 ***
## PBLK -0.110581 0.014623 -7.562 7.21e-14 ***
## P65UP 0.030057 0.035549 0.846 0.3980
## METRO -0.007764 0.002998 -2.590 0.0097 **
## PHSPLUS -0.243284 0.012841 -18.946 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.04261 on 1378 degrees of freedom
## Multiple R-squared: 0.7832, Adjusted R-squared: 0.7819
## F-statistic: 622.3 on 8 and 1378 DF, p-value: < 2.2e-16#Le modèle de régression linéaire multiple permet d’expliquer la proportion d’enfants vivant dans la pauvreté (PPOV) à partir des variables socio-économiques considérées. Le modèle est globalement très significatif (p-value < 2.2e-16) et présente un bon pouvoir explicatif avec un R2 d’environ 0.78. Les variables PUNEM (taux de chômage), PFHH (proportion de femmes chefs de famille) et PEXTR (emploi dans l’industrie extractive) ont un effet positif et significatif sur la pauvreté infantile : une augmentation de ces variables est associée à une hausse de la proportion d’enfants pauvres. La variable PHSP présente également un effet positif significatif. À l’inverse, la variable PHSPLUS (niveau d’éducation) a un effet négatif important, indiquant qu’un niveau d’éducation plus élevé est associé à une diminution de la pauvreté infantile. La variable METRO a également un effet négatif, suggérant que les comtés métropolitains présentent en moyenne une pauvreté infantile plus faible. En revanche, la variable P65UP n’est pas significative dans le modèle.
#Question 7 : #D’après les résultats de la question 5, la variable PBLK (proportion d’afro-américains) est positivement corrélée avec la variable PPOV (corrélation ≈ 0.41), ce qui suggère qu’une augmentation de cette proportion est associée à une augmentation de la pauvreté infantile. Cependant, dans le modèle de régression linéaire multiple, le coefficient estimé pour PBLK est négatif et significatif. Ce résultat apparaît donc incohérent avec la corrélation observée précédemment. Cela signifie qu’une fois les autres variables explicatives contrôlées dans le modèle, l’effet estimé de PBLK devient opposé à celui suggéré par la corrélation simple. Cette incohérence suggère la présence d’une multicolinéarité entre les variables explicatives, en particulier entre PBLK et PFHH, ce qui peut déstabiliser l’estimation des coefficients dans la régression.
vars_exp <- c("PHSP","PFHH","PUNEM","PEXTR",
"PBLK","P65UP","METRO","PHSPLUS")
cor_exp <- cor(st_drop_geometry(south)[,vars_exp],
use = "complete.obs")
round(cor_exp,3)## PHSP PFHH PUNEM PEXTR PBLK P65UP METRO PHSPLUS
## PHSP 1.000 -0.177 0.116 0.402 -0.240 -0.054 -0.062 0.004
## PFHH -0.177 1.000 0.516 -0.206 0.874 -0.055 -0.080 -0.086
## PUNEM 0.116 0.516 1.000 0.064 0.398 -0.086 -0.258 -0.306
## PEXTR 0.402 -0.206 0.064 1.000 -0.153 0.231 -0.411 -0.289
## PBLK -0.240 0.874 0.398 -0.153 1.000 -0.185 -0.014 -0.054
## P65UP -0.054 -0.055 -0.086 0.231 -0.185 1.000 -0.374 -0.180
## METRO -0.062 -0.080 -0.258 -0.411 -0.014 -0.374 1.000 0.444
## PHSPLUS 0.004 -0.086 -0.306 -0.289 -0.054 -0.180 0.444 1.000vif(model_lm)## PHSP PFHH PUNEM PEXTR PBLK P65UP METRO PHSPLUS
## 1.413508 5.867381 1.710764 1.699723 5.224766 1.393145 1.615034 1.368849#L’étude de la matrice de corrélation entre les variables explicatives met en évidence une corrélation très forte entre PFHH et PBLK (0.874). Cette relation indique que ces deux variables contiennent une information très similaire. Les résultats du VIF confirment cette multicolinéarité : les variables PFHH (VIF ≈ 5.87) et PBLK (VIF ≈ 5.22) présentent des valeurs supérieures au seuil indicatif de 4, ce qui suggère une instabilité potentielle des coefficients associés dans la régression. Cette forte multicolinéarité explique pourquoi le coefficient estimé pour PBLK devient négatif dans le modèle de régression alors que sa corrélation simple avec PPOV est positive. L’effet de PBLK est en réalité en grande partie capté par la variable PFHH, ce qui perturbe l’interprétation de son coefficient. Pour éviter ce problème d’instabilité, la variable PBLK est retirée du modèle pour la suite de l’analyse.
model_lm2 <- lm(PPOV ~ PHSP + PFHH + PUNEM + PEXTR +
P65UP + METRO + PHSPLUS,
data = st_drop_geometry(south))
summary(model_lm2)##
## Call:
## lm(formula = PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP + METRO +
## PHSPLUS, data = st_drop_geometry(south))
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.205612 -0.024960 -0.001279 0.024781 0.262327
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.075849 0.009365 8.099 1.21e-15 ***
## PHSP 0.067341 0.009473 7.109 1.87e-12 ***
## PFHH 0.497657 0.020175 24.667 < 2e-16 ***
## PUNEM 1.465323 0.060647 24.162 < 2e-16 ***
## PEXTR 0.343934 0.024277 14.167 < 2e-16 ***
## P65UP 0.119585 0.034195 3.497 0.000485 ***
## METRO -0.008858 0.003054 -2.900 0.003792 **
## PHSPLUS -0.243540 0.013100 -18.591 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.04347 on 1379 degrees of freedom
## Multiple R-squared: 0.7742, Adjusted R-squared: 0.7731
## F-statistic: 675.5 on 7 and 1379 DF, p-value: < 2.2e-16plot(model_lm2)
lm.morantest(model_lm2, lw_queen)##
## Global Moran I for regression residuals
##
## data:
## model: lm(formula = PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP + METRO
## + PHSPLUS, data = st_drop_geometry(south))
## weights: lw_queen
##
## Moran I statistic standard deviate = 19.154, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Observed Moran I Expectation Variance
## 0.3044980977 -0.0030597865 0.0002578173south$residus <- residuals(model_lm2)
# résidus
res <- south$residus
# classes
brks <- seq(min(res), max(res), length.out = 6)
# couleurs
cols <- heat.colors(5)
# découpage des résidus
res_cat <- cut(res, breaks = brks, include.lowest = TRUE)
# carte
plot(st_geometry(south),
col = cols[res_cat],
border = NA,
main = "Carte des résidus du modèle OLS")
# légende
legend("bottomleft",
legend = paste(round(brks[-6],3), "-", round(brks[-1],3)),
fill = cols,
bty = "n",
cex = 0.8)
library(lmtest)## Warning: package 'lmtest' was built under R version 4.4.1## Loading required package: zoo## Warning: package 'zoo' was built under R version 4.4.1##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericbptest(model_lm2)##
## studentized Breusch-Pagan test
##
## data: model_lm2
## BP = 123.52, df = 7, p-value < 2.2e-16#Question 10 — Choisir entre modèle LAG et SEM
# Tests du multiplicateur de Lagrange
lm.LMtests(model_lm2, lw_queen, test = "all")## Please update scripts to use lm.RStests in place of lm.LMtests##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP + METRO
## + PHSPLUS, data = st_drop_geometry(south))
## test weights: listw
##
## RSerr = 354.08, df = 1, p-value < 2.2e-16
##
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP + METRO
## + PHSPLUS, data = st_drop_geometry(south))
## test weights: listw
##
## RSlag = 300.88, df = 1, p-value < 2.2e-16
##
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP + METRO
## + PHSPLUS, data = st_drop_geometry(south))
## test weights: listw
##
## adjRSerr = 125.34, df = 1, p-value < 2.2e-16
##
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP + METRO
## + PHSPLUS, data = st_drop_geometry(south))
## test weights: listw
##
## adjRSlag = 72.145, df = 1, p-value < 2.2e-16
##
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP + METRO
## + PHSPLUS, data = st_drop_geometry(south))
## test weights: listw
##
## SARMA = 426.22, df = 2, p-value < 2.2e-16#Question 11 — Estimation du modèle choisi
library(spatialreg)## Warning: package 'spatialreg' was built under R version 4.4.1## Loading required package: Matrix##
## Attaching package: 'spatialreg'## The following objects are masked from 'package:spdep':
##
## get.ClusterOption, get.coresOption, get.mcOption,
## get.VerboseOption, get.ZeroPolicyOption, set.ClusterOption,
## set.coresOption, set.mcOption, set.VerboseOption,
## set.ZeroPolicyOptionmod_lag <- lagsarlm(
PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP + METRO + PHSPLUS,
data = st_drop_geometry(south),
listw = lw_queen
)
summary(mod_lag)##
## Call:lagsarlm(formula = PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP +
## METRO + PHSPLUS, data = st_drop_geometry(south), listw = lw_queen)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.23171972 -0.02208745 0.00013235 0.02086364 0.21366822
##
## Type: lag
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.0294822 0.0086790 3.3969 0.0006814
## PHSP 0.0583932 0.0085150 6.8577 7.000e-12
## PFHH 0.4244052 0.0194139 21.8609 < 2.2e-16
## PUNEM 1.0712445 0.0572199 18.7215 < 2.2e-16
## PEXTR 0.2410238 0.0225978 10.6658 < 2.2e-16
## P65UP 0.1358662 0.0305602 4.4459 8.754e-06
## METRO -0.0061801 0.0027304 -2.2634 0.0236102
## PHSPLUS -0.2078276 0.0119431 -17.4015 < 2.2e-16
##
## Rho: 0.34252, LR test value: 272.16, p-value: < 2.22e-16
## Asymptotic standard error: 0.020138
## z-value: 17.009, p-value: < 2.22e-16
## Wald statistic: 289.3, p-value: < 2.22e-16
##
## Log likelihood: 2521.369 for lag model
## ML residual variance (sigma squared): 0.0015088, (sigma: 0.038843)
## Number of observations: 1387
## Number of parameters estimated: 10
## AIC: -5022.7, (AIC for lm: -4752.6)
## LM test for residual autocorrelation
## test value: 50.273, p-value: 1.338e-12mod_sem <- errorsarlm(
PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP + METRO + PHSPLUS,
data = st_drop_geometry(south),
listw = lw_queen
)
summary(mod_sem)##
## Call:errorsarlm(formula = PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP +
## METRO + PHSPLUS, data = st_drop_geometry(south), listw = lw_queen)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.23731347 -0.02157246 -0.00087547 0.01974823 0.19201355
##
## Type: error
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.0732903 0.0102571 7.1454 8.977e-13
## PHSP 0.0850482 0.0153964 5.5239 3.316e-08
## PFHH 0.5947809 0.0220961 26.9180 < 2.2e-16
## PUNEM 0.9373405 0.0600132 15.6189 < 2.2e-16
## PEXTR 0.3130227 0.0264418 11.8382 < 2.2e-16
## P65UP 0.1384429 0.0371617 3.7254 0.000195
## METRO -0.0034557 0.0027368 -1.2627 0.206708
## PHSPLUS -0.2194900 0.0124924 -17.5699 < 2.2e-16
##
## Lambda: 0.64842, LR test value: 345.75, p-value: < 2.22e-16
## Asymptotic standard error: 0.026464
## z-value: 24.501, p-value: < 2.22e-16
## Wald statistic: 600.32, p-value: < 2.22e-16
##
## Log likelihood: 2558.163 for error model
## ML residual variance (sigma squared): 0.0013299, (sigma: 0.036467)
## Number of observations: 1387
## Number of parameters estimated: 10
## AIC: -5096.3, (AIC for lm: -4752.6)#Les modèles spatiaux LAG et SEM permettent de prendre en compte la dépendance spatiale observée dans les données. Dans le modèle LAG, le paramètre spatial ρ est positif et significatif, indiquant une influence des comtés voisins sur la pauvreté infantile. Dans le modèle SEM, le paramètre λ est également fortement significatif, ce qui montre que la dépendance spatiale se situe principalement dans les erreurs. La comparaison des critères d’ajustement montre que le modèle SEM présente un AIC plus faible (-5096.3 contre -5022.7 pour le LAG) et une log-vraisemblance plus élevée. Le modèle SEM apparaît donc comme le plus adapté pour représenter la dépendance spatiale des données.
#Question 12 — Comparaison des modèles
library(spatialreg)
anova(mod_lag, mod_sem, model_lm2)## Model df AIC logLik Test L.Ratio p-value
## mod_lag 1 10 -5022.7 2521.4 1
## mod_sem 2 10 -5096.3 2558.2 1
## model_lm2 3 9 -4752.6 2385.3 2 345.75 0#La comparaison des modèles par le test du rapport de vraisemblance montre que le modèle SEM améliore significativement l’ajustement par rapport au modèle OLS (L.Ratio = 345.75, p-value ≈ 0). De plus, le critère AIC est plus faible pour le modèle SEM (-5096.3) que pour le modèle LAG (-5022.7) et le modèle OLS (-4752.6). Le modèle SEM apparaît donc comme le modèle le plus adapté pour représenter la dépendance spatiale des données.
#Question 13 — Test d’hétéroscédasticité
library(lmtest)
bptest(mod_sem$residuals ~ mod_sem$fitted.values)##
## studentized Breusch-Pagan test
##
## data: mod_sem$residuals ~ mod_sem$fitted.values
## BP = 58.209, df = 1, p-value = 2.357e-14moran.mc(mod_sem$residuals, lw_queen, nsim = 999)##
## Monte-Carlo simulation of Moran I
##
## data: mod_sem$residuals
## weights: lw_queen
## number of simulations + 1: 1000
##
## statistic = -0.040583, observed rank = 9, p-value = 0.991
## alternative hypothesis: greater#Le test de Breusch-Pagan appliqué aux résidus du modèle SEM est significatif (BP = 58.21, p-value < 0.001), ce qui indique la présence d’hétéroscédasticité dans les résidus. En revanche, le test de Moran appliqué aux résidus n’est pas significatif (p-value ≈ 0.996), ce qui montre qu’il ne subsiste plus d’autocorrélation spatiale dans les erreurs. Le modèle SEM permet donc de corriger l’autocorrélation spatiale observée dans le modèle OLS, mais il ne résout pas totalement les problèmes d’hétéroscédasticité des résidus.
#Question 15 — Modèle Durbin (SDM)
mod_sdm <- lagsarlm(
PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP + METRO + PHSPLUS,
data = st_drop_geometry(south),
listw = lw_queen,
type = "mixed"
)
summary(mod_sdm)##
## Call:lagsarlm(formula = PPOV ~ PHSP + PFHH + PUNEM + PEXTR + P65UP +
## METRO + PHSPLUS, data = st_drop_geometry(south), listw = lw_queen,
## type = "mixed")
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.23936170 -0.01973886 -0.00076802 0.01833189 0.17428113
##
## Type: mixed
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.0367957 0.0122529 3.0030 0.0026732
## PHSP 0.0893706 0.0216480 4.1283 3.654e-05
## PFHH 0.6047231 0.0226218 26.7319 < 2.2e-16
## PUNEM 0.8498939 0.0605702 14.0316 < 2.2e-16
## PEXTR 0.2882838 0.0276158 10.4391 < 2.2e-16
## P65UP 0.1189128 0.0401249 2.9636 0.0030410
## METRO -0.0022143 0.0027139 -0.8159 0.4145487
## PHSPLUS -0.2046735 0.0129416 -15.8152 < 2.2e-16
## lag.PHSP -0.0674699 0.0232130 -2.9066 0.0036544
## lag.PFHH -0.4431783 0.0328930 -13.4733 < 2.2e-16
## lag.PUNEM 0.2217863 0.1087181 2.0400 0.0413491
## lag.PEXTR -0.1379620 0.0396387 -3.4805 0.0005005
## lag.P65UP -0.0836930 0.0564869 -1.4816 0.1384375
## lag.METRO -0.0148854 0.0052524 -2.8340 0.0045966
## lag.PHSPLUS 0.0971944 0.0226904 4.2835 1.840e-05
##
## Rho: 0.51861, LR test value: 238.82, p-value: < 2.22e-16
## Asymptotic standard error: 0.030981
## z-value: 16.74, p-value: < 2.22e-16
## Wald statistic: 280.21, p-value: < 2.22e-16
##
## Log likelihood: 2612.412 for mixed model
## ML residual variance (sigma squared): 0.001279, (sigma: 0.035763)
## Number of observations: 1387
## Number of parameters estimated: 17
## AIC: -5190.8, (AIC for lm: -4954)
## LM test for residual autocorrelation
## test value: 30.461, p-value: 3.4063e-08#Le modèle SDM (Spatial Durbin Model) prend en compte à la fois l’autocorrélation spatiale de la variable dépendante et celle des variables explicatives. Le paramètre spatial ρ est positif et fortement significatif (ρ = 0.519, p-value < 2.2e-16), ce qui confirme l’existence d’une dépendance spatiale importante entre les comtés. Les variables PFHH, PUNEM et PEXTR ont un effet positif significatif sur la pauvreté infantile, tandis que PHSPLUS a un effet négatif, ce qui est cohérent avec les résultats des modèles précédents. Certains retards spatiaux des variables explicatives (par exemple lag.PFHH, lag.PEXTR ou lag.PHSPLUS) sont également significatifs, ce qui indique l’existence d’effets de diffusion spatiale entre les comtés. Enfin, le modèle SDM présente un AIC plus faible (-5190.8) que les modèles précédents, ce qui suggère qu’il offre le meilleur ajustement aux données.
#Question 16 — Comparaison SDM / SEM
anova(mod_sem, mod_sdm)## Model df AIC logLik Test L.Ratio p-value
## mod_sem 1 10 -5096.3 2558.2 1
## mod_sdm 2 17 -5190.8 2612.4 2 108.5 0#La comparaison des modèles SEM et SDM par le test du rapport de vraisemblance montre que le modèle SDM améliore significativement l’ajustement par rapport au modèle SEM (L.Ratio = 108.5, p-value ≈ 0). De plus, le modèle SDM présente un AIC plus faible (-5190.8 contre -5096.3) et une log-vraisemblance plus élevée. Cela indique que l’introduction des interactions spatiales des variables explicatives permet de mieux expliquer la structure spatiale de la pauvreté infantile. Le modèle SDM apparaît donc comme le modèle le plus adapté parmi les modèles testés.
#Question 17 — Impacts du modèle SDM
imp <- impacts(mod_sdm, listw = lw_queen, R = 1000)
summary(imp)## Impact measures (mixed, exact):
## Direct Indirect Total
## PHSP dy/dx 0.086739763 -0.04124511 0.04549465
## PFHH dy/dx 0.588585317 -0.25300495 0.33558037
## PUNEM dy/dx 0.932418704 1.29380481 2.22622351
## PEXTR dy/dx 0.289721777 0.02254475 0.31226653
## P65UP dy/dx 0.116169610 -0.04300684 0.07316277
## METRO dy/dx -0.004211458 -0.03131028 -0.03552174
## PHSPLUS dy/dx -0.205788497 -0.01748015 -0.22326865
## ========================================================
## Simulation results ( variance matrix):
## Direct:
##
## Iterations = 1:1000
## Thinning interval = 1
## Number of chains = 1
## Sample size per chain = 1000
##
## 1. Empirical mean and standard deviation for each variable,
## plus standard error of the mean:
##
## Mean SD Naive SE Time-series SE
## PHSP dy/dx 0.086088 0.020960 0.0006628 0.0006628
## PFHH dy/dx 0.587795 0.021847 0.0006909 0.0006909
## PUNEM dy/dx 0.935216 0.059707 0.0018881 0.0018881
## PEXTR dy/dx 0.289914 0.026643 0.0008425 0.0008425
## P65UP dy/dx 0.113978 0.039575 0.0012515 0.0012911
## METRO dy/dx -0.004263 0.002798 0.0000885 0.0000885
## PHSPLUS dy/dx -0.205195 0.012814 0.0004052 0.0004052
##
## 2. Quantiles for each variable:
##
## 2.5% 25% 50% 75% 97.5%
## PHSP dy/dx 0.046572 0.071356 0.085501 0.099994 0.127103
## PFHH dy/dx 0.547147 0.573381 0.587765 0.602189 0.632132
## PUNEM dy/dx 0.817438 0.895294 0.934439 0.974533 1.052489
## PEXTR dy/dx 0.236680 0.271990 0.289923 0.308219 0.340237
## P65UP dy/dx 0.035867 0.086689 0.113026 0.140465 0.191292
## METRO dy/dx -0.009806 -0.006081 -0.004247 -0.002533 0.001347
## PHSPLUS dy/dx -0.229278 -0.214469 -0.205062 -0.196577 -0.180093
##
## ========================================================
## Indirect:
##
## Iterations = 1:1000
## Thinning interval = 1
## Number of chains = 1
## Sample size per chain = 1000
##
## 1. Empirical mean and standard deviation for each variable,
## plus standard error of the mean:
##
## Mean SD Naive SE Time-series SE
## PHSP dy/dx -0.04141 0.027465 0.0008685 0.0008185
## PFHH dy/dx -0.25490 0.048513 0.0015341 0.0014501
## PUNEM dy/dx 1.29841 0.162064 0.0051249 0.0051249
## PEXTR dy/dx 0.02185 0.061121 0.0019328 0.0019328
## P65UP dy/dx -0.04398 0.088447 0.0027969 0.0025768
## METRO dy/dx -0.03141 0.009799 0.0003099 0.0003099
## PHSPLUS dy/dx -0.01614 0.038795 0.0012268 0.0012268
##
## 2. Quantiles for each variable:
##
## 2.5% 25% 50% 75% 97.5%
## PHSP dy/dx -0.09998 -0.05776 -0.04195 -0.024484 0.01540
## PFHH dy/dx -0.35282 -0.28620 -0.25481 -0.222254 -0.16159
## PUNEM dy/dx 0.97661 1.18526 1.29773 1.406908 1.61569
## PEXTR dy/dx -0.09687 -0.02219 0.02379 0.064815 0.14075
## P65UP dy/dx -0.22295 -0.10392 -0.03964 0.012817 0.13272
## METRO dy/dx -0.05135 -0.03824 -0.03140 -0.024622 -0.01311
## PHSPLUS dy/dx -0.09043 -0.04218 -0.01593 0.009536 0.05934
##
## ========================================================
## Total:
##
## Iterations = 1:1000
## Thinning interval = 1
## Number of chains = 1
## Sample size per chain = 1000
##
## 1. Empirical mean and standard deviation for each variable,
## plus standard error of the mean:
##
## Mean SD Naive SE Time-series SE
## PHSP dy/dx 0.04468 0.01993 0.0006301 0.0005980
## PFHH dy/dx 0.33290 0.04704 0.0014877 0.0014152
## PUNEM dy/dx 2.23363 0.16925 0.0053522 0.0053522
## PEXTR dy/dx 0.31176 0.06332 0.0020023 0.0020023
## P65UP dy/dx 0.07000 0.09037 0.0028579 0.0026995
## METRO dy/dx -0.03568 0.01061 0.0003355 0.0003355
## PHSPLUS dy/dx -0.22133 0.04196 0.0013269 0.0013269
##
## 2. Quantiles for each variable:
##
## 2.5% 25% 50% 75% 97.5%
## PHSP dy/dx 0.007092 0.030485 0.04456 0.05877 0.08368
## PFHH dy/dx 0.234195 0.302487 0.33370 0.36429 0.42100
## PUNEM dy/dx 1.904878 2.119391 2.22982 2.34612 2.56369
## PEXTR dy/dx 0.192431 0.267149 0.31588 0.35516 0.42983
## P65UP dy/dx -0.111405 0.009693 0.07087 0.12683 0.24428
## METRO dy/dx -0.056946 -0.043075 -0.03560 -0.02841 -0.01490
## PHSPLUS dy/dx -0.300120 -0.248903 -0.22171 -0.19296 -0.13953#Les impacts du modèle SDM permettent de distinguer les effets directs, indirects (spillovers spatiaux) et totaux des variables explicatives sur la pauvreté infantile. Les résultats montrent que PUNEM (taux de chômage) a l’effet total le plus important et positif, indiquant qu’une augmentation du chômage accroît fortement la pauvreté infantile, y compris dans les régions voisines. Les variables PFHH et PEXTR ont également des effets totaux positifs significatifs, tandis que PHSPLUS (niveau d’éducation) a un effet total négatif, suggérant qu’un niveau d’éducation plus élevé réduit la pauvreté infantile. Les effets indirects montrent l’existence de spillovers spatiaux, notamment pour PUNEM, ce qui signifie que les conditions socio-économiques d’une région peuvent influencer la pauvreté dans les régions voisines. Ces résultats confirment l’importance de prendre en compte la dimension spatiale dans l’analyse de la pauvreté infantile.
#Question 18 — Vérification des résidus
library(lmtest)
bptest(mod_sdm$residuals ~ mod_sdm$fitted.values)##
## studentized Breusch-Pagan test
##
## data: mod_sdm$residuals ~ mod_sdm$fitted.values
## BP = 49.799, df = 1, p-value = 1.703e-12moran.mc(mod_sdm$residuals, lw_queen, nsim = 999)##
## Monte-Carlo simulation of Moran I
##
## data: mod_sdm$residuals
## weights: lw_queen
## number of simulations + 1: 1000
##
## statistic = -0.028557, observed rank = 45, p-value = 0.955
## alternative hypothesis: greater#Le test de Breusch-Pagan appliqué aux résidus du modèle SDM est significatif (BP = 49.80, p-value < 0.001), ce qui indique la présence d’hétéroscédasticité dans les résidus. En revanche, le test de Moran appliqué aux résidus n’est pas significatif (p-value ≈ 0.966), ce qui montre qu’il ne subsiste plus d’autocorrélation spatiale dans les erreurs. Le modèle SDM parvient donc à capter la dépendance spatiale présente dans les données, même si une certaine hétéroscédasticité persiste dans les résidus.