Application 1 - Données spatiales

Préparer son environnement pour ce TD

Les données ont été préparées en amont de ce TD. Un lien de lancement rapide est disponible ci-dessus qui met à disposition un environnement prêt à l’emploi sur le SSPCloud.

Après avoir cliqué sur le bouton, il convient de créer un projet RStudio depuis le dossier appli1 (File > New Project) :

Git est normalement préconfiguré dans ce dossier, vous pourrez donc pusher votre travail sur Github si vous créez un dépôt dessus.

Si la récupération des données a échoué pour une raison x ou y, vous pouvez lancer la récupération des données en copiant ce code dans un terminal

#!/bin/bash

mkdir -p appli1
cd appli1
echo "data/" >> .gitignore
git init
git branch -m main

mc cp s3/projet-formation/nouvelles-sources/data/geoparquet/dvf.parquet data/dvf.parquet
mc cp s3/projet-formation/nouvelles-sources/data/geoparquet/carreaux.parquet data/carreaux.parquet
mc cp s3/projet-formation/nouvelles-sources/data/triangle.geojson data/triangle.geojson
mc cp s3/projet-formation/nouvelles-sources/data/malakoff.geojson data/malakoff.geojson
mc cp s3/projet-formation/nouvelles-sources/data/montrouge.geojson data/montrouge.geojson

En premier lieu, ce TD utilise une source administrative nommée DVF (« Demandes de Valeurs Foncières »).

Les données DVF recensent l’ensemble des ventes de biens fonciers réalisées au cours des dernières années, en métropole et dans les départements et territoires d’outre-mer — sauf à Mayotte et en Alsace-Moselle. Les biens concernés peuvent être bâtis (appartement et maison) ou non bâtis (parcelles et exploitations). Les données sont produites par la Direction générale des finances publiques (DGFip). Elles proviennent des actes enregistrés chez les notaires et des informations contenues dans le cadastre. Cette base a été filtrée de manière à être la plus pédagogique possible pour cette formation.

L’analyse de ces données sera complétée des données Filosofi produites par l’Insee :

Les données spatiales carroyées à 200m, produites par l’Insee à partir du dispositif Filosofi, contentant des informations socio-économiques sur les ménages.

Enfin, nous proposons trois contours géographiques ad hoc :

La commune de Malakoff
La commune de Montrouge
Le “Triangle d’or” de Malakoff (autrement dit, son centre-ville à peu de choses près)

L’objectif de ce TD est d’illustrer la manière dont peuvent être traitées des données spatiales de manière flexible avec duckdb.

Préparation de l’environnement

Les librairies suivantes seront utilisées dans ce TD, vous pouvez d’ores et déjà les charger dans votre environnement.

library(duckdb)
library(glue)
library(dplyr)
library(dbplyr)
library(mapview)

Si celles-ci ne sont pas installées, vous pouvez faire en console un install.packages (voir (note-bp-install?)).

Note

Les installations de packages sont à faire en console mais ne doivent pas être écrites dans le code. Bien que ce ne soit pas l’objet de ce cours, il est utile de suivre les bonnes pratiques recommandées à l’Insee et plus largement dans le monde .

Pour en savoir plus, vous pourrez explorer le portail de formation aux bonnes pratiques.

Nous allons avoir besoin des codes Insee suivants pour notre application :

cog_malakoff <- "92046"
cog_montrouge <- "92049"

Import des données

L’import des contours en se fait assez naturellement grâce à sf.

triangle <- sf::st_read("data/triangle.geojson", quiet=TRUE)
malakoff <- sf::st_read("data/malakoff.geojson", quiet=TRUE)
montrouge <- sf::st_read("data/montrouge.geojson", quiet=TRUE)

En premier lieu, on peut visualiser la ville de Malakoff :

mapview(malakoff) + mapview(triangle, col.regions = "#ffff00")

Et ensuite les contours de Montrouge :

mapview(montrouge)

Préparation de DuckDB

DuckDB est un moteur de base de données analytique en mémoire, optimisé pour les requêtes SQL sur des données volumineuses, particulièrement adapté aux fichiers plats comme Parquet ou CSV, et intégrable dans des langages comme Python, R ou SQL.

En principe, duckdb fonctionne à la manière d’une base de données. Autrement dit, on définit une base de données et effectue des requêtes (SQL ou verbes tidyverse) dessus. Pour créer une base de données, il suffit de faire un read_parquet avec le chemin du fichier.

La base de données se crée tout simplement de la manière suivante :

con <- dbConnect(duckdb::duckdb())
dbExecute(con, "INSTALL spatial;")
dbExecute(con, "LOAD spatial;")

Nous verrons ultérieurement pourquoi nous avons besoin de cette extension spatiale.

Cette connexion duckdb peut être utilisée de plusieurs manières. En premier lieu, par le biais d’une requête SQL. dbGetQuery permet d’avoir le résultat sous forme de dataframe puisque la requête est déléguée à l’utilitaire duckdb qui est embarqué dans les fichiers de la librairie :

out <- dbGetQuery(
  con,
  glue(  
    'SELECT * EXCLUDE (geometry) FROM read_parquet("data/dvf.parquet") LIMIT 5'
  )
)
out

La chaîne d’exécution ressemble ainsi à celle-ci :

Même si DuckDB simplifie l’utilisation du SQL en proposant de nombreux verbes auxquels on est familier en R ou Python, SQL n’est néanmoins pas toujours le langage le plus pratique pour chaîner des opérations nombreuses. Pour ce type de besoin, le tidyverse offre une grammaire riche et cohérente. Il est tout à fait possible d’interfacer une base duckdb au tidyverse. On pourra donc utiliser nos verbes préférés (mutate, filter, etc.) sur un objet duckdb : une phase préliminaire de traduction en SQL sera automatiquement mise en oeuvre :

table_logement <- tbl(con, glue('read_parquet("data/dvf.parquet")'))
table_logement %>% head(5)

Partie 1 : Prix immobiliers à Malakoff et à Montrouge

Dans cette partie, l’objectif est d’extraire de l’informations d’une base de données volumineuse à l’aide de DuckDB. Pour le moment, le caractère spatial des données est mis de côté : on découvre et on traite les données via des requêtes attributaires classiques.

Tentons, par une première série d’exercices, de comparer la médiane des prix des transactions immobilières à Malakoff et à Montrouge.

Dans cette partie, nous allons pouvoir faire nos traitements de données avec SQL et/ou tidyverse. Cela illustre l’une des forces de duckdb, à savoir son excellente intégration avec d’autres écosystèmes dont nous sommes familiers.

Lorsque nous irons sur l’aspect spatial, on devra passer en SQL pur, l’écosystème tidyverse n’étant pas encore finalisé pour le traitement de données spatiales avec duckdb.

Premières requêtes SQL : description des données DVF

Tout d’abord, il convient de se familiariser avec les données. Les requêtes proposées pour l’exercice 1 permettent d’obtenir des informations primordiale de manière très rapide et sans nécessité de charger l’ensemble des données dans la mémoire vive.

Exercice 1

Cet exercice nous fera rentrer progressivement dans les données à partir de quelques requêtes basiques.

Lire les 10 premières lignes des données par l’approche SQL et par l’approche tidyverse.
Afficher les noms des colonnes selon les deux approches.
Regarder les valeurs uniques de la colonne nature_mutation selon les deux approches.
Calculer les bornes min et max des prix des transactions selon ces deux approches.

A la question 1, vous devriez avoir :

A la question 2, la liste des colonnes donnera plutôt

 [1] "id_mutation"                  "date_mutation"               
 [3] "numero_disposition"           "nature_mutation"             
 [5] "valeur_fonciere"              "adresse_numero"              
 [7] "adresse_suffixe"              "adresse_nom_voie"            
 [9] "adresse_code_voie"            "code_postal"                 
[11] "code_commune"                 "nom_commune"                 
[13] "code_departement"             "ancien_code_commune"         
[15] "ancien_nom_commune"           "id_parcelle"                 
[17] "ancien_id_parcelle"           "numero_volume"               
[19] "lot1_numero"                  "lot1_surface_carrez"         
[21] "lot2_numero"                  "lot2_surface_carrez"         
[23] "lot3_numero"                  "lot3_surface_carrez"         
[25] "lot4_numero"                  "lot4_surface_carrez"         
[27] "lot5_numero"                  "lot5_surface_carrez"         
[29] "nombre_lots"                  "code_type_local"             
[31] "type_local"                   "surface_reelle_bati"         
[33] "nombre_pieces_principales"    "code_nature_culture"         
[35] "nature_culture"               "code_nature_culture_speciale"
[37] "nature_culture_speciale"      "surface_terrain"             
[39] "longitude"                    "latitude"                    
[41] "geometry"                     "__index_level_0__"

Que contient le champ nature_mutation ? (il a été filtré au ventes classiques pour simplifié cette application ; les vraies données sont plus riches).

A la question 4, vous devriez obtenir des statistiques similaires à celles-ci :

Nous venons de voir comment faire quelques requêtes basiques sur un geoparquet avec duckdb et l’équivalence entre les approches SQL et tidyverse. La dernière question était déjà une introduction au calcul à la volée de statistiques descriptives, ajoutons quelques statistiques avec ce nouvel exercice.

Exercice 2 : statistiques descriptives sur la dimension attributaire

Ne garder que les seules transactions effectuées à Montrouge ou Malakoff et faire une médiane par communes des montants des transactions

Faire ceci avec SQL et dplyr[^chatGPT]

[^chatGPT] : Vous avez le droit d’utiliser chatGPT ou Claude ou vous IA assistante préférée ! Mais ne prenez pas pour argent comptant ce qu’elle vous propose.

Avec l’approche SQL vous devriez obtenir

  code_commune mediane_valeur_fonciere
1        92046                  375000
2        92049                  394500

On peut se rassurer, on obtient la même chose l’approche dplyr :

On peut en conclure que les biens vendus à Montrouge (dans notre base) ont une médiane un peu plus élevée qu’à Malakoff.

Partie 2 : les prix immobiliers à Malakoff, dans le centre et en-dehors.

À présent, nous souhaitons avoir des informations sur les transactions effectuées dans le « fameux » Triangle d’Or de Malakoff (plus prosaïquement, dans son centre-ville commerçant).

Comme il n’est pas possible de distinguer cette zone par requêtes attributaires, nous proposons de :

Via DuckDB, extraire les transactions de l’ensemble de la commune de Malakoff tout en conservant leur caractère spatial (chaque transaction correspond à un point géographique, avec ses coordonnées xy).
Utiliser localement le package sf pour distinguer spatialement les transactions effectuées à l’intérieur ou à l’extérieur du Triangle d’Or (dont nous fournissons les contours).
Calculer la médiane des prix dans les deux sous-zones.

Format des géométries

On extrait les transactions de Malakoff. Pour information, dans le fichier dvf.parquet, les coordonnées spatiales sont stockées dans un format binaire spécifique (Well-Known Binary - WKB). Ce format est efficace pour le stockage et les calculs, mais n’est pas directement lisible ou interprétable par les humains.

En transformant ces géométries en une représentation texte lisible (Well-Known Text - WKT) avec ST_AsText, on rend les données spatiales faciles à afficher, interpréter ou manipuler dans des contextes qui ne supportent pas directement les formats binaires géospatiaux.

Pour le prochain exercice, nous aurons besoin de la structure de requête suivante pour bien interpréter la dimension géographique :

FROM ...
SELECT
  x, y, ..., ST_AsText(geometry) AS geom_text
WHERE ...

Exercice 3

En vous inspirant du template ci-dessus, créer un dataframe transactions_malakoff qui recense les transactions dans cette charmante bourgade.
A ce niveau, les transactions extraites sont maintenant chargées en mémoire et on les transforme dans un format qui facilite leur manipulation en R via le package sf.

transactions_malakoff <- 
  sf::st_as_sf(transactions_malakoff, wkt = "geom_text", crs = 2154) |>
  rename(geometry=geom_text)

Nous allons créer un masque pour reconnaître les transactions qui sont situées ou non dans le triangle d’or. Utiliser la structure suivante pour créer ce masque :

bool_mask <- transactions_malakoff |> 
  # ... |>
  sf::st_intersects(triangle, sparse = FALSE)

⚠️ il faut tenir compte des projections géographiques avant de faire l’opération d’intersection. Ce code est donc à amender à la marge pour pouvoir faire l’intersection.

Cela donne un vecteur de booléen, on peut donc identifier les transactions dans le triangle d’or ou en dehors à partir de celui-ci.

Ci-dessous le dataframe brut extrait via Duckdb (réponse 1).

Ci-dessous, le dataframe transformé en objet sf et prêt pour les opérations spatiales (réponse 2) :

Une fois les données prêtes, on intersecte les points avec le triangle représentant le centre-ville de Malakoff (question 3)

      [,1]
[1,] FALSE
[2,] FALSE
[3,] FALSE
[4,] FALSE
[5,] FALSE
[6,] FALSE

On peut ensuite facilement créer nos deux espaces de Malakoff :

in_triangle <- transactions_malakoff[bool_mask,]
out_triangle <- transactions_malakoff[!bool_mask,]

Une fois que chaque transaction est identifiée comme étant à l’intérieur ou à l’extérieur du Triangle, le calcul de la médiane des prix est immédiat.

median_in <- median(in_triangle$valeur_fonciere)
median_out <- median(out_triangle$valeur_fonciere)

print(glue("Médiane des prix dans le Triangle d'Or de Malakoff : ", median_in))

Médiane des prix dans le Triangle d'Or de Malakoff : 377000

print(glue("Médiane des prix dans le reste de Malakoff : ", median_out))

Médiane des prix dans le reste de Malakoff : 370075

La médiane des prix est un peu plus élevée dans le Triangle qu’en dehors. On peut aller au-delà et étudier la distribution des transactions. Bien que la taille d’échantillon soit réduite, on a ainsi une idée de la diversité des prix dans cette bucolique commune de Malakoff.

Produire la figure sur la distribution du prix des biens

library(ggplot2)
library(scales)

malakoff_identified <- transactions_malakoff %>%
  mutate(
    region = if_else(as.logical(bool_mask), "Triangle d'or", "Hors triangle d'or")
  ) 

ggplot(
  malakoff_identified,
  aes(y = valeur_fonciere, x = region, fill = region)
) +
  geom_violin() +
  scale_y_continuous(
    trans = "log10",
    labels = comma_format(),
    breaks = scales::trans_breaks("log10", function(x) 10^x)
  ) +
  geom_jitter(height = 0, width = 0.1) +
  labs(y = "Valeur de vente (€)") +
  theme_minimal()

Tout ceci ne nous dit rien de la différence entre les biens dans le triangle et en dehors de celui-ci. Nous n’avons fait aucun contrôle sur les caractéristiques des biens. Nous laissons les curieux explorer la mine d’or qu’est cette base.

Partie 3 : Part de ménages pauvres à Malakoff et à Montrouge

Pour finir, on se place dans le cas où :

On souhaite extraire des informations d’un fichier volumineux (les données carroyées de l’Insee).
Mais il n’est pas possible de filtrer les données par des requêtes attributaires (par exemple, il n’est pas possible de faire code_commune = 92049).

Ainsi, nous allons :

Utiliser les contours géographiques des deux communes
Filtrer les données par intersections géographiques des carreaux et des communes, à l’aide de DuckDB
Faire les calculs localement après l’extraction des carreaux d’intérêt.

Pour commencer, on décrit les données carroyées comme précédemment :

describe_dvf <- dbGetQuery(con, "DESCRIBE SELECT * FROM read_parquet('data/carreaux.parquet')")
describe_dvf

preview <- dbGetQuery(con, "SELECT * FROM read_parquet('data/carreaux.parquet') LIMIT 10")
preview

On va faire une petite transformation de données préliminaire à cet exercice afin que la géométrie de Malakoff soit reconnue par DuckDB.

malakoff_2154 <- sf::st_transform(malakoff, 2154)
malakoff_wkt <- sf::st_as_text(sf::st_geometry(malakoff_2154))

Voici comment faire une requête géographique sur les carreaux de Malakoff

geo_query <- glue("
  FROM read_parquet('data/carreaux.parquet')
  SELECT
      *, ST_AsText(geometry) AS geom_text
  WHERE ST_Intersects(
      geometry,
      ST_GeomFromText('{malakoff_wkt}')
  )
")

carr_malakoff <- dbGetQuery(con, geo_query)

carr_malakoff <-
  carr_malakoff |>
  sf::st_as_sf(wkt = "geom_text", crs = 2154) |>
  select(-geometry) |>
  rename(geometry=geom_text)

On peut les visualiser de la manière suivante

mapview(carr_malakoff) + mapview(sf::st_boundary(malakoff))

Exercice 4 : extraction des carreaux intersectant Malakoff

Réitèrer l’opération pour Montrouge
Calculer la proportion moyenne de ménages pauvre dans l’ensemble des carreaux extraits à partir des deux objets obtenus.

Le masque des carreaux de Montrouge est le suivant :

On obtient, in fine, les statistiques suivantes

Part de ménages pauvres dans les carreaux de Malakoff : 12.01

Part de ménages pauvres dans les carreaux de Montrouge : 11.02

Exercice optionnel

Calculer la même statistique dans et hors du triangle d’or de Malakoff

Conclusion

Nous avons donc réussi à lire des données avec DuckDB et à faire des statistiques dessus. Pourquoi est-ce pertinent de passer par DuckDB ? Car ce package permet de faire ceci de manière très efficace sur de gros volumes de données. Il passe très bien à l’échelle.

A noter que notre démarche est une introduction à ce sujet bien plus large qu’est l’analyse géographique. Notre approche serait améliorable sur plusieurs plans :

rationalisation des requêtes,
pertinence statistique des résultats
réplicabilité du code