mais aussi sur les variables centrées-réduites. La première option a l'avant-tage de l'unité de mesure (%) qui est directement interprétable. La seconde option permet une comparaison directe des profils à partir de variables aux ordres de grandeur très différents. C'est cette seconde option qui est choisie ici.

Dans un premier temps, la moyenne des variables est calculée avec le fonction `aggregate()` sur chacun des quatre groupes créés par la classification. On cherche ensuite à faire un diagramme en bâtons (barplot) pour représenter graphiquement chacun des groupes. Ceci pourrait être fait avec une boucle qui afficherait successivement les graphiques. Il faudrait également modifier les paramètres graphiques pour que chaque graphique s'intègre dans la même sortie. Dans l'exemple suivant, cette représentation est produite avec le package `ggplot2` présenté dans le Chapitre 9. Ce package demande un tableau en format long (cf. Section 2.4.3), transformation effectuée avec la fonction `melt()`.

```r
clusProfile <- aggregate(muniCSP[, 1:5],
                       by = list(muniCSP$CLUSCSP),
                       mean)
```

```r
colnames(clusProfile)[1] <- "CLUSTER"
clusLong <- melt(clusProfile, id vars = "CLUSTER")
```

Le graphique prend comme argument le tableau qui vient d'être transformé, les trois variables sont l'identification des classes (`CLUSTER`), l'identification des variables (variable) et les valeurs prises par les variables pour chaque classe (value). Pour construire une planche de graphiques, il faut inclure dans une même sortie graphique un ensemble de graphiques correspondant à des sous-groupes. C'est la fonction `facetWrap()` qui permet cette mise en forme (cf. Section 9.3.4).

```r
ggplot(clusLong) +
  geom_bar(aes(x = variable, y = value, fill = CLUSTER),
          stat = "identity") +
  scale_fill_grey() +
  facetWrap(~ CLUSTER) +
  coordflip() + themeBW()
```
<figure><img src="image_49.png" /></figure>

Ces profils peuvent être cartographiés<fnref n="1" /> avec les méthodes présentées dans les chapitres 9 et 10. L'analyse conjointe de la carte et du graphique montre un Ouest parisien accompagné de quelques communes des Hauts-de-Seine (cluster 1) avec une forte sur-représentation des cadres et des artisans-commerçants. À l'opposé, le Nord de la zone d'étude, correspondant à la moitié Est de la Seine-Saint-Denis (cluster 4), montre une forte sur-représentation d'ouvriers et d'employés. Entre ces deux extrêmes, les profils 2 et 3 sont comparables aux profils 1 et 4 respectivement, mais avec des tendances moins marquées.

Certains cas particuliers pourraient être examinés par une analyse plus fine des profils de classes les comparant au profil moyen de l'ensemble de l'espace d'étude. Il s'agirait de tester, pour chaque variable, si la moyenne de la classe est significativement différente de la moyenne générale. Tous ces calculs pourraient être faits avec les fonctions présentées dans les chapitres 4 et 5. Il existe aussi une fonction tout-en-un, la fonction `catdes()` du package `FactoMineR`.

<footnote n="1">1. Attention, la carte en noir et blanc n'est pas satisfaisante d'un point de vue sémiologique, il serait plus pertinent d'utiliser ici une palette qualitative. Voir à ce propos le Chapitre 9 et en particulier la Section 9.2.</footnote>
<figure><img src="image_50.png" /><figcaption>FIGURE 7.4 – Typologie de communes selon les CSP</figcaption></figure>

Cet exemple simple de classification montre aussi l'intérêt de coupler
les méthodes, en particulier l'analyse factorielle et la classification. Il ar-
rive en effet souvent que les variables d'intérêt soient corrélées les unes
avec les autres, cas systématique lorsque les variables sont des pourcen-
tages qui somment 100 % pour chaque individu statistique. Ici, les com-
munes avec une forte proportion de cadres et d'artisans-commerçants au-
ront mécaniquement une faible proportion des autres catégories. Dans un
tel cas, il est judicieux d'effectuer d'abord une analyse factorielle sur les
variables d'intérêt, analyse qui élimine les problèmes de colinéarité, puis
de faire une classification sur les axes factoriels.
## 7.3 Couplage analyse factorielle - classification

### 7.3.1 Classification sur une métrique du $\chi^2$

Il s'agit ici de dresser des profils de communes en fonction de l'évolution de leur population entre 1936 et 2008. Cette classification constitue la suite logique de l'analyse factorielle des correspondances présentée dans la Section 6.2. Elle débute donc de la même façon par une analyse factorielle sur un tableau de contingence. Comme dans la Section 6.2, le tableau de contingence est le résultat du croisement de deux variables qualitatives, ici les communes (lignes du tableau) et les années de recensement (colonnes).

Dans une telle analyse menée sur les modalités de variables qualitatives, il est d'usage de travailler non pas sur une distance euclidienne mais sur une distance du $\chi^2$. En effet, celle-ci est moins sensible au découpage en modalités du phénomène observé, c'est-à-dire à la définition de typologies ou nomenclatures (propriété d'équivalence distributionnelle). Pour un tableau de $n$ lignes et $p$ colonnes, la distance entre les lignes $i$ et $i'$ s'écrit de la façon suivante :

$$d^2(i, i') = \sum_{j=1}^{p} \left( \frac{1}{f_{.j}} \right) \left( \frac{f_{ij}}{f_{i.}} - \frac{f_{i'j}}{f_{i'..}} \right)^2$$

À la différence d'une classification effectuée sur un tableau d'individus caractérisés par des variables, la classification des lignes ou des colonnes du tableau de contingence travaille sur des effectifs. Il est donc utile de prendre ce fait en compte en pondérant les lignes par l'effectif. Ce type de classification sur des distances du $\chi^2$ prenant en compte une variable de pondération est implémenté dans les logiciels classiques comme SAS ou SPAD, mais pas dans les fonctions R utilisées jusqu'à présent. C'est pour cette raison que ce calcul est effectué avec la fonction `ward.cluster()` du package `FactoClass`.

La première étape consiste à faire l'analyse factorielle telle que présentée à la Section 6.2.
```
popEvol <- popCom3608[, 3:11]
row.names(popEvol) <- popCom3608$CODGEO
coaPop <- dudi.coa(popEvol,
                  scannf = FALSE,
                  nf = ncol(popEvol))
```

Ensuite, la matrice de distances du $\chi^2$ est produite grâce à la fonction `dist.dudi()` du package ade4. L'argument `amongrow` précise si la distance doit être calculée entre les lignes du tableau ou entre ses colonnes. Cette matrice est enfin utilisée par la fonction `ward.cluster()`. L'argument `peso`<fnref n="1" /> permet de pondérer les lignes par l'effectif marginal, c'est-à-dire la somme des effectifs de chaque ligne. Cette somme est calculée directement par la fonction `apply()` présentée à la Section 3.4. L'objet qui stocke les résultats de la classification est de type `hclust` et le dendrogramme peut être découpé comme précédemment avec la fonction `cutree()`.

```
library(FactoClass)
chiDist <- dist.dudi(coaPop, amongrow = TRUE)
cahPop <- ward.cluster(chiDist,
                      peso = apply(popEvol, 1, sum),
                      plots = FALSE, h.clust = 1)
popEvol$CLUSTER <- paste("CLUS", cutree(cahPop, k = 4))
```

### 7.3.2 Description des classes

Le tableau de contingence représente ici les trajectoires des communes en termes de peuplement. La description des classes proposée dans cette section se concentre sur l'analyse des groupes de communes pour lesquelles les trajectoires de peuplement sont semblables. Trois types de tableaux sont créés à partir de résumés numériques calculés sur les classes de communes : un tableau qui donne la somme des populations communales des classes à chaque pas de temps, un tableau qui rapporte cette

<footnote n="1">1. Parmi les quelques 5 000 packages existant à la date de rédaction de ce manuel, un certain nombre d'entre eux sont peu attentifs aux noms d'arguments, mélangeant des noms d'arguments en anglais avec des noms dans la langue des auteurs. Il est donc possible de trouver un argument `peso` (poids) dans un package implémenté par des hispanophones, ou un argument `poids` dans un package implémenté par des francophones.</footnote>
somme à la population totale de l'espace d'étude, un tableau qui donne la population communale moyenne pour chaque classe et à chaque pas de temps.

La fonction `aggregate()` est d'abord utilisée pour calculer les résultats numériques selon la variable de regroupement **CLUSTER** (cf. Section 2.4.2). Les effectifs sont divisés par 1000 pour une meilleure lisibilité (population exprimée en milliers). Le calcul des populations en pourcentage du total à chaque date est fait avec la fonction `prop.table()` présentée précédemment (cf. Section 5.4.2). Cette fonction demande une transformation du tableau en objet de type `matrix`, puis une transformation du résultat en objet de type `data.frame`. Le calcul des populations moyennes se fait comme le premier avec `aggregate()` mais en utilisant la fonction `mean()` au lieu de `sum()`.

```r
# Somme des populations (absolue)
sumClass <- aggregate(popEvol[, 1:9],
                     by = list(popEvol$CLUSTER),
                     FUN = sum)
sumClass[, 2:10] <- sumClass[, 2:10] / 1000
colnames(sumClass)[1] <- "CLUSTER"

# Somme des populations (relative)
sumClassMat <- as.matrix(sumClass[, 2:10])
propClass <- as.data.frame(prop.table(sumClassMat, 2))
propClass$CLUSTER <- sumClass$CLUSTER

# Moyenne des populations
meanClass <- aggregate(popEvol[, 1:9],
                      by = list(popEvol$CLUSTER),
                      FUN = mean)
meanClass[, 2:10] <- meanClass[, 2:10] / 1000
colnames(meanClass)[1] <- "CLUSTER"
```

Comme dans la section précédente, il est possible de produire des représentations graphiques avec les fonctions génériques ou avec le package `ggplot2`. C'est cette option qui est choisie et les tableaux doivent donc être transformés en format long avec la fonction `melt()`. Ensuite, une variable temporelle est créée à partir des noms de variables (`POP1936`, `POP1954`, etc.) : les quatre derniers caractères sont extraits avec la fonction `substr()` et ils sont transformés en nombres entiers.
```
absSum <- melt(sumClass, id_vars = "CLUSTER")
relSum <- melt(propClass, id_vars = "CLUSTER")
popMean <- melt(meanClass, id_vars = "CLUSTER")
absSum$YEAR <- as_integer(substr(absSum$variable, 4, 7))
relSum$YEAR <- as_integer(substr(relSum$variable, 4, 7))
popMean$YEAR <- as_integer(substr(popMean$variable, 4, 7))
```

Ces tableaux servent à produire les trois graphiques : somme des populations (valeur absolue), somme des populations (valeur relative) et moyenne des populations par classes.

```
ggplot(absSum) +
  geom_line(aes(x = YEAR, y = value, linetype = CLUSTER)) +
  scale_xcontinuous("Année",
    breaks = seq(1930, 2010, 10)) +
  scale_ycontinuous("Population (milliers)",
    breaks = seq(0, 2500, 500)) +
  themeBW()
```

<figure><img src="image_51.png" /></figure>

```
ggplot(relSum) +
  geom_line(aes(x = YEAR, y = value, linetype = CLUSTER)) +
  scale_xcontinuous("Année",
    breaks = seq(1930, 2010, 10)) +
  scale_ycontinuous("Population (%)",
    labels = percent) +
  themeBW()
```
<figure><img src="image_52.png" /></figure>

```r
ggplot(popMean) +
  geom_line(aes(x = YEAR, y = value, linetype = CLUSTER)) +
  scale_x continuous("Année",
    breaks = seq(1930, 2010, 10)) +
  scale_y continuous("Population (milliers)",
    breaks = seq(0, 100, 20)) +
  themeBW()
```

<figure><img src="image_53.png" /></figure>

#
Ces graphiques renvoient une image assez claire des dynamiques de peuplement en distinguant ces quatre classes de communes : la première classe se caractérise par une forte diminution, particulièrement marquée entre les années 1960 et 1980, les trois autres se caractérisent par une augmentation continue, assez faible pour la classe 2, et forte pour les classes 3 et 4. Une représentation cartographique des classes de communes<fnref n="1" /> est utile pour spatialiser ce constat : les quatre classes s'organisent, à peu de choses près, en cercles concentriques autour de Paris.

<figure><img src="image_54.png" /><figcaption>FIGURE 7.5 – Typologie de communes en fonction des dynamiques de peuplement</figcaption></figure>

<footnote n="1">1. Attention, la carte en noir et blanc n'est pas satisfaisante d'un point de vue sémiologique, il serait plus pertinent d'utiliser ici une palette qualitative. Voir à ce propos le Chapitre 9 et en particulier la Section 9.2.</footnote>
Ce chapitre a présenté plusieurs techniques de classification permettant de créer des groupes différenciés d'individus homogènes. La classification ascendante hiérarchique est une technique fréquemment utilisée en analyse spatiale, elle peut être mobilisée dans un approche descriptive, exploratoire et/ou modélisatrice.
# CHAPITRE 8

## Analyse de graphes

**Objectifs** : Ce chapitre présente une initiation à l'analyse de graphes avec R : importer les données, les transformer pour créer des graphes, calculer des mesures sur les graphes et les visualiser. Les graphes sont des objets mathématiques composés de sommets (nœuds, vertices) et de liens (arcs, edges). Selon les disciplines, les graphes peuvent représenter des réseaux sociaux, des réseaux de lieux, des réseaux de transport, etc.

*

**Prérequis** : Connaissance des fonctions R de base, du vocabulaire de l'analyse de graphes et des principales mesures calculées sur les graphes.
**Description des packages utilisés** Les deux principaux packages d'analyse de graphe sont présentés : *statnet* et *igraph*. Le premier est en réalité un métapackage incluant divers packages spécifiques tels *sna*, *network*, *ergm* etc. Cet ensemble est développé par des sociologues, il utilise des objets de type *network*. Le second, *igraph*, est développé par des physiciens, il utilise des objets de type *igraph*.

Ces deux packages sont très utilisés et très bien documentés avec des sites dédiés<fnref n="1" />. Le choix d'utiliser l'un ou l'autre, ou les deux, dépend des besoins et des préférences de l'utilisateur. Toutes les fonctionnalités de base sont implémentées dans les deux packages, en revanche certaines fonctionnalités spécifiques de l'analyse de réseaux sociaux (*social network analysis*) ne sont implémentées que dans *statnet* alors que d'autres, plutôt issues de la physique, ne sont implémentées que dans *igraph*.

Il est déconseillé d'utiliser les deux packages en même temps. En effet, il y a un certaine redondance entre les deux packages qui peut être source d'erreur : certaines fonctions différentes de chacun des deux packages ont des noms identiques (cf. Encart p. 43).

## 8.1 Aperçu et préparation des données

Le fichier *MobResid08.txt* est utilisé ici, il s'agit du fichier de mobilités résidentielles entre communes de la petite couronne produit par l'Insee.

```
mobResid <- read.csv2("data/MobResid08.txt",
    stringsAsFactors = FALSE)
```

Le fichier de départ liste les flux résidentiels entre les communes de Paris et petite couronne entre 2003 et 2008. Il comprend 5 colonnes et 15 178 lignes, ses variables sont décrites à la Section 1.6. Il faut d'ores et déjà noter que le tableau est en format long : une ligne du tableau est un couple de lieux assorti d'un flux. Avec les techniques présentées à la

<footnote n="1">1. Le site http://statnet.csde.washington.edu/index.shtml pour statnet et le site http://igraph.org pour igraph.</footnote>
Section 2.4.3, il serait facile de passer à un format large, qui serait ici une matrice origine-destination.

Dans le langage de l'analyse de graphe, le format long est considéré comme une liste de liens (*edges list*) et le format large comme une matrice d'adjacence (*adjacency matrix*). L'exemple ci-dessous montre trois représentations d'un même graphe : une liste de liens, une matrice d'adjacence et une visualisation du graphe correspondant. Ce graphe de trois nœuds a les mêmes caractéristiques que le graphe des communes de l'agglomération parisienne sur lequel porte ce chapitre : c'est un graphe orienté et valué. Il est orienté parce qu'un flux de *i* vers *j* n'est pas équivalent à un flux de *j* vers *i*, ce qui se traduit par une matrice d'adjacence non symétrique. Il est valué parce que les liens sont assortis d'un attribut numérique qui associe une quantité au lien, par exemple le nombre d'individus qui migrent d'une commune vers une autre.

<figure><img src="image_55.png" /><figcaption>FIGURE 8.1 – Trois représentations d'un graphe</figcaption></figure>

Les *packages* utilisés ici permettent de créer un graphe à partir des deux types de tableaux : liste de liens ou matrice d'adjacence. Cependant, avant de faire cette transformation du tableau en graphe, il convient de remarquer que bon nombre de mesures pourraient être calculées directement sur les tableaux, sans passer par la transformation en graphe. L'exemple ci-dessous montre un tel traitement qui vise à calculer un solde résidentiel relatif, mesure qui met en relation les flux entrant et sortant d'une commune.
Dans un premier temps, les flux intracommunaux sont supprimés pour simplifier le traitement et son interprétation<fnref n="1" />. Seuls trois champs sont conservés et réordonnés : la commune d'origine (DCRAN), la commune de destination (CODGEO) et le flux (NBFLUX).

Dans un deuxième temps, le tableau long (liste de liens) est transformé en tableau large (matrice d'adjacence) avec le package `reshape2` (cf. Section 2.4.3). Ce tableau est un objet de type `data.frame`, il est transformé en objet de type `matrix` et les codes des communes d'origine sont ajoutés en tant que noms de lignes (avec `row.names()`) et non en tant que variable du tableau. Enfin, les flux entrants (in) et sortants (out) sont calculés grâce à la fonction `apply()` (cf. Section 3.4).

```r
library(reshape2)
edgesList <- mobResid[mobResid$CODGEO != mobResid$DCRAN,
                     c("DCRAN", "CODGEO", "NBFLUX")]
```

```r
adjMat <- dcast(data = edgesList,
                 formula = DCRAN ~ CODGEO,
                 value.var = "NBFLUX")
```

```r
adjMat <- as.matrix(adjMat[, -1])
row.names(adjMat) <- colnames(adjMat)
```

```r
flowsOut <- apply(adjMat, 1, sum, na.rm = TRUE)
flowsIn  <- apply(adjMat, 2, sum, na.rm = TRUE)
summary(flowsOut)
```

```
##      Min. 1st Qu.  Median   Mean 3rd Qu.  Max.
##       168    2370    3850    5940    7010   34200
```

Une fois les flux entrants et sortants calculés, ils sont regroupés dans un nouveau tableau. Dans le vocabulaire de l'analyse de graphes, ce qui vient d'être calculé est le degré des nœuds pondéré par le flux (ce même calcul sera fait sur le graphe à la Section 8.3). Ici la commune qui présente le flux sortant maximum envoie 34 200 individus vers les autres communes de l'espace d'étude. Enfin, les soldes résidentiels sont calculés : solde absolu d'abord, qui est la différence entre les flux entrants et les flux sortants ;

<footnote n="1">1. Les flux intracommunaux sont les flux d'une commune vers elle-même, c'est-à-dire le nombre d'individus qui n'ont pas changé de commune de résidence. Dans le vocabulaire de l'analyse de graphes, le lien qui relie un noeud à lui-même est appelé « boucle » (loop).</footnote>
solde relatif ensuite, qui est un rapport entre le solde absolu et la somme des flux.

```R
flowsSum <- data.frame(CODGEO = names(flowsIn),
                       FLOWSOUT = flowsOut,
                       FLOWSIN = flowsIn,
                       stringsAsFactors = FALSE)
```

```R
flowsSum$ABSDELTA <- flowsSum$FLOWSOUT - flowsSum$FLOWSIN
flowsSum$RELDELTA <- with(flowsSum, ABSDELTA /
                         (FLOWSOUT + flowsSum$FLOWSIN))
```

Ces variables qui mesurent les soldes résidentiels absolu et relatif peuvent être cartographiées avec les techniques présentées dans le Chapitre 10. Voici le résultat obtenu<fnref n="1" /> : sans surprise, les soldes négatifs correspondent aux arrondissements parisiens et à quelques communes limítrophes, les soldes positifs correspondent au reste des communes de la petite couronne. Ce résultat ne concerne que le système de flux entre Paris et la petite couronne, pour une analyse plus complète il faudrait élargir l'analyse à l'ensemble de la région Île-de-France, de la France et/ou du monde.

## 8.2 Création et exploration du graphe

### 8.2.1 Importation et exportation de graphes

Dans l'exemple ci-dessous, un graphe est créé à partir d'un tableau de flux entre des communes. La manipulation du tableau et sa transformation en graphe sont intégrées dans un flux de travail entièrement effectué avec le logiciel R. Cependant, il arrive que l'utilisateur récupère des données sous un format spécifique, propre aux logiciels d'analyse de graphes, ou bien qu'il veuille exporter un graphe créé avec R vers un autre logiciel de ce type, Ucinet ou Gephi par exemple.

Les deux packages présentés ici, `statnet` et `igraph`, disposent de fonctions pour importer (`read`) et exporter (`write`) les graphes dans des

<footnote n="1">1. Attention, la carte en noir et blanc n'est pas satisfaisante d'un point de vue sémiologique, il serait plus pertinent d'utiliser ici une palette divergente. Voir à ce propos le Chapitre 9 et en particulier la Section 9.2.</footnote>
<figure><img src="image_56.png" /><figcaption>FIGURE 8.2 – Solde résidentiel des communes de Paris et petite couronne</figcaption></figure>

formats de données spécifiques. Avec `statnet`, ces fonctions prennent le nom du format à importer, par exemple `read.paj()` pour un format Pajek. Avec `igraph`, les fonctions sont `read.graph()` et `write.graph()` et le format est précisé dans les arguments des fonctions (pajek, gml, graphml, etc.).

En termes de calculs, R n'a rien à envier aux logiciels spécialisés tels que Ucinet ou Gephi. La principale raison pour laquelle un utilisateur pourrait vouloir utiliser l'un de ces logiciels réside dans la visualisation statique et/ou interactive du graphe.
### 8.2.2 Création d'un graphe avec R

La liste de liens peut être utilisée comme argument pour créer un graphe. Cet objet sera de type *network* avec le package *statnet* et de type *igraph* avec le package *igraph*. Avec *statnet*, il faut d'abord transformer le tableau de liens en matrice, puis ajouter l'attribut, ici le flux, une fois le graphe créé.

```r
library(statnet)
edgesMat <- as.matrix edgesList[, 1:2])
residNet <- network (edgesMat, directed = TRUE)
set.edge.attribute (residNet,
  attrname = "NBFLUX",
  value = edgesList$NBFLUX)
```

La transformation du tableau en graphe se fait de façon semblable avec *igraph*. La fonction utilisée prend directement un *data.frame* comme argument.

```r
library(igraph)
residGraph <- graph.data.frame (edgesList, directed = TRUE)
```

Dans la suite de ce chapitre, les principales fonctions sont données pour les deux *packages*, mais les traitements sont faits uniquement avec *igraph*. Ils seraient bien sûr réalisables avec *statnet* en utilisant les fonctions indiquées.

### 8.2.3 Calcul de mesures globales

Pour explorer la structure du graphe, les deux *packages* disposent d'un ensemble de fonctions permettant de calculer des mesures dites « globales », c'est-à-dire portant sur le graphe dans son ensemble. Les plus couramment utilisées sont les suivantes :

- l'ordre : nombre de sommets,

- la taille : nombre de liens,

- la densité : nombre de liens présents sur nombre de liens maximum,

- le nombre de composantes connexes,
- le diamètre : longueur du plus long des plus courts chemins,

- la transitivité : proportion de triangles fermés.

Ces mesures sont particulièrement intéressantes pour comparer des graphes entre eux. Le tableau suivant donne la syntaxe à utiliser dans les packages étudiés ici, pour un graphe nommé g.

<table><thead><tr><td>Mesure</td><td>Avec statnet</td><td>Avec igraph</td></tr></thead><tbody><tr><td>Ordre</td><td>g</td><td>g</td></tr><tr><td>Taille</td><td>g</td><td>g</td></tr><tr><td>Densité</td><td>gden (g)</td><td>graph.density (g)</td></tr><tr><td>Comp. connexes</td><td>components (g)</td><td>clusters (g)</td></tr><tr><td>Diamètre</td><td>?*</td><td>diameter (g)</td></tr><tr><td>Diades</td><td>dyad.census (g)</td><td>dyad.census (g)</td></tr><tr><td>Triades</td><td>triad.census (g)</td><td>triad.census (g)</td></tr><tr><td>Transitivité</td><td>gtrans (g)</td><td>transitivity (g)</td></tr></tbody></table>

<caption>TABLE 8.1 – Mesures globales d’un graphe g</caption>

Une rapide exploration du graphe est nécessaire pour se faire une idée de ses caractéristiques très particulières : un nombre de liens très importants (14 935) par rapport au nombre de nœuds (143), ce qui se traduit par une densité très forte (0,74). En effet le nombre de liens existants est proche du nombre de liens maximum, égal à $n(n - 1)$ dans le cas d’un graphe orienté et non planaire. Le graphe ne comporte qu’une seule composante connexe et son diamètre est très faible (2).

```
residGraph
graph.density(residGraph)
clusters(residGraph)
diameter(residGraph)
```

## 8.3 Mesures locales et manipulation du graphe

Les mesures locales sont les mesures qui renseignent sur chaque sommet ou chaque lien du graphe, les plus courantes sont les suivantes :

<footnote>* Il ne semble pas y avoir de fonction pour calculer directement le diamètre avec statnet. Une solution est de calculer les plus courts chemins avec la fonction geodist() et de choisir le plus long des plus courts chemins.</footnote>
- le degré : nombre de liens incidents pour un sommet. Dans un graphe orienté, on distingue le degré entrant et le degré sortant ;

- le degré pondéré : nombre de liens incidents pondéré par le poids de ces liens ;

- la famille de mesures (centralité, accessibilité) qui s'appuient sur le calcul des plus courts chemins.

Les mesures qui s'appuient sur le calcul des plus courts chemins peuvent être calculées aussi bien pour les nœuds que pour les liens. L'utilisateur peut souvent choisir de calculer ces mesures sur le graphe en prenant en compte une pondération des liens, c'est-à-dire un attribut quantitatif attaché aux liens. Les fonctions permettant de calculer ces différents indicateurs sont listées dans le tableau suivant.

<table><caption>TABLE 8.2 – Mesures locales d'un graphe g</caption><thead><tr><td>Mesure</td><td>Avec statnet</td><td>Avec igraph</td></tr></thead><tbody><tr><td>Degré</td><td>degree(g)</td><td>degree(g)</td></tr><tr><td>Proximité</td><td>closeness(g)</td><td>closeness(g)</td></tr><tr><td>Intermédiarité</td><td>betweenness(g)</td><td>betweenness(g)</td></tr><tr><td>Plus court chemin</td><td>geodist(g)</td><td>shortest paths(g)</td></tr></tbody></table>

Le calcul des degrés pondérés et non pondérés rejoint le calcul fait dans la première section de ce chapitre. Le degré pondéré donne une idée de la masse qui est envoyée ou reçue alors que le degré non pondéré donne une idée de la variété des destinations ou des origines. La commune qui a le flux sortant maximum envoie 34 200 individus vers les autres communes de l'espace d'étude, la commune qui a la plus grande variété de destinations envoie des individus vers 138 autres communes distinctes, c'est-à-dire vers la presque totalité de cet espace.
```
degOut <- degree(residGraph, mode = "out")
wtdDegOut <- graph.strength(residGraph,
                           mode = "out",
                           weights = E(residGraph)$NBFLUX)
summary(degOut)
```

```
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 23.0   93.5  111.0  104.0  124.0  138.0
```

```
summary(wtdDegOut)
```

```
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 168   2370  3850  5940  7010  34200
```

Pour manipuler les attributs du graphe, la syntaxe est différente selon le package utilisé, statnet ou igraph. Avec igraph, il est possible de lire et d'écrire des attributs de nœuds et des liens avec les fonctions V() (V pour vertices) et E() (E pour edges), suivies de l'opérateur $, suivi du nom de l'attribut. C'est cette syntaxe qui vient d'être utilisée pour pondérer les liens par la variable de flux : E(residGraph)$NBFLUX.

Ces fonctions peuvent être utilisées pour faire appel à un attribut du graphe (lire l'attribut) mais aussi pour enrichir le graphe de nouveaux attributs (écrire l'attribut dans l'objet). Dans cet exemple, un attribut est ajouté aux nœuds, le degré sortant. En exécutant l'objet ou avec la fonction summary(), on vérifie que l'ajout a été fait : il y a bien un nouvel attribut nommé DEGOUT, qui caractérise les nœuds (v) et qui est numérique (n).

```
V(residGraph)$DEGOUT <- degree(residGraph, mode = "out")
residGraph
```

```
## IGRAPH DN-- 143 14935 --
## + attr: name (v/c), DEGOUT (v/n), NBFLUX (e/n)
```

Si l'utilisateur enrichit le graphe sur lequel il travaille avec de nouveaux attributs et qu'il veut récupérer l'ensemble pour le traiter sous forme de tableau, la fonction get.data.frame() permet cette transformation.
## 8.4 Cliques et communautés

C'est une pratique fréquente de chercher à l'intérieur d'un graphe des sous-groupes cohérents, c'est-à-dire des sous-graphes fortement connexes : des cliques (sociologues), des communautés (physiciens) ou des *clusters* (informaticiens). De nombreuses méthodes existent pour détecter ces sous-graphes et seules les principales sont évoquées.

Une première famille de méthodes s'appuie sur la notion de clique, terme qui désigne un sous-graphe maximal complet comprenant au moins trois sommets. En d'autres termes, une clique est un ensemble de sommets (au moins 3), entre lesquels tous les liens possibles sont présents (graphe complet) et il n'est pas possible d'ajouter un sommet sans que la propriété précédente ne disparaisse (graphe maximal). Cette définition est rarement opérationnelle, elle peut être trop restrictive (cas d'un réseau de connaissances dont les liens sont ténus) ou pas assez restrictive. C'est le cas ici sur le graphe des mobilités résidentielles dont la densité est très forte. Plusieurs paramètres permettent de moduler la définition de la clique pour l'assouplir ou au contraire la restreindre.

Il y a principalement deux façons d'assouplir la définition de la clique : jouer sur la distance et jouer sur les degrés. Au lieu de considérer un ensemble de sommets entre lesquels tous les liens possibles sont présents (distance de 1), il est possible de considérer une distance plus grande : par exemple, chaque sommet doit être relié à tous les autres avec une distance de 2. Ceci revient, dans un réseau social, à considérer que deux individus non reliés directement mais qui connaissent un même individu tiers sont effectivement reliés à une distance 2 (l'ami d'un ami est mon ami). Le type de cliques formées de cette façon sont appelés des *n-cliques*.

En assouplissant le critère du degré, on peut créer des sous-ensembles dans lesquels tous les sommets sont connectés à *k* sommets de la clique (*k-core*) ou encore des sous-ensembles dans lesquels tous les liens possibles moins *k* sont présents (*k-plex*).

Les fonctions disponibles pour calculer ces sous-ensembles sont : `clique.census()` et `kcores()` pour le package `statnet`; `cliques()` et `graph.coreness()` pour le package `igraph`. Avec `statnet`, la fonction `clique.census()` permet d'obtenir la taille
des différentes cliques présentes, leur composition, le nombre de cliques auxquelles appartient chaque sommet et enfin les co-appartenances de sommets entre les différentes cliques. Avec *igraph*, il y a des fonctions qui renvoient uniquement le nombre de cliques, *cliques.number()* par exemple, et des fonctions qui renvoient l'appartenance des sommets aux cliques, *cliques()* par exemple.

Ces approches sont plus adaptées pour des réseaux sociaux que pour le cas étudié ici. En effet, le graphe utilisé est très dense et le poids des liens (i.e. le flux résidentiel) porte une information très importante qui doit nécessairement être prise en compte dans le calcul de sous-ensembles. Il existe un grand nombre de méthodes pour créer des sous-ensembles cohérents qui s'appuient sur la modularité. Ces méthodes sont plutôt utilisées par les physiciens, on les retrouvera donc dans le package *igraph* plutôt que dans le package *statnet*.

Pour un graphe dans lequel on distingue plusieurs communautés, la modularité est forte quand les liens intracommunautés sont forts et les liens intercommunautés sont faibles. Cette mesure s'applique à des graphes non orientés mais elle peut prendre en compte le poids des liens. Elle est définie comme la différence entre la proportion observée de liens intracommunauté et la proportion que l'on observerait dans un graphe aléatoire conservant la distribution de degrés du graphe originel :

$$Q = \frac{1}{2m} \sum_{ij} \left[ A_{ij} - \frac{w_i w_j}{2m} \right] \delta(c_i, c_j)$$

Où *m* est la somme de la matrice de poids, où $A_{ij}$ est le poids du lien entre *i* et *j* et où la fonction *δ* est une fonction qui renvoie 1 si *i* et *j* font partie de la même communauté et 0 sinon. La modularité peut être utilisée comme simple mesure de la qualité d'une partition, c'est-à-dire pour évaluer *a posteriori* une partition obtenue par une méthode quelconque. Cependant, la modularité peut aussi être utilisée comme méthode de partitionnement dans un algorithme qui cherche à la maximiser pour produire une partition cohérente. C'est sur ce principe que fonctionne l'algorithme utilisé, dit « méthode de Louvain » : il produit une partition du graphe à deux niveaux qui correspond à un optimum de modularité.
Le graphe est d'abord transformé en graphe non orienté, les arguments de la fonction (`mode` et `edge.attr.comb`) précisent que s'il y a un lien de *i* vers *j* et un lien de *j* vers *i*, il faut les confondre et faire la somme des attributs des deux liens. Ainsi, si Créteil envoie 1 000 individus vers Melun et Melun 500 individus vers Créteil, le sens et l'asymétrie de la relation seront perdus : on considérera qu'il y a un seul lien non orienté entre ces deux communes représentant un flux de 1 500 individus.

La fonction `multilevel.community()` renvoie une liste d'éléments qui contient les noms des sommets et leur appartenance aux communautés créées. On accède à ces éléments avec l'opérateur `$` et on les regroupe dans un tableau.

```r
undirGraph <- as.undirected(residGraph,
                           mode = "collapse",
                           edge.attr.comb = sum)

mltcomResid <- multilevel.community(
  undirGraph,
  weights = E(undirGraph)$NBFLUX
)

mltClust <- t(rbind(mltcomResid$memberships,
                   mltcomResid$names))

V(undirGraph)$CLUST <- mltClust[, 1]
```

Là encore il sera intéressant de cartographier l'appartenance des communes à ces communautés avec les techniques présentées dans le Chapitre 10. Le résultat est affiché sur la carte suivante<fnref n="1" />.

La fonction propose deux niveaux d'agrégation, c'est le plus agrégé qui est cartographié ici. Cinq sous-ensembles sont distingués qui correspondent presque exactement au découpage départemental. Seul le département 75 (Paris) est divisé en deux : les arrondissements de l'Ouest forment une communauté avec l'ensemble des communes des Hauts-de-Seine alors que les arrondissements de l'Est forment une communauté avec quelques communes limitrophes : Bagnolet, Montreuil ou Pantin.

<footnote n="1">1. Attention, la carte en noir et blanc n'est pas satisfaisante d'un point de vue sémiologique, il serait plus pertinent d'utiliser ici une palette qualitative. Voir à ce propos le Chapitre 9 et en particulier la Section 9.2.</footnote>
<figure><img src="image_57.png" /><figcaption>FIGURE 8.3 – Communautés de flux résidentiels (2008)</figcaption></figure>

Cette division Est-Ouest de Paris se retrouve dans bon nombre de variables sociales, économiques ou politiques, la plus fameuse étant la division du vote PS-UMP.
### 8.5 Visualisations

Visualiser les graphes est l'un des réflexes les plus courants, et peut-être l'un des plus trompeurs, en analyse de graphes. En effet, il faut être attentif aux paramètres de visualisation (épaisseur des liens, taille des nœuds) et aux principes de l'algorithme utilisé pour la visualisation. Les deux packages proposent des outils et options comparables. Le résultat par défaut est plus lisible avec `statnet` qu'avec `igraph` mais les arguments des fonctions permettent d'obtenir des visualisations personnalisées de grande qualité quel que soit le package utilisé.

La taille et la densité du graphe des mobilités résidentielles le rendent difficile à visualiser directement. Pour explorer sa structure, il faudrait travailler sur des sous-graphes, par exemple sur les communautés créées dans la section précédente, ou encore sur des sélections des liens. Il existe plusieurs façons de faire émerger des structures simples de graphes très denses, la plus simple étant de sélectionner pour chaque nœud le flux maximum envoyé.

Dans un premier temps, on récupère à partir de la matrice origine-destination l'index de colonne du flux maximum envoyé. Pour cela on applique la fonction `which.max()` ligne à ligne avec la fonction `apply()` (cf. Section 3.4). Cette fonction renvoie un vecteur d'indexation : par exemple, si la première valeur de ce vecteur (correspondant à la commune 75101) est 11, ceci signifie que le flux maximum est atteint à la colonne 11 de la matrice. Il suffit ensuite d'utiliser cet index pour récupérer les noms de colonnes qui correspondent aux destinations des flux. La colonne 11 correspond par exemple à la commune 75111 : le flux résidentiel maximum en provenance du 1er arrondissement se dirige vers le 11e arrondissement. Ces liens sont réorganisés dans un tableau et la fonction `graph.data.frame()` transforme ce tableau en graphe.
```r
indexMax <- apply(adjMat, 1, which.max)
destMax <- colnames(adjMat)[indexMax]
maxFlows <- data.frame(ORI = names(indexMax),
                      DES = destMax,
                      stringsAsFactors = FALSE)
maxGraph <- graph.data.frame(maxFlows, directed = TRUE)
maxGraph

## IGRAPH DN-- 143 143 --
## + attr: name (v/c)
```

Le graphe ainsi créé peut être visualisé directement avec la fonction générique `plot()`. La sortie par défaut n'est pas satisfaisante, il faut spécifier un certain nombre d'arguments graphiques pour obtenir un bon résultat, en particulier la couleur et l'épaisseur des liens (`edge.color` et `edge.width`), la couleur et la taille des nœuds (ici proportionnelle au degré), la disposition du graphe (`layout`). Un grand nombre d'arguments sont disponibles aussi bien dans `igraph` que dans `statnet`. Si l'utilisateur utilise beaucoup ces visualisations avec des réglages spécifiques, il peut facilement les encapsuler dans une fonction créée à cet effet et ainsi éviter de réécrire tous ces arguments à chaque fois.

```r
plot(maxGraph,
     edge.color = "grey60",
     edge.width = 2,
     edge.curved = F,
     edge arrow.mode = "-",
     vertex.frame.color = "white",
     vertex.color = "grey30",
     vertex.label = V(maxGraph)$name,
     vertex.label.cex = 0.8,
     vertex.label.color = "black",
     vertex.size = 2 * degree(maxGraph),
     layout = layout.fruchterman.reingold(maxGraph))
```
<figure><img src="image_58.png" /></figure>

La sélection proposée dans cet exemple est rudimentaire mais elle laisse déjà apparaître certaines structures : les grands arrondissements de l'Est parisien (11ᵉ, 18ᵉ, 19ᵉ, 20ᵉ) comme destination principale des flux intra-parisiens ; le 15ᵉ arrondissement qui fait le lien entre les arrondissements de l'Ouest et les communes des Hauts-de-Seine (92) ; certains pôles d'at-traction départementaux apparaissent, comme Saint-Denis (93066) pour le département de Seine-Saint-Denis (93) ou Créteil (94028) pour les communes du Val-de-Marne (94). Apparaissent également de petits sous-systèmes avec des communes peu nombreuses, souvent limitrophes, qui envoient leurs flux principaux les unes vers les autres, par exemple le trio Tremblay (93073)-Vaujours (93074)-Villepinte (93078).

L'argument *layout* prend comme argument le résultat d'algorithmes classiques de visualisation, comme celui de Fruchterman et Reingold ou celui de Kamada et Kawai. Ces algorithmes sont implémentés aussi bien
dans *igraph* que dans *statnet*. Il est cependant possible pour l’utili-
sateur de proposer ses propres coordonnées sous forme d’une matrice de *n*
lignes (nombre de nœuds) et 2 colonnes pour les coordonnées *X* et *Y*. En
donnant comme argument les coordonnées géographiques des centroïdes
des communes, il serait par exemple possible d’afficher le graphe dans le
système de projection choisi.
# CHAPITRE 9

## Focus sur la visualisation graphique

**Objectifs** : Ce chapitre est une mise au point sur la visualisation graphique avec R. Il traite de l'exportation des tableaux et des graphiques, il présente la gestion des couleurs et il introduit, en préambule au chapitre sur la cartographie, la manipulation du package ggplot2.

*

**Avertissement** Ce chapitre sur la visualisation et la couleur est destiné à être imprimé en noir et blanc, le code utilisé est reproduit mais les sorties graphiques sont rarement affichées. C'est également cette contrainte éditoriale qui amène, dans l'ensemble du manuel, à ne pas toujours respecter les règles de sémiologie graphique concernant les couleurs (cf. Section 9.2).
**Prérequis** Types d'objet présentés dans le Chapitre 2 ; graphiques présentés dans le Chapitre 4 et le Chapitre 5. Notions de sémiologie graphique.

**Description des packages utilisés** Plusieurs packages sont proposés pour l'exportation des tableaux et la création dynamique de contenus : tables pour exporter des tableaux dans des tableurs, xtable pour exporter les tableaux au format LATEX et knitr pour produire des documents dynamiques.

La représentation graphique demande également une manipulation de palettes de couleurs. Il existe aussi des packages spécialisés dans ce domaine, en particulier RColorBrewer.

Plusieurs packages spécialisés dans la représentation graphique sont disponibles. Les deux principaux packages graphiques généralistes sont lattice et ggplot2, d'autres proposent des fonctions spécifiques, comme le package vcd pour la visualisation de variables qualitatives. C'est ggplot2 qui est présenté ici, pour les raisons détaillées par la suite.

**Description des données utilisées** Les applications présentées dans ce chapitre sont effectuées sur les données historiques contenues dans le package HistData déjà utilisé à la Section 3.4 : les données de « statistique morale » d'André-Michel Guerry et les données sur le choléra de John Snow. Les premières comprennent un ensemble de variables touchant à la « moralité » comme le nombre de suicides, de prostituées ou de désertions dans les départements français vers 1830. Les secondes présentent des données spatialisées sur l'épidémie de choléra à Londres dans les années 1850. Les travaux de Guerry et ceux de Snow sont d'une importance maJEURE dans l'histoire de la statistique et de l'épidémiologie, ce pourquoi leurs jeux de données ont été mis à disposition par Michael Friendly et Stéphane Dray <fnref n="1" />.

<footnote n="1">1. Voir en particulier la page dédiée au package Guerry : http://cran.r-project.org/web/packages/Guerry/index.html.</footnote>
## 9.1 Exportation des tableaux et des graphiques

L'utilisateur peut visualiser des résultats numériques et graphiques dans l'environnement RStudio, mais il a souvent besoin d'exporter ces sorties pour les intégrer dans des documents rédigés. Plusieurs cas sont envisagés dans le tableau qui suit, selon le type de sortie (tableau ou image) et le mode de travail, en particulier le logiciel utilisé pour produire du contenu.

Deux types de flux de travail (*workflow*) peuvent être distingués. Le plus classique consiste à découpler l'écriture du texte et la production des sorties numériques ou graphiques : l'utilisateur produit un graphique avec R, l'exporte dans un format d'image et l'insère dans un logiciel de traitement de texte avec lequel il écrit le contenu. L'autre façon de procéder, de plus en plus développée, s'inscrit dans la vague de la *reproducible research*. Elle consiste à combiner l'écriture des contenus textuels et des traitements numériques et graphiques pour générer un document tout-en-un. Ce flux de travail est souvent qualifié de *dynamic report generation*.

<figure><img src="image_59.png" /><figcaption>FIGURE 9.1 – Exportation de tableaux et graphiques</figcaption></figure>

Si la sortie est un tableau et que le document est rédigé sur un traitement de texte de type WYSIWYG (*What You See Is What You Get*) tel que LibreOffice Writer ou MS Word, il faut nécessairement passer par une exportation dans un tableur avec les fonctions `write.table()` présentées dans le premier chapitre. Le *package tables*, avec ses fonction `tabular()` et `write.table.tabular()`, facilite ce travail.

Si la sortie est un tableau et que le document est rédigé avec un langage à balises – LATEX, markdown ou HTML – le package `xtable` propose une
fonction éponyme `xtable()` qui transforme le tableau en code directement utilisable dans le document rédigé. Dans cet exemple, un tableau de dénombrement des départements selon la variable **Region** est exporté au format LATEX :

```r
library(xtable)
library(HistData)
data(Guerry)
xtable(table(Guerry$Region))
```

Si la sortie est une image (graphique ou carte), elle peut être exportée via l'interface graphique de RStudio. Dans l'onglet **Plots** plusieurs options d'exportation des graphiques sont proposés : il est possible de fixer la taille de l'image et de choisir parmi plusieurs formats courants, par exemple `svg`, `png` ou `pdf`. Pour plus de précision dans les réglages et/ou pour intégrer l'exportation dans le script, il faut utiliser fonctions graphiques qui prennent le nom du format à produire : `png()`, `svg()`, `jpeg()`, etc. Ces fonctions initialisent le dispositif graphique, il doit en-suite être fermé par la fonction `dev.off()`. L'exemple suivant exporte un fichier `png` qui contient l'histogramme de la variable **Crime_pers**.

```r
png(filename = "chemin/fichier.png")
hist(Guerry$Crime_pers)
dev.off()
```

Dans la fonction `png()`, l'utilisateur peut fixer une taille (arguments `width` et `height`) dans une unité de longueur (`units`) comme les pouces ou les centimètres, il peut également préciser la résolution en ppi (`pixels per inch`, `pixels par pouce`) pour obtenir une taille en termes de nombre de pixels. L'utilisateur peut aussi fixer directement la taille en pixels en utilisant ces trois arguments. Enfin, l'argument `pointsize` est une mesure relative de la taille des caractères en fonction de la résolution.

Dans ce flux de travail, la meilleure façon d'obtenir des images prêtes à l'emploi est de les paramétrer le mieux possible avec R, de les exporter au format `svg` et éventuellement de les retoucher sur un logiciel de dessin vectoriel comme Inkscape.
Pour produire un document qui fourmille de traitements et de graphiques, l'arme ultime est la combinaison de l'écriture et des traitements dans un même document. R et RStudio disposent d'outils très pratiques pour cela grâce au package *knitr*<fnref n="1" />. Dans l'interface de RStudio, le menu **File > New file** propose plusieurs solutions : **R Sweave** pour une combinaison avec LATEX (ce manuel est écrit de la sorte), **R Markdown**, **R HTML** et **R Presentation** pour produire des documents écrits en markdown et en HTML. À la date d'écriture de ce manuel, la version de RStudio la plus récente permet d'écrire un document en markdown et de compiler le texte et le code pour obtenir une sortie en pdf, en HTML et en docx (Microsoft Word).

## 9.2 Gestion des couleurs

Sans rentrer dans le détail du codage des couleurs pour l'impression et la visualisation numérique, il faut simplement mentionner que plusieurs conventions existent : le code RGB (*red*, *green*, *blue*), RVB en français) qui indique le mélange des trois couleurs primaires pour obtenir une couleur donnée ; le code CMYK (*cyan*, *magenta*, *yellow*, *key*, CMJN en français) issu du monde de l'imprimerie ; le code hexadécimal issu des standards des navigateurs Internet.

Avec R, tous ces codes peuvent être utilisés en cas de contrainte éditoriale forte, mais il s'agit rarement de la méthode la plus simple. Le plus courant est d'assigner des couleurs prédéfinies. Pour obtenir l'ensemble des couleurs disponibles, il faut exécuter la fonction `colors()`. Cette fonction renvoie une liste de 657 couleurs prédéfinies<fnref n="2" />.

Les couleurs peuvent d'abord être désignées par leur nom :

```r
colors()
hist(Guerry$Crime_pers, col = "red", border = "green")
hist(Guerry$Crime_pers,
       col = "springgreen",
       border = "darksalmon")
```

<footnote n="1">1. Voir le site dédié : http://yihui.name/knitr.</footnote>

<footnote n="2">2. Voir à ce propos la page web d'Earl F. Glynn : http://research.stowers-institute.org/efg/R.</footnote>
Les couleurs peuvent aussi être désignées par leur numéro dans le vecteur de 657 couleurs renvoyé par la fonction `colors()`. Ainsi, les deux lignes suivantes sont équivalentes :

```r
hist(Guerry$Crime_pers,
     col = colors()[26],
     border = colors()[133])

hist(Guerry$Crime_pers, col = "blue", border = "firebrick")
```

Les couleurs peuvent être désignées par leur code hexadécimal :

```r
hist(Guerry$Crime_pers,
     col = "#FF0000",
     border = "#FFAA00")
```

Enfin, R dispose de palettes de couleurs prédéfinies, par exemple `rainbow()` ou `heatcolors()`. Ces fonctions demandent comme argument le nombre de couleurs à extraire de la palette et renvoient le code hexadécimal correspondant :

```r
rbPal <- rainbow(n = 5)
boxplot(Guerry$Crime_pers ~ Guerry$Region, col = rbPal)
```

Certains packages proposent des palettes de couleurs supplémentaires, le plus connu étant `RColorBrewer` qui est une implémentation des palettes de couleurs créés par Cynthia Brewer<fnref n="1" />. Ce package contient des palettes avec des variations de couleurs et de valeurs adaptées à trois grands types de variables statistiques : des variables continues, des variables continues avec une valeur seuil, des variables qualitatives.

Pour les premières, la variation statistique est traduite visuellement par une montée en valeur sur une même couleur, du vert clair au vert foncé par exemple. Pour les deuxièmes, la variation statistique est traduite par une montée en valeur sur deux couleurs distinctes divisées autour d'un seuil. Cette palette dite « divergente » est efficace pour représenter par exemple un solde qui prend des valeurs négatives et positives, les valeurs négatives

<footnote n="1">1. Voir le site dédié : http://colorbrewer2.org.</footnote>
seront traduites dans une couleur froide (bleu par exemple) et les posi-
tives dans une couleur chaude (rouge par exemple). Enfin, les variables
qualitatives sont représentées par une variation de couleur.

Toute carte statistique devrait suivre ces règles sémiologiques simples
qui établissent des liens cohérents entre la variation statistique d'un ca-
ractère et la variation visuelle représentée sur la carte. Ces règles sont
difficiles à respecter lorsqu'on ne dispose que de niveaux de gris, c'est
pourquoi plusieurs cartes produites dans ce manuel sont signalées en note
de bas de page comme ne respectant pas les règles qui viennent d'être
énoncées.

Voici le code pour afficher l'ensemble des palettes disponibles, pour
afficher seulement un certain type de palette (divergente dans ce cas) et
pour créer une palette de *n* couleurs (6 dans ce cas).

```r
library(RColorBrewer)
display.brewer.all()
display.brewer.all(type = "div")
pastelPal <- brewer.pal(n = 6, name = "Pastel2")
```

Le *package* RColorBrewer est utilisé par défaut pour les représentations
graphiques produites avec *ggplot2*, *package* présenté dans la section
qui suit.

9.3 Introduction à ggplot2

9.3.1 Intérêt vis-à-vis des fonctions graphiques classiques

C'est *ggplot2* qui est présenté ici pour plusieurs raisons. D'abord,
il s'agit du *package* graphique le plus polyvalent et le plus utilisé, c'est
d'ailleurs l'un des *packages* les plus téléchargés sur l'ensemble des *pa-
ckages* disponibles sur le CRAN. Ensuite, ce *package* est conçu comme
l'implémentation d'une « grammaire graphique » : il s'appuie donc sur
une base théorique robuste qui systématise le passage des attributs sta-
tistiques d'un tableau aux attributs esthétiques d'une représentation gra-
phique. Enfin, il est développé par Hadley Wickham, l'un des dévelop-
peurs les plus productifs de la communauté R : il a vocation à s'étendre et
à communiquer avec d'autres packages, en particulier dans le domaine de la visualisation interactive et de la cartographie.

Le package ggplot2 est décrit dans la plupart des manuels de visualisation avec R, il est présenté en détail dans un manuel dédié, il dispose enfin d'un site web dédié<fnref n="1" />. L'avantage de ce package est qu'il permet une manipulation fine des attributs esthétiques. Les fonctions graphiques de base, utilisées dans les chapitres précédents, sont simples et rapides: boxplot() et hist() permettent par exemple de faire des boîtes à moustaches et des histogrammes, de modifier quelques paramètres comme la couleur ou la taille. Cependant, leur domaine d'utilisation est limité à ces quelques paramétrages.

Avec ggplot2, il est possible de produire tous types de graphiques, pas seulement les grands classiques. On peut produire de nouveaux types de représentations graphiques ou encore reproduire des graphiques uniques : celui de Minard par exemple (cf. Figure 3.4)<fnref n="2" />, certains graphiques oubliés proposés par Tufte comme le slopegraph par exemple.

### 9.3.2 Aperçu de la syntaxe

La polyvalence de ggplot2 est possible parce qu'il permet un travail fin sur les éléments suivants, constitutifs d'une représentation graphique :

- variable visuelle (mappings) : variation esthétique correspondant à la variation statistique, en particulier la taille et la couleur des éléments;

- géométrie (geoms) : implantation géométrique des objets graphiques, principalement des points, des lignes et des polygones;

- échelle (scales) : variation dans l'espace des valeurs numériques mise à l'échelle dans l'espace des variables visuelles.

- système de coordonnées (coord) : projection des valeurs sur le plan, par exemple avec un système de coordonnées cartésien, polaire ou un système de projection cartographique;

<footnote n="1">1. http://docs GGplot2.org.</footnote>

<footnote n="2">2. Voir à ce propos la page web de Michael Friendly: http://www.datavis.ca/gALLERY/re-minard.php.</footnote>
— transformation statistique (stats) : élément optionnel, il est utile par exemple pour représenter graphiquement des résumés numériques (fréquences absolues ou relatives);

— construction d'une planche (facet) : élément optionnel, il est utile pour construire une planche constituée de plusieurs graphiques correspondant à des sous-groupes de données.

Pour utiliser `ggplot2`, il faut d'abord se familiariser avec une syntaxe et un fonctionnement particuliers. D'abord, les données statistiques à graphier doivent nécessairement être contenues dans un objet de type `data.frame`. Ensuite, à la différence de la plupart des fonctions graphiques de base, `ggplot2` peut créer un objet qui contient les données et les paramétrages graphiques. Pour échanger un graphique classique produit avec `plot()` par exemple, il faut mettre à disposition à la fois le code de la représentation graphique et les données nécessaires. L'objet `ggplot2`, en revanche, est autonome, il peut être stocké et échangé.

Pour créer un tel objet, le package propose une fonction de base `ggplot()` à laquelle on ajoute les différents éléments décrits ci-dessus à l'aide de l'opérateur `+`. Voici un exemple fictif qui produirait un nuage de points. Les données sont stockées dans un `data.frame` nommé `mydata`, ces données sont représentées par des points dont les coordonnées sont fixées par deux variables quantitatives `x` et `y`. Enfin, la taille et la couleur des points peuvent varier en fonctions d'autres variables présentes dans le tableau.

```r
# CAS 1
ggplot(mydata) +
  geom_point(aes(x = , y = , size = , color = )) +
  scale_size(...)
  scale_color(...)
```

```r
# CAS 2
ggplot(mydata) +
  geom_point(aes(x = , y = ), size = , color = )
```

La différence entre le premier et le second cas est le statut de la taille et de la couleur. Dans le premier, les arguments sont inclus dans la fonction `aes()` (aesthetics mapping). La taille et la couleur sont donc considérées comme des variables visuelles : elles vont traduire une variation statistique
par une variation visuelle et cette variation sera mise à l'échelle – échelle de taille et de couleur – par les fonctions de type `scale()`.

Dans le second cas, l'utilisateur fixe une couleur et une taille pour l'ensemble des points : il ne s'agit plus de variables visuelles mais de paramètres graphiques d'ensemble. Dans ce cas, les fonctions de type `scale()` ne sont pas utiles.

### 9.3.3 Variables visuelles

Cet exemple produit un nuage de points mettant en relation la variable **Prostitutes**, qui le nombre de prostituées à Paris par département de naissance, et la variable **Distance**, qui indique la distance du département à Paris.

```r
library(ggplot2)
distProst <- ggplot(Guerry) +
  geom_point(aes(x = Distance, y = Prostitutes),
            size = 3, colour = "black")
distProst
```

<figure><img src="image_60.png" /></figure>

```r
class(distProst)
## [1] "gg" "ggplot"
```
L'objet créé est de type *ggplot*, on peut accéder à son contenu, en particulier aux données, avec l'opérateur `$` comme pour les objets de type *list* et *data.frame*. La taille et la couleur sont des paramètres d'ensemble et non des variables visuelles. Le thème par défaut présente un fond gris clair carroyé de lignes blanches, les lignes principales sont fixées avec un pas de 200 pour l'axe des *x* et un pas de 1 000 pour l'axe des *y*, le titre des axes est le nom de variables représentées. Tous ces paramètres peuvent bien sûr être modifiés.

D'un point de vue thématique, la variable **Prostitutes** indique le nombre de prostituées à Paris par département de naissance. Ce nombre décroit avec la distance à Paris selon une fonction exponentielle négative, résul-tat cohérent avec les modèles d'interaction spatiale développés dans le domaine de la géographie et de la démographie des migrations.

L'exemple qui suit est un nuage de points mettant en relation les crimes contre la personne et les crimes contre la propriété. Ces deux variables sont exprimées en population par crime, c'est-à-dire l'inverse de tous les taux calculés de nos jours, par exemple un taux de criminalité exprimé en crimes par population. Le thème par défaut peut être substitué par un thème plus sobre : chaque élément du thème peut être paramétré, la cou-leur du cadre, l'épaisseur des lignes, etc. mais il est plus facile d'utiliser des thèmes tout faits, par exemple le thème `themeBW()`.

La taille et la couleur sont utilisées comme variables visuelles : la taille traduit le nombre total de crimes, variable calculée au préalable, et la cou-leur traduit la région (Nord, Sud, Est, Ouest, Centre)<fnref n="1" />. La valeur de la Corse est modifiée et la région Sud lui est assignée.

<footnote n="1">1. Attention, le graphique en noir et blanc n'est pas satisfaisant d'un point de vue sémiologique, il serait plus pertinent d'utiliser ici une palette qualitative avec différentes couleurs. Voir à ce propos la Section 9.2.</footnote>
```
Guerry$TotalCrimePers <- with(Guerry,
                                1000 * Pop1831 / Crime_pers)
```

```
Guerry$TotalCrimeProp <- with(Guerry,
                                1000 * Pop1831 / CrimeProp)
```

```
Guerry$TotalCrime <- with(Guerry,
                          TotalCrimePers + TotalCrimeProp)
```

```
Guerry$Region[86] <- "S"
```

```
ggplot(Guerry) +
  geom_point(aes(x = Crime_pers,
                y = CrimeProp,
                size = TotalCrime,
                colour = Region)) +
  scale_sizecontinuous(range = c(3, 15)) +
  scale_color_grey() +
  themeBW()
```

<figure><img src="image_61.png" /></figure>

summary(Guerry$TotalCrime)

```
## Min. 1st Qu.  Median  Mean 3rd Qu.  Max.
```

```
## 20.3    49.2    69.2    84.3    93.7   751.0
```

Le graphique fait apparaître de fortes disparités entre ces deux types de crimes ainsi qu'une relation presque inexistante entre les deux. Les deux variables visuelles renseignent sur le nombre absolu de crimes commis
(taille), le département de Seine en tête avec 751 crimes, et sur les grandes régions correspondant aux départements (couleur).

### 9.3.4 Construction d'une planche

Il est souvent utile de faire construire des planches de graphiques ou de cartes, c'est-à-dire d'inclure dans une même sortie graphique un ensemble de cartes ou de graphiques correspondant à des sous-groupes de données. Ce sont les fonctions `facet_grid` et `facetWrap` qui permettent cela, fonctions dites de « maillage » (latticing). Dans cet exemple, on cherche à produire un nuage de points similaire au précédent, mais en différenciant chaque région.

```r
ggplot(Guerry) +
  geom_point(aes(x = Crime_pers, y = CrimeProp),
            size = 2,
            colour = "black") +
  facetWrap(~Region) +
  themeBW()
```

<figure><img src="image_62.png" /></figure>

```r
tapply(Guerry$TotalCrime, Guerry$Region, max)
```

```
## C    E    N    S    W
```

```
## 82.48 153.69 750.62 131.45 108.56
```
riable de regroupement qui indique quels points doivent être reliés pour
former les polygones.

```r
ggplot(geomEx) +
  geom_polygon(aes(x = coordx,
                 y = coordy,
                 group = polyg),
  fill = "red") +
  theme_bw()
```

Ce travail est maintenant effectué sur de véritables données spatialisées, celles de John Snow qui indiquent la localisation des morts de chloréra recensés en 1854 autour de Broad Street. Le tableau **Snow.deaths** est composé de 578 points avec leurs coordonnées et le tableau **Snow.streets** donne la localisation des rues. Ces rues sont des segments composés de points : les champs **x** et **y** donnent leurs coordonnées, le champ **street** indique quels points il faut relier pour former les segments. Dans ce cas, les données nécessaires pour représenter la carte sont contenues dans deux objets différents. Ces objets sont donnés comme argument dans la fonction `geom_` et non dès la fonction `ggplot()`.

La fonction `coord_equal()` est utile dans la production de cartes :
elle force la sortie graphique à respecter un ratio fixe de correspondance
entre l'échelle des **x** et l'échelle des **y**. Par défaut, ce ratio est de 1, ce
qui force à respecter un ratio de 1 pour 1. C'est le ratio qui sera toujours
utilisé en cartographie, mais d'autres ratios seraient paramétrables.

```r
data(Snow.deaths)
data(Snow.streets)

snowMap <- ggplot() +
  geom_path(data = Snow.streets,
            aes(x = x, y = y, group = street)) +
  geom_point(data = Snow.deaths,
            aes(x = x, y = y), size = 2) +
  coord_equal() +
  themebw()

snowMap
```
<figure><img src="image_63.png" /></figure>

Les trois géométries qui viennent d’être présentées, `geom_point()`, `geom_path()` et `geom`Polygon(), constituent la base permettant de produire des cartes avec ` ggplot2`. Elles seront utiles pour cartographier les trois grands types d’implantation d’objets géographiques : implantations ponctuelle, linéaire et zonale.

À l’issue de ce chapitre, l’utilisateur dispose de tous les éléments pour produire des sorties graphiques de qualité, pour bien gérer les palettes de couleurs et pour exporter ces sorties dans différents formats. Ces bases seront utiles pour aborder le chapitre de cartographie, tout en gardant à l’esprit que les règles de sémiologie graphique ne peuvent pas toujours être respectées en ce qui concerne les palettes de couleur vue la contrainte éditoriale d’édition en noir et blanc.
# CHAPITRE 10

## Introduction aux objets spatiaux et à la cartographie

**Objectifs** : Ce chapitre présente les principales techniques de cartographie d'objets vectoriels (points et polygones) et d'images raster. Il mobilise les fonctions graphiques de base mais aussi le package ggplot2 introduit dans le chapitre précédent.

*

**Avertissement** : Compte tenu de la contrainte d'impression en noir et blanc, le code utilisé est reproduit mais les sorties graphiques ne sont pas toujours affichées. C'est également cette contrainte qui amène, dans l'ensemble du manuel, à ne pas toujours respecter les règles de sémiologie graphique concernant les couleurs (cf. Section 9.2).
**Prérequis** Notions de cartographie thématique : discrétisation et sémiologie graphique. Connaissance des formats d'objets spatiaux.

**Description des packages utilisés** L'importation et la manipulation de données spatiales dans R s'appuient sur un certain nombre de packages spécialisés et de types d'objets spatiaux. Les fonctions et les objets spatiaux de base sont définis dans le package *sp* (*spatial*). Tous les autres packages spatiaux dépendent du package *sp*. Ce dernier, ainsi qu'un certain nombre d'autres packages plus spécialisés, sont développés par Roger Bivand, Edzer Pebesma et Virgilio Gómez-Rubio qui sont également les auteurs du manuel de référence en analyse spatiale avec R.

Avant de commencer cette section, il faut donc installer et charger les packages suivants :

- sp : la base de tous les autres packages spatiaux, dans lequel sont définis les types d'objets spatiaux ;

- rgdal : implémentation dans R de la bibliothèque GDAL (Geo-spatial Data Abstraction Library). Ce package permet d'importer de nombreux formats de données spatiales, raster et vectorielles ;

- raster : fonctions pour la manipulation des fichiers raster ;

- rgeos : interface avec la bibliothèque Geometry Engine-Open Source (GEOS) qui permet de manipuler la géométrie des objets spatiaux ;

- mapproj : fonctions de conversion des coordonnées géographiques en coordonnées projetées ;

- fields : ensemble de fonctions dédiées aux objets spatiaux, ce package est utilisé ici pour calculer des matrices de distance entre des semis de points ;

- classInt : permet de discrétiser des variables continues. Il est utilisé ici pour faire des cartes choroplèthes ;

- ggplot2 : package utilisé pour les représentations graphiques, il est aussi utile pour faire des représentations cartographiques (cf. Chapitre 9) ;

- scales : package complémentaire de ggplot2 ;

- RColorBrewer : implémentation dans R de Color Brewer, un ensemble de palettes de couleurs ;
- *animation*: permet de faire des représentations graphiques ou cartographiques animées;

- *reshape2*: contient un ensemble de fonctions permettant de re-structurer l'information d'un tableau (cf Section 2.4.3).

**Description des données utilisées** La distinction entre format vecteur et format *raster* est essentielle pour le traitement informatique d'images<fnref n="1" />. Le format *raster*, ou * bitmap*, représente l'image par une matrice de pixels dont les valeurs se traduisent par des couleurs. Les photographies en sont l'exemple le plus courant. Le format vecteur représente des objets géométriques par leurs caractéristiques (points, segments, polygones) et leurs coordonnées.

Dans un usage géographique, le fichier *raster* ou vectoriel est plus qu'une image, il contient une information géoréférencée. Dans le cas de données *raster*, cette information est contenue dans les valeurs de la matrice de pixels; dans le cas de données vectorielles, l'information intéressante peut être la géométrie des objets elle-même (la structure radio-concentrique d'un ensemble de lignes) mais il s'agit souvent des attributs attachés à ces objets (la population d'une commune).

<figure><img src="image_64.png" /><figcaption>FIGURE 10.1 – Formats vecteur et *raster*</figcaption></figure>

Plusieurs types de fichiers sont mobilisés dans cette section pour montrer les manipulations courantes dans un travail de cartographie : introduire des données externes stockées dans un tableau et avoir à gérer différents formats de données spatiales. Ici seuls deux formats courants sont importés, un format ESRI et un format MapInfo, mais la fonction

<footnote n="1">1. Le terme *image* est employé ici dans son acception du langage courant. Dans le domaine de la géomatique, *image* est synonyme de *raster*.</footnote>
`readOGR()` permet d'importer tous les formats de la bibliothèque OGR (PostGIS, kml, Idrisi, etc.). Les fichiers nécessaires sont les suivants :

- `popCom3608.csv`: fichier des populations dans les communes de Paris et la petite couronne de 1936 à 2008;

- `parispc_com.shp`: limites communales de Paris et petite couronne, données vectorielles (polygones) au format ArcGIS constitué de quatre fichiers avec les extensions `.shp`, `.dbf`, `.shx` et `.prj`;

- `parispc_bnd.shp`: contour de l'espace d'étude (Paris et petite couronne), données vectorielles (polygones) au format ArcGIS constitué de quatre fichiers avec les extensions `.shp`, `.dbf`, `.shx` et `.prj`;

- `parispc_hop.tab`: établissements hospitaliers de Paris et petite couronne, données vectorielles (points) au format MapInfo constitué de deux fichiers avec les extensions `.mif`, `.mid`;

- `Paris_dens.tif`: grille de population, données raster au format `.tif`.

## 10.1 R pour la cartographie

L'utilisation de R pour la cartographie n'est pas un choix évident. L'utilisateur qui souhaite faire une simple carte choroplèthe aura plus vite fait de la produire avec un logiciel de cartographie classique. L'utilisation de R pour faire de la cartographie se justifie si on envisage l'ensemble du flux de travail (*workflow*), c'est-à-dire l'intégration de la cartographie dans une chaîne de traitements.

On peut envisager (au moins) quatre situations dans lesquelles R s'avèrera avantageux. D'abord, pour obtenir rapidement des sorties massives : les logiciels de cartographie et de SIG sont des logiciels à interface graphique (même si certains permettent d'automatiser les traitements) et il devient rapidement coûteux de produire manuellement des ensembles de nombreuses cartes.

Ensuite, il sera intéressant d'utiliser R pour la cartographie dans le cadre d'un flux de travail unifié, c'est-à-dire pour intégrer la cartographie dans un flux de travail « classique » (analyse de données classiques, i.e. non spatiales) entièrement fait avec R et se passer ainsi d'importations et d'exportations parfois plus nombreuses que prévues, et qui sont source d'erreurs. En effet, il arrive fréquemment de faire des calculs dans un logiciel
de traitement de données, puis d'exporter les résultats pour les cartographier sur un autre logiciel, puis de refaire les calculs suite à une erreur ou une modification, puis de réexporter les résultats à cartographier, etc.

Unifier la chaîne de traitement est d'autant plus utile dans un flux de travail qui mobilise des fonctions d'analyse spatiale implémentées sous R. En effet, R dispose de plusieurs packages d'analyse spatiale et d'analyse de graphes (dont certaines sont abordées dans ce manuel) et il est très pratique de pouvoir manipuler des tableaux de données classiques, des objets spatiaux et/ou des graphes dans un flux de travail unifié.

Enfin, dans certains cas, l'utilisateur cherche à combiner les traitements et l'écriture des contenus : il s'agit d'intégrer des sorties numériques, graphiques et cartographiques dans des supports qui sont également écrits avec R (cf. 9.1). Dans ce cas, la cartographie, comme le reste des traitements, s'intègrera dans la production d'ensemble.

## 10.2 Prise en main des données spatiales

Certains packages spatiaux renvoient des messages d'avertissement lorsqu'on les charge, en particulier des messages sur les versions des bibliothèques externes utilisées : c'est le cas de rgdal qui, lors du chargement, donne une précision de ce type sur les bibliothèques GDAL et PROJ.4¹.

```
library(sp)
library(rgdal)
library(rgeos)
library(mapproj)
library(maptools)
library(raster)
```

### 10.2.1 Construction d'un objet spatial

Dans cette première section, un objet spatial est créé à partir de zéro, il s'agira ici de deux sections de route (géométrie linéaire) accompagnées

<footnote n="1">1. Pour installer le package rgdal sur les systèmes Debian ou dérivés, il faut installer au préalable les bibliothèques de développement GDAL et PROJ.4, soit les paquets libgdal-dev et libproj-dev.</footnote>
de deux attributs, la vitesse sur la section et le nombre de voies. Il est très rare qu'un utilisateur ait à créer lui-même des objets spatiaux, la plupart du temps il récupère les données produites par des organismes tels que l'IGN. Il est cependant important de comprendre comment un tel objet est construit.

Un objet spatial vectoriel est la combinaison d'une composante géométrique et d'une composante sémantique. La composante géométrique comporte des points avec des identifiants et des coordonnées. Ces points peuvent être accompagnés d'un attribut de regroupement s'ils forment des géométries plus compliquées, comme des lignes ou des polygones (cf. Section 9.3.5). La composante sémantique est un tableau relié aux géométries par un identifiant et qui contient des attributs correspondant aux objets géométriques.

L'exemple qui suit reproduit une structure de données fréquente pour les données spatiales linéaires telles que les routes ou les voies ferrées. La route est composée d'un ensemble de sections avec des attributs variables (type de route, vitesse, nombre de voies, etc.), chaque section est composée d'un ensemble de lignes, chaque ligne est composée d'un ensemble de points.

<figure><img src="image_65.png" /><figcaption>FIGURE 10.2 – Construction d'un objet linéaire</figcaption></figure>

Dans un premier temps, trois tableaux sont créés qui contiennent des points avec des coordonnées.
```
coordPoints1 <- data.frame(id = c("P1", "P2", "P3"),
                          coordx = c(1, 1, 2),
                          coordy = c(1, 2, 2))
```

```
coordPoints2 <- data.frame(id = c("P3", "P4", "P5"),
                          coordx = c(2, 3, 3),
                          coordy = c(2, 2, 1))
```

```
coordPoints3 <- data.frame(id = c("P5", "P6", "P7"),
                          coordx = c(3, 4, 5),
                          coordy = c(1, 1, 1))
```

Dans un deuxième temps, les trois tableaux sont utilisés pour créer trois lignes avec la fonction `Line()` du package `sp`. Ces lignes ne contiennent pour le moment rien d'autre qu'un ensemble de points avec leurs coordonnées.

```
line1 <- LineLine coordPoints1[, 2:3])
line2 <- Line coordPoints2[, 2:3])
line3 <- Line coordPoints3[, 2:3])
```

Dans un troisième temps, les trois lignes sont regroupées dans des sections avec la fonction `Lines()` et un identifiant leur est attribué. Les deux objets `section1` et `section2` seront les formes géométriques (features) de l'objet spatial. La géométrie finale de l'objet spatial est créée avec la fonction `SpatialLines()` qui contient toutes les formes géométriques ainsi qu'un système de projection (ce point est détaillé dans la section suivante).

```
section1 <- Lines(slinelist = list(line1, line2),
                  ID = "Sec1")

section2 <- Lines(slinelist = list(line3),
                  ID = "Sec2")

geoLines <- SpatialLines(list PUBLIC ("+init=epsg:2154"))
```

Enfin, une table attributaire est créée avec deux champs, la vitesse et le nombre de voies. Ce tableau a deux lignes correspondant aux deux formes
géométriques. Le lien se fait entre l'identifiant des géométries et les noms de ligne du tableau. La fonction `SpatialLinesDataFrame()` crée l'objet spatial final qui combine la composante géométrique et la composante sémantique.

```r
attrTable <- data.frame(SPEED = c(30, 50),
                        LANES = c(2, 3))

row.names(attrTable) <- c("Sec1", "Sec2")

geoDataLines <- SpatialLinesDataFrame(sl = geoLines,
                                     data = attrTable,
                                     match.ID = TRUE)
```

L'objet spatial sans composante sémantique (`geoLines`) peut déjà être manipulé et cartographié. L'objet spatial final contient une information sémantique (`geoDataLines`), il peut être cartographié en utilisant les variables de la table attributaire, en traduisant la variation des attributs numériques par des variations visuelles, de taille ou de couleur par exemple.

```r
plot(geoDataLines,
     lwd = geoDataLines@data$SPEED,
     axes = TRUE)
```

Dans cet exemple, un objet spatial de type vectoriel et linéaire a été construit en assemblant des objets géométriques avec des attributs sémantiques. Dans un flux de travail réel, il est plus fréquent de travailler avec des objets spatiaux déjà construits sur lesquels il faut pouvoir extraire et manipuler des éléments géométriques et/ou sémantiques.

Dans les chapitres précédents, on a vu que certains types d'objets, `data.frame` ou `list` par exemple, sont composés de plusieurs éléments (`slots`) auxquels on accède avec l'opérateur `$`. Les objets spatiaux de type vectoriel comportent deux niveaux d'éléments : on accède au premier avec l'opérateur `@` et au second avec l'opérateur `$`.

Les principaux éléments de premier niveau sont : `@data` qui contient la table attributaire, `@lines` (`@points` ou `@polygons` selon le type de géométrie) qui contient la géométrie des objets et `@proj4string` qui contient l'indication du système de projection. Ensuite, il est possible avec
l'opérateur `$` d'accéder à un champ de la table attributaire de la même façon qu'à un champ de tableau classique. Lors de l'utilisation des opérateurs `@` et `$`, ne pas oublier la touche `Tab` qui complète automatiquement les noms des objets (cf. Section 1.3).

```
geoDataLines@data$TYPE <- "Highway"
attrTable <- geoDataLines@data
```

Pour récupérer et manipuler la composante géométrique, plusieurs fonctions sont disponibles, en particulier `bbox()` qui renvoie l'étendue spatiale des données (`bounding box`) et `coordinates()` qui renvoie les coordonnées des objets géométriques. Si l'objet spatial est ponctuel, la fonction renvoie les coordonnées de ces points ; si l'objet est linéaire, elle renvoie les coordonnées des points qui forment les lignes ; si l'objet est de type polygonal, elle renvoie les centroïdes des polygones.

```
bboxGeolines <- bbox(geoDataLines)
coordPoints <- coordinates(geoDataLines)
```

Cette courte introduction a montré comment construire et manipuler un objet spatial à partir d'un exemple très simple. Les sections qui suivent poursuivent l'étude des communes de Paris et de petite couronne. Deux types de géométries y sont manipulées : des points et des polygones.

### 10.2.2 Importation et prise en main des objets spatiaux

La fonction `readOGR()` du package `rgdal` est utilisée pour charger les données **Communes** (format ArcGIS) et **Hopitaux** (format MapInfo). Cette fonction importe les données et les stocke dans un objet de type *spatial*, selon leur géométrie : *SpatialPolygonsDataFrame*, *SpatialLinesDataFrame* ou *SpatialPointsDataFrame*.

La fonction `readOGR()` peut lire des formats contenant plusieurs couches, c'est pourquoi il ne suffit pas de préciser le nom du fichier (`dsn`), mais également le nom de la couche (`layer`). Dans le cas des fichiers ESRI (`shp`) et MapInfo (`mif`), ils ne contiennent qu'une seule couche (du même nom que le fichier) qu'il faut préciser avec l'argument `layer`.
L'argument (*encoding*) permet de préciser l'encodage, de la même façon que dans l'importation de simples tableaux (cf. Section 2.2.3), et l'argument *StringsAsFactors* sert à spécifier que les champs alphanumériques ne doivent pas être transformés en champs de type *factor*. Pour importer l'image *raster*, la fonction `readGDAL()` du package *raster* est utilisée.

```r
muniBound <- readOGR(dsn = "data/parispc_com.shp",
  layer = "parispc_com",
  encoding = "utf8",
  stringsAsFactors = FALSE)

regBoundary <- readOGR(dsn = "data/parispc_bnd.shp",
  layer = "parispc_bnd",
  encoding = "utf8",
  stringsAsFactors = FALSE)

publicHosp <- readOGR(dsn = "data/parispc_hop.MIF",
  layer = "parispc_hop",
  encoding = "latin1",
  stringsAsFactors = FALSE)

densGrid <- readGDAL("data/paris_dens.tif")
```

Lors de l'importation des fichiers, certaines informations sont renvoyées : le driver utilisé (ESRI et MapInfo dans le cas présent), le nombre d'objets (*features*) et de variables (*fields*) pour les fichiers vectoriels, le type de géométries (*points* et *polygones* dans cet exemple), la taille de la matrice de pixels dans le cas du fichier *raster* (392*380 dans ce cas).

Il est intéressant d'examiner le type et la structure de ces données spatiales grâce aux fonctions *class()* et *str()*. La fonction *str()* renvoie les attributs, le type de données spatiales, elle précise l'étendue des données (*bbox*) et le système de coordonnées (*proj4string*).

```
class(muniBound)
class(publicHosp)

str(muniBound)
str(publicHosp)
```
Une fois les données importées, elles peuvent être visualisées avec la fonction `plot()`. Cette fonction générique s'utilise de la même façon pour les données spatiales que pour les données classiques : l'argument `col` précise la couleur de remplissage, l'argument `border` la couleur du trait, l'argument `pch` le type de point, et finalement l'argument `add` pour superposer chaque nouvelle couche à la précédente.

```r
plot(muniBound,
      col = "grey",
      border = "white",
      axes = TRUE)

plot(publicHosp, pch = 20, col = "black", add = TRUE)
```

<figure><img src="image_66.png" /></figure>

Pour visualiser le raster, on commence par créer une palette de couleurs (cf. Section 9.2) qui est donnée comme argument de la fonction `image()` du package `raster`. Pour travailler directement sur les values
leurs des pixels, il est possible de transformer l'objet raster (*SpatialGrid-Frame*) en une matrice de pixels avec la fonction `as.matrix()`.

```r
colPalHeat <- c("white", rev(heatcolors(12)))
image(densGrid, col = colPalHeat)
plot(muniBound, add = TRUE)
pixGrid <- as.matrix(densGrid)
```

La conception de l'échelle et de la légende est détaillée dans la suite du chapitre.

### 10.2.3 Le système de projection

Jusqu'à présent, les données spatiales ont été importées et visualisées sans se soucier de leur système de projection qui est pourtant essentiel lorsqu'on manipule des données spatiales. Le système de projection est la transformation mathématique qui assure le passage d'une portion d'espace terrestre en trois dimensions à une représentation en deux dimensions. Cette transformation se traduit nécessairement par une déformation, ce qui a un impact en termes de visualisation (cartographie) et en termes de calculs de surfaces ou de distances (statistique spatiale).

Il existe un ensemble de normes d'identification des systèmes de projection. Le code qui identifie le système de projection suit un standard international (SRID - *Spatial Reference System Identifier*) et le standard le plus répandu est le code EPSG (*European Petroleum Survey Group*). Pour connaître le système de projection des objets spatiaux, la fonction `proj4string()` est utilisée, elle renvoie une information contenue dans les formats de fichiers classiques. Par exemple, pour un format ESRI shape, la projection est donnée dans le fichier d'extension .prj.

```
proj4string(muniBound)
```

Le résultat obtenu est une chaîne de caractères contenant plusieurs informations, principalement le système de projection et l'étendue des données (*bounding box*). Ici, le système est de type *Lambert Conformal Conic*.
(lcc) et l'étendue est donnée par les coordonnées géographiques de deux points définissant la « boîte ». Il est utile de trouver le code qui identifie cette projection dans la liste de codes EPSG.

Certains sites web proposent une telle fonctionnalité<fnref n="1" /> en récupérant l'information contenue dans le fichier d'extension .prj du fichier ESRI (shapefile) d'origine. Ceci dit, l'information et la modification du système de projection peuvent se faire intégralement avec R. La projection de l'objet muniBound semble mal référencée et sera convertie en « RGF93-Lambert93 » qui est le système de référence proposé par l'IGN<fnref n="2" />.

Pour repérer le code correspondant à cette projection, il faut utiliser la fonction `make_EPSG()` du package `rgdal` qui renvoie un `data.frame` contenant l'ensemble des 4211 systèmes de projection avec leurs identifiants. Une sélection est effectuée sur ce tableau grâce à la fonction `grep()`, qui recherche les expressions régulières, c'est-à-dire les chaînes de caractères répétés (motifs). En recherchant tous les systèmes qui contiennent la chaîne de caractères « RGF » (Référentiel Géodésique Français), on obtient un tableau de 14 lignes qui contient le code EPSG correspondant au système voulu : il s'agit du code **2154**.

```r
projCodes <- make_EPSG()
selecRgf <- projCodes[grep("RGF", projCodes$note),]
```

Il s'agit maintenant de définir le système de projection des données spatiales, limites communales et hôpitaux, grâce à la fonction `spTransform()` qui ajoute au système initial le nouveau code EPSG. Avec la fonction `proj4string()` utilisée plus haut, on vérifie que le code EPSG a bien été modifié.

```r
muniBound <- spTransform(x = muniBound,
                          CRSobj = CRS("+init=epsg:2154"))
```

```r
publicHosp <- spTransform(x = publicHosp,
                          CRSobj = CRS("+init=epsg:2154"))
```

<footnote n="1">1. Voir par exemple le site http://prj2epsg.org.</footnote>

<footnote n="2">2. Le site de l'IGN est riche en informations sur la géodésie : http://geodesie.ign.fr/index.php?page=geodesie.</footnote>
### 10.2.4 Exportations

Au cours du travail cartographique, on peut souhaiter exporter deux types d'objets : des résultats cartographiques (des images) ou bien des objets spatiaux. Les dispositifs graphiques pour exporter les résultats sous forme d'images sont présentés dans la Section 9.1.

Pour exporter des objets spatiaux, on peut utiliser la fonction `writeOGR()`, pendant de la fonction `readOGR()` utilisée précédemment, et créer des fichiers dans les formats les plus courants, ESRI ou MapInfo par exemple. Le site de la bibliothèque GDAL/OGR<fnref n="1" /> recense tous les formats disponibles. Dans l'exemple qui suit, le fichier de communes est exporté au format ESRI (.shp) dans un dossier intitulé « MyFolder ».

```R
writeOGR(muniBound,
         dsn = "MyFolder",
         layer = "muniBound",
         driver = "ESRI Shapefile")
```

Si le dossier n'existe pas, il est créé au moment de l'exportation. Si seul le nom du dossier est donné sans préciser le chemin qui y mène, ce dossier sera créé dans le répertoire de travail (cf. Section 2.2.1).

## 10.3 Cartographie des objets ponctuels

Le fichier des hôpitaux est doté d'une variable quantitative de stock intitulée **Capacite**. Il s'agit de la capacité totale de l'hôpital mesurée en nombre de lits. Cette variable peut être traduite par une variable visuelle de taille pour produire une carte en cercles proportionnels.

La fonction `plot()` est utilisée en faisant varier la taille des symboles (cex, character expansion factor) selon la capacité des hôpitaux. Par défaut, la variation de la valeur numérique est traduite par une variation proportionnelle du rayon du cercle. Ce paramétrage est corrigé : pour une visualisation appropriée, la propotionnalité doit affecter l'aire du cercle et non son rayon.

<footnote n="1">1. http://www.gdal.org.</footnote>
Il faut également ajouter une légende ainsi qu'une échelle. La légende est créée puis localisée sur le graphique soit par une localisation relative, soit en précisant les coordonnées *x* et *y*. Dans le premier cas, on indique la localisation par des arguments prédéfinis (*left*, *right*, *topleft*, etc.). Dans le second cas, on peut s'aider de la fonction `locator()` qui permet de récupérer les coordonnées d'un point sur le graphique en cliquant dessus.

Finalement, on utilise la fonction `legend()` en veillant à reprendre le style exact des symboles des hôpitaux choisis précédemment, notamment la forme, la couleur et la taille. Les arguments `pch`, `col` et `pt.cex` sont utilisées pour spécifier respectivement le type de figuré, la couleur et la taille des points, et l'argument `bty` (*box type*) permet de dessiner ou non un encadré de la légende.

```r
plot(muniBound, border = "darkgrey")

plot(publicHosp,
       pch = 16,
       cex = 0.2 * sqrt(publicHosp$Capacite / pi),
       col = "black", add = TRUE)

title("Capacité des hôpitaux à Paris-Petite couronne")

legendHosp <- c(100, 400, 700, 1000)

legend("bottomleft",
       legend = legendHosp,
       pch = 16, col = "black",
       pt.cex = 0.2 * sqrt(legendHosp / pi),
       bty = "n", title = "Nombre de lits")
```
Capacité des hôpitaux à Paris-Petite couronne

<figure><img src="image_67.png" /></figure>

Pour ajouter une échelle, il y a deux solutions, qui demandent toutes deux de spécifier manuellement la localisation du trait, la longueur du trait, l'épaisseur et le texte. Dans le cas ci-dessous sont utilisées les fonctions `arrows()` pour le trait de l'échelle, puis `text()` pour le texte de l'échelle. Il est également possible de placer une échelle avec les fonctions `locator()` et `SpatialPolygonsRescale()`. La fonction `locator()` est interactive, elle demande de cliquer sur la carte pour obtenir les coordonnées du point où l'on souhaite placer l'échelle. Puis, la fonction `SpatialPolygonsRescale()` permet de placer l'échelle à partir de ce point.

Finalement, pour visualiser correctement les cercles proportionnels qui se chevauchent, il est nécessaire de spécifier une couleur de bordure qui les différencie (ici le blanc) et de trier les hôpitaux selon leur taille pour que les cercles les plus gros se retrouvent en arrière-plan.
Attention, l’argument *col* désigne deux choses différentes selon que l’on représente un objet zonal ou un objet ponctuel. Pour les objets zonaux, *col* désigne la couleur de remplissage des polygones et *border* la couleur du trait. Pour les objets ponctuels, *col* désigne la couleur du symbole, donc du trait, et l’argument *bg* (*background*) la couleur du fond, donc le remplissage du symbole ponctuel. Le type de point choisi (*pch*) est le numéro 21, parce qu’on cherche à distinguer la couleur de remplissage de la couleur de bordure. Il faut donc un type de symbole avec bordure, ce qui n’est pas le cas des types de symboles 15 à 20. On peut obtenir la liste des symboles ponctuels dans l’aide (`?pch`).

Le fonctionnement de la fonction `arrows()` nécessite une explication. Il s'agit de situer la droite représentant l'échelle par rapport aux marges du graphiques. La fonction `par()` permet de manipuler les paramètres graphiques généraux, l'argument `usr` est un vecteur de quatre valeurs : $x1, x2, y1, y2$ qui désignent les bornes de la sortie graphique. La fonction `arrows()` demande également quatre valeurs, les coordonnées $x$ et $y$ du premier point de la droite et les coordonnées $x$ et $y$ du second point de la droite. Ici le premier point de la droite a pour coordonnées $x1 + 1000$ unités (dans l'unité de mesure de l'objet spatial) et $y1 + 1000$ unités (désigne le coin en bas à gauche de la sortie graphique). Le second point de la droite a pour coordonnées $x1+10000$ unités (soit 10 km) et $y1+1000$ (bien sûr, $y$ reste constant puisque la droite est horizontale). La fonction `text()` fonctionne de la même façon, en situant le point à partir duquel commence le texte et en précisant le contenu de ce texte (10 km).
```R
publicHosp <- publicHosp[order(publicHosp$Capacite,
 decreasing = TRUE), ]
```

```R
plot(muniBound, border = "darkgrey")
```

```R
plot(publicHosp,
     pch = 21,
     cex = 0.2 * sqrt(publicHosp$Capacite / pi),
     col = "white",
     bg = "black",
     add = TRUE)
```

```R
title("Capacité des hôpitaux à Paris-Petite couronne")
```

```R
legendHosp <- c(100, 400, 700, 1000)
```

```R
legend("topleft",
       legend = legendHosp,
       pch = 21,
       col = "white",
       pt.bg = "black",
       pt.cex = 0.2 * sqrt(legendHosp / pi),
       bty = "n",
       title = "Nombre de lits")
```

```R
arrows(par() $ usr[1] + 1000,
        par() $ usr[3] + 1000,
        par() $ usr[1] + 10000,
        par() $ usr[3] + 1000,
        lwd = 2, code = 3,
        angle = 90, length = 0.05)
```

```R
text(par() $ usr[1] + 5050,
      par() $ usr[3] + 2700,
      "10 km", cex = 1.2)
```
Capacité des hôpitaux à Paris-Petite couronne

Nombre de lits

<figure><img src="image_68.png" /></figure>

<figure><img src="image_69.png" /></figure>

Ce code semble exagérément long pour produire une simple carte de cercles proportionnels. Sa longueur vient du fait que chacun des éléments est paramétré de façon peu automatisée. Cependant, le code peut être encapsulé dans une fonction beaucoup plus simple pour un usage plus rapide. Le même résultat aurait aussi pu être produit avec le package ggplot2 présenté dans le chapitre précédent.

10.4 Cartographie des objets zonaux

À la date de préparation de ce manuel, il existe principalement deux façons de produire des cartes avec R : avec la fonction générique plot() qui permet de représenter directement des objets spatiaux, ou bien avec la fonction ggplot() du package ggplot2 très utile pour tous types de représentations graphiques. Ces deux méthodes sont expliquées ici ainsi
que leurs avantages et leurs inconvénients respectifs. L'exemple d'application est la production d'une carte choroplèthe de la densité de population par commune en 2008.

Pour chacune des deux fonctions `plot()` et `ggplot()`, deux critères de discrétisation sont testés : en quartiles et selon l'algorithme de Jenks. Dans ce dernier cas, le package `ClassInt`, est utilisé, comme dans le Chapitre 4.

L'algorithme de Jenks est très utilisé pour produire des cartes choroplèthes parce qu'il a été conçu dans ce but et qu'il est implémenté dans la plupart des logiciels de cartographie. Cependant, d'autres algorithmes existent pour discrétiser une variable continue selon des seuils dits « natuels » et ils sont souvent plus efficaces. Certains sont implémentés dans le package `ClassInt`, comme celui de Fisher ou les algorithmes de classification multivariée (k-means, CAH).

### 10.4.1 Cartographie avec la fonction plot()

On commence par créer des palettes de couleurs à partir du package `RColorBrewer`. Ce n'est pas la seule façon de construire une palette, mais ce package présente l'avantage de prédéfinir des palettes adaptées à différents types de représentations (cf. Section 9.2). La palette créée est une palette continue avec cinq couleurs de même ton et d'intensité variable (montée de valeurs).

```r
library(RColorBrewer)
greyPal <- brewer.pal(n = 5, name = "Greys")
```

Ensuite, une jointure est faite entre les données attributaires de l'objet spatial et les données externes contenues dans le tableau `popCom3608`. Il y a deux précautions à prendre pour faire la jointure, l'une concerne la variable de jointure, l'autre concerne la fonction de jointure. La variable de jointure, ici le code de la commune, doit être du même type dans les deux tableaux. La fonction `class()` permet de s'en assurer. Concernant la fonction de jointure, la fonction `merge()` présentée dans la Section 2.4.1 peut être source d'erreur. Il faut en effet s'assurer que la table attributaire n'est modifiée ni dans l'ordre des lignes ni dans le nombre de lignes. Si, au
cours de la jointure, la table attributaire bouge (l'ordre des lignes change par exemple) la correspondance entre la composante sémantique (i.e. la table attributaire) de la couche à cartographier et sa composante géométrique n'est plus assurée.

La fonction match() garantit cette intégrité de la jointure. Les lignes qui suivent se lisent de la façon suivante : la table attributaire finale est la combinaison de la table attributaire initiale et du tableau externe (popCom3608). Les lignes du tableau externe sont ajoutées quand il y a correspondance entre les identifiants. L'ajout conserve l'ordre et le nombre de lignes du premier argument de la fonction match(), à savoir la table attributaire de l'objet spatial.

```sql
str(muniBound@data$INSEE_COM)
muniBound$INSEE_COM <- as.character(muniBound$INSEE_COM)

muniBound@data <- data.frame(
  muniBound@data,
  popCom3608[match(muniBound$INSEE_COM,
    popCom$CODGEO), ]
)
```

Ce code nécessaire à la jointure est difficilement lisible et il sera fastidieux de le ré-écrire à chaque jointure. Il serait intéressant de l'encapsuler dans une fonction (cf. Section 3.3) pour pouvoir faire des jointures plus facilement. Cette fonction de jointure sur la table attributaire est nommée AttribJoin(), ses arguments sont le tableau externe (df), l'objet spatial (spdf) et les variables de jointure correspondantes (df.field et spdf.field).
```javascript
AttribJoin <- function(df, spdf, df.field, spdf.field) {
  if(is.factor(spdf@data[, spdf.field]) == TRUE) {
    spdf@data[, spdf.field] <- as.character(
      spdf@data[, spdf.field]
    )
  }
  if(is.factor(df[, df.field]) == TRUE) {
    df[, df.field] <- as.character(df[, df.field])
  }
  spdf@data <- data.frame(
    spdf@data,
    df[match(spdf@data[, spdf.field], df[, df.field]), ]
  )
  return(spdf)
}
```

La fonction de jointure commence par tester le type des variables d'entrée : s'il s'agit de facteurs, ils sont transformés en champs alphanumériques. La fonction de jointure est la fonction `match()` présentée plus haut. La fonction qui vient d'être créée est simple et générique, elle peut être utilisée sur n'importe quel type de données spatiales.

```r
muniBound <- AttribJoin(df = popCom3608,
                        spdf = muniBound,
                        df.field = "CODGEO",
                        spdf.field = "INSEE_COM")
```

Une fois la jointure effectuée, on peut travailler directement sur les données attributaires, les discrétiser et les représenter. Un vecteur de seuils (*breaks*) est créé avec la fonction `quantile()`. Celle-ci découpe la variable selon les probabilités données en argument avec la fonction `seq()`, qui correspondent ici à des quintiles.

La discrétisation et l'assignation de la palette de couleurs se font dans la même étape, avec la fonction `cut()`. Cette fonction renvoie un objet de type `factor`, qui désigne chaque catégorie par un entier (ici de 1 à 5) et un label (étiquette de valeur) qui prend ici la couleur de la palette créée précédemment. L'argument `break` est utilisé pour préciser les seuils, l'argument `labels` pour préciser les étiquettes des valeurs, les arguments `include.lowest` et `right` pour créer des intervalles fermés à gauche et ouverts à droite de type `[valeur1; valeur2[.]`.
Finalement la fonction as.character() transforme le facteur en un vecteur alphanumérique qui contient la couleur correspondant à chaque catégorie. Le nouveau champ créé (DENSQ5) dans la table attributaire contiendra donc le code couleur utilisé par la suite avec la fonction plot(). La légende est construite de la même façon, la différence étant dans la fonction levels() qui renvoie la liste de modalités de la variable. Les modalités sont ici les intervalles de discrétisation.

```r
quintBrks <- quantile(muniBound$DENSITE,
                      prosbs = seq(0, 1, 0.2),
                      names = TRUE)

muniBound$DENSQ5 <- as.character(cut(muniBound$DENSITE,
                                     breaks = quintBrks,
                                     labels = greyPal,
                                     include.lowest = TRUE,
                                     right = FALSE))

legendQ5 <- as.character(levels(cut(muniBound$DENSITE,
                                   breaks = quintBrks,
                                   include.lowest = TRUE,
                                   right = FALSE)))
```

La représentation cartographique proprement dite se fait comme dans l'exemple de carte en cercles proportionnels de la section précédente. La variable créée plus haut contient des étiquettes de couleurs : en précisant cette variable dans l'argument col, elle devient une variable visuelle.

```r
plot(muniBound, col = muniBound$DENSQ5, border = "black")

legend("bottomleft",
        legend = legendQ5,
        bty = "n",
        fill = greyPal,
        cex = 0.8,
        title = "Densité (hab/km2)")

title("Densité de population à Paris et en petite couronne")
```
Densité de population à Paris et en petite couronne

<figure><img src="image_70.png" /></figure>

La même manipulation peut être effectuée avec un autre algorithme de
discrétisation, celui de Jenks par exemple. Il suffira de substituer le vecteur
de seuils (breaks) en quintiles par un autre vecteur de seuils.

jenksDiscret <- classIntervals(var = muniBound$DENSITE,
                             n = 5,
                             style = "jenks")
jenksBrks <- jenksDiscret$brks

10.4.2 Cartographie avec la fonction ggplot()

La fonction ggplot() du package ggplot2 permet de personnaliser
les aspects de la représentation plus finement que la fonction plot().
Elle a aussi certains inconvénients : elle ne prend pas comme argument un
objet de type spatial, il faut donc passer par une transformation préalable
expliquée par la suite.

L'usage de `ggplot()` pour la cartographie se justifie quand le résultat
graphique attendu est trop complexe pour pouvoir être fait avec `plot()`
et quand les données spatiales ne sont pas trop lourdes. En effet, la trans-
formation des données et la fonction de représentation sont plus lentes
que la fonction `plot()`. Pour le présent jeu de données, composé de 143
communes et arrondissements, cela ne pose pas de problème, mais il est
par exemple impossible de faire, sur un ordinateur de bureau, des cartes
des 36 000 communes françaises avec `ggplot()`.

En premier lieu, la variable de densité est discrétisée selon la même
méthode que dans la section précédente (quintiles). En revanche, la va-
riable est discrétisée dans le tableau externe (`popCom3608`) et non dans
le tableau de données attributaires de l'objet spatial.

```r
popCom3608$DENSQ5 <- cut(
  popCom3608$DENSITE,
  breaks = quintBrks,
  include.lowest = TRUE,
  right = FALSE,
  labels = c("Q1", "Q2", "Q3", "Q4", "Q5")
)
```

Pour visualiser les communes avec `ggplot()`, il faut transformer l'objet
spatial (*SpatialDataframe*) en un objet que `ggplot()` puisse lire,
à savoir un *data.frame* classique. La fonction `fortify()` du package
`ggplot2` effectue ce travail : pour un objet à géométrie polygonale,
comme ici le fond communal, l'objet résultant est un tableau contenant les
sommets des polygones avec leurs coordonnées (long, lat) ainsi qu'une
variable de regroupement qui identifie le polygone avec les sommets le
constituant. La représentation graphique dessinera les polygones en regroupant ces sommets, de la même façon que dans l'exemple présenté au
chapitre précédent (cf. Section 9.3.5).

C'est à cause de cette transformation du format spatial au format « for-
tifié » qu'il faut nécessairement discrétiser la variable au préalable et non
sur l'objet fortifié. En effet, dans cet objet un polygone est représenté par
*n* lignes du tableau, *n* étant le nombre de sommets le constituant. Si la
discrétisation était faite sur ce tableau, la valeur de la variable associée à chaque polygone serait pondérée par le nombre de sommets, ce qui fausserait le résultat.

Dans un deuxième temps, une jointure est effectuée entre l'objet fortifié et le tableau externe contenant la variable discrétisée. À la différence de la jointure de la section précédente, la fonction merge() peut être utilisée ici, à condition de préciser les arguments all et sort. Cette jointure est l'équivalent d'un left join dans le langage SQL.

```
fortMuni <- fortify(muniBound, region = "INSEE_COM")

fortMuni <- merge(x = fortMuni,
                 y = popCom3608,
                 by.x = "id",
                 by.y = "CODGEO",
                 all.x = TRUE,
                 sort = FALSE)
```

La fonction `geomPolygon()` sert à graphier des objets zonaux : l'argument `group` désigne l'identifiant des polygones et l'argument `fill` désigne la couleur de remplissage, contenue dans le champ `DENSQ5` créé précédemment. Finalement on modifie l'objet graphique (`mapDens`) avec la fonction `coord_equal()` qui assure la conservation de la même échelle `x` et `y` et avec la fonction `scale_fill_brewer()` qui permet de donner un titre à la légende et de désigner une palette de couleurs.

```
mapDens <- ggplot(fortMuni) +
  geomPolygon(aes(x = long,
                 y = lat,
                 group = group,
                 fill = DENSQ5),
  colour = "black") +
  scale_fill_manual(name = "Densité 2008",
                    values = greyPal) +
  coord_equal()
```

L'objet `mapDens` peut maintenant être visualisé, dans ce cas avec le thème défini par défaut `theme_bw()`.
<figure><img src="image_71.png" /></figure>

Il est également possible de définir un thème vide, sans grille, sans axes, sans titres d'axes, etc.

```r
theme_update(axis.ticks = element-blank(),
  axis.text.x = element-blank(),
  axis.title.x = element-blank(),
  axis.text.y = element-blank(),
  axis.title.y = element-blank(),
  panel.grid.minor = element-blank(),
  panel.grid.major = element-blank(),
  panel-background = element-rect(),
  legend(position = "bottom")
)
```
## 10.5 Cartes *take away*

Il faut maintenant revenir aux constats posés en début de chapitre (cf. Section 10.1). Si l'on met dans la balance la longueur du code à écrire pour produire une simple carte, l'utilisation de R n'est pas toujours avantageuse vis-à-vis de logiciels à interface graphique. R est avantageux si on replace la sortie cartographique dans l'ensemble du flux de travail, et en particulier lorsqu'il faut produire des cartes à la chaîne.

Il reste cependant possible de produire rapidement une carte prête à l'emploi, prête à être intégrée dans un document par exemple. Le package *rCarto* présenté ci-dessous<fnref n="1" /> travaille avec un fichier spatial au format ESRI (.shp) et un tableau (*data*, *frame*) existant dans l'espace de travail.

Deux fonctions permettent de reproduire les cartes présentées dans ce chapitre, `mapCircles()` pour la carte en cercles proportionnels et `mapChoropleth()` pour la carte choroplèthe. Voici en exemple qui reproduit la carte de la densité de population discrétisée en quintiles.

```r
library(rCarto)
mapChoropleth(
  shpFile = "data/parispc_com.shp", shpId = "INSEE_COM",
  df = popCom3608, dfId = "CODGEO", var = "DENSITE",
  nclass = 5, palCol = "Greys", style = "quantile",
  posLeg = "topleft", lgdRnd = 0, legend = "Hab./km2",
  title = "Densité de population (2008)",
  author = "Groupe ElementR",
  sources = "Insee, Recensement de la population 2008")
```

<footnote n="1">1. *Package* développé par Timothée Giraud de l'UMS RІATE.</footnote>
<figure><img src="image_72.png" /></figure>

## 10.6 Cartographie à la chaîne

L'un des intérêts majeurs de faire de la cartographie avec R est de pouvoir faire des traitements automatisant la production de cartes, ce qui est coûteux à faire manuellement sur un logiciel à interface graphique. Il s'agit ici de comparer les densités à plusieurs époques (de 1936 à 2008), ce qui pose deux questions : la question du format de sortie et la question de la méthode de discrétisation. Concernant la première question, il y a au moins trois façons d'intégrer dans une seule sortie un ensemble de plusieurs cartes :
- afficher toutes les cartes sur une même page,

- afficher toutes les cartes sur plusieurs pages d'un même document,

- intégrer toutes les cartes dans une animation qui affiche chaque carte successivement.

L'exemple suivant montre la marche à suivre pour la première option. La deuxième option se fait avec les fonctions `pdf()`, `png()` ou autres suivant le format souhaité (cf. Section 9.1). Ces fonctions prendraient comme argument une boucle qui cartographie les données à chaque date du recensement.

La troisième option nécessite l'utilisation de *packages* spécifiques, par exemple `animation` qui permet d'intégrer toutes les cartes dans un format animé tel que `HTML` (fonction `saveHTML()`), `GIF` (fonction `saveGIF()`) ou autre.

Concernant la discrétisation, il faut choisir une méthode qui découpe la série selon des valeurs de centralité et de dispersion (moyenne, écart-type, quantiles, etc.) et non une méthode comme celle de Jenks qui ne s'appuie pas sur des résumés numériques fixes. On peut dès lors discrétiser chaque variable de densité de 1936 à 2008 séparément, ou bien discrétiser sur l'ensemble des valeurs de 1936 à 2008. C'est cette seconde option qui est montrée ici.

Concernant la fonction de représentation graphique, cette planche sera construite avec `ggplot()`, comme dans l'exemple présenté à la Section 9.3.4.

Dans un premier temps, les densités de population de 1936 à 2008 sont calculées. Ensuite, ces variables sont discrétisées ensemble. Il n'est pas nécessaire de joindre toutes ces données dans un seul champ pour calculer les quantiles : la fonction `quantile()` peut être appliquée à l'ensemble du tableau, à condition de le transformer en objet de type `matrix`, qui est un vecteur bidimensionnel.

```
evolDens <- popCom3608[, 3:11] / popCom3608$SURF
evolQuint <- quantile(as.matrix(evolDens),
                     names = TRUE,
                     seq = seq(0, 1, 0.2))
```
Une fois créé le vecteur de seuils, chacune des colonnes du tableau est successivement discrétisée grâce à la fonction `cut()` appliquée colonne à colonne avec la fonction `apply()` (cf. Section 3.4). Le résultat final est transformé en tableau et une jointure est effectuée entre ce tableau et l'objet fortifié.

```r
discrDens <- apply(evolDens, 2, cut,
  breaks = evolQuint,
  include.lowest = TRUE,
  right = FALSE,
  labels = c("Q1", "Q2", "Q3", "Q4", "Q5"))

discrDens <- as.data.frame(discrDens)
discrDens$CODGEO <- popCom3608$CODGEO

colnames(discrDens) <- c(
  "DENS1936", "DENS1954", "DENS1962", "DENS1968", "DENS1975",
  "DENS1982", "DENS1990", "DENS1999", "DENS2008", "CODGEO")

fortMuni <- merge(x = fortMuni,
  y = discrDens,
  by.x = "id",
  by.y = "CODGEO",
  all.x = TRUE,
  sort = FALSE)
```

L'objet `fortMuni` est un tableau « large » contenant neuf champs de densité, un champ pour chaque année du recensement. Pour construire une planche en utilisant une variable désignant des sous-groupes de données (ici des années), il faudrait en faire une tableau « long » qui n'aurait que deux champs à cartographier : un champ contenant une variable qualitative indiquant l'année du recensement et un champ contenant la variable à représenter (ici des classes de densités). L'objet est transformé en format long avec la fonction `melt()` (cf. Section 2.4.3).

```r
fortMuniLong <- melt(
  data = fortMuni,
  id_vars = colnames(fortMuni)[1:7],
  measure_vars = colnames(fortMuni)[22:30]
)
```
Cet objet peut être cartographié, en utilisant l'argument `facetWrap` pour intégrer l'ensemble des cartes dans une seule sortie. Cette fonction constuire une planche de *n* cartes, une pour chaque modalité de la variable donnée comme argument (`facets`).

```r
mapEvolDens <- ggplot(fortMuniLong) +
  geom_polygon(aes(long, lat, group = group, fill = value),
    colour = "black", size = 0.3)
```

```r
mapEvolDens +
  facetWrap(facets = ~variable) +
  scale_fill_manual(name = "Densité (quintiles)",
    values = greyPal) +
  coord_equal()
```

DENS1936

<figure><img src="image_73.png" /></figure>

DENS1954

<figure><img src="image_74.png" /></figure>

DENS1962

<figure><img src="image_75.png" /></figure>

DENS1968

<figure><img src="image_76.png" /></figure>

DENS1975

<figure><img src="image_77.png" /></figure>

DENS1982

<figure><img src="image_78.png" /></figure>

DENS1990

<figure><img src="image_79.png" /></figure>

DENS1999

<figure><img src="image_80.png" /></figure>

DENS2008

<figure><img src="image_81.png" /></figure>

Densité (quintiles) Q1 Q2 Q3 Q4 Q5
Un rapide commentaire thématique s'impose : la population de l'agglomération parisienne, comme c'est le cas dans la plupart des grandes métropoles européennes, a tendu à s'étaler au cours du 20e. D'autres résultats concordants ont été obtenu dans les chapitres précédents : le cœur de l'agglomération, c'est-à-dire Paris *intramuros*, s'est dépeuplé durant la seconde moitié du 20e ; parallèlement, le peuplement de la proche banlieue s'accélérait.

Au terme de ce chapitre, l'utilisateur sait construire et manipuler des
objets spatiaux, il sait travailler sur les tables attributaires et faire des
jointures de données externes, il dispose enfin des outils pour construire
différents types de cartes thématiques. Ces bases seront nécessaires pour
aborder le chapitre suivant sur l'analyse spatiale.
<figure><img src="image_82.png" /></figure>
# CHAPITRE 11

## Initiation aux statistiques spatiales

**Objectifs** : Ce chapitre propose une initiation à des analyses fondées sur les localisations dans l'espace. Elles utilisent des calculs de statistiques spatiales utilisés dans l'analyse des localisations et organisations dans l'espace, mettant en relation proximités géographiques et ressemblances statistiques. Ces calculs peuvent également être utilisés à des fins d'exploration des relations dans l'espace.

*

**Prérequis** Manipulation des objets spatiaux telle que présentée dans le Chapitre 9 et le Chapitre 10. Analyse des semis de points, notions d'auto-corrélation spatiale.
**Description des packages utilisés** Les packages spécifiques pour la manipulation d'objets spatiaux sont présentés dans le chapitre précédent. Deux packages suplémentaires sont utilisés dans ce chapitre : spdep (spatial dependence) pour les statistiques spatiales et ICSNP (Invariant Coordinate System - Non parametric) pour le calcul du point médian.

**Description des données utilisées** Les données utilisées dans ce chapitre sont les mêmes que celles du chapitre précédent : le semis d'hôpitaux, les limites communales et les données socio-économiques attachées à ces communes.

## 11.1 Point moyen et point médian

Le chapitre 4 a présenté les principales mesures de centralité et de dispersion qui résument la distribution d'une variable quantitative. Parmi les mesures de centralité, la médiane et la moyenne sont les plus fréquemment utilisées. La médiane d'une variable est la valeur correspondant au deuxième quartile, elle découpe la distribution en deux sous-ensembles d'effectifs égaux. La moyenne d'une variable est la somme des valeurs divisée par l'effectif.

### 11.1.1 Centralité dans un espace à une dimension

La moyenne, comme la médiane, est une mesure de centralité parce que la somme des écarts entre cette valeur et toutes les valeurs prises par la variable est minimale. Plus précisément, la moyenne est la valeur qui minimise la somme du carré des écarts (1), la médiane est la valeur qui minimise la somme de la valeur absolue des écarts (2).

$$ (1) \quad \arg \min \sum_{i=1}^{n} (x_i - \bar{x})^2 $$

$$ (2) \quad \arg \min \sum_{i=1}^{n} |x_i - \bar{x}| $$

Ces propriétés de la moyenne et de la médiane peuvent être testées sur une variable quantitative quelconque. Dans cet exemple, la variable utilisée est la capacité des hôpitaux (nombre de lits) de l'espace d'étude.
Dans un premier temps, cette variable est extraite de la table attributaire de l'objet spatial. L'histogramme, la moyenne et la médiane de la variable sont calculées. La distribution est dissymétrique à gauche, la moyenne est supérieure à la médiane.

```r
hospCapacity <- publicHosp$Capacite
meanCapacity <- mean(hospCapacity)
medCapacity <- median(hospCapacity)
```

Dans un deuxième temps, une séquence de valeurs est créée qui encadre les valeurs centrales, il s'agit d'une suite d'entiers entre 0 et 1 000 avec un pas de 5. Ces 201 valeurs seront prises successivement comme valeur de référence : la somme des écarts entre chaque valeur de référence (x) et toutes les valeurs prises par la variable ($x_i$) sera calculée (écarts absolus et écarts au carré) avec une fonction de type `apply()` (cf. Section 3.4). L'ensemble est enfin regroupé dans un tableau, c'est-à-dire les séquences de valeurs avec un pas de 5, les 201 sommes de la valeur absolue des écarts, les 201 sommes du carré des écarts et le type correspondant à ces calculs (*absolute* et *square*).

```r
seqVec <- seq(0, 1000, 5)

absDif <- sapply(
  seqVec,
  function(x, xi = hospCapacity) sum(abs(xi - x))
)

squDif <- sapply(
  seqVec,
  function(x, xi = hospCapacity) sum((xi - x) ^ 2)
)

difCurves <- data.frame(
  SUMDIF = c(absDif, squDif),
  SEQVAL = c(seqVec, seqVec),
  TYPE = c(rep("ABSOLUTE", length(absDif)),
           rep("SQUARE", length(squDif)))
)
```

Enfin, la variable est représentée graphiquement avec le package `ggplot2` (cf. Section 9.3) : la somme des écarts est représentée par une
courbe, les valeurs de la moyenne et de la médiane sont représentées par des lignes verticales. On vérifie ici les propriétés de ces deux mesures de centralité qui correspondent au minimum de chacune des deux fonctions.

```r
ggplot() +
  geom_path(data = difCurves,
           aes(x = SEQVAL, y = SUMDIF, linetype = TYPE)) +
  geom_vline(xintercept = medCapacity, linetype = 1) +
  geom_vline(xintercept = meanCapacity, linetype = 2) +
  facetWrap(~ TYPE, scales = "free") +
  xlab("Valeurs prises par la mesure de référence") +
  ylab("Somme des écarts (absolus ou au carré)") +
  themeBW()
```

<figure><img src="image_83.png" /></figure>

```r
medCapacity

## [1] 395.5

seqVec[which.min(absDif)]

## [1] 395
```

Le raisonnement peut être appliqué d'un point de vue spatial dans un espace à une seule dimension. Le modèle de Hotelling en est un bon exemple : l'espace du modèle est une plage linéaire sans largeur où des
vacanciers sont répartis de façon aléatoire<fnref n="1" />. Sur cette plage, des vendeurs de glace cherchent à se rapprocher de leur clientèle potentielle. Le gla-cier qui veut maximiser ses opportunités de vente cherchera à se placer au point qui minimise la somme des distances à tous les vacanciers. Ce point est le point médian, il est facile à trouver dans un espace à une dimension mais il est plus difficile à approcher dans un espace à deux dimensions.

Dans un espace à deux dimensions, la moyenne et la médiane restent des mesures de centralité très utilisées pour décrire un semis de points. Le point moyen, ou centre de gravité, ou barycentre, est le point qui minimise le carré des distances à tous les autres points. Le point médian, quant à lui, minimise la somme des distances à parcourir. Ce point est la solution de problèmes simples de localisation, il est connu sous le nom de « point de Weber », du nom de l'économiste Alfred Weber (frère de Max Weber) qui l'a défini dans un travail sur les localisations industrielles. D'autre points centraux peuvent être utilisés, en particulier le point minimax qui minimise la distance au point le plus éloigné.

Dans la suite de cette section, le point moyen et le point médian sont calculés sur les données de Paris et de la petite couronne : le point médian sur le semis des hôpitaux et le point moyen sur le semis des communes. Le premier calcul répond à la question d'un hypocondriaque extrémiste qui veut se situer au plus près de l'ensemble des hôpitaux. Le second calcul permettra de déterminer le centre géométrique de l'espace d'étude d'abord, son centre démographique ensuite, en prenant en compte la pop-pulation des communes, enfin il permettra de voir le déplacement de ce centre au cours du temps.

### 11.1.2 Point médian dans un espace à deux dimensions

Le point médian est d'abord approché par une méthode à la fois didactique et rudimentaire. Dans un premier temps les coordonnées des hôpitaux sont extraites de l'objet spatial grâce à la fonction *coordinates()*. Ensuite, l'étendue de l'espace d'étude (*bounding box*) est récupéré avec la fonction *bbox()*. Ces deux éléments serviront de base au calcul. La dernière ligne de code utilise la fonction

<footnote n="1">1. Il existe une description similaire de ce modèle avec des magasins à placer dans une rue, ce pourquoi il a parfois été qualifié de *linear city model*.</footnote>
`fortify()` pour extraire la géométrie du polygone (contour de l'espace d'étude) et pouvoir la représenter avec le package `ggplot2` (cf. Section 9.3). Cette dernière étape n'est utile que pour la représentation cartographique, elle n'a pas de rôle à jouer dans le calcul du point médian.

```r
coordHosp <- as.data.frame(coordinates(publicHosp))
colnames(coordHosp) <- c("COORDX", "COORDY")
bboxRegion <- bbox(regBoundary)
fortBoundary <- fortify(regBoundary)
```

Il s'agit ensuite de parcourir l'étendue de l'espace d'étude en s'arrêtant régulièrement pour prendre un point de référence et calculer la somme des distances entre ce point et l'ensemble des hôpitaux. Pour cela, une grille est créée, c'est-à-dire un semis de points régulier qui couvre l'étendue de l'espace d'étude. Pour créer cette grille, on parcourt l'étendue des x et l'étendue des y avec une séquence régulière (pas de 1 000) avec le fonction `seq()`, puis ces deux séquences sont combinées avec la fonction `expand.grid()`.

```r
seqCoordX <- seqierung(bboxRegion[1, 1]),
           ceiling(bboxRegion[1, 2]),
           by = 1000)

seqCoordY <- seqierung(bboxRegion[2, 1]),
           ceiling(bboxRegion[2, 2]),
           by = 1000)

longGrid <- expand.grid(COORDX = seqCoordX,
                       COORDY = seqCoordY)
```

Finalement, une matrice de distances est calculée entre le semis de points des hôpitaux et le semis régulier (grille). C'est la fonction `rdist()` du package `fields` qui effectue ce travail. L'utilisateur pourrait envisager d'écrire lui-même une fonction qui calcule la matrice de distances, cependant une fonction de ce type est rapidement gourmande en calcul et il vaut mieux utiliser une fonction, telle que `rdist()`, qui est programmée en Fortran et qui sera bien plus rapide qu'une fonction programmée en R (cf. Section 3.3.1).
```r
library fields)
matDist <- rdist (coordHosp, longGrid)
sumDist <- apply (matDist, 2, sum)
longGrid$ABSDIST <- sumDist
```

Le résultat est cartographié avec `ggplot2`. Quatre couches sont superposées : la première est une grille (raster) dont chaque pixel représente la somme des distances entre ce pixel et l'ensemble des points qui forment le semis des hôpitaux, la deuxième est une croix correspondant au point pour lequel cette somme est minimale (approximation du point médian), la troisième est le semis des hôpitaux, la quatrième est le contour de l'espace d'étude.

```r
ggplot() +
  geom_raster(data = longGrid,
             aes(x = COORDX, y = COORDY, fill=ABSDIST)) +
  geom_point(data = longGrid[which.min(longGrid$ABSDIST),],
            aes(x = COORDX, y = COORDY),
            colour = "black", shape = 3, size = 10) +
  geom_point(data = coordHosp,
            aes(x = COORDX, y = COORDY),
            colour = "black", size = 2) +
  geomPolygon(data = fortBoundary,
             aes(x = long, y = lat, group = group),
             fill = NA, colour = "black") +
  scale_fillGradientn(colours = c("white", "black")) +
  coord_equal() +
  themeBW()
```
<figure><img src="image_84.png" /></figure>

Il serait possible avec l'algorithmme développé ci-dessus de donner une approximation plus précise du point médian en fixant une résolution plus fine à la grille mais les calculs deviendraient rapidement très lourds. Il existe bien sûr des méthodes bien plus élégantes et plus efficientes pour calculer le point médian, l'algorithmme le plus classique étant implémenté dans plusieurs packages de R, par exemple dans le package ICSNP.

```r
library(ICSNP)
spatial.median(coordHosp)
```

```r
## COORDX COORDY
```

```r
## 652242 6861734
```
### 11.1.3 Point moyen dans un espace à deux dimensions

L'une des caractéristiques de la moyenne vis-à-vis de la médiane est d'être plus sensible aux variations de valeur, que ce soit dans un espace à une ou à deux dimensions. Cette caractéristique peut s'avérer intéressante lorsqu'on cherche à mettre en évidence des évolutions temporelles et le point moyen est souvent utilisé de la sorte.

Il s'agit dans l'exemple suivant de calculer le point moyen sur le se-
mis des communes, du centroïde des communes plus précisément, mais
en tenant compte du poids de chaque commune en termes de population.
Ce qui est calculé est donc un point moyen pondéré par la population,
et ce point moyen sera calculé pour chacune de dates auxquelles on dis-
pose d'une variable de population entre 1936 et 2008. Les coordonnées
de ce point moyen ou centre de gravité ($x_g$ et $y_g$) se calculent en faisant
la moyenne pondérée ($w_i$) des coordonnées des points ($x_i$ et $y_i$)

$$x_g = \frac{\sum_{i=1}^{n} w_i x_i}{\sum_{i=1}^{n} w_i}$$

$$y_g = \frac{\sum_{i=1}^{n} w_i y_i}{\sum_{i=1}^{n} w_i}$$

Les coordonnées des centroïdes des communes sont extraites de l'objet
spatial avec la fonction `coordinates()` et le code des communes est
récupéré dans la table attributaire de cet objet. Ensuite, une jointure est
faite avec le tableau **popCom3608** qui contient les variables de population
à différentes dates.

```R
coordMuni <- data.frame(
  X = coordinates(muniBound)[, 1],
  Y = coordinates(muniBound)[, 2],
  CODGEO = muniBound@data$INSEE_COM,
  stringsAsFactors = FALSE
)
```

```R
coordPop <- merge (
  coordMuni, popCom3608,
  by = "CODGEO",
  sort = TRUE
)
```

Pour calculer automatiquement les points moyens à chaque date, une
boucle est mise en place. Les fonctions de type `apply()` sont souvent
meilleures que les boucles pour appliquer une fonction successivement à
un ensemble (cf. Section 3.4). Elles sont meilleures à la fois en termes
```r
gravityEvol <- data.frame(YEAR = colnames(coordPop)[5:13],
                         XG = xgEvol,
                         YG = ygEvol)
```

```r
plot/gravityEvol$XG, gravityEvol$YG,
    type = "b",
    pch = 20,
    asp = 1,
    xlim = c(652400, 654000))
```

```r
text(x = gravityEvol$XG,
     y = gravityEvol$YG,
     labels = gravityEvol$YEAR,
     cex = 0.8, pos = c(rep(4, 4), rep(3, 2), 4, 1, 3))
```

<figure><img src="image_85.png" /></figure>

Ce graphique peut être brièvement commenté : durant la première moitié du 20ᵉ siècle, l'Ouest de l'agglomération s'est peuplé plus rapidement
que l'Est, tendance qui s'est inversée durant la seconde moitié du siècle. Le centre de gravité pondéré par la population s'est donc déplacé vers l'Est durant la période 1936-2008, résultat cohérent avec les analyses faites à la Section 6.1. Ni le contour de l'espace d'étude, ni les limites municipales n'ont été ajoutées à cette carte parce que l'étendue de ce déplacement reste limitée : entre 1936 et 2008, le centre de gravité traverse le 4ᵉ arrondissement d'Est en Ouest.

## 11.2 Mesurer les espacements dans un semis de points

### 11.2.1 Notion de plus proche voisin et mesures associées

Lorsque l'on dispose d'un semis de points, la description de la structure de répartition du semis se fonde sur des mesures de l'espacement entre les points, résumés à l'ensemble du semis. Dans ces approches, la distance joue un rôle fondamental. Il y a plusieurs façons de mesurer l'éloignement entre deux points, c'est ici la façon la plus simple qui est utilisée : la distance euclidienne à vol d'oiseau. Pour certaines questions, il peut être nécessaire d'utiliser une distance plus appropriée.

À la base de la description des espacements entre points du semis, il faut calculer une matrice de distance entre tous les points, comme cela a été fait dans la section précédente. Il est courant de travailler non sur l'ensemble de la matrice mais sur la notion de distance au plus proche voisin. Cette distance peut aussi être généralisée à des ordres supérieurs comme par exemple la distance moyenne aux 3 ou aux 5 plus proches voisins.

La distance au plus proche voisin peut être utilisée dans un cadre décritif pour décrire des contextes locaux. Elle peut aussi être utilisée dans une approche modélisatrice, par exemple comme variable explicative dans une régression multiple. Elle peut enfin être utilisée pour caractériser le semis de points en le comparant à une distance théorique obtenue par un processus ponctuel aléatoire, un processus de Poisson par exemple.

Le calcul de la distance moyenne au plus proche voisin peut être traité de deux façons :
— Calculer une matrice de distance entre tous les points ou distancier (cf. Section 3.1.2). Cette méthode a l'avantage de créer un objet complet à partir duquel un grand nombre d'indicateurs peuvent être calculés (distance au plus proche voisin, distance moyenne à tous les voisins). Elle a cependant l'inconvénient de créer un objet surdimensionné par rapport à l'utilisation que l'on peut en faire, objet qui peut être très lourd si le semis comporte un très grand nombre de points.

— Utiliser une suite de fonctions implémentées dans le package spdep spécialisé dans les statistiques spatiales. L'inconvénient de cette méthode est qu'il faut manipuler des types d'objets spécifiques avec des fonctions qui souvent ne servent qu'à passer d'un type d'objet à un autre.

### 11.2.2 Calcul d'une matrice de distance

Le distancier est calculé à partir des coordonnées contenues dans le tableau de données en utilisant la fonction `dist()`. Plusieurs distances sont proposées (distance de Manhattan, de Canberra, etc.), c'est ici la distance euclidienne qui est choisie. Utilisée comme suit, elle permet de construire le triangle inférieur hors diagonale d'une matrice de distance. Elle contient dans ce cas $n(n-1)/2$ éléments. Il faut ensuite transformer cet objet en un objet matrice $n×n$ en utilisant la fonction `as.matrix()`.

```r
distHosp <- distcasoHosp, method = "euclidean")
matDistHosp <- as.matrix(distHosp)
```

On effectue ensuite le calcul de la distance au plus proche voisin pour chaque hôpital en ayant eu soin au préalable d'éliminer du calcul les éléments de la diagonale (en leur donnant la valeur NA), sans quoi chaque hôpital serait son plus proche voisin à une distance de 0.

```r
matDistHosp[matDistHosp == 0] <- NA
minDist <- apply號 (matDistHosp, 1, min, na.rm = TRUE)
mean (minDist)
```

```r
## [1] 1392
```
La valeur moyenne de la distance au plus proche voisin est donc de 1 392 mètres. Indépendamment de la méthode des plus proches voisins, on peut s'intéresser à différentes statistiques des distances de chaque hôpital aux autres, en effectuant des statistiques univariées sur chacune des distributions. Les lignes qui suivent permettent de créer un tableau à n lignes et 7 colonnes pour les différents indicateurs créés par la fonction `summary()`. Ce tableau peut être transposé et trié pour en faciliter la lecture.

```
summaryDist <- apply matDistHosp, 1, summary, na.rm = TRUE)
head(t(summaryDist))
```

<table><thead><tr><th></th><th>Min.</th><th>1st Qu.</th><th>Median</th><th>Mean</th><th>3rd Qu.</th><th>Max.</th><th>NA's</th></tr></thead><tbody><tr><td>## 1</td><td>1180</td><td>4190</td><td>8060</td><td>8000</td><td>11300</td><td>16900</td><td>1</td></tr><tr><td>## 2</td><td>922</td><td>6020</td><td>8040</td><td>9130</td><td>11800</td><td>21400</td><td>1</td></tr><tr><td>## 3</td><td>324</td><td>4560</td><td>7510</td><td>8090</td><td>11200</td><td>18700</td><td>1</td></tr><tr><td>## 4</td><td>4130</td><td>13900</td><td>17500</td><td>16900</td><td>21300</td><td>26400</td><td>1</td></tr><tr><td>## 5</td><td>1450</td><td>5120</td><td>8110</td><td>9350</td><td>12700</td><td>19100</td><td>1</td></tr><tr><td>## 6</td><td>2020</td><td>5890</td><td>9440</td><td>10300</td><td>14100</td><td>18000</td><td>1</td></tr></tbody></table>

### 11.2.3 Fonctions spécifiques du package spdep

Le package `spdep` contient un ensemble de fonctions dédiées aux statistiques spatiales. Ces fonctions créent des objets spécifiques, par exemple des listes de voisins, qui doivent souvent être transformées pour alimenter d'autres fonctions du packages, par exemple pour représenter graphiquement les points et leurs voisins.

Dans un premier temps, la fonction `knearneigh()` permet de calculer de manière générique la liste des k plus proches (near) voisins (neighbours) d'un ensemble. Lorsque cette fonction est appliquée à un objet spatial, la fonction récupère directement l'ensemble des informations nécessaires aux calculs et travaille sur le vecteur `$coords$. La proximité est ici entendue au sens de la distance euclidienne. Si on spécifie l'argument `longlat`, la fonction peut aussi faire le calcul sur des coordonnées géographiques (longitudes-latitudes) en degrés décimaux. Dans cet exemple, on commence par créer une liste de plus proches voisins d'ordre 5 (argument k).
```r
library(spdep)
listNearNei <- knearneigh(publicHosp@coords, k = 5)
class(listNearNei)
```

L'objet qui vient d'être créé est une liste de type `knn` dont le premier élément est le vecteur des plus proches voisins. Ce schéma résume les types d'objets manipulés avec `spdep` et les fonctions qui assurent leur transformation.

<figure><img src="image_86.png" /><figcaption>FIGURE 11.1 – Objets et fonctions du package spdep</figcaption></figure>

Le premier élément de cet objet `knn` est, pour chaque hôpital, la liste des cinq hôpitaux les plus proches. À partir de cet objet, la distance est calculée entre chaque hôpital et ses cinq plus proches voisins grâce à la fonction `nbdists()`. Ce calcul ne peut être fait directement sur l'objet `knn` et il faut le transformer en objet de type `nb`. La fonction `nbdists()` demande également de préciser les coordonnées des hôpitaux contenues dans l'objet spatial (`@coords`).

```r
nearNei <- knn2nb(listNearNei,
                  row.names = publicHosp$CodHop)
distNei <- nbdists(nearNei, publicHosp@coords)
class(distNei)
```

La structure de liste renvoyée par la fonction `nbdists()` est simplifiée avec la fonction `unlist()` qui crée un vecteur numérique regroupant toutes les distances. Ce vecteur donne la distance entre tous les couples d'hôpitaux voisins à un ordre 5. Cette variable peut être explorée avec les méthodes présentées dans le Chapitre 4.
```
distNei <- unlist(distNei)
```

```
summary(distNei)
```

```
## Min. 1st Qu. Median Mean 3rd Qu. Max.
```

```
## 0    1550   2450   2700   3890   9580
```

```
hist(distNei, breaks = 20,
```

```
main = "Distance aux 5 plus proches voisins",
```

```
col = "black",
```

```
border = "white",
```

```
xlab = "Distance",
```

```
ylab = "Fréquence")
```

Distance aux 5 plus proches voisins

<figure><img src="image_87.png" /></figure>

Il est parfois intéressant de fixer une distance seuil et non un nombre fixe de voisins. Par exemple, en définissant un seuil de 2 km, il est possible de créer, pour chacun des points, une liste de voisins situés à une distance inférieure à deux kilomètres. La fonction `dnearneigh()` crée cette liste, à partir des coordonnées des points et des seuils de distance (arguments `d1` et `d2`), ici de 0 à 2000 mètres.

```
neiTwoKm <- dnearneigh(publicHosp@coords,
                       d1 = 0,
                       d2 = 2000,
                       row.names=publicHosp$CodHop)
```

Cette fonction renvoie directement un objet de type `nb`.
### 11.2.4 Visualisations des voisinages

Il est utile de compléter les mesures qui viennent d’être calculées avec des représentations cartographiques. Les objets de types nb, qui sont des listes de voisins, peuvent être représentés sous la forme de graphes et ajoutés à la cartographie du semis de points. L’exemple ci-dessous montre les réseaux de voisins à une distance inférieure à deux kilomètres.

```r
plot(muniBound, axes = TRUE,
     col = "grey", border = "white")
plot(publicHosp, col = "black", pch = 20, add = TRUE)
plot(neiTwoKm, publicHosp@coords,
     pch = 20, col = "black", add = TRUE)
```

<figure><img src="image_88.png" /></figure>
# 11.3 Autocorrélation spatiale et ressemblances locales

## 11.3.1 La notion d'autocorrélation spatiale

Les méthodes de statistique spatiale qui viennent d'être introduites s'appliquent à des semis de points. Il existe également toute une gamme de méthodes qui s'appliquent à des phénomènes observés dans des zones, par exemple dans des maillages territoriaux. L'autocorrélation spatiale en fait partie, elle sera travaillée ici sur le maillage communal.

L'autocorrélation spatiale consiste, pour un phénomène donné, à mesurer l'intensité de la relation entre la proximité des lieux et leur degré de ressemblance au regard de ce phénomène. Le fondement de ces mesures est la mise en rapport d'une variabilité locale avec la variabilité globale.

Les premiers développements de ces mesures concernaient des voisinages définis par la contiguïté sur un maillage. La variance locale peut alors s'écrire comme une variance globale pondérée par un poids $w_{ij}$ affecté à chacun des couples $(i, j)$ : $w_{ij} = 1$ si $i$ et $j$ sont des zones contiguës, $w_{ij} = 0$ sinon. Dans ce cas $w$ est une matrice de contiguïté.

Ce principe peut être étendu à la notion de voisinage au sens large, en considérant comme voisines des unités spatiales distantes de moins de 5km, de moins d'une heure ou appartenant à la même région administrative. Dans ce cas, $w_{ij}$ sera aussi une variable binaire et $w$ une matrice de contiguïté. Il est enfin possible de pondérer les couples par l'inverse de la distance (ou du temps), auquel cas $w_{ij}$ sera une variable continue et $w$ une matrice de poids.

Les deux indices d'autocorrélation les plus anciens et les plus utilisés sont l'indice I de Moran (1950) et l'indice C de Geary (1954).

$$I = \frac{N}{\sum_{i} \sum_{j} w_{ij}} \frac{\sum_{i} \sum_{j} w_{ij} (X_{i} - \bar{X})(X_{j} - \bar{X})}{\sum_{i} (X_{i} - \bar{X})^{2}}$$

$$C = \frac{(N - 1)}{\sum_{i} \sum_{j} w_{ij}} \frac{\sum_{i} \sum_{j} w_{ij} (X_{j} - \bar{X})}{\sum_{i} (X_{i} - \bar{X})^{2}}$$
Où N est le nombre d'individus, $X_i$ et $X_j$ sont les valeurs de la variable X en i et j, $\bar{X}$ est la moyenne de la variable X et $w_{ij}$ est une pondération correspondant à la définition du voisinage.

Des développements ont ensuite été proposés pour travailler sur les contributions locales des zones à cette mesure de l'autocorrélation spatiale. Il s'agit des indicateurs locaux d'association spatiale (LISA - local indicators of spatial autocorrelation). L'indice local de Moran correspond ainsi à la contribution de chaque unité à l'autocorrélation globale.

Il s'agit, pour une unité spatiale i, du produit entre l'écart à la moyenne de cette unité spatiale ($z_i$) et la somme des écarts des unités voisines ($z_j$) pondérés par la matrice de contiguïté ou de poids ($w_{ij}$).

$$I_i = z_i \sum_j w_{ij} z_j$$

Cet indicateur permet d'identifier les situations localement homogènes ($I_i$ faible), et les situations localement hétérogènes ($I_i$ élevé). Selon le phénomène étudié, on cherchera à isoler l'une et/ou l'autre de ces situations.

### 11.3.2 Répartition des cadres et des ouvriers

L'exemple suivant cherche à caractériser la répartition spatiale des cadres et des ouvriers dans l'espace d'étude.

La fonction moran.test() du package spdep travaille directement sur les objets de type nb qui donnent la liste des unités voisines. Pour travailler sur un objet de ce type, il faut au préalable transformer l'objet spatial contenant le découpage communal (SpatialPolygonsDataFrame) en objet de type nb.

C'est la fonction poly2nb() qui effectue ce travail. Plusieurs arguments doivent être précisés : l'argument queen (déplacement de la reine au jeu d'échecs) indique qu'un seul sommet commun entre deux polygones est suffisant pour les considérer comme voisins ou qu'il faut un segment commun. L'argument snap permet de fixer la distance en des sous de laquelle deux points d'un polygone sont considérés comme étant
superposés. En effet, il arrive souvent de travailler avec des géométries
de mauvaise qualité (erreurs de topologie), où les sommets formant des
segments ne coïncident pas exactement.

Dans un premier temps, les variables du tableau externe socEco9907 sont jointes aux données attributaires de l'objet spatial. Cette manipulation est faite avec la fonction AttribJoin() créée dans le chapitre précédent (cf. Section 10.4.1). Ensuite, l'objet nb est créé. Enfin l'indice de Moran est calculé sur les deux variables mentionnées : la proportion de cadres et la proportion d'ouvriers dans les communes. La fonction moran.test() prend comme argument la variable à l'étude et la matrice de contiguïté ou de poids. Celle-ci est créée directement avec la fonction nb2listw().

```r
muniBound <- AttribJoin(df = socEco9907,
                          spdf = muniBound,
                          df.field = "CODGEO",
                          spdf.field = "INSEE_COM")

nbCom <- poly2nb(pl = muniBound,
                 row.names = muniBound$INSEE_COM,
                 snap = 50,
                 queen = TRUE)

moran.test(x=muniBound$PCAD07,
           listw = nb2listw nbCom))

moran.test(x=muniBound$POUV07,
           listw = nb2listw nbCom))
```

La fonction renvoie la valeur de l'autocorrélation spatiale (0,75 pour les cadres et 0,69 pour les ouvriers), la valeur de l'écart standardisé (14,59 et 13,43) ainsi que la *p-value*. La proportion de cadres et la proportion d'ouvriers présentent une autocorrélation positive, synonyme dans ce cas de ségrégation spatiale.

Les valeurs de l'indice de Moran peuvent être représentées graphiquement et cartographiquement. Ici la relation est représentée par un nuage de points (*Moran scatterplot*) croisant les valeurs en *i* avec la moyenne des valeurs du voisinage de *i*.
La fonction `moran.plot()` demande, comme la fonction `moran.test()`, une matrice de poids. Celle-ci est produite, comme précédemment, par la fonction `nb2listw()`. L'argument `style` est précisé même s'il s'agit de sa valeur par défaut : il indique que les poids doivent être standardisés par lignes. C'est cet argument qui permet de lire les valeurs sur l'axe y du graphique comme une moyenne pondérée par la distance de la variable à l'étude.

```
moran.plot(x = muniBound$PCAD07,
         listw = nb2listw nbCom, style = "W"),
         xlab = "% de cadres en 2007",
         ylab = "% de cadres dans le voisinage")
```

<figure><img src="image_89.png" /></figure>
```
moran plot(x = muniBound$POUV07,
          listw = nb2listw nbCom, style = "W"),
          xlab = "% d'ouvriers en 2007",
          ylab = "% d'ouvriers dans le voisinage")
```

<figure><img src="image_90.png" /></figure>

Ce type de graphique est très riche, il permet de visualiser l'autocorré-
lation, de classer les unités spatiales en termes de profils et de détecter
les cas exceptionnels. Les droites horizontale et verticale sont centrées sur
les moyennes de chacune des deux variables : ces droites dessinent des
cadrants. Là encore, il sera intéressant de cartographier ces profils <fnref n="1" />.

Pour la proportion de cadres par exemple, ces cadrants guideront la lec-
ture en distinguant les communes comportant peu de cadres et peu de

<footnote n="1">1. Voir par exemple la communication « R et espace » (Commenges, Cura, Mathian) présentée au *Semin-R* du Museum National d'Histoire Naturelle, disponible à l'adresse suivante : http://rug.mnhn.fr/semin-r.</footnote>
cadres dans le voisinage, les communes comportant beaucoup de cadres et beaucoup de cadres dans le voisinage. Les cadrants opposés sont particulièrement intéressants puisqu'ils montrent les discontinuités : les communes comportant beaucoup de cadres situées dans un voisinage qui en comportent peu (Saint-Maur-des-Fossés par exemple, code 94068) et celles qui comportent très peu de cadres dans un voisinage qui en comporte beaucoup (Bagneux par exemple, code 92007).

Ce chapitre a présenté des techniques simples d'analyse spatiale, techniques permettant de caractériser la répartition spatiale d'objets et de phénomènes. Les mesures de centralité, point moyen et point médian, peuvent servir à résumer la distribution d'un semis de points dans une approche descriptive, elles peuvent être utilisées pour saisir des évolutions du semis de point, elles peuvent enfin être utilisées dans un cadre opérationnel (localisation optimale pour un hôpital). Les fonctions du *package spdep* sont très intéressantes pour analyser plus finement les phénomènes spatiaux, par exemple la répartition des semis de points ou l'homogénéité d'un phénomène observé dans des unités spatiales.
<figure><img src="image_91.png" /></figure>
# Annexes
<figure><img src="image_92.png" /></figure>
# ANNEXE A

## Bibliographie thématique

Pour approfondir les méthodes abordées dans le manuel, voici une bibliographie indicative classée par thèmes. Cette liste reste volontairement très brève au regard de la très vaste littérature existante sur l'analyse de données et le logiciel R.

### Programmation

BAATH R. (2012) "The state of naming conventions in R", *The R Journal*, Vol. 4/2, pp.74-75.

CHAMBERS J.M. (2008) *Software for data analysis. Programming with R*, Springer.

GENOLINI CH. (nd) *Petit manuel de S4. Programmation orientée objet sous R*, http://christophe.genolini.free.fr/webTutorial/index.php.

MUENCHEN R. (2008) *R for SAS and SPSS users*, Berlin, Springer.
MUENCHEN R., HILBE J.M. (2012) *R for Stata users*, Berlin, Springer.

WICKHAM H. (2014), *Advanced R*, http://adv-r.had.co.nz.

WICKHAM H. (2011), “The split-apply-combine strategy for data analysis”, *Journal of Statistical Software*, vol.40, num.1, pp. 1-29.

## Analyse statistique

CORNILLON P.-A., MATZNER-LØBER E. (2011) *Régression avec R*, Paris, Springer.

DUMOLARD P., DUBUS N., CHARLEUX L. (2003) *Les statistiques en géographie*, Presses Universitaires de Rennes.

HUSSON F., LÊ S., PAGÈS J. (2009) *Analyse de données avec R*, Presses Universitaires de Rennes.

LEBART L., MORINEAU A., PIRON M. (1995) *Statistique exploratoire multidimensionnelle*, Dunod.

SANDERS L. (1989) *L'analyse statistique des données en géographie*, Montpellier, Alidade-RECLUS.

## Analyse de graphes

BARTHELEMY M. (2011) “Spatial Networks”, *Physics Reports*, vol.499, pp.1-101, http://arxiv.org/abs/1010.0302.

GROUPE FMR (2010-2013) *Synthèses*, http://halshs.archives-ouvertes.fr/FMR.

HANNEMAN, R.A., RIDDLE M. (2005) *Introduction to social network methods*, Riverside, University of California, http://faculty.ucr.edu/~hanneman.

NEWMAN M. (2010) *Networks : an introduction*, Oxford, Oxford University Press.

WASSERMAN S., FAUST K. (1994), *Social network analysis. Methods and applications*, Cambridge, Cambridge University Press.
Visualisation graphique et cartographique

BÉGUIN M., PUMAIN D. (2000) *La représentation des données géographiques*, Paris, Armand Colin.

CHEN C., HÄRDLE W., UNWIN A. (dir.) (2008) *Handbook of data visualization*, Berlin, Springer.

TUFTE E. R. (1983) *The visual display of quantitative information*, Cheshire, Graphics Press.

WICKHAM H. (2011), “Graphical criticism: some historical notes”, *Journal of Computational and Graphical Statistics*, vol.22, num.1, pp. 38-44,
http://vita.had.co.nz/articles.html.

ZANIN CH., TRÉMÉLO M.-L. (2003) *Savoir faire une carte. Aide à la conception et à la réalisation d'une carte thématique univariée*, Paris, Belin.

Analyse spatiale

ANSELIN L. (1995) “Local indicators of spatial association”, *Geographical analysis*, 27/2, pp.93-115.

BIVAND R.S., PEBESMA E.J., Gómez-Rubio V. (2008) *Applied spatial data analysis with R*, Springer.

PUMAIN D., SAINT-JULIEN TH. (2010) *L'analyse spatiale. Les localisations*, tome 1, 2e éd., Paris, Armand Colin.

PUMAIN D., SAINT-JULIEN TH. (2010) *L'analyse spatiale. Les interactions*, tome 2, 2e éd., Paris, Armand Colin.

ZANINETTI J.-M. (2005) *Statistique spatiale : méthodes et applications géomatiques*, Paris, Hermès Science Publications.
<figure><img src="image_93.png" /></figure>
# ANNEXE B

## Liste des packages utilisés

Un certain nombre *packages* sont installés et chargés tout au long de la progression du manuel. En voici la liste exhaustive accompagnée d'une brève description. Lorsqu'un *package* est utilisé dans le manuel pour une utilisation unique et très spécifique, seule cette information est donnée, c'est le cas par exemple du package *ICSNP* :

- ade4 : analyse géométrique des données multivariées

- classInt : discrétisation des variables continues

- cluster : méthodes de classification

- dplyr : manipulation des données

- FactoClass : analyse géométrique des données multivariées

- fields : outils pour le traitement de données spatiales

- ggplot2 : représentations graphiques et cartographiques

- HistData : jeux de données historiques

- Hmisc : mélanges de fonctions courantes

- ICSNP : calcul du point médian

- igraph : analyse de graphes
— mapproj : gestion des systèmes de projection des données spatiales

— maptools : outils de cartographie

— plyr : manipulation des données

— raster : traitement et affichage de données raster

— rCarto : outils de cartographie

— RColorBrewer : création de palettes de couleur

— reshape2 : transposition individus-variables

— rgdal : manipulation des objets spatiaux

— rgeos : manipulation des objets spatiaux

— scales : format des axes des graphiques

— sp : manipulation des objets spatiaux

— spdep : statistiques spatiales

— sqldf : exécution de requêtes SQL
<table><tbody><tr><td>2.2.1 Le répertoire de travail</td><td>25</td></tr><tr><td>2.2.2 Les formats de données</td><td>25</td></tr><tr><td>2.2.3 Codage des caractères et des séparateurs</td><td>25</td></tr><tr><td>2.2.4 Importer des données en format texte</td><td>26</td></tr><tr><td>2.2.5 Exporter des données en format texte</td><td>26</td></tr><tr><td>2.2.6 Importation et exportation depuis/vers d'autres formats</td><td>27</td></tr><tr><td>2.3 Recoder et trier</td><td>27</td></tr><tr><td>2.3.1 Sélectionner et recoder</td><td>27</td></tr><tr><td>2.3.2 Trier</td><td>30</td></tr><tr><td>2.4 Manipulation avancée</td><td>32</td></tr><tr><td>2.4.1 Superposition, concaténation, jointure</td><td>32</td></tr><tr><td>2.4.2 Agrégations et traitements par blocs</td><td>34</td></tr><tr><td>2.4.3 Transposition variables-observations</td><td>37</td></tr><tr><td>2.5 Packages polyvalents</td><td>41</td></tr><tr><td><b>3 Introduction à la programmation</b></td><td><b>45</b></td></tr><tr><td>3.1 Structures itératives</td><td>46</td></tr><tr><td>3.1.1 Exemple 1 : simple boucle for</td><td>47</td></tr><tr><td>3.1.2 Exemple 2 : double boucle for</td><td>48</td></tr><tr><td>3.2 La structure conditionnelle if...else</td><td>51</td></tr><tr><td>3.3 Les fonctions</td><td>52</td></tr><tr><td>3.3.1 Premier aperçu sur les fonctions</td><td>52</td></tr><tr><td>3.3.2 Discrétisation automatique</td><td>55</td></tr><tr><td>3.3.3 Calcul de l'équilibre de Wardrop</td><td>57</td></tr><tr><td>3.3.4 Calcul des vols de Syracuse</td><td>61</td></tr><tr><td>3.4 L'application de fonctions sur des ensembles</td><td>65</td></tr><tr><td><b>4 Analyse univariée</b></td><td><b>71</b></td></tr><tr><td>4.1 Calculs simples et recodages</td><td>72</td></tr><tr><td>4.2 Résumés statistiques</td><td>74</td></tr><tr><td>4.3 Représentation graphique des distributions statistiques</td><td>76</td></tr><tr><td><b>5 Analyse bivariée</b></td><td><b>83</b></td></tr><tr><td>5.1 Préambule</td><td>84</td></tr><tr><td>5.1.1 Types de variables</td><td>84</td></tr><tr><td>5.2 Relation entre deux variables quantitatives</td><td>84</td></tr><tr><td>5.2.1 Représentation graphique de la relation</td><td>85</td></tr><tr><td>5.2.2 Calcul et significativité des corrélations</td><td>88</td></tr><tr><td>5.2.3 Régression linéaire simple</td><td>90</td></tr><tr><td>5.3 Relation entre une variable quantitative et une variable qualitative</td><td>94</td></tr><tr><td>5.3.1 Représentation graphique de la relation</td><td>94</td></tr><tr><td>5.3.2 Analyse de la variance à un facteur</td><td>96</td></tr><tr><td>5.4 Relation entre deux variables qualitatives</td><td>99</td></tr></tbody></table>
<table><tbody><tr><td>5.4.1 Représentation graphique de la relation</td><td>99</td></tr><tr><td>5.4.2 Analyse de la répartition observée</td><td>101</td></tr><tr><td>5.4.3 Analyse des écarts à l'indépendance</td><td>103</td></tr><tr><td><b>6 Analyses factorielles</b></td><td><b>105</b></td></tr><tr><td>6.1 Analyse en composantes principales</td><td>106</td></tr><tr><td>6.1.1 Réaliser l'analyse</td><td>106</td></tr><tr><td>6.1.2 Examiner et interpréter les résultats</td><td>107</td></tr><tr><td>6.1.3 Analyse des contributions</td><td>111</td></tr><tr><td>6.1.4 Analyse des qualités de représentation</td><td>113</td></tr><tr><td>6.2 Analyse factorielle des correspondances</td><td>114</td></tr><tr><td>6.2.1 Réaliser l'analyse</td><td>115</td></tr><tr><td>6.2.2 Interpréter les résultats</td><td>116</td></tr><tr><td><b>7 Méthodes de classification</b></td><td><b>121</b></td></tr><tr><td>7.1 Rappels théoriques</td><td>122</td></tr><tr><td>7.2 Classification ascendante hiérarchique</td><td>126</td></tr><tr><td>7.2.1 Coupure de l'arbre et description des classes</td><td>129</td></tr><tr><td>7.3 Couplage analyse factorielle - classification</td><td>134</td></tr><tr><td>7.3.1 Classification sur une métrique du χ²</td><td>134</td></tr><tr><td>7.3.2 Description des classes</td><td>135</td></tr><tr><td><b>8 Analyse de graphes</b></td><td><b>141</b></td></tr><tr><td>8.1 Aperçu et préparation des données</td><td>142</td></tr><tr><td>8.2 Création et exploration du graphe</td><td>145</td></tr><tr><td>8.2.1 Importation et exportation de graphes</td><td>145</td></tr><tr><td>8.2.2 Création d'un graphe avec R</td><td>147</td></tr><tr><td>8.2.3 Calcul de mesures globales</td><td>147</td></tr><tr><td>8.3 Mesures locales et manipulation du graphe</td><td>148</td></tr><tr><td>8.4 Cliques et communautés</td><td>151</td></tr><tr><td>8.5 Visualisations</td><td>155</td></tr><tr><td><b>9 Focus sur la visualisation graphique</b></td><td><b>159</b></td></tr><tr><td>9.1 Exportation des tableaux et des graphiques</td><td>161</td></tr><tr><td>9.2 Gestion des couleurs</td><td>163</td></tr><tr><td>9.3 Introduction à ggplot2</td><td>165</td></tr><tr><td>9.3.1 Intérêt vis-à-vis des fonctions graphiques classiques</td><td>165</td></tr><tr><td>9.3.2 Aperçu de la syntaxe</td><td>166</td></tr><tr><td>9.3.3 Variables visuelles</td><td>168</td></tr><tr><td>9.3.4 Construction d'une planche</td><td>171</td></tr><tr><td>9.3.5 Variables de regroupement</td><td>172</td></tr></tbody></table>
<table><tbody><tr><td><strong>10 Introduction aux objets spatiaux et à la cartographie</strong></td><td><strong>175</strong></td></tr><tr><td>10.1 R pour la cartographie</td><td>178</td></tr><tr><td>10.2 Prise en main des données spatiales</td><td>179</td></tr><tr><td>10.2.1Construction d'un objet spatial</td><td>179</td></tr><tr><td>10.2.2Importation et prise en main des objets spatiaux</td><td>183</td></tr><tr><td>10.2.3Le système de projection</td><td>186</td></tr><tr><td>10.2.4Exportations</td><td>188</td></tr><tr><td>10.3 Cartographie des objets ponctuels</td><td>188</td></tr><tr><td>10.4 Cartographie des objets zonaux</td><td>193</td></tr><tr><td>10.4.1Cartographie avec la fonction plot()</td><td>194</td></tr><tr><td>10.4.2Cartographie avec la fonction ggplot()</td><td>198</td></tr><tr><td>10.5 Cartes take away</td><td>202</td></tr><tr><td>10.6 Cartographie à la chaîne</td><td>203</td></tr><tr><td><strong>11 Initiation aux statistiques spatiales</strong></td><td><strong>209</strong></td></tr><tr><td>11.1 Point moyen et point médian</td><td>210</td></tr><tr><td>11.1.1Centralité dans un espace à une dimension</td><td>210</td></tr><tr><td>11.1.2Point médian dans un espace à deux dimensions</td><td>213</td></tr><tr><td>11.1.3Point moyen dans un espace à deux dimensions</td><td>217</td></tr><tr><td>11.2 Mesurer les espacements dans un semis de points</td><td>220</td></tr><tr><td>11.2.1Notion de plus proche voisin et mesures associées</td><td>220</td></tr><tr><td>11.2.2Calcul d'une matrice de distance</td><td>221</td></tr><tr><td>11.2.3Fonctions spécifiques du package spdep</td><td>222</td></tr><tr><td>11.2.4Visualisations des voisinages</td><td>225</td></tr><tr><td>11.3 Autocorrélation spatiale et ressemblances locales</td><td>226</td></tr><tr><td>11.3.1La notion d'autocorrélation spatiale</td><td>226</td></tr><tr><td>11.3.2Répartition des cadres et des ouvriers</td><td>227</td></tr></tbody></table>

Annexes 235

A Bibliographie thématique 235

B Liste des packages utilisés 239