Cours complet pour apprendre R avec une pratique pour l'analyse de données sociologiques

L'extension rgrs fournit un exemple d'objet de ce type, nommé lyon, et qui contient le contour des neuf arrondissements de cette commune. On peut le charger dans R de la manière suivante :

R> library(rgrs)
R> data(lyon)

Nous avons désormais à notre disposition un objet nommé lyon que nous pouvons tout de suite représenter graphiquement à l'aide de la fonction plot. Le résultat est indiqué figure 8.1.

R> sp::plot(lyon)

Si l'on veut étudier la structure de l'objet lyon, par exemple en effectuant un str(lyon), on se rend vite compte de la complexité de ce type d'objets. En fait lyon est lui-même composé de plusieurs « sous-objets » (slots) accessibles avec l'opérateur @. Nous décrirons les trois principaux :

• lyon@data est un tableau de données dont chaque ligne correspond à un des polygones (c'est-à-dire ici à un arrondissement de Lyon) et qui lui associe un certain nombre de données (identifiant, nom de l'arrondissement, etc.) ;
• lyon@polygons est une liste de polygones définissant les contours de chaque arrondissement ;
• proj4string est un objet décrivant la projection géographique utilisée pour les données spatiales.

On n'accédera en général jamais directement à ces sous-objets, à part éventuellement data, que l'on peut manipuler comme n'importe quel tableau de données :

R> str(lyon@data)
'data.frame':       9 obs. of 2 variables:
$ DepCom : Factor w/ 301 levels "69001","69002",..: 293 294 295 296 297 298 299 300 301
$ Nom_Com: chr  "LYON 1ER" "LYON 2E" "LYON 3E" "LYON 4E" ...
R> lyon@data$Nom_Com
[1] "LYON 1ER"  "LYON 2E"  "LYON 3E"  "LYON 4E"  "LYON 5E"  "LYON 6E"
[7] "LYON 7E"   "LYON 8E"  "LYON 9E"

La conversion et l'import de données de type spatial sont des opérations relativement complexes, notamment du fait qu'elles mettent en jeu des notions de projection pas toujours faciles à comprendre et à maîtriser pour des non-spécialistes(30).

R propose cependant de nombreux outils pour importer des données de différents formats, notamment via les extensions maptools et rgdal.

On n'entrera pas ici dans le détail de ces opérations (que nous ne maîtrisons guère). Mais voici cependant, à titre indicatif, la marche à suivre pour importer dans R des données de l'IGN telles que fournies par le Centre Maurice Halbwachs (contours Iris de différents départements, par exemple).

Les données étant fournies au format MapInfo, la première opération consiste à les convertir en format ESRI Shapefile. Sous Linux cela se fait très facilement grâce à ogr2ogr avec une commande du type :

ogr2ogr -f "ESRI Shapefile" 69_iris.shp 69_iris.mid

Ceci devrait vous générer trois fichiers portant le même nom, mais avec les extensions shp, shx et dbf. L'import dans R peut alors s'effectuer de la manière suivante :

library(maptools)
library(rgdal)
library(foreign)

## Projection d'origine
proj.string <- "+init=epsg:27572 +proj=lcc +lat_1=45.90287723937
         +lat_2=47.69712276063 +lat_0=46.8 +lon_0=2.337229104484
         +x_0=600000 +y_0=2200000 +units=m +pm=greenwich"

## Projection d'arrivée
proj.string.geo <- "+proj=longlat +datum=WGS84"

## Import du fichier
rhone.iris <- readShapePoly("69_iris.shp", proj4string=CRS(proj.string))

## Transformation de la projection
rhone.iris <- spTransform(rhone.iris, CRS(proj.string.geo))

## Conversion des données en Unicode
rhone.iris$Nom_Com <- iconv(rhone.iris$Nom_Com, from="latin1", to="utf8")
rhone.iris$Nom_Iris <- iconv(rhone.iris$Nom_Iris, from="latin1", to="utf8")

## Sauvegarde
save(rhone.iris, file="rhone_iris.rda")

On peut ensuite charger le contenu du fichier rhone_iris.rda à l'aide de la fonction load dans un autre script.

Nous disposons donc, d'un côté, d'un objet spatial représentant les arrondissements de Lyon, et de l'autre d'un tableau de données contenant un extrait du recensement de 1999 pour les communes du Rhône. Si on regarde un peu plus attentivement la structure de ces deux objets :

R> str(lyon@data)
'data.frame':       9 obs. of 2 variables:
$ DepCom : Factor w/ 301 levels "69001","69002",..: 293 294 295 296 297 298 299 300 301
$ Nom_Com: chr  "LYON 1ER" "LYON 2E" "LYON 3E" "LYON 4E" ...
R> head(rp99$code, 20)
 [1] 69001 69002 69003 69004 69005 69006 69007 69008 69009 69010 69012
[12] 69013 69014 69015 69016 69017 69018 69019 69020 69021

on se rend compte que les deux objets peuvent être « joints » grâce à un champ contenant le code INSEE de chaque commune ou arrondissement. Ce champ se nomme DepCom pour l'objet lyon, et code pour l'objet rp99.

Nous avons dès lors tout ce qu'il nous faut pour afficher une première carte, en l'occurrence le taux de chômage par arrondissement en 1999. On utilise pour cela la fonction carte.prop de rgrs. Le code et le résultat sont indiqués figure 8.2.

R> carte.prop(lyon, rp99, "tx.chom", sp.key = "DepCom", data.key = "code")

La fonction carte.prop admet les arguments suivants :

• le premier argument est l'objet de type spatial contenant les données cartographiques, ici lyon ;
• le second argument est le tableau de données contenant les variables à cartographier, ici rp99 ;
• le troisième argument est le nom de la variable à représenter, ici "tx.chom" ;
• l'argument sp.key correspond au nom du champ de jointure dans l'objet spatial ;
• l'argument data.key correspond au nom du champ de jointure dans le tableau de données.

Le nombre de classes de valeurs et les limites de ces classes sont calculés automatiquement. On peut cependant spécifier soit un nombre de classes à l'aide de l'argument nbcuts(31), soit les limites de ces classes avec l'argument at, comme indiqué figures 8.3 et 8.4.

R> carte.prop(lyon, rp99, "tx.chom", sp.key = "DepCom", data.key = "code", nbcuts = 3)

R> carte.prop(lyon, rp99, "tx.chom", sp.key = "DepCom", data.key = "code",
+     at = c(10, 10.5, 11, 11.5, 12, 15))

De nombreuses options sont disponibles pour personnaliser l'affichage de la carte. On pourra citer :

• main titre de la carte ;
• sub sous-titre de la carte ;
• posleg position de la légende, spécifiée sous forme d'une chaîne de caractères : bottomright, topleft, left, center… ;
• diverg si TRUE, indique que la carte représente à la fois des valeurs négatives et positives ;
• palette.pos palette de couleur utilisée pour les valeurs positives ;
• palette.neg palette de couleur utilisée pour les valeurs négatives ;
• palette palette de couleur spécifiée manuellement.

Enfin, tout argument supplémentaire est transmis à la fonction spplot pour l'affichage de la carte.

Les palettes utilisées en argument de palette.pos et palette.neg sont celles définies par l'extension RColorBrewer, elle-même issue du projet Colorbrewer : http://www.colorbrewer.org.

Le site du projet propose notamment un outil qui permet de visualiser et de choisir une palette de manière interactive : http://www.personal.psu.edu/cab38/ColorBrewer/ColorBrewer.html.

Les noms de palettes passés en argument de palette.pos et palette.neg sont les mêmes que ceux utilisés sur ce site.

Un exemple d'utilisation de ces différents paramètres est donné figure 8.5.

R> carte.prop(lyon, rp99, "tx.chom", sp.key = "DepCom", data.key = "code",
+     main = "Taux de chômage1999", sub = "Source : INSEE, RP 1999",
+     palette.pos = "RdPu", posleg = "topright")

Contrairement à la représentation de proportions, qui s'effectue en affectant une couleur à chaque polygone, les effectifs ou les populations sont en général représentés par un symbole de taille variable. Ce type de carte peut être obtenu avec la fonction carte.eff.

Son utilisation est très semblable à celle de carte.prop pour ses arguments principaux. Un exemple utilisant les options par défaut et représentant la population active de chaque arrondissement est donné figure 8.6.

R> carte.eff(lyon, rp99, "pop.act", sp.key = "DepCom", data.key = "code")

Les options nbcuts et at sont disponibles de la même manière que pour carte.prop, mais elles n'agissent ici que sur la présentation de la légende.

D'autres options de personnalisation sont également disponibles. On retrouve les options main, sub et posleg déjà décrites pour carte.prop, ainsi que les options suivantes :

• col.bg couleur des symboles (rouge par défaut) ;
• col.border couleur de la bordure des symboles (blanc par défaut) ;
• cex facteur d'agrandissement des symboles ;
• pch symbole utilisé ;
• plot.polygons si FALSE, on n'affiche que les symboles et pas les polygones.

Un exemple d'utilisation de ces différents paramètres est donné figure 8.7.

R> carte.eff(lyon, rp99, "pop.act", sp.key = "DepCom", data.key = "code",
+     main = "Population active en 1999", sub = "Source : INSEE, RP 1999",
+     pch = 23, cex = 10, col.bg = "blue", col.border = "yellow",
+     posleg = "topright")

La représentation d'une variable qualitative (typiquement le résultat d'une classification) se fait généralement de la même manière que la représentation d'une proportion, mais en utilisant une palette de couleurs contrastées et n'induisant pas de « hiérarchie » entre les objets représentés. Ce type de carte peut être obtenu avec la fonction carte.qual.

Là encore, son utilisation est très semblable à celle de carte.prop et carte.eff. Un exemple utilisant les options par défaut et représentant la population active de chaque arrondissement est donné figure 8.8.

R> rp99$qual <- sample(c("A", "B", "C", "D", "E"), nrow(rp99), replace = TRUE)
R> carte.qual(lyon, rp99, "qual", sp.key = "DepCom", data.key = "code")

Le tableau de données rp99 ne contenant pas de variable de type qualitative, on a simulé les données de manière tout à fait aléatoire.

Là encore, plusieurs options sont disponibles pour personnaliser l'affichage de la carte. Ce sont en fait les mêmes options que celles de carte.prop (main, sub, posleg, palette…), ainsi que l'option palette.qual, qui permet de sélectionner une palette de l'extension RColorBrewer. Un exemple de personnalisation est donné figure 8.9.

R> carte.qual(lyon, rp99, "qual", sp.key = "DepCom", data.key = "code",
+     main = "Catégories d'arrondissements", sub = "Source : aléatoire",
+     posleg = "bottomright", palette.qual = "Set2")

Il peut être intéressant d'ajouter une bordure à une carte soit pour des raisons esthétiques (par exemple, pour délimiter son contour), soit pour mettre en valeur une zone particulière.

Imaginons, par exemple que nous souhaitons mettre en valeur le premier arrondissement de Lyon sur notre carte. Un objet de type spatial, malgré sa complexité, supporte l'indexation d'une manière semblable aux tableaux de données. On peut donc sélectionner un sous-ensemble de notre carte des arrondissements de Lyon de la manière suivante :

R> lyon1 <- lyon[lyon@data$Nom_Com == "LYON 1ER", ]

Cette commande sélectionne le polygone dont le nom est LYON 1ER et le place dans un nouvel objet. Il est alors très simple de rajouter une bordure autour de cet arrondissement en utilisant la fonction plot, comme indiqué figure 8.10.

R> carte.prop(lyon, rp99, "tx.chom", sp.key = "DepCom", data.key = "code",
+     main = "Taux de chômage 1999")
R> lyon1 <- lyon[lyon@data$Nom_Com == "LYON 1ER", ]
R> plot(lyon1, lwd = 5, border = "red", add = TRUE)

Si on souhaite délimiter le contour de notre carte avec une bordure, on doit effectuer une opération supplémentaire, qui consiste à « fusionner » l'ensemble des polygones de notre carte pour n'en garder que le contour. Cette opération peut se faire dans R à l'aide de la fonction unionSpatialPolygons de l'extension maptools. Celle-ci permet de fusionner les polygones d'une carte en fonction des valeurs d'une variable : les polygones ayant la même valeur sont alors regroupés pour n'en former plus qu'un.

Dans notre cas, nous voulons fusionner tous les polygones, nous pouvons donc créer une variable artificielle (ici nommée fusion) contenant la même valeur pour toutes nos zones, puis appliquer la fonction unionSpatialPolygons, ce qui donne :

R> library(maptools)
R> fusion <- rep(1, nrow(lyon@data))
R> lyon.contour <- unionSpatialPolygons(lyon, fusion)

On peut ensuite utiliser l'objet ainsi calculé pour ajouter une bordure globale à notre carte, comme dans la figure 8.11.

R> carte.prop(lyon, rp99, "tx.chom", sp.key = "DepCom", data.key = "code",
+     main = "Taux de chômage 1999")
R> plot(lyon1, lwd = 5, border = "red", add = TRUE)
R> plot(lyon.contour, lwd = 3, border = "black", add = TRUE)

L'ajout de labels est souvent indispensable pour augmenter la lisibilité d'une carte et permettre le repérage des zonages géographiques représentés. La fonction carte.labels est faite pour cela.

Elle accepte comme principaux arguments :

• le nom de l'objet spatial ;
• un vecteur de chaînes de caractères contenant les labels.

On peut spécifier les coordonnées de placement des labels via l'argument coords. Si cet argument vaut NULL (ce qui est le cas par défaut), la position des labels est automatiquement calculée en fonction du polygone auquel il appartient (ce qui n'exclut cependant pas les chevauchements).

La figure 8.12 montre comment ajouter les noms des arrondissements à notre carte de Lyon.

R> carte.prop(lyon, rp99, "tx.chom", sp.key = "DepCom", data.key = "code",
+     main = "Taux de chômage 1999")
R> carte.labels(lyon, lyon@data$Nom_Com)

Plusieurs options sont disponibles pour personnaliser l'affichage des labels, notamment :

• cex facteur d'agrandissement ;
• font style de police de caractères (gras par défaut). Voir la page d'aide de par pour plus de détails ;
• col couleur du texte ;
• outline, outline.decal, outline.col permettent d'ajouter une « bordure » autour des labels.

La figure 8.13 montre comment ajouter la valeur du taux de chômage à notre carte, en leur ajoutant une bordure blanche.

R> carte.prop(lyon, rp99, "tx.chom", sp.key = "DepCom", data.key = "code",
+     main = "Taux de chômage 1999")
R> lyon.tx.chom <- round(rp99$tx.chom[rp99$code %in% lyon@data$DepCom], 1)
R> carte.labels(lyon, lyon.tx.chom, outline = TRUE)