[R] Combiner R et PostGIS, quelle bonne idée !

R_postgis

Pourquoi Combiner R et PostGIS ?

Combiner R et PostGIS permet d’utiliser dans un environnement de travail intégré regroupant à la fois la puissance de R (gestion des données, traitements statistiques, …) et la puissance de PostGIS (géotraitements, requêtes spatiales, …) sans rupture de charge. Cette alliance permet notamment de : traiter, analyser et cartographier les données (R) ; bénéficier des dernières évolutions de PostGIS 2 (gestion topologique, indexation, géotraitements, …) ; structurer les données de façon efficace (PostGIS) ; favoriser le recherche reproductible (R & PostGIS). Et tout cela avec les avantages du libre et du gratuit (nombreuses mises à jour, grosse communauté d’utilisateurs). Bref,  Cette combinaison permet d’offrir un environnement de développement complet et puissant adapté à la très grande majorité des traitements utiles en géographie, en statistique et géomatique. Que demander de plus…

Prérequis

Pour commencer, il faut faire quelques installations :

installer R (http://www.r-project.org/)
installer et configurer PostGIS (http://postgis.net/install/ & http://postgis.net/docs/postgis_installation.html)

Comment transférer les données de R à PostGIS et inversement ?

Théoriquement, il est possible d’importer et d’exporter des spatialdataframes dans PostGis  avec le package rgdal. Néanmoins, après de nombreux essais infructueux,  une solution basée sur l’import/export des géométries au format WKT s’avère nettement  plus efficace. La meilleure solution trouvée sur le web est proposée sur le site de Philipp Hunziker. Sur ce site, l’auteur propose deux fonctions permettant d’importer ou d’exporter les géométries et les tables attributaires.  Les fonctions proposées sur ce site ne gérant ni les projections ni l’import de polygones uniques, je propose ci-dessous une version un petit peu mise mise à jour.

download

Une fois le script Rpostgis.R téléchargé et sourcé, la fonction dbWriteSpatial() permet de transférer un spatialdataframe vers PostGIS.
La fonction dbReadSpatial() permet de de transférer un une table PostGIS vers un spatialdataframe

Exemple d’utilisation :

# Paramètres de connexion
con <- dbConnect(dbDriver("PostgreSQL"), dbname="test")

# Shapfile -> R -> PostGIS
df <- readOGR(dsn ="geom/",layer = "polys")
dbWriteSpatial(con, df, schemaname="public", tablename="newtable", replace=T,srid=3035)

# PostGIS -> R
spdf <- dbReadSpatial(con, schemaname="public", tablename="newtable", geomcol="geom")
proj4string(spdf)
head(spdf@data)
plot(spdf)

Exécuter une requête PostGIS dans R

L’exécution de requêtes PostGIS dans R s’effectue avec le package RPostgreSQL

# Appel du package RPostgreSQL
require(RPostgreSQL)
# Définition des paramètres de connexion
con <- dbConnect(dbDriver("PostgreSQL"), dbname="test")
# Requête (e.g zone tampon)
q<-"CREATE TABLE nom-de-ma-nouvelle-table-postgis- AS SELECT id,ST_Buffer(geom,200000) AS geom FROM nom-de-ma-table-postgis"
# Envoi de la requête vers PostGIS
res<-dbSendQuery(con,q)
# Libération des ressources pour effectuer éventuellement une autre requête
dbClearResult(res)

Un exemple de traitement simple (e.g. Polygons -> polylines)

img1 

source("sources/Rpostgis.R")

# 1. Connexion
con <- dbConnect(dbDriver("PostgreSQL"), dbname="test")

# 2. Transfert d'une couche sur les pays européens dans PostGIS
df <- readOGR(dsn ="data/",layer = "NUTS0_RG_20M_2006")
dbWriteSpatial(con, df, schemaname="public", tablename="pays", replace=T,srid=3035)

# 3. Exécution d'un requête (Polygones -> polylignes)
q<-"CREATE TABLE bound AS SELECT id,ST_Boundary(geom) AS geom FROM pays"
res<-dbSendQuery(con,q)
dbClearResult(res)

# 4. Récupération du résultat dans R
df2 <- dbReadSpatial(con, schemaname="public", tablename="bound", geomcol="geom")

# 5.  Affichage
par(mfrow=c(1,2))
plot(df,col="red")
plot(df2,col="red")
par(mfrow=c(1,1))

# 6. Nettoyage de la base PostGIS
res1<-dbSendQuery(con,"SELECT DropGeometryTable('bound')")
dbClearResult(res1)
res2<-dbSendQuery(con,"SELECT DropGeometryTable('pays')")
dbClearResult(res2)

process 

Conclusion

Si l’exemple montré ici est simple (il aurait pu être effectué intégralement dans R sans PostGIS), il montre néanmoins un aperçu des possibilités offertes par l’articulation R<->PostGIS ainsi que la logique d’une telle utilisation.  Dans un prochain billet, je montrerai comment utiliser cet environnement de développement pour réaliser une généralisation cartographique qui conserve la topologie.


Nicolas Lambert

Ingénieur de recherche CNRS en sciences de l'information géographique. Membre de l'UMS RIATE et du réseau MIGREUROP / CNRS research engineer in geographical information sciences. Member of UMS RIATE and the MIGREUROP network.

Vous aimerez aussi...

1 réponse

  1. Rapenne dit :

    Bonjour,

    Merci beaucoup pour ce code.
    J’ai toutefois du faire des modifications aussi bien dans dbWriteSpatial que dans dbReadSpatial pour que ces fonctions fonctionnent.
    Dans dbWriteSpatial, j’ai modifié la ligne 43 : query5 <- paste("SELECT UpdateGeometrySRID('",tablename,"','geom',",srid,")",sep="")
    par query5 <- paste("SELECT UpdateGeometrySRID('",schemaname,"','",tablename,"','geom',",srid,")",sep="")

    Dans dbReadSpatial j'ai modifié schema.table = paste(schemaname, ".", tablename, sep="") en schema.table = paste0("\"",schemaname, "\".\"", tablename,"\"") ainsi que query1 <- paste("SELECT ST_SRID(geom) FROM",tablename,"LIMIT 1")
    en query1 <- paste("SELECT ST_SRID(\"",tablename,"\".geom) FROM \"",schemaname,"\".\"",tablename,"\" LIMIT 1",sep="")

    Cordialement

    Serge

Laisser un commentaire

Votre adresse e-mail ne sera pas publiée. Les champs obligatoires sont indiqués avec *

Ce site utilise Akismet pour réduire les indésirables. En savoir plus sur comment les données de vos commentaires sont utilisées.

Rechercher dans OpenEdition Search

Vous allez être redirigé vers OpenEdition Search