# Introduction à la géomatique pour le statisticien :
quelques concepts et outils innovants de gestion,
traitement et diffusion de l'information spatiale

François Smécurbe, Elise Coudin

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

## To cite this version:

François Smécurbe, Elise Coudin. Introduction à la géomatique pour le statisticien : quelques concepts et outils innovants de gestion, traitement et diffusion de l'information spatiale. 2022. hal-05376241

HAL Id: hal-05376241

https://insee.hal.science/hal-05376241v1

Preprint submitted on 21 Nov 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

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

Distributed under a Creative Commons Attribution 4.0 International License
# Introduction à la géomatique pour le statisticien : quelques concepts et outils innovants de gestion, traitement et diffusion de l'information spatiale

## Documents de travail

N° 2022-01 – Février 2022

<figure><img src="image_3.png" /></figure>
<figure><img src="image_4.png" /></figure>
INSTITUT NATIONAL DE LA STATISTIQUE ET DES ÉTUDES ÉCONOMIQUES

Série des documents de travail « Méthodologie Statistique »
de la Direction de la Méthodologie et de la Coordination Statistique et Internationale

M 2022/01

Introduction à la géomatique pour le statisticien :
quelques concepts et outils innovants de gestion, traitement et
diffusion de l'information spatiale

François SEMECURBE

Élise COUDIN

Insee

Février 2022

Direction de la méthodologie et de la coordination statistique et internationale

Département des Méthodes Statistiques -Timbre L001 -

88 Avenue Verdier - CS 70058 - 92541 Montrouge Cedex - France -

Tél. : 33 (1) 87 69 55 00 - E-mail : -DG75-L001@insee.fr - Site Web Insee : http://www.insee.fr

Ces documents de travail ne reflètent pas la position de l'Insee et n'engagent que leurs auteurs.

Working papers do not reflect the position of INSEE but only their author's views.
<figure><img src="image_5.png" /></figure>
# Introduction à la géomatique pour le statisticien :
quelques concepts et outils innovants de gestion,
traitement et diffusion de l'information spatiale*

François Sémécurbe‡ et Elise Coudin§

Janvier 2022

## Résumé

Ce document vise une présentation simple à l'attention des statisticiens des outils géomatiques récents qui permettent de stocker, traiter et diffuser l'information spatiale. Les logiciels comme R ou Python intègrent désormais les caractéristiques géographiques rendant plus accessibles leur traitement. Pour autant, le foisonnement des technologies, les possibilités offertes par le web, les technologies web, sont autant d'obstacles à dépasser pour celles et ceux souhaitant réaliser des cartographies thématiques percutantes. Ce document propose une présentation unifiée des concepts géomatiques, avec des extraits de code en R, Python et PostGIS. Il se concentre sur les données vectorielles et décrit les traitements et manipulations classiques à connaître pour construire une statistique spatiale. Il aborde aussi les outils et les flux permettant une diffusion dynamique (cartes interactives) de l'information spatiale. Il discute enfin le rôle de la spatialisation dans la représentation des données statistiques.

**Mots clés** : géomatique, cartographie dynamique, cartographie thématique, statistiques spatiales, données spatiales, systèmes d'information géographique.

<footnote>*Ce document s'appuie sur et complète les travaux du groupe de travail Great sur les nouveaux outils cartographiques qui s'est tenu en 2019. Les auteurs remercient Julie Djiriguian, Yves-Laurent Benichou, Sébastien Calvet, Sébastien Novella, Fabrice Garnes, David Levy, Vianney Costemalle, Sonia Oujia, Marc Branchu, Julien Jamme, Arlindo Dos Santos, Aurélien D'Isanto, Frédéric Comte, Christophe Yon, Lionel Cacheux, Olivier Levitt pour leur contribution aux travaux de ce groupe. Julie Djiriguian et François Sémécurbe ont poursuivi ces travaux sous la forme d'un site web de prise en main de ces outils, dont certaines parties de ce document s'inspirent très largement. Les auteurs remercient à ce titre une deuxième fois Julie Djiriguian. Ce document a bénéficié des relectures d'Aurélien Fortin, Clément Guillot, Benoit Hurpeau, Patrick Sillard que les auteurs remercient. Elise Coudin tient à préciser qu'elle a tenu un rôle de relectrice active plus qu'autre chose et que le mérite revient au premier auteur.</footnote>

<footnote>‡ Insee et ThèMA - Besançon</footnote>

<footnote>§ Insee et CREST</footnote>
# Introduction to geomatics for statisticians: some innovative concepts and tools for storing, processing and disseminating spatial information

## Abstract

This document is intended as a simple presentation for statisticians of recent geomatics tools, which allow the storage, processing and dissemination of spatial information. Softwares like R or Python now integrate directly geographical characteristics making their processing more accessible. However, the proliferation of technologies, the possibilities offered by the web, web technologies, are all obstacles to overcome for those wishing to produce powerful thematic maps. This document offers a unified presentation of geomatics concepts, with code snippets in R, Python and PostGIS. It focuses on vector data and describes the classic treatments and manipulations to know to build a spatial statistic. It also addresses the tools and flows allowing dynamic dissemination (interactive maps) of spatial information. Finally, it discusses the role of spatialization in the representation of statistical data.

**Keywords:** geomatics, dynamic cartography, thematic cartography, spatial statistics, spatial data, geographic information systems.
# Table des matières

<table><tr><td>Introduction</td><td>8</td></tr><tr><td>1 Les données spatiales</td><td>12</td></tr><tr><td>1.1 La géographie</td><td>12</td></tr><tr><td>1.2 Stockage des données spatiales</td><td>13</td></tr><tr><td>1.2.1 Fichiers cartographiques</td><td>14</td></tr><tr><td>1.2.2 PostGIS, une base de données spatiales</td><td>16</td></tr><tr><td>1.2.3 Système de projection</td><td>19</td></tr><tr><td>1.3 Exemple : création d'une couche carroyée à partir de données ponctuelles</td><td>21</td></tr><tr><td>2 Traitement des données spatiales</td><td>24</td></tr><tr><td>2.1 Fonctions et transformations géométriques élémentaires</td><td>24</td></tr><tr><td>2.1.1 Transformation de données ponctuelles en données surfaciques</td><td>24</td></tr><tr><td>2.1.2 Simplification de données surfaciques ou transformation de données surfaciques en données ponctuelles</td><td>27</td></tr><tr><td>2.2 Transformations complexes mettant en relation plusieurs éléments d'une couche ou de plusieurs couches</td><td>29</td></tr><tr><td>2.2.1 Agrégation des données spatiales</td><td>29</td></tr><tr><td>2.2.2 Intersection et jointure des données spatiales</td><td>30</td></tr><tr><td>2.3 Exemples de manipulation de données spatiales</td><td>31</td></tr><tr><td>2.3.1 Définir l'emprise spatiale des zones bâties</td><td>31</td></tr><tr><td>2.3.2 Déterminer l'accessibilité à un équipement</td><td>33</td></tr><tr><td>2.3.3 Agréger les données carroyées à l'échelle des communes</td><td>38</td></tr><tr><td>3 Diffusion dynamique de l'information spatiale</td><td>40</td></tr><tr><td>3.1 Les flux géographiques</td><td>40</td></tr><tr><td>3.1.1 Les flux WMS</td><td>40</td></tr><tr><td>3.1.2 Les flux WFS</td><td>42</td></tr><tr><td>3.1.3 Exemple d'application dynamique</td><td>44</td></tr><tr><td>3.2 GeoServer : diffuser des flux géographiques</td><td>46</td></tr><tr><td>3.3 Pour aller plus loin : les tuiles vectorielles</td><td>47</td></tr><tr><td>3.3.1 Le principe des tuiles vectorielles</td><td>47</td></tr><tr><td>3.3.2 Mise en oeuvre avec Tegola</td><td>49</td></tr><tr><td>Discussion : le territoire des statisticiens</td><td>53</td></tr></table>
## Table des figures

<table><tr><td>1</td><td>Exemples de cartes thématiques</td><td>8</td></tr><tr><td>2</td><td>Répartition des terrasses éphémères parisiennes en 2021 - données ponctuelles ou surfaciques</td><td>13</td></tr><tr><td>3</td><td>Les trois principales projections utilisées par les statisticiens : Lambert 93, LAEA (projection européenne) et WGS 84 (utilisée par le système de positionnement par satellite GPS)</td><td>20</td></tr><tr><td>4</td><td>Associer une surface à un point : exemple avec les terrasses éphémères</td><td>25</td></tr><tr><td>5</td><td>Empresses géométriques minimales de la ville de Paris</td><td>26</td></tr><tr><td>6</td><td>Simplification d'un fond de carte : simplifier le contour des arrondissements parisiens</td><td>28</td></tr><tr><td>7</td><td>Les étapes de construction des taches bâties</td><td>33</td></tr><tr><td>8</td><td>Localisation des terrasses éphémères disposant de plus de 100 terrasses voisines dans un rayon de 200m</td><td>35</td></tr><tr><td>9</td><td>Carreaux à moins de 500m d'un parc</td><td>36</td></tr><tr><td>10</td><td>Contenu du flux IGN WMS</td><td>41</td></tr><tr><td>11</td><td>Bâtiments autour de la Direction Générale de l'Insee à Montrouge suivant différentes définitions de l'image (paramètres width et height)</td><td>42</td></tr><tr><td>12</td><td>Résultat de la requête getfeature précédente</td><td>43</td></tr><tr><td>13</td><td>Interface de GeoServer</td><td>47</td></tr><tr><td>14</td><td>Indexation des tuiles en fonction de l'échelle et de leur localisation. Source : Documentation QGIS</td><td>48</td></tr><tr><td>15</td><td>Application de prévisualisation des tuiles vectorielles de Tegola.</td><td>52</td></tr><tr><td>16</td><td>"Land doesn't vote, people do."</td><td>55</td></tr><tr><td>17</td><td>Deux visions d'un même phénomène : répartition des ménages pauvres dans l'agglomération parisienne. Source : Filosofi 2015</td><td>56</td></tr><tr><td>18</td><td>MAUP : effet de zonage</td><td>57</td></tr><tr><td>19</td><td>MAUP : effet d'échelle.</td><td>58</td></tr><tr><td>20</td><td>Are larger cities greener or smoggier ? Relation entre l'émission de CO2 et la taille des villes aux Etats-Unis.</td><td>59</td></tr></table>

## Liste des tableaux

<table><tr><td>1</td><td>Exemple de table spatiale : Position des terrasses éphémères parisiennes en 2021</td><td>13</td></tr><tr><td>2</td><td>Précision géographique en fonction de l'échelle d'observation. Source : Mapbox</td><td>48</td></tr></table>
## Liste des codes

<table><tr><td>1</td><td>Lire et écrire des données spatiales en R</td><td>15</td></tr><tr><td>2</td><td>Lire et écrire des données spatiales en Python</td><td>15</td></tr><tr><td>3</td><td>Transformer des coordonnées géographiques en données spatiales avec R</td><td>15</td></tr><tr><td>4</td><td>Transformer des coordonnées géographiques en données spatiales avec Python</td><td>15</td></tr><tr><td>5</td><td>Rajouter l'extension PostGIS à une base PostgreSQL depuis R</td><td>16</td></tr><tr><td>6</td><td>Rajouter l'extension PostGIS à une base PostgreSQL depuis Python</td><td>17</td></tr><tr><td>7</td><td>Envoyer et récupérer une table sous PostGIS depuis R</td><td>18</td></tr><tr><td>8</td><td>Envoyer et récupérer une table sous PostGIS depuis Python</td><td>18</td></tr><tr><td>9</td><td>Reprojection d'un fond départemental en R</td><td>21</td></tr><tr><td>10</td><td>Reprojection d'un fond départemental en Python</td><td>21</td></tr><tr><td>11</td><td>Reprojection d'un fond départemental en PostGIS</td><td>21</td></tr><tr><td>12</td><td>Transformation de données ponctuelles en données surfaciques en R</td><td>22</td></tr><tr><td>13</td><td>Transformation de données ponctuelles en données surfaciques en Python</td><td>22</td></tr><tr><td>14</td><td>Transformation données ponctuelles en données surfaciques en PostGIS</td><td>23</td></tr><tr><td>15</td><td>Diagramme de Voronoi des terrasses éphémères en R</td><td>25</td></tr><tr><td>16</td><td>Diagramme de Voronoi des terrasses éphémères en Python</td><td>25</td></tr><tr><td>17</td><td>Diagramme de Voronoi des terrasses éphémères en PostGIS</td><td>26</td></tr><tr><td>18</td><td>Polygones d'emprise spatiale minimale en R</td><td>26</td></tr><tr><td>19</td><td>Polygones d'emprise spatiale minimale en Python</td><td>27</td></tr><tr><td>20</td><td>Polygones d'emprise spatiale minimale en PostGIS</td><td>27</td></tr><tr><td>21</td><td>Compression du fond de carte des arrondissements parisiens à l'aide de mapshaper</td><td>28</td></tr><tr><td>22</td><td>Agrégation d'une couche géographique en R</td><td>29</td></tr><tr><td>23</td><td>Agrégation d'une couche géographique en Python</td><td>29</td></tr><tr><td>24</td><td>Agrégation d'une couche géographique en SQL</td><td>30</td></tr><tr><td>25</td><td>Intersection spatiale en R</td><td>30</td></tr><tr><td>26</td><td>Intersection spatiale en Python</td><td>30</td></tr><tr><td>27</td><td>Intersection spatiale en PostGIS</td><td>31</td></tr><tr><td>28</td><td>Délimitation des taches en R</td><td>32</td></tr><tr><td>29</td><td>Délimitation des taches en Python</td><td>32</td></tr><tr><td>30</td><td>Délimitation des taches en PostGIS</td><td>33</td></tr><tr><td>31</td><td>Requête spatiale reposant sur les structures k-d tree en R</td><td>34</td></tr><tr><td>32</td><td>Requête spatiale reposant sur les structures k-d tree en Python</td><td>34</td></tr><tr><td>33</td><td>Déterminer les carreaux à moins de 500 m d'un parc ou d'un jardin en R</td><td>35</td></tr></table>
<table>
  <tbody>
    <tr>
      <td>34</td>
      <td>Déterminer les carreaux à moins de 500 m d'un parc ou d'un jardin en Python</td>
      <td>36</td>
    </tr>
    <tr>
      <td>35</td>
      <td>Application de la DBSCAN à la détection des agrégats de terrasses en R</td>
      <td>37</td>
    </tr>
    <tr>
      <td>36</td>
      <td>Application de la DBSCAN à la détection des agrégats de terrasses en Python</td>
      <td>38</td>
    </tr>
    <tr>
      <td>37</td>
      <td>Calcul d'une population communale à partir de données carroyées au prorata de la surface commune</td>
      <td>39</td>
    </tr>
    <tr>
      <td>38</td>
      <td>Requête WMS pour récupérer une image</td>
      <td>41</td>
    </tr>
    <tr>
      <td>39</td>
      <td>Téléchargement d'une image à partir d'une requête WMS en R</td>
      <td>41</td>
    </tr>
    <tr>
      <td>40</td>
      <td>Téléchargement d'une image à partir d'une requête WMS en Python</td>
      <td>42</td>
    </tr>
    <tr>
      <td>41</td>
      <td>Requête WFS pour récupérer les données vectorielles</td>
      <td>43</td>
    </tr>
    <tr>
      <td>42</td>
      <td>Traitement d'une requête WFS en R</td>
      <td>43</td>
    </tr>
    <tr>
      <td>43</td>
      <td>Traitement d'une requête WFS en Python</td>
      <td>44</td>
    </tr>
    <tr>
      <td>44</td>
      <td>Cartographie dynamique d'un flux WFS à l'aide de la librairie Open-Layer</td>
      <td>44</td>
    </tr>
  </tbody>
</table>
Introduction

Les **cartes thématiques**, appelées également **cartes statistiques** ou **quanti-tatives**, sont des médiations symboliques qui visent à restituer graphiquement l'organisation spatiale des phénomènes sociaux (Palsky, 1990). Leur origine est étroitement liée à celle de la statistique. La première carte thématique moderne est attribuée au statisticien Charles Dupin. Publiée en 1826, elle décrit la répartition de l'instruction populaire en France (fig 1a).

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

(a) Une première carte thématique : instruction populaire en France, Charles Dupin, 1826

Lecture : nombre d'élèves masculins scolarisés en primaire rapporté au nombre d'habitants dans le département selon dix tranches en aplat de couleur. Ce rapport est faible dans les zones les plus foncées.

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

(b) Part du bâti concentré en 2005, Himpens, Poulhes, et Sémécurbe (2021)

Lecture : en 2005, dans la majorité des territoires du quart nord-est de la France, plus de 85 % du bâti est concentré dans des noyaux villageois.

<caption>FIGURE 1 – Exemples de cartes thématiques</caption>

Cette carte marque une rupture avec les cartes topographiques. Pour la première fois, une représentation graphique vise à représenter les variations spatiales d'une grandeur abstraite, en l'occurrence une statistique sociale<fnref n="1" />. Ceci est rendu

<footnote n="1">1. La carte figurative de l'instruction populaire sépare la France en deux : dans le nord-est la France « éclairée » s'oppose à la France « obscure » dans le reste du pays selon les commentaires de Dupin lors de la présentation de sa carte au Conservatoire des Arts et Métiers (Palsky, 1996). Le message est clair et a nourri tout le long du 19e les représentations politiques et sociales, ce</footnote>
possible grâce à ce qui s'appellera à la fin du 20ᵉ siècle, la **sémiologie cartographique**. Plus les teintes sont claires et plus la part de l'instruction populaire est importante dans les départements. Contrairement, aux cartes topographiques dont la lecture est plus ou moins immédiate, la compréhension des cartes thématiques nécessite d'acquérir certains réflexes perceptifs<sup>2</sup>. En contrepartie de cet effort interprétatif, la cartographie thématique produit une vision utile à l'action publique et à l'appropriation par les citoyens des enjeux sociaux. En particulier, la comparaison de cartes permet de mettre en relation des phénomènes sociaux a priori indépendants : tels que l'instruction populaire et la concentration de l'habitat<sup>3</sup> (fig 1b).

Autrefois réalisées à la main, les cartographies thématiques actuelles nécessitent l'acquisition de savoir-faire **géomatiques**. La **géomatique**, qui provient de la contraction des termes « géographie » et « informatique », désigne l'ensemble des outils, technologies et méthodes permettant de représenter, d'analyser et d'intégrer des données géographiques<sup>4</sup>. Elle fait appel aux sciences, aux technologies de mesure de la terre ainsi qu'aux technologies de l'information pour faciliter l'acquisition, le traitement et la diffusion de ces données. Elle couvre les bases de données spatiales, les systèmes d'information géographique (SIG)<sup>5</sup>, l'imagerie spatiale, mais aussi les technologies du web, la géolocalisation<sup>6</sup>...

La profusion actuelle de l'information spatiale n'a pas encore favorisé sa démo-
cratisation. Au contraire, le foisonnement des technologies et la nécessité d'acquérir
des compétences plurielles en analyse spatiale, en gestion des réseaux et en dévelop-
pement web représentent des obstacles à dépasser. Si des solutions simples existent
pour diffuser de l'information spatiale à petite échelle, telle que Rshiny, la diffusion
massive d'une information nécessite des infrastructures et des connaissances spéci-
fiques. C'est à ce besoin que ce document cherche à répondre. Ainsi, il ne constitue
pas un guide en matière de cartographie spatiale, de sémiologie spatiale, ou encore
d'analyses statistiques spatiales. Pour ce faire, le lecteur pourra se référer au très

<footnote>qui fait de la carte un outil au statut particulier dont la lecture a un impact supérieur à celle d'un simple tableau. En tant qu'objet symbolique la carte peut être porteuse de présupposés voire de préjugés, ainsi le recours aux couleurs sombres peut refléter « l'obscurantisme ».</footnote>

<footnote n="2">2. Palsky (1991) reprend Playfair (1805) qui écrit « the reader will find, five minutes attention to the principle on which they are constructed, a saving of much labor and time; but without that trifling attention, he may as well look at a blank sheet of paper. »</footnote>

<footnote n="3">3. voir Le Bras et Todd (2013) pour une explication à cette corrélation spatiale.</footnote>

<footnote n="4">4. voir le site de l'Institut National de Géographie : https://www sig-geomatique.fr/ sig-geomatique.html</footnote>

<footnote n="5">5. Un système d'information géographique ou SIG (en anglais, geographic information system ou GIS) est un système d'information conçu pour recueillir, stocker, traiter, analyser, gérer et présenter tous les types de données spatiales et géographiques Wikipédia (2021c).</footnote>

<footnote n="6">6. Pour une revue des outils de géolocalisation disponibles sur le web ou l'Insee et d'autres exemples cartographiques voir https://juliedjidji.github.io/memocarto/geocodage.html</footnote>
complet manuel d'analyse spatiale, Loonis et Bellefon (2018), pour les méthodes d'analyse, à Lambert et Zanin (2016), ou encore au portail de la cartographie proposé par le pôle de service de l'action régionale sur l'analyse territoriale de l'Insee. Ce document s'inscrit en complément comme une initiation, une explicitation, une présentation, à destination d'un public de statisticiens des outils et les méthodes actuels et innovants de la géomatique, en se concentrant sur le cas des données spatiales vectorielles<sup>7</sup>, c'est-à-dire définies en contours des territoires.

Le traitement des données spatiales a été bouleversé par une prise en compte intégrée de la géographie dans les langages de programmation. Autrefois, la prise en compte de la géographie sous R reposait sur le *package* sp qui traitait séparément les données statistiques et les données spatiales (Giraud et Pecout, 2021). La manipulation simultanée de l'information statistique et spatiale était délicate et lente. La géométrie, c'est-à-dire l'information purement spatiale, est à présent une information comme une autre dans les tables de données ce qui simplifie considérablement les traitements. Cette approche a été rendue possible par le recours systématique aux librairies spatiales GEOS, GDAL et PROJ dans les logiciels Python avec le *package* geopandas (extension spatiale de la librairie de traitement des tables de données pandas), et pour R avec le *package* sf compatible avec la suite tidyverse et PostgreSQL et l'extension spatiale PostGIS. Ces solutions logicielles étant libres et basées sur les mêmes librairies garantissent l'uniformité des traitements cartographiques. Ce document essaye tant que faire se peut d'illustrer cette proximité en proposant des exemples de codes dans les trois langages R, Python et SQL.

La suite de ce document se compose de quatre parties. La première partie sur les **données spatiales** décrit la façon dont l'espace est pris en compte dans les systèmes d'information géographique. La deuxième est consacrée aux **traitements et manipulations sur les données spatiales**. Elle décrit les opérations couramment utilisées pour produire une connaissance statistique à l'aide de données géographiques et de statistiques élémentaires. La troisième s'intéresse à la **diffusion dynamique de l'information spatiale**. Elle explore les outils disponibles pour mettre à disposition du grand public des cartes thématiques interactives. La dernière partie de ce document plus conclusive discute des mésusages de la cartographie thématique. La dimension graphique voire esthétique des cartes thématiques en fait des outils puissants pour « faire parler les chiffres » et raconter une histoire et ce, quelle que soit la qualité des analyses sous-jacentes. Il est donc toujours primordial pour les statisticiens de se questionner sur le rôle de la spatialisation et celui de la représentation des données statistiques dans la production d'une connaissance géographique.

<footnote n="7">7. Ce document n'aborde pas les images raster qui recourent en effet à des outils différents : reconnaissance d'images...</footnote>
L'ensemble des outils présentés dans ce document peuvent être exécutés en local sur un ordinateur personnel et dans ce cas, sauf manipulation spécifique, ils sont accessibles uniquement à son utilisateur. Il est également possible de les déployer sur des serveurs *adhoc* ou à l'aide de conteneurs Docker, ceci afin de disposer d'une puissance de traitement plus importante et de mettre à disposition d'un large public l'information spatiale.
# 1 Les données spatiales

Les données spatiales désignent les informations liées à des objets ou éléments présentant une composante spatiale. Elles recouvrent deux types de données :

**Vectoriel** Ce type de données repose sur l'usage de figures géométriques telles que des points, des lignes et des polygones définis à l'aide de listes de coordonnées pour représenter les objets spatiaux, tels que la répartition des équipements ou de la population dans l'espace ou les contours des départements et des régions.

**Raster** Ce type utilise des images, soit des matrices de nombres, pour décrire les différents constituants de l'espace. Les images aériennes mais aussi les plans d'occupation du sol sont généralement diffusées en format raster.

Ce document se concentre sur les **données spatiales vectorielles**. Les bases de données spatiales vectorielles sont optimisées pour stocker et interroger des données reliées à des objets référencés géographiquement, comme des points, les lignes et des polygones, communément appelés **géométrie** de la base. Les données numériques et textuelles de la base spatiale sont dénommées **données attributaires**.

## 1.1 La géographie

Dans les logiciels R, Python et les bases de données PostGIS, l'information spatiale une fois chargée en mémoire se présente sous la forme de tableaux de données classiques (*dataframe*) auxquels a été rajoutée la colonne de géométrie (*geometry*), laquelle décrit la représentation spatiale des objets. Cette géométrie suit une norme commune, le standard ISO 19125 *simple feature access*, appelé plus simplement Simple Features et défini conjointement par l' *Open Geospatial Consortium* (OGC) et l' *International Organization for Standardization* (ISO). Cette norme permet de retranscrire en format numérique la géométrie des objets et des entités, ce qui garantit la reproductibilité des traitements cartographiques entre les trois logiciels.

La grammaire de la norme Simple Features est décrite sur la page Wikipédia associée. On décrit ci-dessous ses représentations les plus courantes. Les points et les polygones sont les représentations des objets les plus utilisées dans les cartes thématiques. Le point sert à représenter l'emplacement des équipements ou des activités dont la surface est négligeable par rapport aux phénomènes d'intérêt. Il est codé à l'aide du mot clef *POINT*(X₁, Y₁). Un polygone est codé à l'aide du mot clef *POLYGON*(X₁, Y₁, X₂, Y₂...Xₙ, Y₁) suivi d'une liste de coordonnées<fnref n="8" />. On l'utilise pour décrire l'emprise spatiale d'une commune ou d'un carreau appartenant à un carroyage par exemple. Les fleuves et les routes qui servent à agrémenter une

<footnote n="8">8. la première coordonnée listée correspond à la dernière, de façon à clore le polygone.</footnote>
<table><thead><tr><td>Nom du commerce</td><td>geometry</td></tr></thead><tbody><tr><td>Kifolie</td><td>POINT (654406.186 6863626.372)</td></tr><tr><td>Le bistrot d'André</td><td>POINT (646934.770 6859933.893)</td></tr><tr><td>Les frérots</td><td>POINT (649919.421 6864956.475)</td></tr><tr><td>PRET A MANGER</td><td>POINT (650048.122 6864069.131)</td></tr><tr><td>la Traversée</td><td>POINT (652155.733 6865449.423)</td></tr></tbody></table>

Source : https://opendata.paris.fr/explore/dataset/terrasses-ephemeres

<caption>TABLEAU 1 – Exemple de table spatiale : Position des terrasses éphémères parisiennes en 2021</caption>

carte sont généralement représentés sous la forme de structure linéaire codée à l'aide du mot clef LINE(X₁, Y₁, X₂, Y₂...Xₙ, Yₙ). Des structures plus complexes, tels que des objets composés de plusieurs parties (territoires avec des exclaves), peuvent être décrites à l'aide des mots clefs MULTIPOINT, MULTILINE et MULTIPOLYGON.

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

Source : https://opendata.paris.fr/explore/dataset/terrasses-ephemeres

<caption>FIGURE 2 – Répartition des terrasses éphémères parisiennes en 2021 - données ponctuelles ou surfaciques</caption>

## 1.2 Stockage des données spatiales

Les données spatiales sont stockées soit dans des fichiers cartographiques soit dans des bases de données. Les fichiers cartographiques sont appelés **fonds de cartes** en géomatique. Le type de stockage conditionne l'accès aux données et les
traitements géomatiques possibles ensuite de réaliser. Les fichiers cartographiques sont efficaces pour diffuser l'information spatiale mais ils s'avèrent inefficaces pour traiter et manipuler de très gros volumes de données. Les bases de données (Post-GIS par exemple) répondent à ce besoin mais nécessitent en contrepartie de maîtriser un langage de requêtage (le langage SQL pour PostGIS) et de savoir se connecter à une base de données.

### 1.2.1 Fichiers cartographiques

Une multitude de formats de fichiers géographiques existent. Le plus connu et le plus utilisé jusqu'à présent est le format *shapefile* associé au logiciel ArcGIS de la société ESRI. Dans ce format, un fond de carte est composé de cinq fichiers : le `.shp` qui est la porte d'entrée sur les données, un fichier `.dbf` qui contient les données attributaires (les données associées aux objets contenus dans le *shapefile*) <fnref n="9" /> et d'autres fichiers qui gèrent la projection et la géométrie. Tous les fichiers doivent être localisés au même endroit pour que le fond puisse être ouvert correctement.

D'autres formats tendent progressivement à s'imposer tel que celui proposé par la communauté QGIS, le Geopackage, qui a l'avantage de pouvoir intégrer plusieurs bases/fonds de carte en même temps. Les données sont contenues dans une base de données SQLite, l'extension est `*.gpkg`. Le format geojson est très utilisé par les librairies de cartographie dynamique telles que Leaflet, mais est particulièrement volumineux. geojson est une extension spatiale de json, les données ne sont pas compressées et sont organisées selon une logique de dictionnaire, ce qui permet de le lire facilement, comme en témoigne le fichier ci-dessous contenant deux points et deux variables (*lat* pour la latitude et *lon* pour la longitude) <fnref n="10" />.

```json
{
  "type": "FeatureCollection",
  "crs": {"type":"name", "properties": {"name":"urn:ogc:def:crs:EPSG::4328"}}},
  "features": [
    {
      "type": "Feature", "properties": { "lon": 2.416, "lat": 48.848 },
      "geometry": { "type": "Point", "coordinates": [ 2.416, 48.848 ] }
    },
    {
      "type": "Feature", "properties": { "lon": 2.361, "lat": 48.888 },
      "geometry": { "type": "Point", "coordinates": [ 2.361, 48.888 ] }
    }
  ]
}
```

<footnote n="9">9. L'ouverture de ce fichier dans Python ou R pour permet d'accéder aux données attributaires sans avoir à charger les données spatiales.</footnote>

<footnote n="10">10. A la ligne 3, apparaît `crs:EPSG::4328` qui détermine la projection du fond de carte. Les projections sont abordées dans la partie 1.2.3.</footnote>
10 }

Les outils modernes sont généralement capables de lire tous les formats disponibles. Une fois ouvert sous la forme d'une table spatiale, le fond peut être ensuite enregistré dans un autre format.

En R, les commandes `st_read` et `st_write` servent à ouvrir et sauvegarder un fond spatial.

```r
1 library sf)
2 url = 'https://france-geojson GregoireDavid.fr/repo/regions.
geojson'
3 region = st_read(url)
4 st_write(region, 'region.shp')
```

### Code 1 – Lire et écrire des données spatiales en R

En Python, l'ouverture et l'enregistrement d'un fond de carte s'effectue à l'aide des commandes `read_file` et `to_file`.

```python
1 import geopandas as gpd
2 url = 'https://france-geojson GregoireDavid.fr/repo/regions.
geojson'
3 region = gpd.read_file(url)
4 region.to_file('region.shp')
```

### Code 2 – Lire et écrire des données spatiales en Python

Très souvent, les données ponctuelles, telles que l'emplacement des terrasses éphémères sont diffusées sous la forme d'un fichier csv contenant deux colonnes de coordonnées. L'instruction `st_as_sf` et `gpd points_from_xy` permettent de transformer ces fichiers en fond de carte respectivement sous R et Python. Attention, ces instructions ne fonctionnent pas en présence de valeurs manquantes.

```r
1 library (sf)
2 library (tidyverse)
3
4 terrasses = read_delim('terrazz es-eph emeres.csv', lim em=':')
5 sterrasses = st_as_sf(terrazz es, co ords = c('Coord X', 'Coord Y'),
                       crs=2154)
```

### Code 3 – Transformer des coordonnées géographiques en données spatiales avec R

```r
1 import geopandas as gpd
2 import pandas as pd
3
4 terrasses = pd.read_csv('terrazz es-eph emeres.csv', sep=";")
5 geom = gpd points_from_xy(terrazz es['Coord X'], terrasse['Coord Y'])
6 sterrasses = gpd.GeoDataFrame(terrazz es, geometry=geom, crs=2154)
```

### Code 4 – Transformer des coordonnées géographiques en données spatiales avec Python
### 1.2.2 PostGIS, une base de données spatiales

PostGIS<fnref n="11" /> est l'extension spatiale de PostgreSQL (comme Oracle spatial est celle d'Oracle) (PostGIS, 2021). Une table spatiale sous PostGIS est analogue à un tableau de données spatiales sous R ou python : une colonne de géométrie décrit à l'aide de la norme Simple Features l'emprise spatiale des objets de la table. PostGIS permet de manipuler des objets géographiques grâce à ces géométries (points, lignes, polygones...), la géométrie étant alors un type à part entière et non une simple succession de coordonnées. Ce type est alors manipulable au même titre que les autres, avec un grand nombre d'opérations possibles. Par rapport aux fichiers cartographiques, l'apport de PostGIS est de pouvoir stocker et traiter de grandes quantités d'informations. Ceci est rendu possible par l'indexation spatiale basée sur deux principes :

- chaque objet est indexé par sa bounding box, soit les 4 couples de coordonnées qui définissent le plus petit rectangle qui le contient;

- ensuite, un arbre spatial permet de classer et regrouper les objets selon leur proximité et rompre avec une lecture séquentielle des objets contenus dans une table pour prendre en compte leur positionnement relatif.

Le recours à l'indexation spatiale permet d'éviter dans les requêtes spatiales, telles que les jointures spatiales, des calculs inutiles entre des objets très éloignés sans intersection.

L'accès aux données peut se faire à l'aide d'un client SQL (PGAdmin voir encadré, DBeaver...), à travers un logiciel tiers comme QGIS, en ligne de commande (en utilisant les utilitaires de la bibliothèque GDAL), ou à l'aide d'un serveur cartographique (geoserver...). En particulier, deux approches permettent de communiquer avec PostGIS depuis R ou Python. La première consiste à lancer des requêtes SQL qui s'exécutent sous PostGIS. Les résultats restent sous PostGIS, il n'y a pas de retour sous R ou Python. Cette approche permet de configurer ou de créer de nouvelles vues ou tables dans la base de données.

En R, l'accès à une base de données s'effectue à l'aide de la librairie *DBI*. De plus amples informations sur cette librairie et le requêtage d'une base de données depuis R sont disponibles sur le site utilitR.

1. library(DBI)

2. library(RPostgres)

3. conn <- DBI::dbConnect(drv = RPostgres::Postgres())

4. host = "mypostgis",

<footnote n="11">11. Voir http://postgis.net le site de référence avec la documentation de la dernière version http://postgis.net/docs/manual-dev/postgis-fr.html voir aussi un tutoriel pour débuter avec PostGIS https://www.sigterritoires.fr/index.php/debuter-avec-postgresql-postgis-introduction-a-pgadmin4/ Le site de référence</footnote>
```
6
7
8
9
10
11
```

Code 5 – Rajouter l'extension PostGIS à une base PostgreSQL depuis R

En Python, la communication avec une base de données, telle que Postgres, est dévolue à la librairie `psycogp2`.

```
1
2
3
4
5
6
7
8
9
10
11
```

Code 6 – Rajouter l'extension PostGIS à une base PostgreSQL depuis Python

La connexion à une base de données PostGIS nécessite de spécifier cinq informations :

- Le *port* logiciel permet, sur un ordinateur donné, de distinguer différents interlocuteurs. Ces interlocuteurs sont des programmes informatiques qui, selon les cas, écoutent ou émettent des informations sur ces ports. Un port est distingué par son numéro (Wikipédia, 2022). Par défaut le port de PostGIS est 5432;

- l'hôte (*host*) est le nom du serveur;

- *database* est le nom de la base de données;

- *user* est le nom de l'utilisateur. Une même base peut avoir plusieurs utilisateurs différents;

- Enfin *password* est le mot de passe associé à l'utilisateur.

La deuxième approche vise à envoyer ou récupérer de l'information entre une base de données PostGIS et R ou Python. Les données pouvant être très volumineuses, les transferts peuvent être très longs, de sorte, qu'il est préférable de concentrer les calculs sur la base de données de PostGIS et de récupérer les bases finales en fin de traitement.

En R, les commandes *sf*, *st_read* et *st_write* servent également à communiquer des informations avec PostGIS. Dans l'exemple ci-dessous, une table spatiale *sf* de
points aléatoires dans l'emprise spatiale de Paris est créée puis exportée dans
PostGIS. Cette dernière table est enfin requêtée sous R.

1 library(DBI)
2 library(RPostgres)
3 libraryzf)
4
5 # Creation d'une table de points sf et export vers PostGIS (dans 1
'emprise spatiale de Paris)
6 lon = runif(20, 2.25, 2.42)
7 lat = runif(20, 48.8, 48.9)
8 df = data.frame(lon, lat)
9 df = st_as_sf(df, coordinates = c("lon", "lat"),
                       crs = 4326, agr = "constant")
10
11 st_write(df, dsn = conn, layer = "points",
          delete layer = TRUE)
13
14 # Import de la table points sous sf
15 sql = 'SELECT * FROM points'
16 points <- st_read(con, query = sql)

Code 7 – Envoyer et récupérer une table sous PostGIS depuis R

En Python, la librairie `sqlalchemy` permet à la librairie `geopandas` de communiquer avec les bases de données spatiales PostGIS. En revanche, contrairement à R, les transferts s'effectuent à l'aide de deux nouvelles commandes `to_postgis` et `read_postgis`.

1 import pandas as pd
2 import geopandas as gpd
3 from sqlalchemylib create_engine
4
5 db_connection_url = " postgresql://myusername:mypassword@myhost
:5432/mydatabase"
6 engine = create_engine(db_connection_url)
7
8 # Creation d'une table de points Geopandas et export vers PostGIS
(dans l'emprise spatiale de Paris)
9 df = pd. Frame({ 'lon' : np.randomcas uniform(2.25, 2.42,20),
'lat' : np.randomcas uniform(48.8, 48.9, 20)})
10
df = gpd.GeoDataFrame(df,
  geometry = gpd. points_from_xy(df. lon, df. lat),
  crs = ' epsg : 4328')
13
df.to_postgis(' points ', engine, if_exists='replace')
15
16 # Import de la table points sous Python
17 sql = 'SELECT * FROM points'
18 points = gpd.read_postgis/sql, engine, geom_col='geometry')

Code 8 – Envoyer et récupérer une table sous PostGIS depuis Python
de France métropolitaine. Elle a vocation à déformer le moins possible les entités spatiales contenues en France Métropolitaine. La distance entre Lille et Marseille calculée avec la formule de Pythagore sur les coordonnées en Lambert 93 est de 833.812 km. La distance géodésique, tenant compte de la courbure de la Terre, est de 834.005 km. L'approximation de Pythagore est très bonne en France métropolitaine. C'est d'ailleurs cette distance qui est utilisée en statistique spatiale. En revanche, une distance vol d'oiseau entre Lille et l'Australie génère une erreur de plus de 1500 km. La projection LAEA est le système en vigueur en Europe, utilisé entre autres par Eurostat pour diffuser de l'information spatiale. Ce système vise à déformer le moins possible la géométrie des pays européens.

<figure><img src="image_9.png" /><figcaption>FIGURE 3 – Les trois principales projections utilisées par les statisticiens : Lambert 93, LAEA (projection européenne) et WGS 84 (utilisée par le système de positionnement par satellite GPS)</figcaption></figure>

En fonction de la projection retenue, les représentations des territoires de la France métropolitaine peuvent être très différentes (fig 3). En particulier, la carte en WGS 84 se différencie fortement des deux autres. La projection WGS 84 est très utilisée pour la cartographie dynamique et les systèmes de positionnement par satellite. Contrairement aux projections Lambert 93 et LAEA, la projection
WGS 84 n'est pas une projection en tant que telle, mais une modélisation du
globe terrestre. Les coordonnées sont des angles qui permettent de situer tout
point sur la surface de la terre à partir d'une référence. Par exemple, le centre
de Lille a pour coordonnées (3.063, 50.636) et celui de Marseille (5.369, 43.296).
La projection WGS84 ne permet pas de calculer directement des distances et est
donc inadaptée à l'analyse spatiale. Une erreur courante est d'intervvertir les deux
angles.

Une projection est identifiée à l'aide d'un numéro EPSG (European Petroleum Survey Group). Par exemple, 2154 pour le lambert 93, 3035 pour le LAEA (projection européenne) et 4326 pour le WGS 84. On passe d'une projection à une autre, en effectuant une opération de reprojection. Cette opération peut être très longue en fonction de la taille de la couche cartographique. Il est donc recommandé d'adopter une projection unique dès le début d'un projet cartographique.

En R, l'opération de reprojection s'effectue à l'aide de la commande `st_transform` :

1 library sf)

2 departement = readtras('depfrancemetro_2020.geojson')

3 departement = st_tnansform(departement, 2154)

Code 9 – Reprojection d'un fond départemental en R

En Python, cette opération s'obtient avec la commande .to_crs :

1 import geopandas as gpd

2 departement = gpd.read_file('depfrancemetro_2020.geojson')

3 departement = departement.to_crs("EPSG:2154")

Code 10 – Reprojection d'un fond départemental en Python

En PostGIS, cette opération s'effectue à l'aide de la commande `ST_TRANSFORM` :

1 ALTER TABLE departement

2 ALTER COLUMN geometry TYPE geometry (MultiPolygon, 2154)

3 USING ST_Transform (geometry, 2154)

Code 11 – Reprojection d'un fond départemental en PostGIS

1.3 Exemple : création d'une couche carroyée à partir de données ponctuelles

La transformation d'une donnée ponctuelle en une donnée surfacique ne nécessite pas toujours l'emploi d'opérations géomatiques spécifiques. La transformation du fond de carte des terrasses éphémères en des données carroyées de 200m de côté
proposée ici repose sur l'application de fonctions mathématiques élémentaires sur les coordonnées géographiques<fnref n="13" />.

La première étape consiste à déterminer les coordonnées du coin sud-est du carreau d'appartenance de chaque terrasse en recourant à l'aide de la fonction valeur entière : $X_{car} = \lfloor floor(X/s) \times s \rfloor$ et $Y_{car} = \lfloor floor(Y/s) \times s \rfloor$. Où *size* est la taille désirée des carreaux. Le regroupement des données selon ces deux nouvelles coordonnées permet de créer la partie attributaire de la future table de données carroyées. Enfin, on déduit la géométrie des carreaux à partir des coordonnées de leur coin Sud-Est :

Polygon($X_{car}, Y_{car}, X_{car}, Y_{car}+size, X_{car}+size, Y_{car}+size, X_{car}+size, Y_{car}, X_{car}, Y_{car}$)

```r
library sf)
library (tidyverse)
size = 200
terrasses = read_delim('terrasses-ephemeres.csv', delimiter=';')
terrasses['x'] = floor(terrasses$'Coord X' / size) * size
terrasses['y'] = floor(terrasses$'Coord Y' / size) * size
terrasses$densite = 1 / (size / 1000)^2
lterrasses = terrasses %>%
  group_by(x, y) %>%
  summarize(densite = sum(densite)
            , .groups = 'drop')
lterrasses$geometry <- sprintf("POLYGON ((%f %f, %f %f, %f %f, %f %f, %f %f))",
                 lterrasses$x, lterrasses$y,
                 lterrasses$x+size, lterrasses$y,
                 lterrasses$x+size, lterrasses$y+
                 size,
                 lterrasses$x, lterrasses$y+size,
                 lterrasses$x, lterrasses$y)
lterrasses <- st_as_SF(lterrasses, wkt = "geometry", crs = 2154)
plot(lterrasses)
```

Code 12 – Transformation de données ponctuelles en données surfaciques en R

```r
import geopandas as gpd
import pandas as pd
import numpy as np
from shapelygeometry import Polygon
```

<footnote n="13">13. Cette transformation fonctionne uniquement si les coordonnées sont exprimées dans une projection métrique, voir la sous-section précédente pour plus de renseignements sur les projections.</footnote>
```
terrasses = pd.read_csv('terrasses-ephemeres.csv', sep=";")
size = 200

terrasses['x'] = np.floor(terrasses['Coord X'] / 200) * 200
terrasses['y'] = np.floor(terrasses['Coord Y'] / 200) * 200
terrasses['densite'] = 1 / (size / 1000) **2
lterrasses = terrasses.groupby(['x', 'y'],
                              as_index=False)['densite'].sum()

geometry = [Polygon([(
    (x + size, y),
    (x + size, y + size),
    (x, y + size)])
             for x,y in zip(lterrasses.x,lterrasses.y)]
lterrasses = gpd.GeoDataFrame(lterrasses, geometry=geometry, crs=
                              'EPSG:2154')
```

Code 13 – Transformation de données ponctuelles en données surfaciques en Python

```
CREATE table lterrasses as
select temp.x, temp.y, sum(temp.densite) as "densite",
ST_Polygon(CONCAT('LINESTRING', temp.x, ' ', temp.y, ' ',
              temp.x + 200, ' ', temp.y, ' ',
              temp.x + 200, ' ', temp.y + 200, ' ',
              temp.x, ' ', temp.y + 200, ' ',
              temp.x, ' ', temp.y, ' ))), 2154 ) as geometry
from
(select floor(t."Coord X" / 200) * 200 as "x",
    floor(t."Coord Y" / 200) * 200 as "y",
    1 / (200.00/1000.00) as "densite"
from terrasses as t) as "temp"
group by x,y;
```

Code 14 – Transformation données ponctuelles en données surfaciques en PostGIS

Cette approche simple fonctionne uniquement parce que les carreaux sont agencés de façon régulière dans l'espace. Pour des objets plus complexes il est indispensable de recourir à des traitements géomatiques.
# 2 Traitement des données spatiales

Un même phénomène tel que l'organisation des terrasses éphémères peut être représenté de diverses manières. L'enjeu pour le statisticien est de savoir transformer la représentation géométrique d'un phénomène en fonction des objectifs statistiques et esthétiques visés. Ces transformations peuvent être simples et ne pas nécessiter de fonctions spatiales complexes comme vu dans l'exemple de la figure 2b) mais plus généralement, les transformations sont le résultat d'un enchaînement de traitements géomatiques plus complexes.

De nombreuses fonctions géomatiques existent, seules les plus utiles par les statisticiens sont discutées dans cette partie. Pour plus de renseignements, le lecteur peut consulter l'aide mémoire en R des fonctions géomatiques (Pebesma, 2021), les fonctions géomatiques en Python et la liste des fonctions géomatiques en PostGIS.

## 2.1 Fonctions et transformations géométriques élémentaires

Les fonctions et transformations élémentaires visent à transformer la nature ou le type de géométrie d'un fond de carte spatiale sans avoir à recourir à une information spatiale tierce. Elles sont classées dans cette section en fonction de la géométrie initiale de la donnée et de la géométrie visée.

### 2.1.1 Transformation de données ponctuelles en données surfaciques

Les zones tampons (*buffers*) et les diagrammes de Voronoi sont utilisés pour associer à chaque point d'une couche ponctuelle une surface.

**Zones tampons (buffers)** Les zones tampons associent à chaque point de la couche considérée un disque centré sur celui-ci de rayon défini par l'utilisateur (voir figure 4 a). Les zones tampons servent souvent à déterminer des zones d'accessibilité selon une logique vol d'oiseau autour d'un équipement, d'un point d'intérêt ou à définir l'environnement d'un ménage ou d'un individu. Sous R et PostGIS, les zones tampons sont obtenues à l'aide de la fonction *st_buffer* et sous Python avec *buffer*. Un exemple d'application des zones tampons est présenté dans la section 2.3.1.

**Diagramme de Voronoi** Le diagramme de Voronoi est un pavage de l'espace qui associe à chaque point de la couche considérée l'ensemble des points du plan plus proches de ce point que d'aucun autre (voir figure 4 b). Cette modélisation de l'espace permet de définir des zones de chalandise sans recouvrement, contrairement aux espaces définis à l'aide d'une zone tampon.
<figure><img src="image_10.png" /><figcaption>FIGURE 4 – Associer une surface à un point : exemple avec les terrasses éphémères</figcaption></figure>

En R, le calcul des diagrammes de Voronoi se réalise à l'aide de la fonction `st_voronoi` disponible directement dans la librairie sf.

```r
1 library sf
2
3 terrasses = st_read('sterrasses.gpkg')
4
5 voronoi <- terrasses %>%
6   st_union() %>%
7   st_voronoi() %>%
8   st_collection_extract()
```

Code 15 – Diagramme de Voronoi des terrasses éphémères en R

En Python, le calcul des diagrammes de Voronoi s'effectue à l'aide du package `geovoronoi`. La première étape consiste à récupérer les coordonnées $(x, y)$ à partir de la colonne géométrie.

```python
1 import geopandas as gpd
2 import numpy as np
3 from geovoronoi import voronoiRegionsFromCoordinates
4
5 terrasses = gpd.read_file('sterrasses.gpkg')
6
7 coord = np.vstack([terrasses_geometry.x.T, terrasses_geometry.y.T
8                      ).T
8 pol, _ = voronoiRegionsFromCoordinatescroix, terrasses.unary_union.
9 convex_hull)
9 Voronoi = gpd.GeoSeries(list pol.values())
```

Code 16 – Diagramme de Voronoi des terrasses éphémères en Python
En PostGIS,

1 create table voronoi_terrasses as

2 SELECT

3 ST_VoronoiPolygons(st_union(geometry)) As "geometry"

4 FROM terrasses;

<caption>Code 17 – Diagramme de Voronoi des terrasses éphémères en PostGIS</caption>

**Emprise géométrique minimale** Déterminer l'emprise géométrique minimale d'un ensemble de points, c'est-à-dire déterminer la surface englobant cet ensemble de points, est une opération courante en géomatique. Cette opération est utile par exemple pour affecter une surface à un ensemble de points afin par exemple de calculer une densité. Parmi l'infinité de méthodes permettant de trouver une telle surface, la recherche du polygone convexe minimal (enveloppe convexe ou *convex hull*) englobant tous les points est l'une des méthodes les plus utilisées. Une autre solution simple consiste à déterminer le plus petit rectangle qui contient tous les points et dont les côtés sont parallèles aux axes de la projection (rectangle d'emprise minimale ou *envelope*), voir figure 5 pour une illustration concernant la ville de Paris.

<figure><img src="image_11.png" /><figcaption>(a) Polygone convexe d'emprise minimale (b) Rectangle d'emprise minimale (*envelope*)
(*convex hull*)

FIGURE 5 – Emprises géométriques minimales de la ville de Paris</figcaption></figure>

En R, le polygone convexe minimal (*convex hull*) s'obtient en regroupant tous les points en un unique objet à l'aide de la fonction `st_union` puis en appliquant la commande `stconvex hull`. Il n'existe pas de fonction pour calculer directement le rectangle d'emprise minimale (*envelope*), celle-ci s'obtient en combinant la commande `st_bbox` qui retourne les coordonnées de l'enveloppe et la commande `st_as_sfc` qui permet de les transformer en un polygone.

1 library sf
```
2
3 terrasses = st_read('sterrasses.gpkg')
4 convex_hull = terrasses %>% st_union() %>% stoj convex_hull()
5 envelope = terrasses %>% st_union() %>% stbbbox() %>% stas sfc
```

### Code 18 – Polygones d'emprise spatiale minimale en R

En Python, le polygone convexe minimal et le rectangle d'emprise minimale, s'obtiennent en regroupant tous les points en un unique objet à l'aide de la fonction unary_union puis en appliquant la fonction convex_hull ou la fonction envelope.

```
1 import geopandas as gpd
2 terrasses = gpd.read_file('sterrasses.gpkg')
3 convex_hull = gpd.GeoSeries(terrasses.unary_union(convex_hull, crs='EPSG:2154'))
4 envelope = gpd.GeoSeries(terrasses.unary_union(envelope, crs='EPSG:2154'))
```

### Code 19 – Polygones d'emprise spatiale minimale en Python

En PostGIS, le polygone convexe minimal et le rectangle d'emprise minimale s'obtiennent à l'aide des fonctions st_j convex hull et st envelope.

```
1 create table convex_hull as
2 SELECT
3 ST_j convex hull(st_union(geometry)) As "geometry"
4 FROM terrasses;
5
6 create table envelope as
7 SELECT
8 ST_j envelope(st_union(geometry)) As "geometry"
9 FROM terrasses;
```

### Code 20 – Polygones d'emprise spatiale minimale en PostGIS

## 2.1.2 Simplification de données surfaciques ou transformation de données surfaciques en données ponctuelles

**Centres de gravité ou centroids** Les centres de gravité sont utilisés pour résumer l'emprise spatiale d'une couche de polygones par des points. En R sf, et PostGIS, ils s'obtiennent à l'aide de la fonction st_j centroid et en Python, la fonction centroid. Ce passage par une donnée ponctuelle est parfois nécessaire pour obtenir une transformation d'une donnée surfacique en un autre type de donnée surfacique. Par exemple, un diagramme de Voronoi sur les centres de gravité d'une couche de polygones permet d'obtenir une version simplifiée de cette couche plus adaptée à une diffusion web (voir figure 6).
**Simplifier la géométrie d'une couche avec mapshaper.** La simplification de la géométrie d'une couche de polygones peut être également obtenue à l'aide de la solution *open source mapshaper* <fnref n="14" />. **mapshaper** a l'avantage de simplifier la géométrie sans altérer son organisation globale : les contigüités entre les objets sont conservées. L'installation de **mapshaper** est détaillée sur la page github as-sociée. Un exemple d'utilisation en ligne de commande est proposée ci-dessous :

1. curl https://minio.lab.sspcloud.fr/h529p3/diffusion/armf_bdtopo Dep_75_2021.zip --output armf.zip

2. unzip armf.zip

3. mapshaper armf_bdtopo Dep_75_2021.shp -simplify dp 20% -o format=shapefile armf_simplify.shp

**Code 21 – Compression du fond de carte des arrondissements parisiens à l'aide de mapshaper**

Le package Rmapshaper permet en R d'accéder à mapshaper sans avoir à écrire des lignes de commande dans un terminal. https://mapshaper.org/ permet d'accéder directement aux fonctions de mapshaper en ligne dans une interface graphique, les transformations sont effectuées en local et les fichiers ne sont pas envoyés sur le serveur de l'application.

<figure><img src="image_12.png" /><figcaption>**FIGURE 6 – Simplification d'un fond de carte : simplifier le contour des arrondissements parisiens**</figcaption></figure>

**Autres modifications de couches surfaciques** Il est aussi possible de créer des zones tampons, des polygones d'emprise minimale autour des couches surfaciques, en appliquant les fonctions décrites dans la partie 2.1.1. Un exemple d'application des zones tampons avec un rayon négatif est proposé dans la partie 2.3.1. Dans ce cas, les tampons servent à éroder les polygones.

<footnote n="14">14. https://mapshaper.org/</footnote>
## 2.2 Transformations complexes mettant en relation plusieurs éléments d'une couche ou de plusieurs couches

Les transformations présentées dans cette partie visent à produire une nouvelle information spatiale et statistique à partir d'une ou de plusieurs couches. Elles peuvent être bien plus longues à calculer que les transformations précédentes et peuvent nécessiter beaucoup de mémoire. Le recours à PostGIS permet généralement de réaliser des calculs impossibles à terminer dans une session R ou Python.

### 2.2.1 Agrégation des données spatiales

L'agrégation spatiale vise à regrouper les unités spatiales entre elles pour former des unités spatiales plus larges. Cette transformation est utile pour déduire à partir d'une couche spatiale disposant d'une granularité fine une couche à la granularité plus large. Par exemple, il arrive souvent de vouloir construire à partir d'un fond départemental un fond régional, ou à partir d'un fond communal un fond des Établissement public de coopération intercommunale (EPCI).

Pour ce faire, deux opérations sont réalisées conjointement : la géométrie des entités appartenant à une même macro-zone est fusionnée, et les données statistiques associées sont agrégées.

En R, avec le package *dplyr*, l'agrégation des données spatiales repose sur la fonction *group_by* et *st_union* comme celle sur données classiques :

```r
library sf)
library(dplyr)

departement = read('dep_francemetro_2020.geojson')

region = departement %>
  group_by(reg) %>
  summarize(geometry = st_union(geometry),
            surf = sum(surf))
```

Code 22 – Agrégation d'une couche géographique en R

En Python, l'agrégation s'obtient avec la commande *dissolve*. Le paramètre *aggfunc* détermine la fonction d'agrégation à appliquer aux données numériques associées.

```python
import geopandas as gpd
departement = gpd.read_file('dep_francemetro_2020.geojson')
region = departement.dissolve(by='reg', aggfunc='sum', as_index=False)
```

Code 23 – Agrégation d'une couche géographique en Python
En SQL, tout comme en R, l'agrégation s'obtient à l'aide de la commande `group by` et de la fonction `st_union` qui regroupent les objets géographiques appartenant à un même groupe.

```sql
1 create table region as
2 SELECT ST_Union(geometry) as "geometry", sum(surf) as "surf"
3 FROM departement
4 GROUP BY reg;
```

Code 24 – Agrégation d'une couche géographique en SQL

### 2.2.2 Intersection et jointure des données spatiales

Une jointure spatiale vise à joindre deux couches spatiales selon la position relative de leurs objets respectifs. Dans une grande majorité des cas, c'est leur intersection qui est déterminée : la jointure relie alors deux objets s'ils ont au moins un point du plan en commun. La jointure spatiale est un outil important pour le statisticien car elle permet de produire les statistiques, comme dénombrer un phénomène selon un découpage spatial. Les exemples ci-dessous comptabilisent les terrasses éphémères contenues dans l'emprise spatiale des carreaux de la source Filosofi. Une fois la jointure réalisée, une agrégation permet de produire les statistiques au carreau.

En R, la jointure s'obtient à l'aide de la fonction `st_join`, par défaut selon l'intersection.

```r
1 library sf)
2 library (tidyverse)
3
4 carreaux = read_SF('Filosofi2015_carreaux_200m_Paris.gpkg')
5 terrasses = read_SF('terrasses.gpkg')
6
7 temp = st_join(carreaux, terrasses, left=F)
8 #left=F indique une jointure droite
9 temp = temp %>
10 st_drop_geometry %>
11 group_by(IdINSPIRE) %>
12 summarize(densite=n())
```

Code 25 – Intersection spatiale en R

En Python, la jointure s'obtient à l'aide de la fonction `sjoin`, par défaut également selon l'intersection.

```python
1 import geopandas as gpd
2 import pandas as pd
3
4 carreaux = gpd.read_file('Filosofi2015_carreaux_200m_Paris.gpkg',
 driver="GPKG")
```
```
5 terrasses = gpd.read_csv('terrasses-ephemeres.csv', sep=";")
```

```
7 temp = gpd.sjoin(carreaux, terrasses)
```

```
8 temp['densite'] = 1
```

```
9 temp = temp.groupby(['IdINSPIRE'], as_index=False)['densite'].sum()
```

## Code 26 – Intersection spatiale en Python

En PostGIS, la jointure spatiale et l'agrégation sont réalisées au sein d'une même requête.

```
1 create table carreaux_terrasses as
```

```
2 SELECT
```

```
3 carreaux."IdINSPIRE",
```

```
4 carreaux_geometry,
```

```
5 count(*) / power(200.0/1000.0,2) as "densite"
```

```
6 FROM carreaux
```

```
7 JOIN terrasses
```

```
8 ON ST_intersects(carreaux_geometry, terrasses_geometry)
```

```
9 group by carreaux."IdINSPIRE", carreaux_geometry
```

## Code 27 – Intersection spatiale en PostGIS

Avec *intersects*, deux objets sont reliés dès qu'ils partagent au moins un point en commun. D'autres types de jointure, tels que *contains* ou *within* peuvent être mobilisés. Avec *contains* les objets sont reliés si aucun point du deuxième objet se trouve en dehors du premier. *within* est identique à *contains* sauf que les deux objets sont permutés dans l'appel de la fonction<fnref n="15" /> Dans le cas d'une jointure d'un fond de carte surfacique avec un fond ponctuel, *intersects* et *contains* et *within* donnent des résultats similaires.

## 2.3 Exemples de manipulation de données spatiales

Cette partie illustre les manipulations de fonds de carte sur trois exemples classiques en statistique publique : la détermination d'une emprise spatiale, ici de zones bâties, un calcul d'accessibilité, ici aux terrasses éphémères parisiennes supposées ponctuelles, ou aux parcs et jardins dont l'emprise spatiale ne se résume pas à un point, et enfin l'agrégation de données carroyées à l'échelle communale.

### 2.3.1 Définir l'emprise spatiale des zones bâties

Définir les limites morphologiques d'une ville est une opération complexe qui repose sur une série d'abstractions géographiques<fnref n="16" />. La principale étant de pro-

<footnote n="15">15. Des explications détaillées peuvent être trouvées dans Entin (2019).</footnote>

<footnote n="16">16. Les bases de données de la BD Topo disponibles sur le site de l'IGN (ign.fr) contiennent l'ensemble des bâtiments de France et permettent de s'exercer à la constitution des zones bâties.</footnote>
duire une continuité à partir de bâtiments à l'emprise spatiale finie et isolés les uns des autres. Pour définir les unités urbaines par exemple, l'Insee retient un seuil de distance de 200 m pour délimiter les zones agglomérées (tâches bâties). Deux bâtiments situés à moins de 200 m appartiennent à la même zone. Les zones agglomérées finales sont obtenues en appliquant cette règle de proche en proche aux agrégats déjà constitués. En géomatique, la production de telles taches bâties repose sur un enchaînement d'opérations élémentaires :

1. Calcul des zones tampons de rayon 100 m autour de chaque bâtiment. Le résultat est une liste de `POLYGON` correspondant à chaque tampon;

2. Regroupement des zones tampons qui se chevauchent pour former une unique zone soit un `MULTIPOLYGON`;

3. Application d'un tampon négatif de 100 m sur la zone précédente pour éroder celle-ci et obtenir une tache dont les frontières se rapprochent des bâtiments;

4. Séparation des différentes sous-parties composant la tache érodée pour obtenir une liste de `POLYGON` représentant chaque tache indépendante.

En R, l'union des tampons (étape 2) est réalisée avec la commande `st_union`. La transformation (étape 4) d'un `MULTIPOLYGON` en une liste de `POLYGON` s'obtient à l'aide de l'instruction `stCast` avec le paramètre `'POLYGON'`.

```r
library sf)
batiment = read.sf("BATIMENT_extrait.geojson")
batiment = st_transform(batiment, 2154)
taches_urbaines = st_buffer(batiment, 100) %>%
  st_union() %>%
  st_buffer(-100) %>%
  stCast('POLYGON')
```

### Code 28 – Délimitation des taches en R

En Python, la commande `unary_union` regroupe les tampons (étape 2) et la fonction ` explode ` réalise la désagrégation d'un `MULTIPOLYGON` en `POLYGON` (étape 4).

```python
import geopandas as gpd
batiment = gpd.read_file('BATIMENT_extrait.geojson').to_crs(2154)
taches_urbaines = batiment_geometry.buffer(100).unary_union.buffer(-100)
taches_urbaines = gpd.GeoSeries(taches_urbaines). explode()
```

### Code 29 – Délimitation des taches en Python

En PostGIS, la désagrégation d'un `MULTIPOLYGON` en `POLYGON` (étape 4) s'effectue à l'aide de l'instruction `st_split`.
<figure><img src="image_13.png" /><figcaption>(a) Fond de carte des bâtiments</figcaption></figure>

<caption>(b) Étape 1 : Calcul des zones tampons de rayon 100 m</caption>

<caption>(c) Étape 2 : Regroupement des zones (d) Étapes 3 et 4 : Application d'un tampon négatif et séparation des différentes sous-parties</caption>

<caption>FIGURE 7 – Les étapes de construction des taches bâties</caption>

1. CREATE table tache_urbaine as
2. select st_split(st_buffer(ST_Union(st_buffer(geometry, 100))), -100)) as "geometry"
3. from batiment;

Code 30 – Délimitation des taches en PostGIS

### 2.3.2 Déterminer l'accessibilité à un équipement

Un deuxième exemple classique consiste à déterminer l'accessibilité d'une population à des équipements. C'est une opération apparemment simple mais qui se révèle difficile à mener lorsque les volumétries des objets considérés sont importantes. Naïvement, l'accessibilité d'une population carroyée dont on a résumé la spatialité au point centroïde de chaque carreau à des équipements eux aussi réduits à un point consiste à calculer deux à deux tous les couples de distances et conserver uniquement ceux dont la distance est inférieure à un certain seuil. Lorsque les ensembles de carreaux et d'équipements sont grands, le produit cartésien entre les
deux ensembles peut être gigantesque et ces calculs impossibles à réaliser dans leur globalité. Pour autant, le calcul d'une grande partie des couples n'a pas d'intérêt, et ce problème peut largement être simplifié. En effet, si le seuil de distance est de 200 m, l'appréciation des distances entre les carreaux de Paris et des équipements à Marseille n'est pas pertinente.

Les structures en arbre, telles que les **k-d trees** ont été développées pour répondre à ce type de question (Bentley, 1975). La structure en arbre permet de stocker des points et de requêter plus rapidement qu'en parcourant le tableau de points. Elle évite de calculer explicitement les distances pour des couples de points éloignés, comme entre Paris et Marseille pour un seuil de distance de 200 m par exemple. Les structures en arbre permettent d'appliquer des algorithmes gourmands en mémoire et en temps de calcul sur des ensembles de données de très grandes tailles comme les algorithmes de *clustering* spatial qui relèvent du *Machine Learning* notamment la *DBSCAN* (voir encadré). Ci-dessous, on cherche à déterminer pour chaque terrasse éphémère ses terrasses voisines dans un rayon de 200 m (voir Figure 8 pour le résultat).

En R, la structure **k-d tree** est implémentée dans la librairie *dbscan*, au sein de laquelle La fonction *frNN* recherche les points voisins dans un rayon donné (*nearest neighbors*).

1 library(dbscan)

2 library(tidyverse)

3

4 terrasses = read_delim('https://opendata.paris.fr/explore/dataset/terasses-ephemeres/download/?format=csv&timezone=Europe/Berlin&lang=fr&use labels for header=true&csv separator=%3B', delimiter =';')

5 voisins = frNN(terasses[,c('Coord X', 'Coord Y')], 200, terrasses[,c('Coord X', 'Coord Y')])

6 voisins = lapply(voisinsid[1:100], length))

7 # la fonction frNN renvoie une liste de points voisins pour chaque point. L'aggregation s'effectue dans un deuxième temps.

Code 31 – Requête spatiale reposant sur les structures **k-d tree** en R

En Python, la librairie *scikit-learn* contient une implémentation de la structure **k-d tree**.

1 import pandas as pd

2 from sklearn.neighbors import KDTree

3

4 #terasses = pd.read_csv(’=Europe/Berlin&lang=fr& use labels for header=true&csv separator=%3B’, sep=";")

5 terrasses = terrasses[terasses['Coord X'] >0]

6 tree = KDTree(terasses[[ 'x', 'y']])

7
```
8 voisins = tree.query_range(terrasses[[ 'x', 'y']], 200,
                             count_only=True)
9 #count_only permet de calculer directement le nombre de voisins
```

Code 32 – Requête spatiale reposant sur les structures k-d tree en Python

<figure><img src="image_14.png" /><figcaption>FIGURE 8 – Localisation des terrasses éphémères disposant de plus de 100 terrasses voisines dans un rayon de 200m</figcaption></figure>

Le recours à une approche géomatique permet des calculs d'accessibilité plus complexes, par exemple l'accessibilité à des équipements dont l'emprise spatiale ne se résume pas à un point. Des zones tampons dont le rayon correspond au seuil de distance permettent alors de déterminer l'emprise spatiale d'une accessibilité à vol d'oiseau. Les codes ci-dessous déterminent la population à moins de 500 m d'un jardin dans la ville de Paris. Les parcs et jardins proviennent d'une base en open data proposée par la ville de Paris et accessible ici (fig 9).

```
1 library sf)
2 population = st_read('Filosofi2015_carreaux_200m_metropole.shp')
3
4 parc = st_read('jardin.geojson')
5 parc = st_transform(parc, 2154)
6 parc = parc[as.numeric(st_area(parc)) > 10000,]
7 buffer = st_buffer(parc, 500)
8
```
```r
9 res_intersection = st_join(population, buffer, left=FALSE)
10 plot(res_intersection$geometry)
```

Code 33 – Déterminer les carreaux à moins de 500 m d'un parc ou d'un jardin en R

```r
1 import geopandas as gpd
```

```r
2 population = gpd.read_file('Filosofi2015_carreaux_200m_metropole.shp')
```

```r
3 parc = gpd.read_file('jardin.geojson').to_crs('EPSG:2154')
```

```r
4 parc = parc[parc.area>10000]
```

```r
5 buffer = parc.copy()
```

```r
6 buffer_geometry = buffer.buffer(500)
```

```r
8 res_intersection = gpd.sjoin(population, buffer, how="inner", op='intersects')
```

Code 34 – Déterminer les carreaux à moins de 500 m d'un parc ou d'un jardin en Python

<figure><img src="image_15.png" /><figcaption>(a) Les parcs et leur tampon de 500 m</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>(b) En vert les carreaux qui intersectent un tampon et en rouge les carreaux isolés</figcaption></figure>

FIGURE 9 – Carreaux à moins de 500m d'un parc
Un exemple de clustering spatial : DBSCAN

DBSCAN pour density-based spatial clustering of applications with noise est un algorithme de clustering proposé par Ester, Kriegel, Sander, et Xu (1996) qui vise à repérer les agrégats spatiaux de points pour former des clusters soient les zones de forte densité de points. Avant cet algorithme, ces zones étaient délimitées en calculant par exemple une densité à l'aide d'une estimation non paramétrique à noyau. Les clusters étaient obtenus en conservant comme agglomérées les zones de densité supérieure à un seuil défini par l'utilisateur. La DBSCAN évite de recourir à une estimation non paramétrique, ce qui en fait un algorithme très rapide. Les implantations sous R (package DBSCAN), Python (librairie scikit-learn) et PostGIS (fonction ST Clus terDBSCAN) repose sur l'exploitation des arbres K-D trees. Les deux paramètres eps et Minpts déterminent la taille des clusters. eps est analogue au rayon de lissage (bandwidth) de l'estimation de densité non paramétrique et Minpts correspond au seuil de troncature des densités. De façon simplifiée, l'algorithme se décompose en trois étapes :

— pour chaque point, son nombre de points voisins dans un rayon *eps* est calculé;

— seuls les points (points cœurs) dont le nombre de voisins est supérieur à *Minpts* sont conservés pour la suite de l'analyse;

— un graphe non orienté est constitué en reliant tous les points cœurs situés à une distance inférieure à *eps*. Les agrégats sont les composantes connexes de ce graphe (tous les points d'une même composante peuvent être reliés par un ou plusieurs chemins).

La DBSCAN a été utilisée dans des travaux de l'Insee pour déterminer l'em-
prise des centres villes en recherchant les polarités commerçantes (Baëhr et
Courthial, 2020) ou encore l'emprise des villages dans le cadre de la mesure
du mitage occasionné par les nouvelles constructions (Himpens, Poulhes, et
Sémécurbe, 2021).

Un exemple d'application de la DBSCAN est proposé ci-dessous sur les ter-
rasses éphémères. Les paramètres de la DBSCAN sont ceux de la partie
2.3.2, eps = 200 et Minpts = 100.

1 library(dbscan)

2 agregats = dbscan(terrasses[,c('Coord X', 'Coord Y')], 200,
50)

3 agregats$cluster # renvoie l'agregat d'appartenance de chaque
point. 0 caracterise les points isoles.

Code 35 – Application de la DBSCAN à la détection des agrégats de terrasses en R
1 from sklearn.cluster import DBSCAN

2 agregats = DBSCAN eps=200, min_sample s=100).fit(sterrasses[[ 'x','y']])

3 agregats.labels_ # renvoie l'agregats d'appartenance. -1
 caracterise les points isoles

Code 36 – Application de la DBSCAN à la détection des agrégats de terrasses en Python

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

Les agrégats contiennent plus de points que ceux obtenus avec la recherche
des points voisins (fig 8), ceci s'explique par le fait que l'algorithme intègre,
par défaut, dans les agrégats les points à une distance inférieur à eps. Ce
comportement peut être modifié sous R avec l'option borderPoints.

2.3.3 Agréger les données carroyées à l'échelle des communes

L'agrégation de données ponctuelles dans des données surfaciques ne posent
pas de difficulté particulière car chaque point appartient à un unique carreau. Il
n'en est pas de même lorsque l'on recherche à déduire d'une information statistique
contenue dans une couche surfacique (A) dans une autre couche surfacique (B).
La stratégie consistant à réaliser l'intersection spatiale puis l'agrégation n'est pas
satisfaisante. Une unité spatiale de la couche (A) pouvant chevaucher plusieurs
unités spatiales de la couche (B), les statistiques agrégées sur la couche (B) sont
susceptibles de contenir de nombreux doubles comptes. Afin de tenir compte des doubles comptes, il est nécessaire de désagréger l'information au prorata de la surface commune. Dans ce cas, la population communale totale estimée correspond exactement dans ce cas à la population carroyée totale. Cette opération augmente le temps de calcul car il est nécessaire de calculer explicitement les intersections entre les communes et les carreaux et est difficilement réalisable sous Python ou R. Elle est réalisable avec PostGIS.

```sql
1 select "INSEE_COM",
2     sum("Men"*interarea/area) "Men",
3     sum("Ind"*interarea/area) "Ind"
4 from (
5     select c."INSEE_COM", f."Men", f."Ind",
6         ST_Area(ST_Intersection(c Official geometry, f Official geometry)) "interarea"
7     ,
8         ST_Area(f Official geometry) "area"
9     from commune as c inner join filosofi as f on
10    ST_Intersects(c Official geometry, f Official geometry)) as cf
    group by cf."INSEE_COM"
```

Code 37 – Calcul d'une population communale à partir de données carroyées au prorata de la surface commune
# 3 Diffusion dynamique de l'information spatiale

La cartographie dynamique est le processus automatisé de construction, de représentation et de manipulation de cartes géographiques. Elle permet de créer des cartes interactives. Elle repose sur l'exploitation de flux mettant à disposition l'information spatiale en fonction des actions de l'utilisateur. Ces flux sont hébergés sur des serveurs cartographiques accessibles à l'aide de requêtes HTML. Les librairies Javascript de cartographies dynamiques telles que Leaflet ou Openlayer permettent de représenter directement ces flux. Il en est de même des logiciels dédiés tels que QGIS. Sous R et Python, l'exploitation des résultats de requêtes HTML permet de charger en mémoire sous forme de fond de carte les données contenues dans les flux.

Dans cette partie seront présentés les principaux flux géographiques et Geoserver l'un des serveurs de flux géographiques le plus utilisé actuellement.

## 3.1 Les flux géographiques

Parmi les nombreux types de flux existants, deux se distinguent par leur popularité : les flux WMS (Web Map Service) et les flux WFS (Web Feature Service).

Les flux WMS envoient des images géoréférencées. Ces images ne contiennent pas d'informations en dehors des coordonnées où elles doivent s'afficher. Ce type de flux est adapté à la diffusion d'orthophotographies (images aériennes ou satellites de la surface terrestre) mais également de couches vectorielles, tel que le contour des communes par exemple. Dans ce cas, l'apparence des communes est déterminée par le serveur cartographique et ne peut être modifiée par l'utilisateur.

Les flux WFS renvoient des données vectorielles dans différents formats (json, csv, shapefile, gml). Les informations géométriques récupérées peuvent s'accompagner de nombreuses informations complémentaires. L'apparence des objets récupérés est modifiable par l'utilisateur.

### 3.1.1 Les flux WMS

Un flux WMS est une adresse url qui se termine par généralement par WMS. L'IGN propose par exemple deux adresses pour répondre à la directive Inspire d'accès aux données géographiques<fnref n="17" /> : https://wxs.ign.fr/inspire/inspire/r/WMS pour accéder aux orthophotographies) et https://wxs.ign.fr/inspire/

<footnote n="17">17. https://geoservices.ign.fr/services-web-inspire</footnote>
inspire/v/WMS pour récupérer sous forme d'images les bâtiments et les routes de France.

L'accès aux données s'obtient en complétant les adresses à l'aide d'un opérateur (à définir avec le mot clé *request*) et des paramètres. Le premier opérateur à connaître est *GetCapabilities* qui permet de décrire les données disponibles dans le flux. Le résultat de la requête https://wxs.ign.fr/inspire/inspire/v/WMS?request=GetCapabilities est un fichier XML (cf. fig 10).

<figure><img src="image_18.png" /><figcaption>FIGURE 10 – Contenu du flux IGN WMS</figcaption></figure>

L'opérateur *GetMap* permet ensuite de récupérer les images d'un flux WMS. La requête suivante permet de récupérer la couche des bâtiments autour du bâtiment White, emplacement de la direction générale de l'Insee :

1. https://wxs.ign.fr/inspire/inspire/v/WMS?layer=BU Building&

2. exception=text/xml&format=image/geotiff&

3. service=WMS&version=1.3.0&request=GetMap&styles=&crs=EPSG:4326&

4. bbox=48.8155,2.3072,48.8173,2.3099&width=256&height=256

Code 38 – Requête WMS pour récupérer une image

Ceci permet d'obtenir les images .png de la figure 11.

En R, la librairie *ows4R* permet d'accéder directement aux flux WFS et WMS sans avoir à écrire une requête et d'enregistrer l'image.

1. library(tiff)

2. url = 'https://wxs.ign.fr/inspire/inspire/v/WMS?LAYERS=BU Building &
<figure><img src="image_19.png" /><figcaption>FIGURE 11 – Bâtiments autour de la Direction Générale de l’Insee à Montrouge suivant différentes définitions de l’image (paramètres *width* et *height*)</figcaption></figure>

```python
3 EXCEPTIONS=text/xml&FORMAT=image/geotiff&SERVICE=WMS&
4 VERSION=1.3.0&REQUEST=GetMap&STYLES=&CRS=EPSG:4326&
5 BBOX=48.8155,2.3072,48.8173,2.3099&WIDTH=256&HEIGHT=256'
6 download.file(url, 'white.tif')
7 white = readTIFF('white.tif')
```

Code 39 – Téléchargement d'une image à partir d'une requête WMS en R

Sous Python, le téléchargement des images issues d'un flux WMS peut s'effectuer à l'aide de la librairie ow (également valable pour les flux WFS). Tout comme sous R, il est également possible de procéder sans librairie *adhoc* en enregistrant directement les résultats d'une requête sous la forme d'une image.

```python
1 import urllib.request
2
3 url = 'https://wxs.ign.fr/inspire/inspire/v/WMS?LAYERS=BU Building &'
4 'EXCEPTIONS=text/xml&FORMAT=image/geotiff&'
5 'SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&STYLES=&CRS=EPSG:4326&'
6 'BBOX=48.8155,2.3072,48.8173,2.3099&WIDTH=256&HEIGHT=256'
7 urllib.request.urlopen(url, 'white.tif')
8 I = plt-imread('white.tif')
```

Code 40 – Téléchargement d'une image à partir d'une requête WMS en Python

### 3.1.2 Les flux WFS

Tout comme les flux WMS, un flux WFS est une adresse HTML qui se termine par WFS et qui doit être également complétée par des opérateurs et des paramètres. L'opérateur *GetCapabilities* reste valable pour décrire les données disponibles dans le flux. En revanche, l'accès aux données s'effectue à l'aide de
5

6

bbox=48.8155,2.3072,48.8173,2.3099'
download.file(url, 'white.geojson')

## Code 42 – Traitement d'une requête WFS en R

Sous Python, tout comme pour les flux WMS, la commande `urllib.request` peut être utilisée pour enregistrer en local les résultats d'une requête WFS :

1

import urllib.request

2

url = 'https://wxs.ign.fr/beta/geoportail/WFS?SERVICE=WFS&'\ 
'VERSION=2.0.0&REQUEST getfeature&'\ 
'typename=CADASTRALPARCELS.PARCELLAIRE_EXPRESS:batiment&'\ 
'outputformat=application/json&'\ 
bbox=48.8155,2.3072,48.8173,2.3099'

5

6

7

urllib.request.urlretrieve(url, 'batiment.geojson')

## Code 43 – Traitement d'une requête WFS en Python

### 3.1.3 Exemple d'application dynamique

Le code HTML suivant vise à illustrer l'intégration d'un flux WFS dans une carte dynamique. La librairie de fonctions JavaScript utilisée est OpenLayer<fnref n="19" /> qui permet une intégration native de ce type de flux contrairement à la libraire Leaflet<fnref n="20" />. De façon classique l'application se décompose en deux parties. Dans la première partie, les instructions servent à charger la librairie OpenLayer et à définir le style de la page. La deuxième partie écrite en JavaScript intègre le flux WFS dans le composant de cartographie dynamique d'OpenLayer.

1

<!DOCTYPE html>

2

<html>

3

<head>

4

<title>Test WFS</title>

5

<script src="https://cdn.jsdelivr.net/gh/openlayers/openlayers.github.io@master/en/v6.6.1/build/ol.js"></script>

6

<link rel="stylesheet" href="https://cdn.jsdelivr.net/gh/openlayers/openlayers.github.io@master/en/v6.6.1/css/ol.css"/>

7

<style>

8

.map {

9

width: 100%;

10

height:400px;

11

}

12

</style>

13

</head>

14

<body>

15

<div id="map" class="map"></div>

16

<script>

<footnote n="19">19. https://openlayers.org/</footnote>

<footnote n="20">20. En revanche, OpenLayer et Leaflet sont à même de traiter de façon native les flux WMS</footnote>
```javascript
var vector = new ol layer Vector({
  //data sources
  source: new ol source Vector({
    format: new ol format GeoJSON(),
    url: 'https://wxs.ign.fr/beta/geoportail/WFS?
    SERVICE=WFS&'+
    'VERSION=2.0.0&REQUEST getfeature&'+
    '=CADASTRALPARCELS.PARCELLAIRE_EXPRESSION:batiment&'+
    'outputformat=application/json&bbox
    =48.8155,2.3072,48.8173,2.3099'
  }),
  //layer style
  style: function(feature, resolution) {
    return new ol style Style({
      stroke: new ol style Stroke({
        color: 'blue'
      })
    });
  }
));

var map = new ol Map({
  layers: [new ol layer Tile({source: new ol source OSM
    ()}), vector],
    target: 'map',
    view: new ol View({center: [2.308,48.816], maxZoom:
      19,
      zoom: 16, projection: 'EPSG:4326'})
});
</script>
</body>
</html>
```

Code 44 – Cartographie dynamique d'un flux WFS à l'aide de la librairie OpenLayer
Pour s'initier en ligne à **Leaflet**, le lecteur peut se référer à la mini-application web proposée par Fabrice Garnes dont le code est accessible à cette adresse https://codepen.io/fabcg/pen/wvWGQdW. Cette application, dont le code (html, css, et JavaScript) est éditable, utilise Leaflet pour afficher des informations cartographiques dans un navigateur et récupère des fonds de carte sur les serveurs d'OpenStreetMap. Une couche image est récupérée à l'aide d'une requête WMS ainsi qu'une couche vectorielle à l'aide d'une requête WFS. Elle mobilise aussi le géocodeur ª de l'IGN basé sur la Base Adresse Nationale. Il est également possible de dupliquer le code sur un projet personnel par exemple, sous codepen.io ou jsfiddle.net en créant un compte ce qui permet d'enregistrer les modifications.

a. L'objectif d'un géocodeur est d'associer à une adresse une position spatiale

## 3.2 GeoServer : diffuser des flux géographiques

Un serveur cartographique transforme une source de données spatiales, par exemple, un fichier *shapefile* ou une table spatiale contenue dans une base de données PostGIS en un flux interrogeable (WMS ou WFS) à l'aide d'une requête HTML. De nombreuses applications informatiques ont été développées pour réaliser cette tâche. Parmi celles-ci, GeoServer est l'une des plus plébiscitées par les cartographes du fait de sa polyvalence et de son interface graphique qui permet une administration simple des données sources et des flux en sortie.

La définition d'un flux s'obtient à l'aide de trois étapes :

1. Tout d'abord, il est nécessaire de définir un espace de travail. Un espace de travail permet d'organiser et de regrouper des flux similaires. Cet espace de travail se retrouve dans les requêtes WFS et WMS et correspond à la première partie des couches contenues dans les paramètres *layers* et *type-name* des requêtes WMS et WFS. Par exemple, dans la requête WFS sur le bâti, CADASTRALPARCELS.PARCELLAIRE_EXPRESS est l'espace de travail qui contient la couche bâtiment.

2. Ensuite, un entrepôt de données, à savoir une source de données qui contient des données *raster* ou vecteurs, (*workspace*), doit être rattaché à l'espace de travail. Cet entrepôt peut être une base PostGIS ou un répertoire de travail contenant des fichiers *shapefile* par exemple.

3. La dernière étape consiste à sélectionner la couche de données spatiales contenue dans l'entrepôt de données qui doit être diffusé. Pour le format de sortie WMS, un style doit être défini et rajouté à la couche.
<figure><img src="image_20.png" /></figure>

GeoServer

<figure><img src="image_21.png" /><figcaption>FIGURE 13 – Interface de GeoServer</figcaption></figure>

Le menu Prévisualisation de la couche permet de tester les flux directement mis en ligne sans avoir à écrire une requête HTML.

L'interface de GeoServer permet d'apprendre ses principaux concepts et fonctionnalités. En revanche, le déploiement de nombreux flux y est fastidieux à réaliser à la main. GeoServer est accessible à l'aide d'une API pour automatiser la mise à disposition des flux. La librairie R geosapi permet d'accéder à cette API sans recourir à des lignes de commande. Il n'existe pas d'équivalent de cette librairie en Python. Les requêtes doivent être écrites explicitement et envoyées au GeoServer directement à l'aide de la librairie *request*. Un exemple est proposé dans ce post de blog.

## 3.3 Pour aller plus loin : les tuiles vectorielles

### 3.3.1 Le principe des tuiles vectorielles

Les flux WFS ne sont pas adaptés à la diffusion de couches géographiques vectorielles volumineuses. En effet, quelle que soit l'échelle d'affichage retenue, l'intégralité de la couche est téléchargée en local, ceci a pour conséquence de ralentir considérablement la réactivité des cartes dynamiques. Les outils de cartographie en ligne tels que Google Maps ou OpenStreetMap reposent sur une technologie adaptée aux données géographiques massives, à savoir les tuiles spatiales<fnref n="21" />. Les tuiles

<footnote n="21">21. L'emprise des tuiles s'observe dans Google Maps lorsque l'on change le zoom : l'information géographique s'affiche par bloc carré. Chaque bloc est une tuile.</footnote>
ont été développées pour représenter des cartes évolutives en fonction de l'échelle cartographique et de l'emprise spatiale choisies. Cette technologie affiche uniquement les détails pertinents pour une échelle d'utilisation donnée : plus l'échelle d'observation est fine, plus les fonds de carte contenus dans les tuiles sont précis.

La figure 14 et le tableau 2 présentent l'organisation d'une tuile vectorielle. L'espace géographique est divisé en tuiles indexées par trois paramètres : l'échelle d'observation, la longitude et la latitude. Les tuiles sont caractérisées par les deux propriétés suivantes :

- pour une échelle donnée, les tuiles forment une partition du plan géographique;

- les tuiles sont organisées selon un principe d'emboîtement d'échelle : les tuiles les plus fines sont strictement contenues dans les tuiles les plus grossières. Plus précisément, chaque tuile se décompose en quatre sous-tuiles à l'échelle plus fine suivante.

<figure><img src="image_22.png" /><figcaption>FIGURE 14 – Indexation des tuiles en fonction de l'échelle et de leur localisation. Source : Documentation QGIS</figcaption></figure>

<table><caption>TABLEAU 2 – Précision géographique en fonction de l'échelle d'observation. Source : Mapbox</caption><thead><tr><td>Niveau d'échelle</td><td>Information visible</td></tr></thead><tbody><tr><td>0</td><td>La terre</td></tr><tr><td>3</td><td>Un continent</td></tr><tr><td>4</td><td>Une grande île</td></tr><tr><td>6</td><td>Une grande rivière</td></tr><tr><td>10</td><td>Une autoroute</td></tr><tr><td>15</td><td>Un bâtiment</td></tr></tbody></table>

En pratique, seule l'information contenue dans les tuiles intersectant l'emprise spatiale de la zone géographique d'observation est téléchargée ce qui réduit le temps d'affichage de la carte 22 (le choix de la zone d'observation définie implicitement l'échelle d'observation). L'accès à des données provenant de tuiles vectorielles dans une application web cartographique s'obtient à l'aide d'une adresse web du type nom_domaine/z/x/y, où z est l'échelle (le zoom), x et y l'emplacement de la tuile.

La société Mapbox, spécialisée en cartographie en ligne, propose le format de

<footnote n="22">22. Sous Google Maps, l'échelle d'observation est le dernier chiffre de l'adresse url. L'adresse www.google.fr/maps/@48.8565926,2.3533182,7z produit une carte centrée sur Paris et couvrant le grand bassin parisien. L'échelle de restitution est égale à 7. En augmentant celle-ci, par exemple en passant de 7 à 12, l'emprise spatiale se réduit et la précision de l'information géographique augmente. La carte recouvre alors la commune de Paris et sa proche périphérie</footnote>
13

14

15

16

17

18

[[providers layers]]
name = "departement"
geometry_fieldname = "geometry"
id_fieldname = "dep"
sql = "SELECT ST_AsBinary(ST_MULTI(geometry)) AS geometry, dep
FROM departement WHERE geometry && !BBOX!"
"

19

20

21

22

23

24

[[providers layers]]
name = "commune"
geometry_fieldname = "geometry"
id_fieldname = "commune"
sql = "SELECT ST_AsBinary(st_multi(geometry)) AS geometry,
commune FROM commune WHERE geometry && !BBOX!"
"

25

26

27

28

29

30

[[providers layers]]
name = "iris"
geometry_fieldname = "geometry"
id_fieldname = "iris"
sql = "SELECT ST_AsBinary(st_multi(geometry)) AS geometry, iris
FROM iris WHERE geometry && !BBOX!"
"

31

32

33

34

[[maps]]
name = "paris"
center = [2.3626, 48.8744, 7.0] # carte centr e sur Paris

35

36

37

38

39

[[maps layers]]
provider layer = "defaultdb.departement"
min zoom = 6
max zoom = 8

40

41

42

43

44

[[maps layers]]
provider layer = "defaultdb.commune"
min zoom = 9
max zoom = 12

45

46

47

48

49

[[maps layers]]
provider layer = "defaultdb.iris"
min zoom = 13
max zoom = 18

[[providers]] sert à établir la connexion avec le serveur PostGIS. [[providers layers]] définit les couches utilisées. En particulier, l'instruction sql est utilisée pour écrire la requête qui permet de sélectionner l'information contenue dans les couches. [[maps]] sert à paramétrer l'échelle et la localisation initiale et [[maps layers]] les échelles d'affichage des couches.

Le lancement du serveur Tegola s'effectue avec la ligne de commande suivante,
le serveur étant lancé en local :
1. `./tegola serve --config=config.toml`

Tegola dispose d'une application de prévisualisation des tuiles dans un navigateur internet à l'adresse `localhost:8080` (fig 15).

L'affichage d'une carte choroplèthe dynamique à l'aide de tuiles vectorielles peut être réalisée directement à l'aide de la librairie `mapbox.gl`<fnref n="26" /> et `OpenLayers`.

Sous Leaflet, il faut rajouter les plugins `Leaflet VECTOR Grid`.

<footnote n="26">26. https://docs.mapbox.com/mapbox-gl-js/example/updating-choropleth/</footnote>
<figure><img src="image_23.png" /><figcaption>(a) Affichage des départements</figcaption></figure>

<figure><img src="image_24.png" /><figcaption>(b) Affichage des communes</figcaption></figure>

<figure><img src="image_25.png" /><figcaption>(c) Affichage des iris</figcaption></figure>

Note : les lignes rouges représentent l'emprise des tuiles. Les numéros dans les tuiles désignent leur indices : l'échelle, la longitude, la latitude.

FIGURE 15 – Application de prévisualisation des tuiles vectorielles de Tegola.
Discussion : le territoire des statisticiens

Construire une statistique territoriale revient à associer à un territoire ou à une portion de l'espace, la mesure statistique d'un phénomène. Cette opération est une abstraction géographique qui ne va pas de soit et mérite d'être questionnée. Dans cette dernière partie, en guise de conclusion, on rappelle les questions que les statisticiens doivent se poser avant de faire une carte ainsi que les risques d'erreurs classiques ou de biais perceptifs propres à la cartographie thématique à garder en mémoire dès qu'on utilise des données spatiales.

Quels messages passent mieux par une carte ? Toute quantité numérique peut être représentée sous forme d'une carte. En pratique, une carte n'a d'intérêt que si elle illustre la dimension spatiale d'un phénomène socio-économique et qu'elle aide à sa compréhension. Une carte doit être structurée pour que le message qu'elle diffuse soit lisible (par exemple opposition nord-sud comme sur la carte de Dupin). En géographie, cette structuration de l'espace peut se formaliser sous la forme de « lois » géographiques. Par exemple, pour décrire et expliquer les concentrations et dispersions des implantations humaines, Tannier (2017) propose cinq lois :

— Friction de la distance : des éléments proches ont davantage de chance d'être en relation (pour des processus) ou d'être similaires (états) que des éléments éloignés car la friction de la distance atténue l'influence mutuelle des phénomènes au fur et à mesure.

— Rôle de l'appartenance à un même milieu : des éléments d'un même milieu ont davantage de chances d'être en relation (ou de se ressembler) que des éléments de milieux différents. Un même milieu peut être autant une même maille administrative (commune, région...), une même nation, les mêmes conditions climatiques, une même langue, un même groupe social...

— Influence de la taille des entités considérées : des entités de grande taille, ayant une large emprise spatiale ou contenant une quantité importante d'individus, d'activités, d'information..., ont davantage d'influence sur les éléments de leur environnement que des entités de petite taille.

— Dépendance d'échelles : pour un type donné d'implantations humaines, l'existence d'un processus de concentration considérant un certain voisinage implique un processus de dispersion considérant un voisinage plus large (ou plus retreint).

— Dépendance dans le temps : chaque phénomène influence sa propre variation dans le temps.

Ces lois sont utiles car elles permettent de soutenir les réflexions sur le rôle de l'espace dans la production des comportements sociaux. Les garder en mémoire est utile pour les statisticiens pour questionner les messages transmis par les cartes.
## Le « piège territorial »
Les territoires (nationaux, régionaux, communaux, locaux ou les carreaux d'un carroyage), appelés également unités spatiales, occupent une place centrale dans l'action publique. Outils de gouvernance, ils sont utilisés pour adapter les politiques publiques aux contextes locaux. Les cartes thématiques réalisées par les statisticiens participent à l'appropriation par les acteurs publics de la connaissance de leurs territoires. Il n'en reste pas moins que le territoire du statisticien n'est pas celui de l'acteur public. John Agnew dénomme « piège territorial » l'erreur consistant à confondre les deux, voir Agnew et Dufoix (2014). Le piège territorial selon John Agnew s'incarne dans trois préjugés géographiques :

- la réification des espaces territoriaux de l'État en unités figées d'espace souverain sûr;

- opposition entre l'intérieur et l'extérieur;

- supposer que l'État territorial préexiste à la société et qu'il en est le conteur alors même que les populations se déplacent dans l'espace et ont plusieurs lieux d'encrage.

Les cartes thématiques reposant sur des découpages territoriaux sont susceptibles de participer à la reproduction de ce piège territorial. La cartographie des résultats de l'élection présidentielle américaine de 2016 illustre la façon dont les cartes peuvent façonner différents discours territoriaux, voir Figure 3.3.2. La première carte du vote majoritaire par comté est une analyse en aplat de couleurs avec uniquement deux teintes séparées par la borne 50% (carte de taux ou carte de ratio). Elle suggère une victoire écrasante du candidat républicain au détriment de la candidate démocrate. Les deux autres cartes nuacent ce constat : les territoires peu denses ont voté majoritairement républicain mais contribuent peu en valeur absolue à la géographie du vote. Ce constat a fait dire à l'auteur de ces cartes, Karim Douieb, que « les territoires ne votent pas. Les gens, oui ».

Faut-il pour autant bannir la carte de ratio des analyses statistiques ? Dans le cas des États-Unis, celle-ci est porteuse d'une information pertinente, car le système électoral privilégie une base territoriale à travers le système des grands électeurs contrairement à un suffrage majoritaire qui ne tient pas compte de l'inégale répartition spatiale des voix. La carte de ratio permet de comprendre les résultats de l'élection, la victoire du candidat républicain. En revanche, les cartes en volume relativisent cette victoire et permettent d'affirmer que l'élection ne fut pas plébiscite pour le vainqueur<fnref n="27" />.

En France, la majorité des personnes défavorisées sont concentrées en volume dans les centres urbains (fig 17). Ces sous-populations sont donc inégalement représentées dans l'espace. Très présentes dans les banlieues des grands pôles, elles peuvent aussi être ségrégées. Les politiques de mixité sociale et d'actions sur les

<footnote n="27">27. Ceci est d'autant plus vrai que la candidate démocrate avait plus de voix en volume que le candidat républicain.</footnote>
<figure><img src="image_26.png" /></figure>

Note : trois manières de visualiser les résultats de l'élection américaine de 2016 (en rouge le vote républicain et en bleu le vote démocrate) : en haut, le gagnant par comté, au milieu rapporté à la population plutôt qu'à la surface du comté et en bas la répartition des votes rapportée à la taille de la population par comté.

Source : Karim Douïeb, sur Twitter.

<caption>FIGURE 16 – “Land doesn’t vote, people do.”</caption>
quartiers de la politique de la ville (dont le découpage correspond à des zones de concentration de ménages en situation de pauvreté) ont vocation à agir sur les effets de la ségrégation sur les populations concernées. Le piège territorial dans ce cas est de considérer que seuls les territoires de concentration décrivent l'organisation des populations défavorisées. Les cartes de volume nuacent ce constat et sont utiles pour les politiques de droit commun qui visent à dimensionner par exemple des équipements dans les centres-villes à destination des populations défavorisées.

<figure><img src="image_27.png" /><figcaption>(a) Part des ménages pauvres parmi les ménages
(b) Densité de ménages pauvres par km²</figcaption></figure>

FIGURE 17 – Deux visions d’un même phénomène : répartition des ménages pauvres dans l’agglomération parisienne. Source : Filosofi 2015

Opposer les représentations en volume et en part n'est pas pertinent, les deux apportent une connaissance différente sur un même phénomène social. Les outils de visualisation dynamique modernes permettent de changer les points de vue en passant d'une représentation à l'autre et ainsi de dépasser le piège territorial.

### L'erreur écologique ou Modifiable Areal Unit Problem en géographie

Même si les données ont une précision métrique, la production d'un savoir géographique nécessite de les relier entre elles à l'aide de règles spatiales, autrement dit de spatialiser les données. La spatialisation consiste à associer à chaque objet une portion, une quantité d'espace bien défini. Cette mise en relation spatiale s'apparente à un point de vue que l'on se donne afin de rendre intelligible les distributions spatiales. Ce point de vue opère comme un filtre qui met en valeur certaines propriétés spatiales et en délaisse d'autres. Il en résulte que les analyses géographiques ont donc toujours une part d'arbitraire et sont fonction de la spatialisation choisie.

Les données spatiales sont ainsi sensibles aux effets d'agrégation appelés
erreur écologique en sciences sociales (Courgeau, 2000) et *Modifiable Areal Unit Problem* en géographie (Openshaw et Pumain, 1981).

L'erreur écologique consiste à inférer des relations statistiques à l'échelle des individus à partir de données agrégées, (Wikipédia, 2021a).<fnref n="28" />

Openshaw (1983) a formalisé sous le concept de *Modifiable Areal Unit Problem* (MAUP) les effets d'agrégation propres à la géographie. Le MAUP se décompose en deux effets interdépendants : l'**effet de zonage** et l'effet d'échelle. Le premier est lié à la forme des unités spatiales élémentaires<fnref n="29" /> tandis que le deuxième désigne la dépendance des analyses géographiques et statistiques à l'échelle d'observation.

<figure><img src="image_28.png" /><figcaption>FIGURE 18 – MAUP : effet de zonage</figcaption></figure>

Dans les analyses statistiques, on peut se prémunir en partie des effets de zonage

<footnote n="28">28. L'analyse de Durkheim sur la prévalence des suicides selon la confession religieuse à l'échelle des provinces de la Prusse est souvent citée comme un exemple manifeste d'erreur écologique, voir Courgeau (2000). Durkeim extrapole un taux de suicide pour les protestants quatre fois supérieur à celui observé sur données individuelles en s'appuyant sur la relation presque linéaire qu'il existe entre taux de suicide et part des protestants dans les provinces de Prusse. Or cette relation reflète surtout des variations de conditions de vie (corrélées positivement avec la part de protestants et explicatives du taux de suicide). Il s'agit là d'un effet de structure, appelé aussi effet de composition, voir Courgeau (2000).</footnote>

<footnote n="29">29. Le redécoupage des circonscriptions électorales par les élus est un exemple classique d'effet de zonage (fig 18) . Aux États-Unis cette pratique porte le nom de *gerrymandering* en référence à Elbridge Gerry qui fut accusé d'avoir redécoupé la circonscription d'un comté afin de favoriser son parti (Wikipédia, 2021b).</footnote>
<figure><img src="image_29.png" /></figure>

Note : densités de population à 4 échelles différentes. Source : Recensement de la population 2017 et Filosofi 2015.

<caption>FIGURE 19 – MAUP : effet d'échelle.</caption>

en recourant à des découpages réguliers tels que des carroyages ou bien encore
en réalisant des lissages spatiaux. En revanche, on ne peut quasiment rien faire
quant aux effets d'échelle - le choix du rayon dans la définition de la concentration
d'un ensemble est un effet d'échelle - ils sont inhérents à la prise en compte de
l'espace dans les analyses statistiques. La figure 19 témoigne des effets d'échelle
sur un indicateur densité. Elle montre comment la variabilité spatiale des densités
augmente avec la finesse de l'échelle d'analyse.

Les effets d'échelle du MAUP ne sont pas uniquement cartographiques et
concernent également les modélisations statistiques. Louf et Barthelemy (2014)
explorent le lien existant entre la population des villes et la production de CO₂ (fig
??). Ils constatent qu'en fonction de la définition territoriale des villes américaines
choisie, les villes les plus vertueuses ne sont plus les mêmes. Avec le découpage en
aire urbaine (définition fonctionnelle), plus les villes sont petites et plus elles sont
vertes alors que la relation s'inverse avec le découpage en unité urbaine (défini-
tion morphologique). Ils concluent qu'en l'absence d'un modèle sous-jacent, il est
dangereux d'interpréter les résultats d'une analyse de données.

Jusqu'à récemment, le MAUP était considéré comme une limite à la connais-
<figure><img src="image_30.png" /></figure>

En rouge, analyse avec le découpage en unité urbaine et en vert avec le découpage en aire urbaine. Source : Louf et Barthelemy (2014)

<caption>FIGURE 20 – Are larger cities greener or smoggier ? Relation entre l'émission de CO₂ et la taille des villes aux Etats-Unis.</caption>

sance géographique. Dans le meilleur des cas, les statisticiens et les géographes le résolvaient en utilisant l'échelle d'action des acteurs publics (mais est-elle la plus pertinente ?) ou en recherchant une échelle caractéristique au phénomène étudié (mais sur quel critère la déterminer?). A présent, le MAUP est pensé comme une source d'information sur l'organisation spatiale des distributions étudiées. Plus précisément, la variation des indicateurs à travers les échelles dévoile le caractère multiscalaire des distributions spatiales étudiées (Grasland et Madelin, 2006). Pour comprendre l'organisation d'une distribution spatiale, on doit nécessairement l'observer à différentes échelles. L'usage de tuiles vectorielles permet dans les cartes dynamiques de faire comprendre le rôle des échelles sur l'observation des phénomènes par les utilisateurs.
Références

AGNEW, J., ET S. DUFOIX (2014) : "Le piège territorial," *Raisons politiques*, (2), 23-51.

BAËHR, A., ET M. COURTHIAL (2020) : "Le commerce de centre-ville recule dans les villes moyennes," *Insee Flash Hauts-de-France*, (91).

BENTLEY, J. L. (1975) : "Multidimensional Binary Search Trees Used for Associative Searching," *Commun. ACM*, 18(9), 509-517.

COURGEAU, D. (2000) : "Réflexions sur la causalité en sciences sociales," *Revue des politiques sociales et familiales*, 60(1), 49-60.

ENTIN, M. (2019) : "Contains/Covers/Intersects/Within?," https://mentin.medium.com/which-predicate-cb608b470471, Medium Post. Accessed : 2022-01-19.

ESTER, M., H.-P. KRIEGEL, J. SANDER, ET X. XU (1996) : "A Density-Based Algorithm for Discovering Clusters in Large Spatial Databases with Noise," in *Proceedings of the Second International Conference on Knowledge Discovery and Data Mining*, KDD'96, p. 226-231. AAAI Press.

GIRAUD, T., ET H. PECOUT (2021) : "Cartographie avec R," https://rcarto.github.io/carto_avec_r/, Accessed : 2021-08-16.

GRASLAND, C., ET M. MADELIN (2006) : "Espon 3.4. 3 The Modifiable Areas Unit Problem. Final Report, 254 p,".

HIMPENS, S., M. POULHES, ET F. SÉMÉCURBE (2021) : "Bâti dispersé, bâti concentré, des disparités territoriales persistantes," *Insee Analyses*, (63).

LAMBERT, N., ET C. ZANIN (2016) : *Manuel de cartographie : principes, méthodes, applications*. Armand Colin.

LE BRAS, H., ET E. TODD (2013) : *Le mystère français*. Seuil Paris.

LOONIS, V., ET M.-P. BELLEFON (2018) : *Manuel d'analyse spatiale. Théorie et mise en œuvre pratique avec R* Insee et Eurostat.

LOUF, R., ET M. BARTHELEMY (2014) : "Scaling : lost in the smog," *Environment and Planning B : Planning and Design*, 41(5), 767-769.
O'BRIEN, R. (2020) : "How To Load Your Pandas DataFrame To Your Database 10x Faster," https://towardsdatascience.com/upload-your-pandas-dataframe-to-your-database-10x-faster-eb6dc6609ddf, Accessed : 2021-08-16.

OPENSHAW, S. (1983) : "The modifiable areal unit problem (Vol. 38)," *Norwich : Geo Books*.

OPENSHAW, S., ET D. PUMAIN (1981) : "Le problème de l'agrégation spatiale en géographie," *L'Espace géographique*, 10(1), 15-24.

PALSKY, G. (1990) : "La cartographie thématique en France : recherches sur ses origines et son évolution jusqu'à la fin du 19e siècle," Ph.D. thesis, Paris 1.

——— (1991) : "La cartographie statistique de la population au XIXe siècle," *Espace Populations Sociétés*, 9(3), 451-458.

——— (1996) : *Des chiffres et des cartes : naissance et développement de la cartographie quantitative française au XIXe siècle*. Comité des travaux historiques et scientifiques Paris.

PEBESMA, E. (2021) : "Simple Features for R," https://r-spatial.github.io/sf/, Accessed : 2021-08-16.

PLAYFAIR, W. (1805) : *An inquiry into the permanent causes of the decline and fall of powerful and wealthy nations*. Books4x Company.

POSTGIS (2021) : "Postgis," https://postgis.net/, Accessed : 2021-08-16.

TANNIER, C. (2017) : "Analyse et simulation de la concentration et de la dispersion des implantations humaines de l'échelle micro-locale à l'échelle régionale - Modèles multi-échelles et trans-échelles," Habilitation à diriger des recherches, Université Bourgogne Franche-Comté.

WIKIPÉDIA (2021a) : "Erreur écologique — Wikipédia, l'encyclopédie libre," En ligne ; Page disponible le 24-août-2021.

——— (2021b) : "Gerrymandering — Wikipédia, l'encyclopédie libre," [En ligne ; Page disponible le 30-novembre-2021].

——— (2021c) : "Système d'information géographique — Wikipédia, l'encyclopédie libre,".

——— (2022) : "Port (logiciel) — Wikipédia, l'encyclopédie libre," [En ligne ; Page disponible le 1-janvier-2022].
Série des Documents de Travail
« Méthodologie Statistique »

9601 : Une méthode synthétique, robuste et efficace pour réaliser des estimations locales de population.

G. DECAUDIN, J.-C. LABAT

9602 : Estimation de la précision d'un solde dans les enquêtes de conjoncture auprès des entreprises.

N. CARON, P. RAVALET, O. SAUTORY

9603 : La procédure FREQ de SAS - Tests d'indépendance et mesures d'association dans un tableau de contingence.

J. CONFAIS, Y. GRELET, M. LE GUEN

9604 : Les principales techniques de correction de la non-réponse et les modèles associés.

N. CARON

9605 : L'estimation du taux d'évolution des dépenses d'équipement dans l'enquête de conjoncture : analyse et voies d'amélioration.

P. RAVALET

9606 : L'économétrie et l'étude des comportements. Présentation et mise en œuvre de modèles de régression qualitatifs. Les modèles univariés à résidus logistiques ou normaux (LOGIT, PROBIT).

S. LOLLIVIER, M. MARPSAT, D. VERGER

9607 : Enquêtes régionales sur les déplacements des ménages : l'expérience de Rhône-Alpes.

N. CARON, D. LE BLANC

9701 : Une bonne petite enquête vaut-elle mieux qu'un mauvais recensement ?

J.-C. DEVILLE

9702 : Modèles univariés et modèles de durée sur données individuelles.

S. LOLLIVIER

9703 : Comparaison de deux estimateurs par le ratio stratifiés et application

aux enquêtes auprès des entreprises.

N. CARON, J.-C. DEVILLE

9704 : La faisabilité d'une enquête auprès des ménages.

1. au mois d'août.

2. à un rythme hebdomadaire

C. LAGARENNE, C. THIESSET

9705 : Méthodologie de l'enquête sur les déplacements dans l'agglomération toulousaine.

P. GIRARD.

9801 : Les logiciels de désaisonnalisation TRAMO & SEATS : philosophie, principes et mise en œuvre sous SAS.

K. ATTAL-TOUBERT, D. LADIRAY

9802 : Estimation de variance pour des statistiques complexes : technique des résidus et de linéarisation.

J.-C. DEVILLE

9803 : Pour essayer d'en finir avec l'individu Kish.

J.-C. DEVILLE

9804 : Une nouvelle (encore une !) méthode de tirage à probabilités inégales.

J.-C. DEVILLE

9805 : Variance et estimation de variance en cas d'erreurs de mesure non corrélées ou de l'intrusion d'un individu Kish.

J.-C. DEVILLE

9806 : Estimation de précision de données issues d'enquêtes : document méthodologique sur le logiciel POULPE.

N. CARON, J.-C. DEVILLE, O. SAUTORY

9807 : Estimation de données régionales à l'aide de techniques d'analyse multidimensionnelle.

K. ATTAL-TOUBERT, O. SAUTORY

9808 : Matrices de mobilité et calcul de la précision associée.

N. CARON, C. CHAMBAZ

9809 : Échantillonnage et stratification : une étude empirique des gains de précision.

J. LE GUENNEC

9810 : Le Kish : les problèmes de réalisation du tirage et de son extrapolation.

C. BERTHIER, N. CARON, B. NEROS

9901 : Perte de précision liée au tirage d'un ou plusieurs individus Kish.

N. CARON

9902 : Estimation de variance en présence de données imputées : un exemple à partir de l'enquête Panel Européen.

N. CARON

0001 : L'économétrie et l'étude des comportements. Présentation et mise en œuvre de modèles de régression qualitatifs. Les modèles univariés à résidus logistiques ou normaux (LOGIT, PROBIT) (version actualisée).

S. LOLLIVIER, M. MARPSAT, D. VERGER

0002 : Modèles structurels et variables explicatives endogènes.
J.-M. ROBIN

0003 : L'enquête 1997-1998 sur le devenir des personnes sorties du RMI - Une présentation de son déroulement.

D. ENEAU, D. GUILLEMOT

0004 : Plus d'amis, plus proches ? Essai de comparaison de deux enquêtes peu comparables.
O. GODECHOT

0005 : Estimation dans les enquêtes répétées : application à l'Enquête Emploi en Continu.

N. CARON, P. RAVALET

0006 : Non-parametric approach to the cost-of-living index.

F. MAGNIEN, J. POUGNARD

0101 : Diverses macros SAS : Analyse exploratoire des données, Analyse des séries temporelles.

D. LADIRAY

0102 : Économétrie linéaire des panels : une introduction.

T. MAGNAC

0201 : Application des méthodes de calages à l'enquête EAE-Commerce.

N. CARON

C 0201 : Comportement face au risque et à l'avenir et accumulation patrimoniale - Bilan d'une expérimentation.

L. ARRONDEL, A. MASSON, D. VERGER

C 0202 : Enquête Méthodologique Information et Vie Quotidienne - Tome 1 : bilan du test 1, novembre 2002.

J.-A. VALLET, G. BONNET, J.-C. EMIN, J. LEVASSEUR, T. ROCHER, P. VRIGNAUD, X. D'HAULTFOEUILLE, F. MURAT, D. VERGER, P. ZAMORA

0203 : General principles for data editing in business surveys and how to optimise it.

P. RIVIERE

0301 : Les modèles logit polytomiques non ordonnés : théories et applications.

C. AFSA ESSAFI

0401 : Enquête sur le patrimoine des ménages - Synthèse des entretiens monographiques.

V. COHEN, C. DEMMER

0402 : La macro SAS CUBE d'échantillonnage équilibré

S. ROUSSEAU, F. TARDIEU

0501 : Correction de la non-réponse et calage de l'enquêtes Santé 2002

N. CARON, S. ROUSSEAU
**0502** : Correction de la non-réponse par répondération et par imputation
**N. CARON**

**0503** : Introduction à la pratique des indices statistiques - notes de cours
**J-P BERTHIER**

**0601** : La difficile mesure des pratiques dans le domaine du sport et de la culture - bilan d'une opération méthodologique
**C. LANDRE, D. VERGER**

**0801** : Rapport du groupe de réflexion sur la qualité des enquêtes auprès des ménages
**D. VERGER**

**M2013/01** : La régression quantile en pratique
**P. GIVORD, X. D'HAULTFOEUILLE**

**M2014/01** : La microsimulation dynamique : principes généraux et exemples en langage R
**D. BLANCHET**

**M2015/01** : la collecte multimode et le paradigme de l'erreur d'enquête totale
**T. RAZAFINDROVONA**

**M2015/02** : Les méthodes de Pseudo-Panel
**M. GUILLERM**

**M2015/03** : Les méthodes d'estimation de la précision pour les enquêtes ménages de l'Insee tirées dans Octopus
**E. GROS, K. MOUSSALAM**

**M2016/01** : Le modèle Logit Théorie et application.
**C. AFSA**

**M2016/02** : Les méthodes d'estimation de la précision de l'Enquête Emploi en Continu
**E. GROS, K.MOUSSALAM**

**M2016/03** : Exploitation de l'enquête expérimentale Vols, violence et sécurité.
**T. RAZAFINDROVONA**

**M2016/04** : Savoir compter, savoir coder. Bonnes pratiques du statisticien en programmation.
**E. L'HOUR, R. LE SAOUT, B. ROUPPERT**

**M2016/05** : Les modèles multiniveaux
**P. GIVORD**

**M. GUILLERM**

**M2016/06** : Econométrie spatiale : une introduction pratique
**P. GIVORD, R. LE SAOUT**

**M2016/07** : La gestion de la confidentialité pour les données individuelles
**M. BERGEAT**

**M2016/08** : Exploitation de l'enquête expérimentale Logement internet-papier
**T. RAZAFINDROVONA**

**M2017/01** : Exploitation de l'enquête expérimentale Qualité de vie au travail
**T. RAZAFINDROVONA**

**M2018/01** : Estimation avec le score de propension sous R
**S. QUANTIN**

**M2018/02** : Modèles semi-paramétriques de survie en temps continu sous R
**S. QUANTIN**

**M2019/01** : Les méthodes de décomposition appliquées à l'analyse des inégalités
**B. BOUTCHENIK, E. COUDIN**

**S. MAILLARD**

**M2020/01** : L'économétrie en grande dimension
**J. L'HOUR**

**M2021/01** : R Tools for JDemetra+ - Seasonal adjustment made easier
**A. SMYK, A. TCHANG**

**M2021/02** : Le traitement du biais de sélection endogène dans les enquêtes auprès des ménages par modèle de Heckman
**L. CASTELL, P. SILLARD**

**M2021/03** : Conception de questionnaires auto-administrés
**H. KOUMARIANOS, A. SCHREIBER**

**M2022/01** : Introduction à la géomatique pour le statisticien : quelques concepts et outils innovants de gestion, traitement et diffusion de l'information spatiale
**F. SEMECURBE, E. COUDIN**