3. Géolocalisation

Parmi le volume exponentiel de données générées par l’activité humaine, un très grand nombre d’entre elles sont géolocalisées. Ceci signifie, qu’au delà de la donnée elle même, on dispose de l’information (éventuellement partielle) du lieu géographique auquel elle est attachée. On entend par information partielle, une adresse, une région administrative, un pays, etc… A cette information partielle, pour les besoins de l’affichage, on doit faire correspondre des coordonnées géographiques, qui sont le plus souvent exprimées par un couple (latitude, longitude). Pour rappel, la latitude et la longitude sont deux mesures angulaires définissant les coordonnées sphériques d’un point à la surface du globe terrestre.

3.1. Géocoding

Le geocoding est l’action de faire correspondre une chaîne de caractère caractéristique d’un lieu géographique avec le couple (latitude, longitude) évoqué ci dessus.

Nous allons travailler avec la Base Adresse Nationale. Elle est co-produite par la Direction Interministérielle du Numérique (la DINSIC), par l’IGN, par La Poste, par la Direction Générale des Finances Publiques (DGFiP) et par l’association OpenStreetMap France, dans le cadre d’une convention. Elle fait partie du Service Public des Données de référence.

L’API correspondante est disponible à cette adresse : https://geo.api.gouv.fr/adresse.

Par exemple, le code suivant permet de récupérer une chaîne de caractères contenant une structure GeoJSON:

import requests, json
import urllib.parse

api_url = "https://api-adresse.data.gouv.fr/search/?q="
adr = "2, boulevard Blaise Pascal, 93160 Noisy le Grand"
r = requests.get(api_url + urllib.parse.quote(adr))

print(r.content.decode('unicode_escape'))

Plusieurs packages ou fonctions sont utilisés:

• requests permet de faire des requêtes HTTP

• la fonction quote() permet d’encoder l’adresse en une chaîne de caractères valide pour le protocole HTTP en remplaçant les caractères réservés par un code. En particulier, le protocole HTTP interdisant les espaces et la ponctuation dans les requêtes, l’adresse doit être transformée (url encoding) avant transmission effective sur le réseau. Ainsi, la chaîne "2, boulevard Blaise Pascal, 93160 Noisy le Grand" est remplacée par "2%2C%20boulevard%20Blaise%20Pascal%2C%2093160%20Noisy%20le%20Grand"

• la fonction decode() utilisée avec le paramètre unicode_escape permet de récupérer les caractères accentués. Ainsi "\u00cele-de-France" est transformé en "Île-de-France".

La requête

https://api-adresse.data.gouv.fr/search/?q=2%2C+boulevard+Blaise+Pascal%2C+93160+Noisy+le+Grand

renvoie

{
    "type": "FeatureCollection",
    "version": "draft",
    "features": [
        {
            "type": "Feature",
            "geometry": {
                "type": "Point",
                "coordinates": [
                    2.583275,
                    48.838151
                ]
            },
            "properties": {
                "label": "2 Boulevard Blaise Pascal 93160 Noisy-le-Grand",
                "score": 0.9577360238468353,
                "housenumber": "2",
                "id": "93051_0258_00002",
                "type": "housenumber",
                "name": "2 Boulevard Blaise Pascal",
                "postcode": "93160",
                "citycode": "93051",
                "x": 669412.64,
                "y": 6859869.58,
                "city": "Noisy-le-Grand",
                "context": "93, Seine-Saint-Denis, Île-de-France",
                "importance": 0.5350962623151893,
                "street": "Boulevard Blaise Pascal"
            }
        }
    ],
    "attribution": "BAN",
    "licence": "ODbL 1.0",
    "query": "2, boulevard Blaise Pascal, 93160 Noisy le Grand",
    "limit": 5
}

On trouvera plus d’informations sur le format GeoJSON sur la page Wikipedia et se familiariser avec le format en générant ses propres données sur geojson.io.

On peut vérifier la validité de l’opération de geocoding avec Google Maps.

http://maps.google.com/?q=48.838151,2.583275.

3.2. Reverse geocoding

L’opération de reverse geocoding permet de convertir une coordonnée (latitude, longitude) en une chaîne de caractères contenant une adresse ou un lieu.

Le processus est tout à fait similaire au geocoding et nécessite de forger une URL à partir d’un couple (latitude, longitude). Ainsi la requête

https://api-adresse.data.gouv.fr/reverse/?lon=2.583&lat=48.838

fournit une structure JSON contenant le résultat du reverse geocoding:

"label": "2 Boulevard Blaise Pascal 93160 Noisy-le-Grand"

3.2.1. Automatiser le processus

Pour automatiser le processus de (reverse) geocoding, un script Python devra:

• forger l’URL à partir de l’adresse à géocoder, ou des coordonnées à « reverse géocoder »

• encoder cette URL avec urllib.parse.urlencode()

• utiliser le module requests pour récupérer le résultat

• convertir l’objet JSON retourné en un dictionnaire Python avec la fonction json.loads()

• exploiter la structure du dictionnaire obtenu pour récupérer les données intéressantes

3.3. Tracer une carte

Python dispose de plusieurs modules permettant la représentation de données géolocalisées. Celui que nous allons utiliser dans ce cours est Folium qui permet d’utiliser l’écosystème Python pour créer des cartes interactives basées sur la bibliothèque JavaScript Leaflet.js. Ces cartes sont visualisables dans un navigateur web.

Cette bibliothèque est très largement utilisée:

• https://whatismyipaddress.com/

• https://fr.foursquare.com

• https://www.accuweather.com

• https://www.openstreetmap.org

• http://www.ign.com/maps/the-elder-scrolls-5-skyrim/skyrim

• https://earthquake.usgs.gov/earthquakes/eventpage/usp000jkn8#map

• https://www.velib-metropole.fr/map#/

Visitez la gallerie pour avoir une idée plus globale de ce que l’on peut faire. Et de comment on peut le faire…

3.3.1. Afficher une carte

La figure ci dessous est une capture d’écran d’une carte générée avec Folium dont on trouvera la documentation complète ici.

Le code utilisé pour produire cette figure est donné ci dessous:

import folium

coords = (48.8398094,2.5840685)
map = folium.Map(location=coords, tiles='OpenStreetMap', zoom_start=15)
map.save(outfile='map.html')

Comme on le voit, afficher une carte nécessite simplement deux actions :

• instancier la classe Map en passant au constructeur le paramètre location qui définit le centre de la carte. tiles définit la cartographie utilisée et zoom_start le niveau de zoom appliqué.

• appeler la méthode save() sur l’objet créé à l’étape précédente. La carte ainsi générée est visible dans un navigateur web.

Pour plus de renseignements, consulter la documentation de la classe Map.

3.3.2. Positionner un marqueur

Folium permet de positionner un marqueur sur un fond de carte, pour identifier un POI (Point Of Interest) particulier.

Le code utilisé pour produire cette figure est donné ci dessous:

import folium

coords = (48.8398094,2.5840685)
map = folium.Map(location=coords, tiles='OpenStreetMap', zoom_start=15)
coords = [48.8490591,2.577023]
folium.Marker(location=coords, popup = "ESIEE Paris").add_to(map)
map.save(outfile='map.html')

Ajouter un marqueur, c’est créer un objet de type Marker dont au moins l’attribut location est renseigné. Pour plus de renseignements, consulter la documentation de la classe Marker.

3.3.3. Afficher des données sur la carte

Les données sont représentées sur une carte grâce à des symboles de forme, de taille et de couleur variables.

Le fichier meteo2014.json contient les températures relevées en 2014 pour l’ensemble les stations météo du territoire français. Une fois importées dans Python avec la fonction json.load(), les données sont organisées dans des dictionnaires imbriqués au format METEO[year][month][day][hour][code].

Pour ce qui va suivre, on va considérer uniquement un sous ensemble de ces données:

• les stations météo de France métropolitaine (ordonnées alphabétiquement) placées dans la liste STATIONS

• les latitudes de ces stations placées dans la liste ordonnée LATS

• les longitudes de ces stations placées dans la liste ordonnée LONGS

• les températures moyennes du mois de janvier 2014, relevées à 12:00, placées dans la liste ordonnée TEMPS

Les données extraites du fichier JSON sont les suivantes:

STATIONS = ['ABBEVILLE', 'AJACCIO', 'ALENCON', 'BALE-MULHOUSE',
            'BELLE ILE-LE TALUT', 'BORDEAUX-MERIGNAC', 'BOURGES',
            'BREST-GUIPAVAS', 'CAEN-CARPIQUET', 'CAP CEPET',
            'CLERMONT-FD', 'DIJON-LONGVIC', 'EMBRUN', 'GOURDON',
            'LE PUY-LOUDES', 'LILLE-LESQUIN', 'LIMOGES-BELLEGARDE',
            'LYON-ST EXUPERY', 'MARIGNANE', 'MILLAU', 'MONT-DE-MARSAN',
            'MONTELIMAR', 'MONTPELLIER', 'NANCY-OCHEY',
            'NANTES-BOUGUENAIS', 'NICE', 'ORLY', 'PERPIGNAN',
            "PLOUMANAC'H", 'POITIERS-BIARD', 'PTE DE CHASSIRON',
            'PTE DE LA HAGUE', 'REIMS-PRUNAY', 'RENNES-ST JACQUES',
            'ROUEN-BOOS', 'ST GIRONS', 'STRASBOURG-ENTZHEIM',
            'TARBES-OSSUN', 'TOULOUSE-BLAGNAC', 'TOURS', 'TROYES-BARBEREY']

LATS = [50.136, 41.918, 48.4455, 47.614333, 47.294333,
        44.830667, 47.059167, 48.444167, 49.18, 43.079333,
        45.786833, 47.267833, 44.565667, 44.745, 45.0745,
        50.57, 45.861167, 45.7265, 43.437667, 44.1185,
        43.909833, 44.581167, 43.577, 48.581, 47.15,
        43.648833, 48.716833, 42.737167, 48.825833,
        46.593833, 46.046833, 49.725167, 49.209667,
        48.068833, 49.383, 43.005333, 48.5495, 43.188,
        43.621, 47.4445, 48.324667]

LONGS = [1.834, 8.792667, 0.110167, 7.51, -3.218333, -0.691333,
         2.359833, -4.412, -0.456167, 5.940833, 3.149333, 5.088333,
         6.502333, 1.396667, 3.764, 3.0975, 1.175, 5.077833, 5.216,
         3.0195, -0.500167, 4.733, 3.963167, 5.959833, -1.608833,
         7.209, 2.384333, 2.872833, -3.473167, 0.314333, -1.4115,
         -1.939833, 4.155333, -1.734, 1.181667, 1.106833, 7.640333,
         0.0, 1.378833, 0.727333, 4.02]

TEMPS = [7.6, 13.5, 7.6, 6.8, 10.5, 11.5, 8.5, 9.7, 8.6, 11.8, 9.1,
         7.2, 5.7, 9.2, 6.0, 7.2, 7.6, 8.4, 12.0, 6.1, 11.6, 9.6, 11.7,
         6.5, 10.0, 11.7, 8.1, 12.6, 9.9, 9.1, 10.8, 9.5, 7.4, 9.0,
         7.1, 10.3, 6.7, 10.8, 10.6, 8.4, 8.1]

Pour tracer une carte, il faut importer le module folium et créer un objet Map:

import folium
coords = (46.539758, 2.430331)
map = folium.Map(location=coords, tiles='OpenStreetMap', zoom_start=6)

Note

Les coordonnées utilisées ici sont celles du centre de la France.

Pour chaque station on crée un CircleMarker dont on définit:

• la position ;

• le rayon (proportionnel à la température) ;

• la couleur.

Ce CircleMarker est ajouté à la carte avec la méthode add_to().

for i in range(len(STATIONS)):
    folium.CircleMarker(
        location = (LATS[i], LONGS[i]),
        radius = TEMPS[i]*2,
        color = 'crimson',
        fill = True,
        fill_color = 'crimson'
    ).add_to(map)

Le résultat obtenu est conforme aux attentes:

3.3.4. Customiser l’affichage

Pour produire quelque chose de plus visuel, on peut utiliser également la couleur pour coder la température. L’utilisation d’une colormap est appropriée.

Une colormap (ou palette de couleurs) est un ensemble prédéfini de couleurs que l’on peut mapper (linéairement par exemple) avec des valeurs numériques (les données). On utilisera ici la convention classique, d’utiliser le bleu pour les températures basses, et le rouge pour les températures élevées.

import folium, branca

coords = (46.539758, 2.430331)
map = folium.Map(location=coords, tiles='OpenStreetMap', zoom_start=6)

cm = branca.colormap.LinearColormap(['blue', 'red'], vmin=min(TEMPS), vmax=max(TEMPS))
map.add_child(cm) # add this colormap on the display

for lat, lng, size, color in zip(LATS, LONGS, TEMPS, TEMPS):
    folium.CircleMarker(
        location=[lat, lng],
        radius=size,
        color=cm(color),
        fill=True,
        fill_color=cm(color),
        fill_opacity=0.6
    ).add_to(map)

map.save(outfile='map.html')

Le résultat obtenu montre sans surprise que les températures sont plus élevées dans le sud de la France et sur la côte Atlantique:

3.3.5. Ajout d’objets géométriques

Il peut être parfois nécessaire de superposer à la carte des informations géographiques complémentaires. Par exemple, les contours administratifs des entités administratives (communes, cantons, départements).

Le format le plus répandu pour le stockage des informations géographiques est le format shapefile (SHP) développé par l’entreprise ESRI (La spécification complète).

Cependant, ce format présente un certain nombre d’inconvénients. Pour le développement web, le format GeoJSON est utilisé.

Astuce

Comme un grand nombre de ressources géographiques sont encore au format SHP, il pourra être utile d’utiliser le convertisseur en ligne MapShaper pour les convertir au format GeoJSON.

3.3.5.1. Les données

data.gouv.fr fournit le découpage administratif des communes françaises au format GeoJSON. Les données sont disponibles dans le fichier datagouv-communes.geojson (310 MB).

Pour illustrer l’ajout d’objets géométriques, on va travailler avec un sous ensemble des données globales et ne conserver que les données relatives aux contours administratifs des communes d’Ile de France.

L’extraction des données qui concernent l’Ile de France est effectuée à partir du fichier global en utilisant le package geopandas.

Astuce

Les deux méthodes ci dessous donnent de bons résultats pour l’installation de geopandas qui peut s’avérer complexe à cause des dépendances.

Dans un environnement Anaconda (Windows, Mac, Linux)

$ conda install geopandas

Avec une installation Python (Windows) standard

$ python -m pip install pipwin
$ python -m pipwin install gdal
$ python -m pipwin install fiona
$ python -m pip install geopandas

geopandas permet de lire, traiter et écrire les données au format GeoJSON de façon simple et efficace en reprenant le concept de data frames de pandas. En particulier, la fonction pandas.concat() permet d’agréger plusieurs data frames. Le code correspondant:

import geojson, geopandas, pandas

# lecture du fichier global
france = geopandas.read_file("datagouv-communes.geojson")

l = []
# sélection des données d'Ile de France
for dpt in ["75", "77", "78", "91", "92", "93", "94", "95"]:
    dptidf = france[france["code_commune"].str.startswith(dpt)]
    l.append(dptidf)

# construction de la GeoDataFrame correspondante
idf = pandas.concat(l)

# écriture dans un fichier
with open("idf.geojson", "w") as f:
    geojson.dump(idf, f)

Après exécution du code ci dessus, le fichier idf.geojson doit être présent dans le répertoire de travail. On peut également le trouver ici (6.5 MB).

3.3.5.2. L’affichage

A partir des données du fichier idf.geojson, on peut afficher le contour des communes d’Ile de France sur la carte.

import folium

coords = (48.7453229,2.5073644)
map = folium.Map(location=coords, tiles='OpenStreetMap', zoom_start=9)

#style function
sf = lambda x :{'fillColor':'#E88300', 'fillOpacity':0.5, 'color':'#E84000', 'weight':1, 'opacity':1}

folium.GeoJson(
    data="idf.geojson",
    name="idf",
    style_function= sf
).add_to(map)

map.save(outfile='map.html')

Le résultat obtenu:

Note

https://htmlcolorcodes.com/, https://color.adobe.com, http://www.paletton.com sont d’excellentes ressources pour le choix des couleurs définies dans la lambda function sf().

3.3.6. Cartes choroplèthes

Une carte choroplèthe est une carte pour laquelle des surfaces élémentaires (définies par un ou plusieurs fichiers GeoJSON) sont colorées conformément à une donnée (une statistique par exemple).

On peut par exemple produire une carte de la population par commune en Ile de France en croisant le contour administratif des communes avec les données de recensement de l”INSEE. Ces données de recensement sont disponibles dans le fichier insee-pop-communes.csv (1 MB).

Tracer une telle carte requiert plusieurs étapes:

• préparation des données géographiques

• préparation des données numériques

• création d’une instance de Folium.Map

• utilisation de la classe Choropleth à cette instance

3.3.6.1. Préparation des données géographiques

L’observation de la structure du fichier GeoJSON montre que:

• les données sont contenues dans une liste de dictionnaires, un par commune

• pour chaque chaque dictionnaire

  • la clé properties contient le code de la commune. Ce code servira à faire le lien avec le dataset de l’INSEE

  • la clé geometry contient les données géographiques de la commune sous la forme d’un polygone dont les sommets sont définis par les coordonnées GPS.

{   "type":"FeatureCollection",
    "features": [
        {   "type": "Feature",
            "geometry":
                {   "type": "Polygon",
                    "coordinates": [[
                        [2.304785272419358, 48.67849158613953],
                        ...
                        [2.304785272419358, 48.67849158613953]
                    ]]
                },
            "properties":
                {   "code_epci": "200056232",
                    "commune": "Longjumeau",
                    ...
                    "code_commune": "91345",
                    ...
                    "intitule_commune": "91345 Longjumeau"}
                },
        },
        ...
    ]
}

3.3.6.2. Préparation des données numériques

Les données doivent être lues et stockées dans une data frame Pandas:

>>> pop_data = pd.read_csv('insee-pop-communes.csv', sep=';')
>>> pop_data
    DEPCOM                       COM   PMUN  PCAP   PTOT
0      01001  L' Abergement-Clémenciat    776    18    794
1      01002    L' Abergement-de-Varey    248     1    249
2      01004         Ambérieu-en-Bugey  14035   393  14428
3      01005       Ambérieux-en-Dombes   1689    34   1723
4      01006                   Ambléon    111     6    117
...      ...                       ...    ...   ...    ...
34990  97419               Sainte-Rose   6418    79   6497
34991  97420            Sainte-Suzanne  23505   199  23704
34992  97421                   Salazie   7312    75   7387
34993  97422                 Le Tampon  78629  1076  79705
34994  97423         Les Trois-Bassins   7139    95   7234

[34995 rows x 5 columns]

La colonne DEPCOM contient le code INSEE des communes:

>>> dpt_code = pop_data['DEPCOM']
>>> dpt_code
0        01001
1        01002
2        01004
3        01005
4        01006
        ...
34990    97419
34991    97420
34992    97421
34993    97422
34994    97423
Name: DEPCOM, Length: 34995, dtype: object

La création d’un masque booléen permet de construire le sous ensemble recherché:

>>> mask = ( ( dpt_code.str.startswith('75')  )
...         | ( dpt_code.str.startswith('77') )
...         | ( dpt_code.str.startswith('78') )
...         | ( dpt_code.str.startswith('91') )
...         | ( dpt_code.str.startswith('92') )
...         | ( dpt_code.str.startswith('93') )
...         | ( dpt_code.str.startswith('94') )
...         | ( dpt_code.str.startswith('95') ) )
>>> mask
0        False
1        False
2        False
3        False
4        False
        ...
34990    False
34991    False
34992    False
34993    False
34994    False
Name: DEPCOM, Length: 34995, dtype: bool

Il est utilisé pour effectuer le filtrage par ligne:

>>> pop_data = pop_data[mask]
>>> pop_data
    DEPCOM                       COM   PMUN  PCAP   PTOT
29297  75101  Paris 1er Arrondissement  16266   129  16395
29298  75102   Paris 2e Arrondissement  20900   142  21042
29299  75103   Paris 3e Arrondissement  34115   274  34389
29300  75104   Paris 4e Arrondissement  28088   282  28370
29301  75105   Paris 5e Arrondissement  58850   781  59631
...      ...                       ...    ...   ...    ...
34878  95676        Villers-en-Arthies    503    10    513
34879  95678             Villiers-Adam    859    11    870
34880  95680           Villiers-le-Bel  27676   132  27808
34881  95682           Villiers-le-Sec    185     3    188
34882  95690       Wy-dit-Joli-Village    333     7    340

[1287 rows x 5 columns]

3.3.6.3. La carte proprement dite

Après création de la carte, il faut utiliser la classe Choropleth dont les paramètres importants sont:

• geo_data : les données géographiques au format GeoJSON ;

• pop_data : les données numériques au format Pandas ;

• columns : la paire clé/valeur des données numériques ;

• key_on : nom de la variable géographique permettant de faire le lien entre données géographiques et numériques.

coords = (48.7190835,2.4609723)
map = folium.Map(location=coords, tiles='OpenStreetMap', zoom_start=9)

folium.Choropleth(
    geo_data=geo_data,                              # geographical data
    name='choropleth',
    data=pop_data,                                  # numerical data
    columns=['DEPCOM', 'PTOT'],                     # numerical data key/value pair
    key_on='feature.properties.code_commune',       # geographical property used to establish correspondance with numerical data
    fill_color='YlGn',
    fill_opacity=0.7,
    line_opacity=0.2,
    legend_name='Population'
).add_to(map)

map.save(outfile='map.html')

Note

La ressource Python Folium: Create Web Maps From Your Data est un bon complément à ce qui précède.

Avertissement

Le fichier GeoJSON décrit l’agglomération parisienne comme une seule et même entité alors que les données numériques de population découpent la capitale en arrondissements. Lorsque la correspondance ne peut être établie entre les données géographiques et les données numériques, aucune donnée n’est affichée pour le polygone en question.

Deux solutions sont possibles :

• on insère les données géographiques des arrondissements dans le fichier GeoJSON ;

• ou on aggrège les données de population des arrondissements pour constituer une population unique pour Paris dans les données numériques.

3.4. Exercice

Il n’y a pas de fichier squelette associé à cet exercice.

La population totale de la commune est une donnée importante mais elle ne reflète pas complètement l’urbanisation du territoire. La densité de population serait une grandeur plus intéressante dans ce cas. La surface de chacune des communes est une donnée qui n’apparait pas dans les datasets utilisés dans le cours.

On pourrait imaginer la calculer à partir des coordonnées du polygone correspondant. Mais cette opération est complexe (voir l’algorithme du lacet).

On peut également la rechercher dans une source de données différente. Utiliser ce dataset pour produire une carte choroplèthe de la densité de population des communes d’Ile de France.