load("data/Lyon.Rdata")
library(spdep)
## Matrice de contiguïté selon le partage d'un nœud (Queen)
nb_queen <- poly2nb(LyonIris, queen = TRUE)
w_queen <- nb2listw(nb_queen, zero.policy = TRUE, style = "W")
# I de Moran sur la variable lden elon l’hypothèse de la normalité
moran.test(LyonIris$Lden, listw=w_queen, zero.policy = TRUE, randomisation = FALSE)
# I de Moran sur la variable lden selon la normalité selon l’hypothèse de la randomisation
moran.test(LyonIris$Lden, listw=w_queen, zero.policy = TRUE, randomisation = TRUE)
# I de Moran sur la variable lden selon des permutations Monte-Carlo
moran.mc(LyonIris$Lden, listw=w_queen, zero.policy = TRUE, nsim = 999)11 Correction des exercices
11.1 Exercices du chapitre 1
11.1.1 Exercice 1
11.1.2 Exercice 2
library(spdep)
load("data/Lyon.Rdata")
# Modèle MCO
formule <- "Lden ~ Pct0_14 + Pct_65 + Pct_Img + Pct_brevet + NivVieMed"
modele_mco <- lm(formule, data = LyonIris)
# Résultats du modèle
summary(modele_mco)
## Matrice de contiguïté (Rook)
nb_rook <- poly2nb(LyonIris, queen = FALSE)
w_rook <- nb2listw(nb_rook, zero.policy = TRUE, style = "W")
# I de Moran sur les résidus du modèle global (MCO)
lm.morantest(modele_mco, w_rook)
# Cartographie des résidus
LyonIris$MCO.Residus <- modele_mco$residuals
tmap_mode("plot")
tm_shape(LyonIris)+
tm_borders(col="gray25", lwd=.5)+
tm_fill(col="MCO.Residus", n = 6, style = "pretty",
legend.format = list(text.separator = "à",
decimal.mark = ","),
midpoint = 0,
palette = "-RdBu",
title = "MCO") +
tm_layout(frame = FALSE)+
tm_scale_bar(breaks = c(0,5))11.2 Exercices du chapitre 2
11.2.1 Exercice 1
11.2.2 Exercice 2
11.2.3 Exercice 3
11.3 Exercices du chapitre 3
11.3.1 Exercice 1
library(sf)
library(spatialreg)
load("data/Lyon.Rdata")
# Matrice de contiguïté selon le partage d'un segment
nb_rook <- poly2nb(LyonIris, queen=FALSE)
w_rook <- nb2listw(nb_rook, zero.policy = TRUE, style = "W")
# Modèles
formule <- "PM25 ~ Pct0_14+Pct_65+Pct_Img+Pct_brevet+NivVieMed"
# Modèles MOC, SLX, SDM, SDEM, SAR, SEM, Manski
MCO <- lm(formule, data = LyonIris)
SLX <- lmSLX(formule, listw = w_rook, data = LyonIris)
Manski <- sacsarlm(formule, listw=w_rook, data = LyonIris, type="sacmixed")
SDM <- lagsarlm(formule, listw = w_rook, data = LyonIris, type = "mixed")
SDEM <- errorsarlm(formule, listw=w_rook, data = LyonIris, etype = 'emixed')
SAR <- lagsarlm(formule,listw=w_rook, data = LyonIris, type = 'lag')
SEM <- errorsarlm(formule, listw=w_rook, data = LyonIris)11.3.2 Exercice 2
summary(lm.LMtests(model = MCO,
listw = w_rook,
test = c("LMlag","LMerr","RLMlag","RLMerr")))11.3.3 Exercice 3
## Valeurs d'AIC et de BIC
AICs <- AIC(MCO, SLX, SDM, SDEM, SAR, SEM)
BICs <- BIC(MCO, SLX, SDM, SDEM, SAR, SEM)
## Autocorrélation spatiale des résidus
IMoran_MCO <- moran.mc(resid(MCO), w_rook, nsim = 999)
IMoran_SLX <- moran.mc(resid(SLX), w_rook, nsim = 999)
IMoran_SAR <- moran.mc(resid(SAR), w_rook, nsim = 999)
IMoran_SEM <- moran.mc(resid(SEM), w_rook, nsim = 999)
IMoran_SDM <- moran.mc(resid(SDM), w_rook, nsim = 999)
IMoran_SDEM <- moran.mc(resid(SDEM), w_rook, nsim = 999)
MoranI_s <- c(IMoran_MCO$statistic, IMoran_SLX$statistic,
IMoran_SAR$statistic, IMoran_SEM$statistic,
IMoran_SDM$statistic, IMoran_SDEM$statistic)
MoranI_p <- c(IMoran_MCO$p.value, IMoran_SLX$p.value,
IMoran_SAR$p.value, IMoran_SEM$p.value,
IMoran_SDM$p.value, IMoran_SDEM$p.value)
## Tableau
Comparaison <- data.frame(Modele = c("MCO", "SLX", "SAR", "SEM", "SDM", "SDEM"),
AIC = round(AICs$AIC,2),
BIC = round(BICs$BIC,2),
dl = AICs$df,
MoranI = round(MoranI_s,3),
MoranIp = round(MoranI_p,3))
Comparaison11.4 Exercices du chapitre 4
11.4.1 Exercice 1
11.4.2 Exercice 2
11.4.3 Exercice 3
11.5 Exercices du chapitre 5
11.5.1 Exercice 1
11.5.2 Exercice 2
11.5.3 Exercice 3
11.6 Exercices du chapitre 6
11.6.1 Exercice 1
library(sf)
library(mgcv)
library(car)
library(DHARMa)
library(mgcViz)
# Chargement du jeu de données sur l'agglomération Lyonnaise
load("data/Lyon.Rdata")
# Vérification de la multicolinéarité (VIF)
vif(lm(PM25 ~ Lden+Pct_65+Pct_Img+Pct_brevet+NivVieMed, data = LyonIris))
# Construction du modèle GLM gaussien
modele_gam1 <- gam(PM25 ~ Pct0_14+Pct_65+Pct_Img+Pct_brevet+NivVieMed,
family = gaussian, data = LyonIris)
# Résidus simulés du modèle gaussien
gam1_res <- simulateResiduals(modele_gam1, plot = FALSE)
plot(gam1_res)
# Construction du modèle GLM avec une distribution de Student
modele_gam2 <- gam(PM25 ~ Pct0_14+Pct_65+Pct_Img+Pct_brevet+NivVieMed,
family = scat, data = LyonIris)
# Résidus simulés avec une distribution de Student
gam2_res <- simulateResiduals(modele_gam2, plot = FALSE)
plot(gam2_res)
# Comparaison des deux modèles
AIC(modele_gam1, modele_gam2)
# Résultats du modèle gaussien
summary(modele_gam2)
# Matrice de contiguïté selon le partage d'un nœud
queen_nb <- poly2nb(LyonIris, queen = TRUE)
queen_w <- nb2listw(queen_nb, style = 'W', zero.policy = TRUE)
# I de Moran sur les résidus gaussiens
moran.mc(residuals(modele_gam2,
type = 'pearson'),
listw = queen_w, nsim = 999,
zero.policy = TRUE)
# I de Moran sur les résidus simulés
moran.mc(residuals(gam2_res,
type = 'pearson'),
listw = queen_w, nsim = 999,
zero.policy = TRUE)11.6.2 Exercice 2
# Chargement du jeu de données sur l'agglomération Lyonnaise
load("data/Lyon.Rdata")
# Ajout des coordonnées x et y dans la couche sf LyonIris
XY <- st_coordinates(st_centroid(LyonIris))
LyonIris$X <- XY[,1]
LyonIris$Y <- XY[,2]
# Modèle GAM avec la distribution de student et une spline
# sur les coordonnées x et y avec la distribution de Student
modele_gam_spline <- gam(NO2 ~ Pct0_14+Pct_65+Pct_Img+Pct_brevet+NivVieMed+
s(X,Y, bs = 'gp', m = 3),
data = LyonIris, family = scat)
summary(modele_gam_spline)11.7 Exercices du chapitre 7
11.7.1 Exercice 1
11.7.2 Exercice 2
11.7.3 Exercice 3
11.8 Exercices du chapitre 8
11.8.1 Exercice 1
library(sf)
library(spgwr)
load("data/Lyon.Rdata")
# Ajout des coordonnées x et y
xy <- st_coordinates(st_centroid(LyonIris))
LyonIris$X <- xy[,1]
LyonIris$Y <- xy[,2]
# Optimisation du nombre de voisins avec le CV
formule <- "PM25 ~ Pct0_14+Pct_65+Pct_Img+Pct_brevet+NivVieMed"
bwaCV_voisins <- gwr.sel(formule,
data = LyonIris,
method = "cv",
gweight=gwr.bisquare,
adapt=TRUE,
verbose = FALSE,
RMSE = TRUE,
longlat = FALSE,
coords=cbind(LyonIris$X,LyonIris$Y))
# Optimisation du nombre de voisins avec l'AIC
formule <- "PM25 ~ Pct0_14+Pct_65+Pct_Img+Pct_brevet+NivVieMed"
bwaCV_voisins <- gwr.sel(formule,
data = LyonIris,
method = "CV",
gweight=gwr.bisquare,
adapt=TRUE,
verbose = FALSE,
RMSE = TRUE,
longlat = FALSE,
coords=cbind(LyonIris$X,LyonIris$Y))
# Réalisation de la GWR
modele_gwr <- gwr(formule,
data = LyonIris,
adapt=bwaCV_voisins,
gweight=gwr.bisquare,
hatmatrix=TRUE,
se.fit=TRUE,
coords=cbind(LyonIris$X,LyonIris$Y),
longlat=FALSE)
# Affichage des résultats
modele_gwr11.8.2 Exercice 2
# Modèle MCO
modele_global <- lm(formule, data = LyonIris)
# Comparaison des R2
r2_global <- summary(modele_global)$r.squared
rss <- sum((modele_gwr$lm$y - modele_gwr$SDF$pred)^2)
tss <- sum((modele_gwr$lm$y - mean(modele_gwr$SDF$pred))^2)
R2_GWRquasiglobal <- 1 - (rss / tss)
cat("R2 global (LM) = ", round(r2_global, 3),
"\nR2 quasi-global (GWR) : ", round(R2_GWRquasiglobal, 3))
# Comparaison des R2
anova(modele_gwr)
LMZ.F1GWR.test(modele_gwr)
LMZ.F2GWR.test(modele_gwr)
# Les coefficients du modèle GWR varient-ils significativement dans l'espace?
LMZ.F3GWR.test(modele_gwr)11.8.3 Exercice 3
library(tmap)
# Récupération du R carré locaux et valeurs de t locales
LyonIris$GWR.R2 <- modele_gwr$SDF$localR2
vars_indep <- c("Pct0_14", "Pct_65", "Pct_Img", "Pct_brevet", "NivVieMed")
for(e in vars_indep){
# Nom des nouvelles variables
var.coef <- paste0("GWR.", "B_", e)
var.t <- paste0("GWR.", "T_", e)
# Récupération des coefficients pour les variables indépendantes
LyonIris[[var.coef]] <- modele_gwr$SDF[[e]]
# Calcul des valeurs de t pour les variables indépendantes
LyonIris[[var.t]] <- modele_gwr$SDF[[e]] / modele_gwr$SDF[[paste0(e, "_se")]]
}
# Cartographie des R2 locaux
legendes_parametres <- list(text.separator = "à",
decimal.mark = ",",
big.mark = " ")
tm_shape(LyonIris)+
tm_borders(col="gray25", lwd=.5)+
tm_fill(col="GWR.R2",
palette="YlOrBr",
n=5, style="quantile",
legend.format = legendes_parametres,
title = "R2 locaux")+
tm_layout(frame = FALSE)+
tm_scale_bar(breaks=c(0,5))
# Cartographie des valeurs de t locales
## Paramètres pour la légende
legende.parametres <- list(text.separator = "à",
text.less.than = "Moins de",
text.or.more = "et plus",
decimal.mark = ",",
big.mark = " ")
## Cartographie
classes_intervalles = c(-Inf, -3.29, -2.58, -1.96, 1.96, 2.58, 3.29, Inf)
carte1 <- tm_shape(LyonIris)+ tm_borders(col="gray25", lwd=.5)+
tm_fill(col="GWR.T_Pct0_14", palette="-RdBu",
midpoint = NA,
breaks = classes_intervalles,
legend.format = legende.parametres,
title = "Moins de 15 ans (%)")+
tm_layout(frame = FALSE, legend.outside = TRUE)
carte2 <- tm_shape(LyonIris)+ tm_borders(col="gray25", lwd=.5)+
tm_fill(col="GWR.T_Pct_65", palette="-RdBu",
midpoint = NA,
breaks = classes_intervalles,
legend.format = legende.parametres,
title = "65 ans et plus (%)")+
tm_layout(frame = FALSE, legend.outside = TRUE)
carte3 <- tm_shape(LyonIris)+ tm_borders(col="gray25", lwd=.5)+
tm_fill(col="GWR.T_Pct_Img", palette="-RdBu",
midpoint = NA,
breaks = classes_intervalles,
legend.format = legende.parametres,
title = "Immigrants (%)")+
tm_layout(frame = FALSE, legend.outside = TRUE)
carte4 <- tm_shape(LyonIris)+ tm_borders(col="gray25", lwd=.5)+
tm_fill(col="GWR.B_Pct_brevet", palette="-RdBu",
midpoint = NA,
breaks = classes_intervalles,
legend.format = legende.parametres,
title = "Faible scolarité (%)")+
tm_layout(frame = FALSE, legend.outside = TRUE)
carte5 <- tm_shape(LyonIris)+ tm_borders(col="gray25", lwd=.5)+
tm_fill(col="GWR.T_NivVieMed", palette="-RdBu",
midpoint = NA,
breaks = classes_intervalles,
legend.format = legende.parametres,
title = "Niveau de vie (€1000)")+
tm_layout(frame = FALSE, legend.outside = TRUE)+
tm_scale_bar(breaks=c(0,5))
tmap_arrange(carte1, carte2, carte3, carte4, carte5, ncol = 2, nrow=3)