Réalisez des cartes administratives avec GADMTools

Parfois, les hasards font vraiment bien les choses ! Il y a quelques mois, Jean Pierre Decorps, développeur de package R dans le domaine de l’épidémiologie, chez Epiconcept, m’a contacté pour m’informer de la publication prochaine de la nouvelle version de son package GAMTools. Il s’agit d’un package permettant de manipuler, d’assembler et d’exporter des cartes à partir des fichiers shapefile du projet GADM (dont le but est de réaliser des cartes administratives de tous les pays du monde, à différents niveau de sous division). Lorsque j’ai reçu ce mail, je n’avais jamais fait de cartographie avec R, alors ça m’a paru un peu obscure.

Et puis quelques semaines plus tard, un premier client m’a demandé de réaliser une carte sous R. J’avais alors utilisé leaflet, mais surtout ça m’a permis de me familiariser avec les fichiers shapefile (le fichier qui contiennent les coordonnées des contours). Et puis encore quelques semaines plus tard, dans le cadre d’une prestation en épidémiologie, un autre client m’a demandé de l’aider à réaliser une carte des taux d’incidence d’une infection virale, et en m’orientant vers le package GADMTools ! J’ai donc exploré ce package, et je l’ai trouvé vraiment très intéressant et très pratique. Et comme vous le verrez dans la troisième partie de cet article, vous pouvez l‘utiliser dans d’autres domaines que l’épidémiologie !

Les différents niveaux de cartes administratives

Comme expliqué en introduction, le package GADMTools permet de réaliser des cartes administratives de n’importe quel pays ou région du monde très facilement. La première étape consiste à charger le fichier shapefile du pays désiré, en utilisant la fonction `gadm_sf_loadCountries()`. Le niveau de sous-division souhaité (région, département etc…) peut être choisi via l’argument `level`.

Les codes des pays et les différents niveaux de découpage disponibles pour chacun d’eux sont décrits dans ce document .

Cartes administratives globales (level=0)

Par exemple, pour la France :

install.packages("GADMTools")
library(GADMTools)
carte_globale <- gadm_sf_loadCountries("FRA", level=0 )
gadm_plot(carte_globale) 

Le fichier shapefile (le fichier qui contient les coordonnées des contours) de la France (code “FRA”) est chargé dans le working directory, par la fonction `gadm_sf_loadCountries()`.

L’argument level=0 correspond au pays sélectionné dans son ensemble, c’est à dire sans découpage supplémentaire (par région par exemple).

La fonction `gadm_plot()` permet simplement plotter le fichier shapefile.

Cartes administratives des régions (level=1)

Pour obtenir la carte administrative de  la  France avec un découpage par région, il suffit d’indiquer `level=1` :

carte_regions<- gadm_sf_loadCountries("FRA", level=1 )
gadm_plot(carte_regions) 

Très pratique pour faire réviser les régions de France aux enfants !

Remarque : A la place de la fonction gadm_plot(), on peut utiliser la fonction `plot_map()` :

plotmap(carte_regions) 

Cartes administratives des départements (level=2)

Pour obtenir la carte administrative de la France avec un découpage par département, il suffit d’indiquer `level=2` :

carte_dep <- gadm_sf_loadCountries("FRA", level=2)
gadm_plot(carte_dep) 

Le package est tidyverse compatible, si bien qu’on peut utiliser aussi cette syntaxe:

gadm_sf_loadCountries("FRA", level=2 ) %>% 
gadm_plot() 

On peut ajouter :

– un titre, en utilisant l’argument `title` au sein de la fonction `gadm_plot()`,
– une flèche nord avec la fonction `gadm_showNorth()`,
– une échelle avec la fonction `gadm_showScale()`

gadm_sf_loadCountries("FRA", level=2 ) %>% 
gadm_plot(title="Carte des départements français") %>% 
gadm_showNorth("tl") %>% 
gadm_showScale('bl')
 

Et on peut aussi ajouter des points. Pour cela, on crée un data.frame avec les latitudes et longitudes où l’on souhaite placer les points. Par exemple, ici l’Isère et les Alpes de Haute Provence :

longitude <- c(5.93, 6.24)
latitude <- c(44.996, 44.08) 
mydf <- data.frame(longitude,latitude) 

J’ai trouvé les coordonnées gps ici.

France <- gadm_sf_loadCountries("FRA", level=2 ) 
dots(France, mydf, color="red", size = 3) 

On peut également réaliser une carte administrative de plusieurs pays  :

gadm_sf_loadCountries(c("FRA", "ITA"), level=1 ) %>% 
gadm_plot(title="Carte des régions françaises et italiennes") %>% 
gadm_showNorth("tl") %>% 
gadm_showScale('bl') 

Dessiner des portions de cartes administratives

Le package GADMTools permet de faire des subsets de cartes administratives, en utilisant la fonction `gadm_subset()`. Par exemple pour représenter les départements de la région PACA :

PACA <- gadm_sf_loadCountries("FRA", level = 2) %>% 
        gadm_subset(level=1, regions="Provence-Alpes-Côte d'Azur")  
plotmap(PACA) 

Utiliser `level=2` lors du chargement du fichier shapefile, permet d’obtenir un fichier comportant les contours des départements. En revanche, on utilise `level=1` dans la fonction `gadm_subset()` afin de faire un découpage au niveau de la région.

Les labels (étiquettes) des régions, employées par le package GADMTools peuvent être obtenues avec cette syntaxe :

 FRA <- gadm_sf_loadCountries("FRA", level=1 )
    listNames(FRA, level=1)
    ##  [1] "Auvergne-Rhône-Alpes"       "Bourgogne-Franche-Comté"   
    ##  [3] "Bretagne"                   "Centre-Val de Loire"       
    ##  [5] "Corse"                      "Grand Est"                 
    ##  [7] "Hauts-de-France"            "Île-de-France"             
    ##  [9] "Normandie"                  "Nouvelle-Aquitaine"        
    ## [11] "Occitanie"                  "Pays de la Loire"          
    ## [13] "Provence-Alpes-Côte d'Azur" 

Pour obtenir les labels de départements :

FRA <- gadm_sf_loadCountries("FRA", level=2 )
head(listNames(FRA, level=2))
## [1] "Ain"         "Allier"      "Ardèche"     "Cantal"      "Drôme"      
## [6] "Haute-Loire 

Le package permet également d’ajouter un fond de carte openstreet map. Pour plus d’info, consultez la vignette du package.

Réaliser des cartes choroplèthes

Une carte choroplèthe, c’est “une carte thématique où les régions sont colorées ou remplies d’un motif qui montre une mesure statistique, tels la densité de population ou le revenu par habitant. Ce type de carte facilite la comparaison d’une mesure statistique d’une région à l’autre ou montre la variabilité de celle-ci pour une région donnée.” (cf. Wikipedia).

Récemment, sur le ton de la plaisanterie, j’ai dit à des proches que j’allais m’acheter un Bubbly ou une Dax :

Du coup, je me suis dit que ça pourrait être intéressant de représenter, pour chaque département, le nombre de vols de deux roues rapporté à la taille de sa population (un peu comme un taux d’incidence)  au cours du dernier mois (enfin le dernier mois dont les données sont disponibles).

Les données de vols de 2 roues en France

Acquisition des données de vols de 2 roues

Le fichier Excel importé, nommé “tableaux-4001-ts” contient les chiffes d’une centaines de délits recensés entre janvier 2000 et mars 2019. La première feuille contient les données pour l’ensemble de la France, la seconde feuille contient les données de la France métropolitaine, et les feuilles suivantes contiennent les données pour chaque département, individuellement :

Importation des données de délits dans R

Le code ci dessous, permet de lire toutes les feuilles du fichier Excel, et de compiler ensemble les données de tous les départements (en ajoutant le numéro de chaque département). Pour plus détail sur ce code, consulter mon article “Importer automatiquement les feuilles d’un fichier Excel dans R“.

library(readxl)

sheet_names <- excel_sheets(here::here("data","tableaux-4001-ts.xlsx"))
delits_fr <- NULL

for (i in 3:98){
      name <- sheet_names[i]
      data <- read_excel(here::here("data","tableaux-4001-ts.xlsx"), sheet = i)
      data$num_dep=name
      delits_fr<- bind_rows(delits_fr,data)        
    } 

Sélection des données

Il s’agit de conserver :

– les données du dernier mois disponible (ici mars 2019)
– les données concernant les vols de deux roues (correspond à l’Index 36)
– la variable contenant le numéro du département

 vol2roues <- delits_fr %>% 
        filter(Index==36) %>% 
        select("libellé index", "2019_03","num_dep") 

Puis, j’ai renommé la a variable “2019_03” en “nb_vols”

vol2roues <- vol2roues %>% 
        rename(nb_vols='2019_03')

head(vol2roues)
    ## # A tibble: 6 x 3
    ##   `libellé index`                       nb_vols num_dep
    ##                                         
    ## 1 Vols de véhicules motorisés à 2 roues      27 01     
    ## 2 Vols de véhicules motorisés à 2 roues      20 02     
    ## 3 Vols de véhicules motorisés à 2 roues      11 03     
    ## 4 Vols de véhicules motorisés à 2 roues       6 04     
    ## 5 Vols de véhicules motorisés à 2 roues       3 05     
    ## 6 Vols de véhicules motorisés à 2 roues     137 06 

Les données de population

Acquisition des données de population

J’ai simplement copier-coller les colonnes contenant :

  • le numéro du département
  • le nom du département
  • la taille de la population totale

dans un fichier csv, que j’ai nommé “pop”.

Importation des données de population

pop <- read.csv2(here::here("data", "pop.csv"))
str(pop)
    ## 'data.frame':    96 obs. of  3 variables:
    ##  $ num       : Factor w/ 96 levels "01","02","03",..: 1 2 3 4 5 6 7 8 9 10 ...
    ##  $ nom       : Factor w/ 96 levels "Ain","Aisne",..: 1 2 3 4 43 5 6 7 8 9 ...
    ##  $ pop_totale: Factor w/ 96 levels "1 036 153","1 075 649",..: 81 69 53 27 22 4 50 44 25 49 ... 

La variable “pop_totale” est considérée comme un facteur. Pour la convertir en variable numérique, il faut :

  • la passer en chaîne de caractères
  • supprimer les espaces qui sépares les centaines des milliers
  • la passer en format numeric
pop$pop_totale <- as.character(pop$pop_totale)

library(stringr)
pop$pop_totale <- str_replace_all(pop$pop_totale,pattern=" ",repl="")
pop$pop_totale <- as.numeric(pop$pop_totale) 

Assemblage des données de vols et de population

Pour réaliser une carte choroplèthe, mon fichier contenant des données numériques, doit contenir une variable avec le nom des départements, tels qu’ils existent dans les fichiers shapefile importées par le package GADMtools.

Comme expliqué précédemment, ces noms sont accessibles via :

FRA <- gadm_sf_loadCountries("FRA", level=2 )

listNames(FRA, level=2)
    ##  [1] "Ain"                     "Allier"                 
    ##  [3] "Ardèche"                 "Cantal"                 
    ##  [5] "Drôme"                   "Haute-Loire"            
    ##  [7] "Haute-Savoie"            "Isère"     

Puisque ces noms de départements sont dans le fichier de données de population, j’ai choisi d’assembler les deux jeux de données en ajoutant le nombre de vols aux données de population avec une jointure gauche. J’ai nommé ce fichier assemblé “mydata” :

mydata <- pop %>% 
        left_join(vol2roues, by=c("num"="num_dep"))

head(mydata)
    ##   num                     nom pop_totale
    ## 1  01                     Ain     653688
    ## 2  02                   Aisne     528016
    ## 3  03                  Allier     333065
    ## 4  04 Alpes-de-Haute-Provence     161980
    ## 5  05            Hautes-Alpes     141784
    ## 6  06         Alpes-Maritimes    1080899
    ##                           libellé index nb_vols
    ## 1 Vols de véhicules motorisés à 2 roues      27
    ## 2 Vols de véhicules motorisés à 2 roues      20
    ## 3 Vols de véhicules motorisés à 2 roues      11
    ## 4 Vols de véhicules motorisés à 2 roues       6
    ## 5 Vols de véhicules motorisés à 2 roues       3
    ## 6 Vols de véhicules motorisés à 2 roues     137 

Création de la variable "taux de vols"

J’ai commencé par calculer le taux de vols en divisant le nombre de vols par la taille de la population.

mydata <- mydata %>% 
   mutate(Taux_vols=nb_vols/pop_totale)  

Puis dans un second temps, j’ai exprimé ce taux pour 100 000 habitants ne pas être embêtée avec les puissance du format scientifique (e-05).

 mydata <- mydata %>% 
        mutate(TV=100000*Taux_vols) 

head(mydata)
    ##   num                     nom pop_totale
    ## 1  01                     Ain     653688
    ## 2  02                   Aisne     528016
    ## 3  03                  Allier     333065
    ## 4  04 Alpes-de-Haute-Provence     161980
    ## 5  05            Hautes-Alpes     141784
    ## 6  06         Alpes-Maritimes    1080899
    ##                           libellé index nb_vols    Taux_vols        TV
    ## 1 Vols de véhicules motorisés à 2 roues      27 4.130411e-05  4.130411
    ## 2 Vols de véhicules motorisés à 2 roues      20 3.787764e-05  3.787764
    ## 3 Vols de véhicules motorisés à 2 roues      11 3.302659e-05  3.302659
    ## 4 Vols de véhicules motorisés à 2 roues       6 3.704161e-05  3.704161
    ## 5 Vols de véhicules motorisés à 2 roues       3 2.115895e-05  2.115895 

Carte du nombre de vols par habitants en mars 2019

La carte est créée avec les commandes suivantes : 

 France <- gadm_sf_loadCountries("FRA", level=2 )

mydata <- data.frame(mydata)
    choropleth(France, 
               data = mydata, 
               step=4,
               value = "TV", 
               adm.join = "nom",
               breaks = "quantile", 
               palette = "Oranges",
               legend="Taux de vols \n de deux roues\n pour 100 000 hab",
                title="Vols de deux roues en France") 

Vu comme ça c’est un peu chaud les Alpes de Haute Provence ! Mais si on trie les données par ordre décroissant des taux de vols :

 mydata <- mydata %>% 
        arrange(desc(TV)) %>% 
        mutate(place=1:nrow(.))
 

On peut voir que le 04 (Alpes de Haute Provence) est à la 48ème place, avec 3.7 vols pour 100 000 habitants ! Pile poile au milieu !

mydata[mydata$num=="04",]

    ##    num                     nom pop_totale
    ## 48  04 Alpes-de-Haute-Provence     161980
    ##                            libellé.index nb_vols    Taux_vols       TV
    ## 48 Vols de véhicules motorisés à 2 roues       6 3.704161e-05 3.704161
    ##    place
    ## 48    48 

Donc au final, on en s’en sort pas trop mal dans les Alpes de Haute Provence !

Conclusion

Et vous, que pensez vous de ce package GADMTools ?

Est ce que je vous ai convaincu de son utilité ?

Une dernière question : si j’achète un deux roues, je choisis lequel ? La Dax ou le Bubbly ? 😉

Si cet article vous a plu, ou vous a été utile, et si vous le souhaitez, vous pouvez soutenir ce blog en faisant un don sur sa page Tipeee 🙏

15 réponses

  1. Bonjour Madame Claire,

    merci pour tout le travail que vous faites sur ce blog.
    Et particulièrement aujourd’hui pour cet article.

    Par contre, je rencontre un problème lors de l’installation du package GADMTools. En effet, R m’indique que ce package n’est pas disponible pour la version 3.4.3 de R, est-ce normal ?

    Cordialement.

  2. Bonjour,

    Attention, les nombre de vols sont des données quantitatives absolues et doivent être représentées avec des cercles proportionnels, et non avec une plage de couleur. Pour pouvoir créer une carte choroplèthe, il faut rapporter les données au nombre d’habitants, comme vous y pensez d’ailleurs en fin d’article. Sinon, il est tout à fait normal de noter plus de vols à Paris et à Marseille, où la population est plus nombreuse.
    En l’état, cette carte est fausse.

    N’oubliez pas également l’échelle, indispensable à toute bonne carte et les sources date, indispensables à sa lecture critique.

    Quelques ressources :
    Représenter les données cartographiques : http://philcarto.free.fr/fortunel/FichesRepresentationDonnees.pdf, p20
    Représenter des données cartographiques sur R : https://riatelab.github.io/anfdataviz/

    Bien cordialement,

  3. Bonjour Claire,
    Article très intéressant, comme d’ailleurs vos articles antérieurs sur le logiciel R. Enfin, je suis comblé. Je cherchais un package R permettant de realiser des cartes administratives de maniere conviviale. Celui ci je le trouve tout simplement chic! Je voudrais savoir comment faire pour avoir l’abreviation des pays. La France c’est “FRA”, l’Italie c’est “ITA”. Et le Burkina Faso par exemple mon pays?
    Merci encore, c’est sympa de votre part.

  4. Bonjour Claire,

    Vraiment très intéressant surtout avec en plus la possibilité de mettre un fond de cartographie osm. Je me permets de répondre à Relwendé : c’est BFA comme indiqué dans le fichier disponible après le premier paragraphe de “Les différents niveaux de cartes administratives”.

    Merci encore pour vos articles.
    Thierry

  5. C’est bien bien claire cher @Maxime et @Thierry. Merci beaucoup et grand merci à @Claire pour cet article très intéressant.

  6. Bonjour, je trouve cet article super intéressant. J’ai essayé de refaire l’exemple pour ensuite pouvoir l’adapter à mes besoins mais dès les premières lignes de code j’ai un message d’erreur..
    > carte_globale <- gadm_sf_loadCountries("FRA", level=0 )
    trying URL 'https://biogeo.ucdavis.edu/data/gadm3.6/Rsf/gadm36_FRA_0_sf.rds&#039;
    Content type '`=»Îø' length 745650 bytes (728 KB)
    downloaded 728 KB
    Error in starts_with("GID") : could not find function "starts_with"

    Savez-vous d'où provient l'erreur ? Et surtout comment la résoudre ? D'avance merci beaucoup !
    Vincent

    1. Bonjour Vincent,

      non malheureusement ,je ne sais pas d’où vient erreur, ni comment la résoudre.
      Bonne continuation

    2. Perdón escribir en español, no sé Francés, yo lo solucioné instalando y cargando
      library(tidyverse)
      library(sf)
      library(GADMTools)
      library(mapview)

  7. Bonjour,

    Merci pour cet article avec des explications super détaillées. Ça fait plaisir de trouver ce type de contenu d’excellente qualité. Je me permets de vous partager quelques doutes:

    1) Comment vous faites pour changer l’échelle des valeurs dans la légende pour avoir des nombres standard ? J’obtiens des valeurs en log de mon côté.
    2) Quel est le format souhaité pour les vecteurs ? Je n’arrive pas à faire les catégories pour la légende manuellement.
    3) Comment faire si on ne souhaite pas avoir des breaks mais une échelle continue de valeurs ? En profitant que le package est construit sur ggplot j’ai essayé avec les scale_x_continous() et scale_y_continous() mais ça ne marche pas vu que choropleth utilise la fonction cut() en interne.
    4) Le package contient une option pour les territoires d’Outre-Mer ? Un collage avec ggrid (bien que possible) me semble peu pratique.

    Je vous remercie de votre aide,

  8. bonjour,

    j’ai essayé de faire une carte choroplethe mais j’ai l’erreur suivante “Error in classIntervals(.data[, .value], n = .steps, style = breaks) : var is not numeric”. ma variable que je veux présenter est numérique.
    j’ai desintallé et réinstallé le package ainsi que celui de ClassInt mais rien n’y fait.

    Merci de votre aide

    Bonne journée

    Valentin

    voici mon code
    library(GADMTools)
    data<-read_csv("https://www.data.gouv.fr/fr/datasets/r/4acad602-d8b1-4516-bc71-7d5574d5f33e&quot😉
    data2%
    filter(extract_date == today())

    France % as.numeric(data2$libelle_dep)
    replace(x = data2$taux_occupation_sae, “.”, “,”)

    choropleth(France,
    data = data2,
    step=4,
    value = “taux_occupation_sae”,
    adm.join = “libelle_dep”,
    breaks = “quantile”,
    palette = “Reds”,
    legend=”Taux d’occupation\n réanimation (%)”,
    title=”Taux d’occupation réanimation\n par rapport niveau habituel”)

Laisser un commentaire

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

Bonjour !

vous venez souvent ?

Identifiez-vous pour avoir accès à toutes les fontionnalités !

Aide mémoire off'R ;)

Enregistrez vous pour recevoir gratuitement mes fiches “aide mémoire” (ou cheat sheets) qui vous permettront de réaliser facilement les principales analyses biostatistiques avec le logiciel R et pour être informés des mises à jour du site.