GÉO-VISUALISATION AVEC R

Action nationale de formation
MATE SHS

Cette formation de 3 heures 30 porte sur la visualisation de données géographiques sous R. Y sont abordées les traitements SIG de base, la cartographie thématique (cartes en figurés proportionnels, cartes choroplèthes, cartes de typologie, etc) et des cartes reposant sur des techniques plus avancées comme les cartes sur grille ou les cartes de discontinuités.

Sont ici principalement utilisés les packages sf (manipulation de données spatiales) et cartography (cartographie thématique). Chaque exemple propose une chaîne de traitements depuis le chargement des données, leur mise en forme, jusqu’à la mobilisation de méthodes adaptées permettant de répondre à des questions spatiales.

Les exemples proposés sont réalisés à différentes échelles, du local (région Occitanie) au Monde, en passant par l’Europe.

Pourquoi faire de la cartographie avec R ?

“La science est infaillible, mais les savants se trompent toujours.” (Anatole France, 1889). Ce principe peut également être appliqué à la cartographie. En effet, toute carte est issue d’un processus complexe de choix, de selections, d’opérations statistiques ou géomatiques. Certains auteurs énoncent que les cartes sont subjectives (Brunet). D’autres auteurs disent carrément que les cartes mentent (Monmonnier). Quoi qu’il en soit, toute carte résulte de l’acte créateur des choix de son auteur. Lorsque l’on se situe dans une démarche scientifique, ces choix doivent être traçables, partageables et soumis à la discussion scientifique (ce qui est difficilement faisable quand ces cartes sont réalisées dans un environnement “clic-bouton”).

La réalisation des cartes dans le langage R permet de tracer toutes les opérations nécessaires à une réalisation cartographique de qualité. Réaliser des cartes dans ce langage unique permet, en diffusant le code source en même temps que les cartes, de jouer “cartes sur table”. Cela permet de détailler les choix qui ont été faits et s’exposer à la controverse scientifique. Cela permet aussi de travailler à plusieurs sur une carte, en associant des compétences complémentaires (sémiologie graphique, statistique, géomatique, etc.) et de faciliter la mise à jour de documents déjà rélisés (en rééxecutant un code préalablement réalisé, par exemple).

Au final, l’utilisation de R nécessite un effort non négligeable pour ceux qui ne sont pas habitués à l’univers de la programmation informatique. Mais définitivement, l’investissement n’a que des avantages.

Packages & données

Environnement de travail

Pour exécuter les programmes proposés par la formation, installez / mettez à jour R (version 3.4.4 minimum), R Studio. Lancez ensuite R studio et installez les packages au moyen des commandes suivantes (prévoir au moins 10 minutes).

Les packages

Installation des packages utilisés pendant la formation. Ouvrez R-Studio et tapez :

install.packages("sf")
install.packages("cartography")
install.packages("cartogram")
install.packages("readxl")
install.packages("countrycode")

Optionnel

install.packages("rmapshaper")
install.packages("osmdata")
install.packages("reshape2")
install.packages("eurostat")
install.packages("rnaturalearth")

Version des packages utilisés dans ce document

## [1] "sf (version 0.6.3), cartography (version 2.1.2), countrycode (version 1.0.0), readxl (version 1.1.0), rnaturalearth (version 0.1.0), cartogram (version 0.1.0), reshape2 (version 1.4.3), eurostat (version 3.2.9), rmapshaper (version 0.4.0), osmdata (version 0.0.7.1)"

Packages dédiés à visualisation de données géographiques

  • Le package cartography permet de créer et intégrer des cartes thématiques dans sa chaîne de traitements en R. Il permet des représentations cartographiques tels que les cartes en symboles proportionnels, des cartes choroplèthes, des typologies, des cartes de flux ou des cartes de discontinuités. Il offre également des fonctions qui permettent d’améliorer la réalisation de la carte, comme des palettes de couleur, des éléments d’habillage (échelle, flèche du nord, titre, légende…), d’y rattacher des labels ou d’accéder à des APIs cartographiques.
  • Le package cartogram comme son nom l’indique permet de construire des cartogrames continus ou discontinus.

Packages dédiés à la manipulation de données (spatiales ou non)

  • Le package sf est un package qui permet d’importer, gérer et transformer des données géographiques (gestion des projections, opérations SIG).
  • Le package readxl permet d’importer des fichiers Excel sous R.
  • Le package reshape2 permet de transformer et mettre en forme des données sous R.
  • Le package rmapshaper permet de généraliser le contours d’un fond de carte

Packages fournisseurs de données

  • Le package eurostat permet d’importer des données d’Eurostat, fournisseur de référence pour les territoires européens.
  • Le package rnaturalearth permet d’importer les fonds de carte vectoriels de référence pour le Monde (niveaux État ou infra-nationaux) ainsi que des éléments d’habillage (graticules, traits de côte, etc.)
  • Le package osmdata permet d’importer des données d’OpenStreetMap.

Exo 1 - Bienvenue à Sète

Objectifs

  • Espace d’étude : Hérault, communes
  • Prise en main des fonctions spatiales de R

Téléchargement des données

Téléchargez les données ici : Télécharger

Décompressez le dossier et cliquez sur proj.Rproj. Le fichier exo1.R s’ouvre dans R-Studio.

Commandes de base

Vérifiez votre répertoire de travail courant.

getwd()

Consulter le contenu. Regarder ce qu’il y a dans le répertoire “data”

list.files()
list.files("data")

Import d’une couche SIG dans R avec le package sf. Pacakge basé sur GEOS (Geometry Engine Open Source) et GDAL (Geospatial Data Abstraction Library).

library("sf")
## Linking to GEOS 3.6.2, GDAL 2.2.3, proj.4 4.9.3
communes <- st_read(dsn = "data/communes.shp", stringsAsFactors = F)
## Reading layer `communes' from data source `/home/nlambert/ownCloud/ANF2018 - R geoviz/data/communes.shp' using driver `ESRI Shapefile'
## Simple feature collection with 343 features and 5 fields
## geometry type:  POLYGON
## dimension:      XY
## bbox:           xmin: 662675 ymin: 6234874 xmax: 796333 ymax: 6319348
## epsg (SRID):    NA
## proj4string:    +proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +units=m +no_defs

Voir la table atributaire

communes
## Simple feature collection with 343 features and 5 fields
## geometry type:  POLYGON
## dimension:      XY
## bbox:           xmin: 662675 ymin: 6234874 xmax: 796333 ymax: 6319348
## epsg (SRID):    NA
## proj4string:    +proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +units=m +no_defs
## First 10 features:
##    INSEE_COM                       NOM_COMM           STATUT SUPERFICIE
## 1      34029                        BELARGA   Commune simple        415
## 2      34290 SAINT-VINCENT-DE-BARBEYRARGUES   Commune simple        222
## 3      34032                        BEZIERS   Sous-prfecture       9567
## 4      34082                    COMBAILLAUX   Commune simple        902
## 5      34108                     FRONTIGNAN Chef-lieu canton       3984
## 6      34166                      MONTBLANC   Commune simple       2716
## 7      34066                    CAZEVIEILLE   Commune simple       1619
## 8      34296                      SAUSSINES   Commune simple        630
## 9      34269        SAINT-JEAN-DE-MINERVOIS   Commune simple       3277
## 10     34012                     ARGELLIERS   Commune simple       5065
##    NOM_DEPT                       geometry
## 1   HERAULT POLYGON ((742075 6272005, 7...
## 2   HERAULT POLYGON ((770896 6289405, 7...
## 3   HERAULT POLYGON ((718759 6244261, 7...
## 4   HERAULT POLYGON ((761606 6284489, 7...
## 5   HERAULT POLYGON ((758738 6257679, 7...
## 6   HERAULT POLYGON ((728067 6248191, 7...
## 7   HERAULT POLYGON ((762395 6294296, 7...
## 8   HERAULT POLYGON ((784582 6294234, 7...
## 9   HERAULT POLYGON ((686897 6252592, 6...
## 10  HERAULT POLYGON ((756922 6287661, 7...
head(communes,3)
## Simple feature collection with 3 features and 5 fields
## geometry type:  POLYGON
## dimension:      XY
## bbox:           xmin: 710431 ymin: 6244261 xmax: 771300 ymax: 6292043
## epsg (SRID):    NA
## proj4string:    +proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +units=m +no_defs
##   INSEE_COM                       NOM_COMM         STATUT SUPERFICIE
## 1     34029                        BELARGA Commune simple        415
## 2     34290 SAINT-VINCENT-DE-BARBEYRARGUES Commune simple        222
## 3     34032                        BEZIERS Sous-prfecture       9567
##   NOM_DEPT                       geometry
## 1  HERAULT POLYGON ((742075 6272005, 7...
## 2  HERAULT POLYGON ((770896 6289405, 7...
## 3  HERAULT POLYGON ((718759 6244261, 7...
View(communes)

L’instruction plot permet d’afficher la couche. L’instruction st_geometry permet d’accéder à la variable définissant les géométries.

plot(st_geometry(communes))

Afficher la couche avec des parametres graphiques

plot(st_geometry(communes), col="#aec8f2", border="darkblue", lwd=1)

Ajuster les marges

par(mar = c(0.5,0.5,1.5,0.5)) 
plot(st_geometry(communes), col="#aec8f2", border="darkblue", lwd=1)

Introduction au SIG avec R / Manipuler les informations spatiales

Nous cherchons à identifier les communes qui ont une partie de leur territoire situé à moins de 20 km du centre de Sète.

Nous commençons par extraire la commune de Sète du fond de carte communal de référence.

macommune <- "SETE"
monpoly <- communes[communes$NOM_COMM == macommune,]
par(mar = c(0.5,0.5,1.5,0.5)) 
plot(st_geometry(communes), col="#aec8f2", border="darkblue", lwd=1)
plot(st_geometry(monpoly), col="red", border="purple", lwd=1, add=T)

Calculer le centroïde de la commune de Sète.

moncentre <- st_centroid(x = monpoly)
## Warning in st_centroid.sf(x = monpoly): st_centroid assumes attributes are
## constant over geometries of x
par(mar = c(0.5,0.5,1.5,0.5)) 
plot(st_geometry(communes), col="#aec8f2", border="darkblue", lwd=1)
plot(st_geometry(monpoly), col="red", border="purple", lwd=1, add=T)
plot(st_geometry(moncentre), pch=20, col="black", cex=2, add=T)

Calculer une zone tampon

mydist <- 20000
buff <- st_buffer(x = st_geometry(moncentre), dist=mydist)
par(mar = c(0.5,0.5,1.5,0.5)) 
plot(st_geometry(communes), col="#aec8f2", border="darkblue", lwd=1)
plot(st_geometry(monpoly), col="red", border="purple", lwd=1, add=T)
plot(st_geometry(moncentre), pch=20, col="black", cex=2, add=T)
plot(buff, col=NA, border="black", lty=2,add=T)

Récupérer la liste des communes

communes$buff <- st_intersects(st_geometry(communes), st_geometry(buff), sparse = FALSE)
head(communes)
## Simple feature collection with 6 features and 6 fields
## geometry type:  POLYGON
## dimension:      XY
## bbox:           xmin: 710431 ymin: 6244261 xmax: 771300 ymax: 6292043
## epsg (SRID):    NA
## proj4string:    +proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +units=m +no_defs
##   INSEE_COM                       NOM_COMM           STATUT SUPERFICIE
## 1     34029                        BELARGA   Commune simple        415
## 2     34290 SAINT-VINCENT-DE-BARBEYRARGUES   Commune simple        222
## 3     34032                        BEZIERS   Sous-prfecture       9567
## 4     34082                    COMBAILLAUX   Commune simple        902
## 5     34108                     FRONTIGNAN Chef-lieu canton       3984
## 6     34166                      MONTBLANC   Commune simple       2716
##   NOM_DEPT                       geometry  buff
## 1  HERAULT POLYGON ((742075 6272005, 7...  TRUE
## 2  HERAULT POLYGON ((770896 6289405, 7... FALSE
## 3  HERAULT POLYGON ((718759 6244261, 7... FALSE
## 4  HERAULT POLYGON ((761606 6284489, 7... FALSE
## 5  HERAULT POLYGON ((758738 6257679, 7...  TRUE
## 6  HERAULT POLYGON ((728067 6248191, 7... FALSE

Récupérer le nombre de communes situées à moins de 20 km de Sète.

nb <- dim(communes[communes$buff == T,])
nb[1]
## [1] 40

Extraction et affichage des communes

communes20km <- communes[communes$buff == TRUE,]
par(mar = c(0.5,0.5,1.5,0.5)) 
plot(st_geometry(communes), col="#aec8f2", border="darkblue", lwd=1)
plot(st_geometry(monpoly), col="red", border="purple", lwd=1, add=T)
plot(st_geometry(moncentre), pch=20, col="black", cex=2, add=T)
plot(buff, col=NA, border="black", lty=2,add=T)
plot(st_geometry(communes20km), col=NA, lwd=2, border="red", add=T)