3ème séance du cycle de formation à l’analyse de réseaux avec R

par Ryma Hachi · Publié 09/01/2018 · Mis à jour 01/09/2020

La troisième séance du cycle de formation organisé par fmr sur l’analyse de réseaux avec R a eu lieu le 19 décembre 2017, à l’Institut de géographie de Paris. Elle était animée par Hadrien Commenges, maître de conférence à l’université Paris 1. Le contenu de la séance est détaillé ci-dessous. Les données utilisées sont téléchargeables ici.

La séance avait pour objectifs de:

Savoir construire un graphe à partir d’une liste de liens
Connaître les principes de la détection de communautés
Savoir choisir parmi les différents algorithmes de détection
Être capable de concevoir une chaîne de traitement efficace

Présentation des données

Le fichier citibike042017.Rds contient tous les trajets effectués sur le réseau Citibike durant le mois d’avril 2017. Les données brutes sont téléchargeables sur le site http://www.citibikenyc.com/system-data. Elles sont enrichies par :

des variables calculées : distance euclidienne, distance de Manhattan, vitesse (calculée sur la distance de Manhattan)
une variable météorologique : pluie, variable binaire 0-non, 1-oui

Traitements préalables

Charger les packages

library(sf)
library(cartography)
library(igraph)
library(dplyr)

Charger les données

citiBike <- readRDS("DATA/citibike042017.Rds")

Créer un fichier de stations : tableau de points (data.frame) et objet spatial (sf)

# vérifier le nombre total de stations
citiBike$start_station_id %>% unique() %>% length()

## [1] 609

citiBike$end_station_id %>% unique() %>% length()

## [1] 608

c(citiBike$start_station_id, citiBike$end_station_id) %>% unique() %>% length()

## [1] 609

# créer un index de valeurs dupliquées
idUnic <- !duplicated(citiBike$start_station_id)

# filtrer et mettre en forme le tableau de points
tabStations <- citiBike[idUnic, c(1, 3, 4)] %>% 
  rename(id = start_station_id,
         lat = start_station_latitude,
         lon = start_station_longitude) %>% 
  arrange(id)

# transformer le tableau en objet spatial
spStations <- st_as_sf(tabStations,
                       coords = c("lon", "lat"),
                       agr = "constant",
                       crs = 4326, 
                       stringsAsFactors = FALSE)

Cartographier les stations

osmTiles <- getTiles(x = spStations, type = "osm", zoom = 14, crop = FALSE)

## Loading required namespace: raster

tilesLayer(osmTiles)
plot(spStations$geometry, pch = 20, add = TRUE)

Pour l’exploration utiliser leaflet

library(leaflet) leaflet() %>% addTiles() %>% addCircleMarkers(data = spStations)

Construction du graphe

Le tableau dont on dispose (citiBike) peut être considéré comme une liste de liens : chaque station représente un noeud (sommet) du graphe, chaque trajet reliant deux stations représente un arc (arrête). En l’état il s’agit donc d’une liste de liens multiples non pondérés : entre deux stations il y a un nombre n de trajets. Pour pouvoir analyser ce réseau il sera utile d’agréger ces liens multiples : si 10 trajets relient la station A à la station B, on ne considère pas qu’il y a dix arcs non valués mais un seul arc de poids 10.

Pour faire l’agrégation il est possible de créer le graphe non valué qui aura 609 noeuds et 12,4 millions d’arcs puis de collapser les arcs multiple avec la fonction simplify() et l’argument remove.multiple. Cette option est assez mauvaise parce qu’elle demande une opération coûteuse de création d’un graphe à 12,4 millions d’arcs. Il est donc préférable de faire l’agrégation en amont, avec aggregate() ou avec les fonctions du packages dplyr.

tabAgr <- citiBike %>% 
  group_by(start_station_id, end_station_id) %>% 
  summarise(NB = n())

Explorons ce réseau.

summary(tabAgr$NB)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   1.000   1.000   3.000   9.111  10.000 761.000

nrow(tabAgr) / (609 * 608)

## [1] 0.3686587

Il y a donc 136504 arcs pondérés reliant les 609 stations. C’est un réseau très dense : le rapport entre nombre d’arcs observés et nombre d’arcs maximal est de 0,37. Les poids varient entre 1 et 761, avec une distributions très concentrée à gauche : min à 1, Q1 à 1, médiane à 3, Q3 à 10.

Ce tableau peut maintenant être transformé en objet de type igraph.

citiNet <- graph.data.frame(tabAgr, directed = TRUE, vertices = tabStations)

On peut explorer l’objet, mais éviter de le représenter : avec 136504 arcs ce serait long et le résultat serait illisible.

summary(citiNet)

## IGRAPH DN-- 609 136504 -- 
## + attr: name (v/c), lat (v/n), lon (v/n), NB (e/n)

components(citiNet)

## $membership
##   72   79   82   83  116  119  120  127  128  143  144  146  150  151  152 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  153  157  161  164  167  168  173  174  195  212  216  217  223  225  228 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  229  232  236  237  238  239  241  242  243  244  245  247  248  249  251 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  252  253  254  257  258  259  260  261  262  264  265  266  267  268  270 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  274  275  276  278  279  280  281  282  284  289  291  293  295  296  297 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  301  302  303  304  305  306  307  308  309  310  311  312  313  314  315 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  316  317  319  320  321  322  323  324  325  326  327  328  330  331  332 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  334  335  336  337  339  340  341  342  343  344  345  346  347  348  349 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  350  351  353  354  355  356  357  358  359  360  361  362  363  364  365 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  366  368  369  372  373  376  377  379  380  382  383  384  385  386  387 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  388  389  390  391  392  393  394  395  396  397  398  399  400  401  402 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  403  405  406  407  408  409  410  411  412  414  415  416  417  418  419 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  420  421  422  423  426  427  428  430  432  433  434  435  436  437  438 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  439  440  441  442  443  444  445  446  447  448  449  450  453  454  455 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  456  457  458  459  460  461  462  465  466  467  468  469  470  471  472 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  473  474  476  477  478  479  480  481  482  483  484  485  486  487  488 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  490  491  492  494  495  496  497  498  499  500  501  502  503  504  505 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  507  508  509  511  513  514  515  516  517  518  519  520  522  523  524 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  525  526  527  528  529  530  531  532  533  534  536  537  539  540  545 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
##  546 2000 2001 2002 2003 2005 2006 2008 2009 2010 2012 2021 2022 2023 3002 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3016 3041 3042 3043 3044 3046 3047 3048 3049 3050 3052 3053 3054 3055 3056 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3057 3058 3059 3060 3061 3062 3063 3064 3065 3066 3067 3068 3069 3070 3071 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3072 3073 3074 3075 3076 3077 3078 3080 3081 3082 3083 3085 3086 3087 3088 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3090 3091 3092 3093 3094 3095 3096 3100 3101 3102 3103 3105 3106 3107 3108 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3109 3110 3111 3112 3113 3114 3115 3116 3117 3118 3119 3120 3121 3122 3123 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3124 3125 3126 3127 3128 3129 3130 3131 3132 3134 3135 3136 3137 3139 3140 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3141 3142 3143 3144 3145 3146 3147 3148 3150 3151 3152 3155 3156 3157 3158 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3159 3160 3161 3162 3163 3164 3165 3166 3167 3168 3169 3170 3171 3172 3173 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3175 3176 3177 3178 3179 3180 3182 3221 3223 3224 3226 3231 3232 3233 3235 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3236 3238 3241 3242 3243 3244 3245 3249 3254 3255 3256 3258 3259 3260 3263 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3282 3283 3284 3285 3286 3288 3289 3290 3292 3293 3294 3295 3297 3298 3300 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3301 3302 3303 3304 3305 3306 3307 3308 3309 3310 3311 3312 3313 3314 3315 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3316 3317 3318 3319 3320 3321 3322 3323 3324 3325 3326 3327 3328 3329 3330 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3331 3332 3333 3335 3336 3337 3338 3339 3340 3341 3342 3343 3344 3345 3346 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3347 3348 3349 3350 3351 3352 3353 3354 3355 3356 3357 3358 3359 3360 3361 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3362 3363 3364 3365 3366 3367 3368 3370 3371 3373 3374 3375 3376 3377 3378 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3379 3381 3382 3383 3384 3386 3387 3388 3389 3390 3391 3392 3393 3394 3395 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3396 3397 3398 3399 3400 3401 3402 3403 3404 3405 3407 3408 3409 3410 3411 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3412 3413 3414 3415 3416 3417 3418 3419 3420 3421 3422 3423 3424 3425 3427 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3428 3429 3430 3431 3434 3435 3436 3438 3440 3441 3443 3445 3447 3449 3450 
##    1    1    1    1    1    1    1    1    1    1    1    1    1    1    1 
## 3452 3453 3454 3455 3456 3457 3459 3462 3463 
##    1    1    1    1    1    1    1    1    1 
## 
## $csize
## [1] 609
## 
## $no
## [1] 1

diameter(citiNet)

## [1] 4

Ces mesures confirment le diagnostic précédent, une seule composante connexe et un diamètre égal à 4, le réseau est très dense. Quelques mesures sur les degrés des noeuds permettent de se faire une idée plus précise de la structure des flux.

degIn <- strength(citiNet, mode = "in", weights = E(citiNet)$NB)
degOut <- strength(citiNet, mode = "out", weights = E(citiNet)$NB)
hist(degIn, col = "firebrick", border = "white")

hist(degOut, col = "firebrick", border = "white")

plot(degIn, degOut, pch = 20)
abline(a = 0, b = 1, col = "red")

Détection de communautés

Sur un graphe simple, tous les algorithmes de détection de communautés donneront un résultat correct et rapide. Voici un exemple.

g <- make_full_graph(5) %du% make_full_graph(5) %du% make_full_graph(5)
g <- add_edges(g, c(1,6, 1,11, 6, 11))

Avec l’algorithme basé sur l’intermédiarité des liens :

betweenClus <- cluster_edge_betweenness(g, weights = NULL, directed = FALSE)
plot(g, vertex.color = betweenClus$membership)

Avec l’algorithme fast greedy :

greedyClus <- cluster_fast_greedy(g, weights = NULL)
plot(g, vertex.color = greedyClus$membership)

Avec l’algorithme de Louvain :

louvainClus <- cluster_louvain(g, weights = NULL)
plot(g, vertex.color = louvainClus$membership)

La notion principale mobilisée dans la détection de communauté est la modularité (modularity), elle s’écrit :

où est la demi-somme des poids, est le poids de l’arc reliant les noeuds et , et k sont les sommes marginales des noeuds et (i.e. leur degré entrant et sortant) et la fonction δ une fonction binaire qui prend la valeur 1 si et sont dans la même communauté (et 0 sinon).

La modularité peut être calculée a posteriori sur n’importe quel objet de type communities pour caractériser le partitionnement :

modularity(betweenClus)

## [1] 0.5757576

modularity(greedyClus)

## [1] 0.5757576

modularity(louvainClus)

## [1] 0.5757576

Dans un grand nombre d’algorithmes, la modularité n’est pas utilisé comme mesure a posteriori de la qualité du partionnement mais comme quantité à maximiser. C’est ce que font les algorithmes suivants :

fast greedy avec cluster_fast_greedy()
Louvain avec cluster_louvain()
optimum de modularité avec cluster_optimal()
principale valeur propre avec cluster_leading_eigen()

Application aux données citibike

Parmi les algorithmes présentés (il existe d’autres algorithmes qui ne s’appuient pas sur la modularité et qui n’ont pas été présentés) :

celui de l’intermédiarité des liens est clairement hors de propos : il est itératif, donc lent, et surtout la notion d’intermédiarité a peu de sens dans un graphe aussi dense (voir plus haut)
celui de la valeur propre est rarement applicable parce qu’il crée seulement deux clusters
celui de l’optimum de modularité donnerait le meilleur résultat mais il est inapplicable sur un graphe d’une telle taille
les deux autres – fast greedy et Louvain – sont applicables et donneront des résultats semblables, il donnent un approximation rapide de l’optimum de modularité

Attention la plupart des algorithmes demandent un graphe non orienté. On commence donc par transformer le réseau : on considère que s’il y a un lien bidirectionnel entre une station i et une station j il convient de sommer les deux sens au moment de l’agrégation. Bien sûr d’autres fonctions pourraient être appliquées.

citiNetUn <- as.undirected(graph = citiNet, mode = "collapse", edge.attr.comb = sum)
summary(citiNetUn)

## IGRAPH UN-- 609 84105 -- 
## + attr: name (v/c), lat (v/n), lon (v/n), NB (e/n)

On applique maintenant l’algorithme de Louvain et on récupère l’appartenance des noeuds aux communautés. L’intérêt de cette méthode est de produire plusieurs niveaux d’agrégation (multilevel) si le réseau est complexe. Dans ce cas, le dernier niveau est donné par le vecteur membership, l’ensemble des niveaux est donné dans la matrice memberships. Dans l’exemple de Citibike le réseau est simple et un seul niveau d’agrégation est produit, la matrice memberships ne contient qu’une seule ligne :

multCom <- cluster_louvain(graph = citiNetUn, weights = E(citiNetUn)$NB)
vecMember <- multCom$membership
matMember <- multCom$memberships

Le réseau reste difficile à représenter avec ses 84105 arcs, on sélectionne une centaine d’arcs parmi les flux principaux en supprimant tous les arcs dont le poids inférieur à 300.

mainNet <- delete_edges(graph = citiNetUn, edges = E(citiNetUn)[E(citiNetUn)$NB < 300])
summary(mainNet)

## IGRAPH UN-- 609 124 -- 
## + attr: name (v/c), lat (v/n), lon (v/n), NB (e/n)

On cherche à afficher les noeuds aux coordonnées géographiques correspondantes aux stations. Il faut donner dans la fonction plot() une matrice de coordonnées : les noeuds du graphe (qui contiennent comme attribut la longitude et la lattitude), on fabrique une matrice contenant long et lat (dans cet ordre : la longitude est l’axe des x et la latitude l’axe des y), on désigne cette matrice avec l’argument layout =.

tabCoord <- get.data.frame(citiNetUn, what = "vertices")
matCoord <- as.matrix(tabCoord[, c(3, 2)])

plot(mainNet, 
     vertex.color = multCom$membership, 
     vertex.label = NA, 
     vertex.frame.color = NA,
     vertex.size = 2,
     edge.color = "firebrick",
     edge.width = 2,
     layout = matCoord)

Si on veut cartographier le résultat il faut extraire un tableau, éventuellement l’enrichir de mesures calculées sur le graphe (ici les clusters et les degrés des noeuds), puis créer un objet spatial.

tabCoord$CLUS <- multCom$membership
tabCoord$DEGIN <- degIn
tabCoord$DEGOUT <- degOut

spClus <- st_as_sf(tabCoord,
                   coords = c("lon", "lat"),
                   agr = "constant",
                   crs = 4326, 
                   stringsAsFactors = FALSE)

Cet objet peut être représenté :

avec la fonction plot()
avec les fonctions du packagecartography
avec le packageleaflet si on veut explorer et comprendre pourquoi certaines stations du Sud de Manhattan font partie du même cluster que Brooklyn.

Option 1

plot(spClus$geometry)

Option 2

tilesLayer(osmTiles)
propSymbolsLayer(x = spClus, var = "DEGIN", inches = 0.1, col = "black", border = "white")

## Warning in st_centroid.sfc(st_geometry(x), of_largest_polygon =
## of_largest_polygon): st_centroid does not give correct centroids for
## longitude/latitude data

Option 3

Ici on ajoute des tuiles CartoDB (voir la liste des tuiles disponibles sur ce site)

colPal <- c("firebrick", "grey", "springgreen", "blue")
leaflet() %>% 
  addProviderTiles(provider = "CartoDB.DarkMatter") %>% 
  addCircleMarkers(data = spClus, 
                   stroke = FALSE, 
                   fillColor = colPal[spClus$CLUS],
                   fillOpacity = 0.6)

Poursuite de l’analyse

Créer des sous-graphes en filtrant le tableau des déplacements en amont, faire par exemple :

le graphe des déplacements en jour de pluie
le graphe des déplacements des abonnés vs. utilisateurs occasionnels
le graphe des déplacements en heure de pointe du matin

Ryma Hachi

Ryma Hachi est doctorante à l’UMR Géographie-cités depuis octobre 2014, sous la direction de Lena Sanders. Architecte de formation, sa thèse vise à révéler, grâce à la théorie des graphes, les processus de formation et de transformation des réseaux viaires.

OpenEdition vous propose de citer ce billet de la manière suivante :
Ryma Hachi (9 janvier 2018). 3ème séance du cycle de formation à l’analyse de réseaux avec R. groupe fmr (flux, matrices, réseaux). Consulté le 19 septembre 2024 à l’adresse https://doi.org/10.58079/p85d

3ème séance du cycle de formation à l’analyse de réseaux avec R

Présentation des données

Traitements préalables

Construction du graphe

Détection de communautés

Application aux données citibike

Poursuite de l’analyse

Ryma Hachi

Vous aimerez aussi...

Laisser un commentaire Annuler la réponse

Catégories

Articles récents

Groupe Res-Hist – Réseaux et histoire

QUANTI / Sciences sociales

Liens

Catégories

Billets récents

GDR Analyse de réseaux en SHS

RES-HIST

QUANTI / sciences sociales

Liens

3ème séance du cycle de formation à l’analyse de réseaux avec R

Présentation des données

Traitements préalables

Construction du graphe

Détection de communautés

Application aux données citibike

Poursuite de l’analyse

Vous aimerez aussi...

Laisser un commentaire Annuler la réponse

Catégories

Articles récents

Liens

Catégories

Billets récents

Mots-clef

Liens