Étude de la démographie française du XIXe siècle à partir de données collaboratives de généalogie

Arthur Charpentier et Ewen Gallic
freakonometrics.hypotheses.org
egallic.fr
Janvier 2018.

---

1 Aspects méthodologique : préparation des données

Le nombre de données relevant quasiment du milliard d’observations, un ordinateur classique se révèle pour l’heure un peu capricieux lors de l’import. La stratégie adoptée consiste à segmenter les observations par morceaux.

1.1 Découpage en tronçons

1.1.1 Une segmentation par tronçons plus petits

Les données dont nous disposons sont partagées en trois tronçons. Malgré ce découpage, le nombre de lignes reste trop conséquent pour l’ordinateur. Nous décidons d’effectuer un découpage plus important.

Pour importer un tronçons, nous utilisons la fonction fread. Nous lisons les observations par petits morceaux de \(20\times10^6\) lignes. Dans un premier temps, il nous est nécessaire de connaître le nombre total d’observations du tronçons sur lequel nous travaillons. L’instruction suivante nous permet de révéler cette information :

```r library(dplyr) library(stringr) library(ggplot2) library(readr) library(tidyr) library(stringi) library(data.table) library(pbapply)

nrow(fread(fichier, select = 1L)) ```

Le nom des colonnes s’obtient facilement :

col_label <- fread(fichier, nrows = 0)
col_label <- colnames(col_label)

Le premier tronçon contient \(412~406~274\) lignes. Nous pouvons créer une variable indiquant les valeurs de débuts et fins de lignes à importer, de manière à lire uniquement \(20\times10^6\) lignes au maximum à chaque importation.

seq_nbs <- seq(1, 412406274, by = 20*10^6)
seq_nbs <- c(seq_nbs, 412406274)
seq_nbs[1] <- 2 # La première ligne correspond aux noms de colonnes
a_parcourir <- data.frame(debut = seq_nbs[-length(seq_nbs)], fin = seq_nbs[-1]-1)

L’importation de chaque tronçon se fait à l’aide de la fonction fread, comme suit (ici, pour le premier tronçon) :

numero_chunk <- 1 # Tronçon 1
val_deb <- a_parcourir$debut[numero_chunk]
val_fin <- a_parcourir$fin[numero_chunk]
val_skip <- val_deb-1
val_nrows <- val_fin-val_deb+1
chunk <- fread(fichier,
encoding = "Latin-1",
header = FALSE, sep = ";",
fill=TRUE, col.names=col_label, quote='',
stringsAsFactors = FALSE,
colClasses=c(rep("character",19)),
na.strings=c("NA","NaN"," ","","0"),
skip = val_skip,
nrows = val_nrows)

Pour ce qui nous intéresse ici, nous pouvons nous permettre de séparer les observations les unes des autres, à condition de conserver celles provenant d’un même arbre d’utilisateur. Nous regroupons chaque observations en fonction du premier caractère alpha numérique des utilistaurs.

# Creation de chunks en fonction de la premiere lettre du nom du prorietaire
chunk[, prem_lettre:=str_to_lower(str_sub(sourcename, 1, 1))]

# Les valeurs différentes pour le nom d'utilisateur
motifs <- unique(chunk$prem_lettre) %>% sort()

# On conserve uniquement les lettres, et on ajoute une valeur pour "non-lettre" (pour les chiffres par exemple)
motifs <- c(motifs[str_detect(motifs, "[:letter:]")], "[^[:letter:]]")

# Sauvegarde dans des tronçons en fonction de la première lettre du nom d'utilisateur
pb2 <- txtProgressBar(min = 1, max = length(motifs), style = 3, char = "@")
for(ii in 1:length(motifs)){
  prem <- motifs[ii]
  if(ii == length(motifs)) prem <- "0"
  chunk_lettre <- chunk[str_detect(prem_lettre, motifs[ii])]
  chunk_lettre[,prem_lettre:=NULL]
  save(chunk_lettre, file = str_c("../data/Geneanet/raw/chunks/chunk_", sprintf("%02s", no_chunks_geneanet), "_", sprintf("%02s", numero_chunk), "_", prem, ".rda"))
  setTxtProgressBar(pb2, ii)
}

Il suffit alors de créer une petite fonction pour faire tourner sur les trois fichiers fournis par Geneanet, et de sauvegarder le résultat dans un fichier de données.

1.1.2 Regroupement des petits tronçons par noms d’utilisateurs

Comme les données initiales sont réparties dans trois fichiers, les petits tronçons obtenus dans l’étape précédente doivent être regroupés (si jamais des données de l’utilisateur toto se trouvent dans au moins deux fichiers différents).

# Chemin vers les tronçons
path_to_files <- "../data/Geneanet/raw/chunks/"
# Liste des fichiers
fichiers <- list.files(path_to_files, pattern = "\\.rda", full.names = TRUE)

motifs <- list.files(path_to_files, pattern = "\\.rda")
motifs <- str_replace(motifs, "chunk_0(1|2|3)_[[:digit:]]{2}_", "") %>% 
  str_replace("\\.rda", "")
motifs <- unique(motifs) %>% sort()
motifs <- c(motifs[str_detect(motifs, "[:letter:]")], "[^[:letter:]]")

# Créé le dossier d'enregistrements des tronçons, si besoin
if(!dir.exists("../data/Geneanet/raw/chunks_letter/")) dir.create("../data/Geneanet/raw/chunks_letter/", recursive = TRUE)

# motif <- "z"
save_chunks_lettre <- function(motif){
  fichiers_charger <- fichiers[str_detect(fichiers, str_c(motif, ".rda"))]
  # Charger les 
  chunks <- pblapply(fichiers_charger, function(x) {load(x) ; unique(chunk_lettre)})
  chunk_lettre <- rbindlist(chunks)
  chunk_lettre <- unique(chunk_lettre)
  if(motif == "[^[:letter:]]") motif <- "0"
  save(chunk_lettre, file = str_c("../data/Geneanet/raw/chunks_letter/chunk_", motif, ".rda"))
}

# Regroupe les tronçons et sauvegarde le résultat
pblapply(motifs, save_chunks_lettre)

1.1.3 Découpage en 5 fichiers pour chaque lettre d’utilisateurs

Afin d’accélérer les étapes suivantes, les tronçons d’utilisateurs sont découpés en 5 parties contenant plus ou moins le même nombre d’utilisateurs différents.

# Les tronçons
N <- list.files("../data/Geneanet/raw/chunks_letter/", pattern = "*.rda", full.names = TRUE)

# Création du dossier d'enregistrement des petits tronçons
if(!dir.exists("../data/Geneanet/raw/chunks_letter_2/")) dir.create("../data/Geneanet/raw/chunks_letter_2/", recursive = TRUE)

# Fonction pour découper des vecteurs en parties à peu près égales
chunk2 <- function(x,n) split(x, cut(seq_along(x), n, labels = FALSE)) 

#' decouper_chunks
#' @i: indice de position du fichier dans N
decouper_chunks <- function(i){
  # Obtenir la lettre des noms d'utilisateurs du tronçon courant
  lettre <- str_replace(N[i], "../data/Geneanet/raw/chunks_letter//chunk_", "")
  lettre <- str_replace(lettre, "\\.rda", "")
  # Chargement du tronçon
  load(N[i], envir = globalenv())
  # Les noms d'utilisateurs
  sourcenames <- unique(chunk_lettre$sourcename)
  # Partage équitable des noms d'utilisateurs (!= des observations)
  chunk_sourcenames <- chunk2(sourcenames, 5)
  
  #' sauvegarder_chunk
  #' @j: indice de chunk_sourcenames à extraire pour traiter
  #'     les utilisateurs qu'il contient
  sauvegarder_chunk <- function(j){
    sourcenames_chunk <- chunk_sourcenames[[j]]
    chunk_lettre_tmp <- chunk_lettre[sourcename %in% sourcenames_chunk]
    save(chunk_lettre_tmp, file = str_c("../data/Geneanet/raw/chunks_letter_2//chunk_lettre_tmp_",
                                        lettre, "_", str_pad(j, width = 2, pad = "0"), ".rda"))
  }# Fin de sauvegarder_chunk()
  
  # Création de 5 petits tronçons pour le fichier courant
  pblapply(1:5, sauvegarder_chunk)
  
}# Fin de decouper_chunks()

# Découper en 5 petits tronçons les fichiers
pblapply(1:length(N), decouper_chunks)

1.2 Regroupement par département et première simplification

Pour déployer les étapes de traitement des données suivantes sur plusieurs machines en même temps, nous regroupons les individus par département. Dans cette même étape, nous regroupons les informations de chaque individu sur une seule ligne. En effet, jusqu’à présent, une ligne correspond à un événement pour un individu dans l’arbre d’un utilisateur.

Pour commencer, nous consignons dans un tableau les différentes régions françaises présentes dans les données :

# Choix d'un tronçon
load("../data/Geneanet/raw/chunks_letter_2/chunk_lettre_tmp_0_05.rda")
liste_depts <-
  chunk_lettre_tmp %>%
  rename(dept = `sous-region`) %>% 
  select(pays, region, dept) %>%
  filter(pays == "FRA") %>%
  unique() %>%
  arrange(region, dept) %>%
  filter(!is.na(dept)) %>%
  tbl_df()

save(liste_depts, file = "liste_depts.rda")

Nous travaillons par tronçons obtenu dans l’étape précédente. Il est de fait nécessaire de les lister.

N <- list.files("../data/Geneanet/raw/chunks_letter_2/", full.names = TRUE)

Le code qui suit permet de traiter un seul tronçon, dont la position dans N est notée i_fichier. Nous prenons le premier fichier en exemple dans le code qui suit. Un simple boucle sur les indices des éléments de N permet de traiter chaque tronçon.

i_fichier <- 1
fichier <- N[i_fichier]
# La première lettre du nom d'utilisateur
lettre <- str_replace(fichier, "../data/Geneanet/raw/chunks_letter_2//chunk_lettre_tmp_", "") %>% 
  str_replace(., "_..\\.rda", "")
# Le numéro de tronçon
num_chunk <- str_replace(fichier, "../data/Geneanet/raw/chunks_letter_2//chunk_lettre_tmp_", "") %>% 
  str_replace(., "._", "") %>% 
  str_replace(., "\\.rda", "")

1.2.1 Quelques fonctions utiles

1.2.1.1 Lieux des événements renseignés

Trois types d’événements sont présents dans la base : naissance, mariage et décès. Nous créons une fonction pour indiquer les informations relatives à ces événements : géographie et date.

#' endroit_acte
#' 
#' Retourne un tableau de donnees indiquant la geographie de l'enregistrement
#' pour un type d'acte donnee
#' 
#' @x: (data.table) informations sur un individu
#' @acte : (string) type d'acte : "N" pour naissance, "M" pour mariage, ou "D" pour deces
#' x <- nai_per ; acte <- "N"
endroit_acte <- function(x, acte){
  if(nrow(x)>0){
    res <- 
      data.table(lieu = x$lieu,
                 dept = x$dept,
                 region = x$region,
                 pays = x$pays,
                 lat = x$latitude,
                 long = x$longitude,
                 stringsAsFactors = FALSE)
    
    if(acte == "N"){
      date <- x$date_naissance
    }else if(acte == "M"){
      date <- NA
    }else{
      date <- x$date_deces
    }
    
    res$date <- date
  }else{
    res <- 
      data.table(lieu = NA, dept = NA, region = NA, pays = NA, lat = NA, long = NA, date = NA, stringsAsFactors = FALSE)
  }
  
  res <- unique(res)
  names(res) <- str_c(names(res), "_", acte)
  res %>% tbl_df()
}# Fin de endroit_acte()

1.2.1.2 Simplification d’une personne

La fonction simplifier_personne() permet de créer un seul enregistrement par individu, pour un arbre d’utilisateur donné. Ainsi, au lieu d’avoir une ligne par événement, nous n’en conservons plus qu’une. Cela a une incidence sur les différents mariages : nous ne pouvons de fait en conserver qu’un seul.

#' simplifier_personne
#' 
#' Retourne les informations relatives a la naissance, le mariage et le deces
#' d'une personne
#' 
#' @un_id_personne: (string) id unique de la personne
#' @dt_source: (data.table) tableau de donnees dans lequel chercher les
#' informations pour cette personne
#' dt_source <- chunk_user ; un_id_personne <- ids_personnes[1]
simplifier_personne <- function(un_id_personne, dt_source){
  
  enregistrements_personne <- dt_source[id_personne == un_id_personne]
  
  sourcename_per <- enregistrements_personne$sourcename %>% unique()
  nom_per <- enregistrements_personne %>% 
    arrange(desc(str_length(nom))) %>% 
    slice(1) %>% 
    .$nom
  prenoms_per <- enregistrements_personne %>% 
    arrange(desc(str_length(prenoms))) %>% 
    slice(1) %>% 
    .$prenoms
  
  sexe_per <- enregistrements_personne$sexe %>% unique()
  
  id_mere_per <- enregistrements_personne$id_mere %>% unique()
  id_pere_per <- enregistrements_personne$id_pere %>% unique()
  
  # Naissance
  nai_per <- enregistrements_personne[stri_detect_regex(event_type, "N")]
  nai <- endroit_acte(nai_per, "N")
  # S'il y a plusieurs lieux de naissance, choisir celui le plus court
  if(nrow(nai)>1){
    nai <- 
      nai %>% 
      filter(!is.na(lieu_N)) %>% 
      arrange(str_length(lieu_N)) %>% 
      slice(1)
  }
  # Mariage
  mar_per <- enregistrements_personne[stri_detect_regex(event_type, "M")]
  mar <- endroit_acte(mar_per, "M")
  # S'il y a plusieurs dates de mariages : on ne retient que la plus ancienne
  if(nrow(mar)>1){
    mar <- 
      mar %>% 
      arrange(date_M) %>% 
      slice(1)
  }
  # Deces
  dec_per <- enregistrements_personne[stri_detect_regex(event_type, "D")]
  dec <- endroit_acte(dec_per, "D")
  # S'il y a plusieurs lieux de deces, choisir celui le plus court
  if(nrow(dec)>1){
    dec <- 
      dec %>% 
      filter(!is.na(lieu_D)) %>% 
      arrange(str_length(lieu_D)) %>% 
      slice(1)
  }
  
  nai$date_N <- enregistrements_personne$date_naissance %>% unique()
  dec$date_D <- enregistrements_personne$date_deces %>% unique()
  
  data.table(sourcename = sourcename_per, id_personne = un_id_personne,
             nom = nom_per, prenoms = prenoms_per, sexe = sexe_per,
             id_mere = id_mere_per, id_pere = id_pere_per, stringsAsFactors = FALSE) %>% 
    cbind(., nai) %>% 
    cbind(., mar) %>% 
    cbind(., dec) %>% 
    unique()
}# Fin de simplifier_personne()

1.2.1.3 Simplification des individus d’un même arbre

La fonction simplifier_proprietaire() permet de simplifier l’ensemble des individus présents dans l’arbre d’un usager de Geneanet. Elle s’appuie sur la fonction simplifier_personne() précédemment définie.

#' simplifier_proprietaire
#' 
#' Pour un proprietaire, simplifie les informations de tous les individus
#' 
#' @id_user: (string) identifiant d'un proprietaire d'arbre
#' @chunk_lettre: (tbl.df) base de donnees contenant les informations des utilisateurs dans le chunk
#' id_user <- "61p"
simplifier_proprietaire <- function(id_user, chunk_lettre){
  # Se restreinde aux individus de l'arbre du proprietaire
  chunk_user <- chunk_lettre[sourcename == id_user]
  
  # Creation d'un identifiant unique par personne
  chunk_user[, id_personne := str_c(sourcename, ref_locale, sep = "@")]
  
  # Ajout de l'identifiant des parents
  identifiants_mere <- chunk_user[, list(ID_num, id_personne)] %>% 
    rename(ID_num_mere = ID_num, id_mere = id_personne)
  
  identifiants_pere <- chunk_user[, list(ID_num, id_personne)] %>% 
    rename(ID_num_pere = ID_num, id_pere = id_personne)
  
  chunk_user <- identifiants_mere[chunk_user, on = "ID_num_mere"]
  chunk_user <- identifiants_pere[chunk_user, on = "ID_num_pere"]
  
  # Suppression de variables qui vont perturber la simplification
  # chunk_user[, c("ref_locale", "ID_num", "ID_num_pere", "ID_num_mere", "ID_num_conjoint") := NULL]
  chunk_user[, c("ID_num", "ID_num_pere", "ID_num_mere", "ID_num_conjoint") := NULL]
  
  # Pour chaque personne, un enregistrement correspond soit :
  # a une naissance (N), un mariage (M), un deces (D),
  # ou un mixe de ces evenements (e.g., NM pour naissance et mariage)
  # Cette distinction est faite parce que les lieux de N, M ou D peuvent varier
  # On va prendre uniquement la date du premier mariage en compte ici.
  ids_personnes <- chunk_user[!is.na(id_personne)]$id_personne %>% unique()
  
  lapply(ids_personnes, simplifier_personne, dt_source = chunk_user) %>% 
    rbindlist
  
}# simplifier_proprietaire

1.2.2 Simplification pour un tronçon

Nous chargeons en mémoire un tronçon précédemment créé, et nous travaillons sur celui-ci. Nous donnons ici les étapes de traitement pour un tronçon ; il est aisé de créer par la suite une fonction et de l’appliquer à chaque tronçon.

# Chargement du chunk
fichier <- N[1]
load(fichier, envir = .GlobalEnv)
chunk_lettre <- chunk_lettre_tmp
rm(chunk_lettre_tmp)

Nous listons les différents usagers (ou propriétaires d’arbre), et procédons à quelques opérations de nommage de variables, afin de faciliter l’écriture du code.

ids_sourcename <- chunk_lettre$sourcename %>% unique()
nbs_obs_ids_sourcename <- length(ids_sourcename)

chunk_lettre <- chunk_lettre %>% rename(dept = `sous-region`)
chunk_lettre <- chunk_lettre %>% rename(ID_num = increment,
                                        ID_num_mere = numero_mere,
                                        ID_num_pere = numero_pere, 
                                        ID_num_conjoint = numero_conjoint,
                                        prenoms = prenom)

Nous avons pris le parti de suivre les descendants d’individus nés entre 1800 et 1804 en France. La base est filtrée en conséquence.

gen_0 <- chunk_lettre
gen_0[, annee_nai := stri_sub(date_naissance, 1, 4)]
gen_0 <- gen_0[annee_nai %in% seq(1800, 1804)]
gen_0 <- gen_0[stri_detect_regex(event_type, "N")]
# Se restreindre uniquement aux abres dont un des individus est né entre 1800 et 1804
chunk_lettre <- chunk_lettre[sourcename %in% unique(gen_0$sourcename)]

1.2.2.1 Pour un tronçon, focus sur un département

Afin d’alléger les opérations à faire subir aux machines, nous travaillons par départements. Nous proposons ici la méthode pour gérer le cas d’un département. À nouveau, une fonction englobant le code qui suit permet de traiter l’ensemble des départements ultérieurement.

i_departement <- "F01"
departement <- liste_depts$dept[i_departement]

# Création du dossier de sauvegarde des résultats
if(!dir.exists(str_c("../data/individuals/migration/", departement, "/"))) dir.create(str_c("../data/individuals/migration/", departement, "/"), recursive = TRUE)

Nous filtrons la base pour se restreindre aux individus de première génération nés dans le département :

gen_0 <- gen_0[dept %in% departement]
# Noms d'utilisateurs de la sous-partie de la base de données
ids_sourcename <- gen_0$sourcename %>% unique()

Il s’agit à présent de se restreindre aux parents et descendants des individus ainsi sélectionnés.

# Initialisation

# La generation courante
gen_n <- gen_0[, list(sourcename, ID_num, ID_num_mere, ID_num_pere)]
nb_obs_gen_0 <- nrow(gen_n)
conserver <- gen_n[, list(sourcename, ID_num)]

# Les generations restantes (on retire la generation courante)
chunk_lettre <- chunk_lettre[sourcename %in% unique(gen_0$sourcename)]
gen_restants <- chunk_lettre[,list(sourcename, ID_num, ID_num_mere, ID_num_pere)]
gen_restants <- fsetdiff(gen_restants, gen_n, all = FALSE)

# Les parents
parents_gen_0 <-
  gen_n %>% 
  select(sourcename, ID_num_mere) %>% 
  rename(ID_num = ID_num_mere) %>% 
  bind_rows(
    gen_n %>% 
      select(sourcename, ID_num_pere) %>% 
      rename(ID_num = ID_num_pere)
  ) %>% 
  unique()

parents_gen_0_complet <- gen_restants[parents_gen_0, on = c("sourcename", "ID_num")]

# On retire les individus trouvés
# Bémol : cette étape retire les individus nés d'une relation incestueuse
gen_restants <- fsetdiff(gen_restants, parents_gen_0_complet, all = FALSE)


# Boucle
# Parmi ces personnes, lesquelles ont pour parent quelqu'un de la generation precedente
compteur <- 0
while(nrow(gen_n) > 0){
  compteur <- compteur+1
  if(compteur>=15) stop("Trop de tours")
  
  parents_m <- data.table(gen_n %>% select(sourcename, ID_num) %>% rename(ID_num_mere = ID_num))
  parents_p <- data.table(gen_n %>% select(sourcename, ID_num) %>% rename(ID_num_pere = ID_num))
  enfants <- data.table(gen_restants %>% select(sourcename, ID_num_mere, ID_num_pere, ID_num))
  
  
  gen_n <- enfants[parents_m, on = c("sourcename", "ID_num_mere")][!is.na(ID_num)] %>% 
    rbind(
      enfants[parents_p, on = c("sourcename", "ID_num_pere")][!is.na(ID_num)]
    ) %>% 
    unique()
  
  conserver <- rbind(conserver, gen_n %>% select(sourcename, ID_num)) %>% unique()
  # Les generations restantes
  gen_restants <- gen_restants[!gen_n, on = c("sourcename", "ID_num")]
}# Fin du while

# Les individus à conserver
conserver <- rbind(parents_gen_0, conserver)
# La base concernant ces individus
chunk_partiel <- chunk_lettre[conserver, on = c("sourcename", "ID_num")]

Reste alors à procéder à la simplification des lignes pour tous ces individus. Le code s’attache à le faire par utilisateur de Geneanet.

res <- pblapply(ids_sourcename, simplifier_proprietaire, chunk_lettre = chunk_partiel)

Comme cette dernière opération informatique prend beaucoup de temps, il est possible de proposer une version parallélisée :

library(parallel)
# Nombre de clusters
ncl <- detectCores()-1
(cl <- makeCluster(ncl))

invisible(clusterEvalQ(cl, library(dplyr, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(stringr, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(tidyr, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(stringi, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(data.table, warn.conflicts=FALSE, quietly=TRUE)))

# Evnoi des fonctions/données aux clusters
clusterExport(cl, c("simplifier_personne", "endroit_acte"))

res <- pblapply(ids_sourcename, simplifier_proprietaire, cl = cl, chunk_lettre = chunk_partiel)
stopCluster(cl)

Reste à rassembler chaque résultat dans un seul tableau :

res <- 
  res %>% rbindlist

Et enfin de sauvegarder le résultat :

save(res, file = str_c("../data/individuals/migration/", departement,"/chunk_", lettre, "_", num_chunk, ".rda"))

1.3 Rassemblement des arbres et élimination de doublons

Comme chaque usager du site Geneanet peut construire son propre arbre, il existe de nombreux doublons. Nous adoptons une stratégie en plusieurs étapes pour tenter de joindre les arbres entre-eux, tout en évitant les doublons. Du fait de l’important des données manquantes, cette tâche s’avère délicate. Il persiste sûrement au final des doublons, qui ne seront, pour la plupart, pas utilisés dans l’analyse. En effet, s’ils n’ont pas été repérés en tant que doublons, c’est qu’ils sont porteur de très peu d’informations, et seront donc écartés par la suite.

Cette étape dans le processus de nettoyage des données permet non seulement de relier les arbres d’utilisateurs entre eux, mais également de compléter certaines informations qui pourraient manquer dans l’arbre d’un usager mais être présentes dans celui d’un autre.

1.3.1 Quelques fonctions utiles

Nous présentons ici quelques fonctions permettant d’effectuer les regroupements.

1.3.1.1 Valeur la plus probable

La fonction most_probable_value() permet de trouver la valeur la plus probable parmi les différentes proposées. Elle prend en entrée le nom de la variable à laquelle on s’intéresse, le tableau de données contenant les valeurs, et une variable indiquant si des poids dans les observations sont fournis. La valeur finale sera celle dont la fréquence (pondérée le cas échéant) est la plus élevée parmi les propositions.

Cette fonction s’utilise sur des données concernant des individus étant identifiés comme désignant la même personne.

#' most_probable_value
#' Find the most probable value for a variable
#' using the values found within all candidates
#' (excluding NAs)
#' @variable_name: (string) name of the variable
#' @df: (data.table)
#' @weights: (logic) should the simplification consider the weights? (default: FALSE)
#' variable_name <- "date_N"
most_probable_value <- function(df, variable_name, weights = FALSE){
  valeurs <- df[[variable_name]]
  if(!all(is.na(valeurs))){
    if(weights){
      poids <- df[["weight"]]
      res <- data.frame(valeurs, poids, stringsAsFactors = FALSE) %>% 
        group_by(valeurs) %>% 
        summarise(poids = sum(poids)) %>% 
        ungroup()
      
      # Si la variable est une date
      # on va privilegier les informatins relatives aux dates completes
      if(variable_name %in% c("date_N", "date_M", "date_D")){
        tmp <- res %>% 
          mutate(mois = str_sub(valeurs, 5, 6),
                 jour = str_sub(valeurs, 7, 8)) %>% 
          filter(!(mois == "00" | jour == "00"))
        # S'il reste des informations, on se base sur elles
        if(nrow(tmp)>0) res <- tmp
      }
      
      res <- 
        res %>%
        arrange(desc(poids)) %>% 
        slice(1) %>% 
        magrittr::extract2("valeurs")
    }else{# Si pas de poids
      table_freq <- sort(table(valeurs))
      res <- names(table_freq)[1]
    }
  }else{
    res <- NA
  }
  res
}# End of most_probable_value()

1.3.1.2 Valeur la plus probable

La fonction simplifier_personne(), à partir d’individus identifiés comme désignant la même personne, simplifie les valeurs de chaque variable pour n’obtenir qu’une seule observation en sortie. Les valeurs de chaque variable sont obtenues en faisant appel à la fonction most_probable_value() précédemment définie.

# df_corresp <- corresp ; un_groupe_num <- 2 ; weights <- TRUE
# rm(df_corresp, un_groupe_num, groupe_cour, individus_cour, nouvelles_valeurs, weight)
#' @weights: (logic) should the simplification consider the weights? (default: TRUE)
simplifier_personne <- function(df_corresp, un_groupe_num, weights = TRUE){
  groupe_cour <- df_corresp[id_personne_num_groupe %in% un_groupe_num]
  
  individus_cour <- individus_simple[id_personne_num %in% groupe_cour$id_personne_num]
  
  weight <- sum(individus_cour$weight)
  
  if(nrow(individus_cour) == 1){
    nouvelles_valeurs <- individus_cour[, mget(c(variables_names))]
  }else{
    # Les valeurs probables pour les variables d'interet
    nouvelles_valeurs <- 
      variables_names %>% 
      sapply(., most_probable_value, df = individus_cour, weights = weights) %>% 
      t() %>% 
      data.table(stringsAsFactors = FALSE)
  }
  
  
  cbind(id_personne_num = un_groupe_num, nouvelles_valeurs, weight = weight)
}# End of simplifier_personne()

1.3.1.3 Une fonction de minimum

Une fonction permettant de gérer les valeurs manquantes lors de l’utilisation de la fonction min().

my_min <- function(x) ifelse( !all(is.na(x)), min(x, na.rm=T), NA)

1.3.2 Chargement des données

Nous allons retenir uniquement certaines variables pour les individus. Les variables textuelles de lieu sont évincées ici.

# Les vairables qui nous interessent
variables_names <- c("nom", "prenoms",
"sexe",
"lieu_N", "dept_N", "region_N", "pays_N", "lat_N", "long_N", "date_N",
"lieu_M", "dept_M", "region_M", "pays_M", "lat_M", "long_M", "date_M",
"lieu_D", "dept_D", "region_D", "pays_D", "lat_D", "long_D", "date_D",
"prenom")

Il s’agit ensuite de charger les données en mémoire dans la session. Nous repartons des données par département, en se concentrant uniquement sur les arbres dans lesquels se trouvent des individus nés entre 1800 et 1804. Nous ne montrons le traitement des donnée que pour un département. Il est aisé d’enrober le code dans une fonction pour ensuite le déployer sur l’ensemble des départements.

departement <- liste_depts$dept[i_departement]
# L'ensemble des données pour ce département
N <- list.files(str_c("../data/individuals/migration/", departement, "/"), pattern = "*.rda", full.names = TRUE)
# Chargement
individus <- 
  pblapply(N, function(x){load(x) ; res}) %>% 
  rbindlist() %>% 
  unique()

Les observations pour lesquelles les coordonnées de lieux sont 0.00 sont manquantes, donc à interpréter comme non disponibles dans R.

individus <- 
  individus %>% 
  mutate(long_N = ifelse(long_N == "0.00000", yes = NA, no = long_N),
         long_M = ifelse(long_M == "0.00000", yes = NA, no = long_M),
         long_D = ifelse(long_D == "0.00000", yes = NA, no = long_D)) %>% 
  mutate(lat_N = ifelse(lat_N == "0.00000", yes = NA, no = lat_N),
         lat_M = ifelse(lat_M == "0.00000", yes = NA, no = lat_M),
         lat_D = ifelse(lat_D == "0.00000", yes = NA, no = lat_D))

1.3.3 Préparation avant la première simplification

1.3.3.1 Création d’une base d’individus

Nous allons nous appuyer sur des identifiants numériques pour les individus. Nous créons une base individus recensant tous les individus et leur attribuant un identifiant.

# Les individus sont identifies de maniere unique par la colonne
# id_personne. Il s'agit d'une chaine de caracteres qui peut etre assez longue
# On va utiliser un identifiant numerique, ce qui sera plus pratique.
individus <- 
  individus %>% 
  tbl_df() %>% 
  mutate(id_personne_num = row_number()) %>% 
  data.table()

On va creer une variable contenant les prénoms de manière simplifiée, idem pour les noms. Nous retirons les espaces et caractères spéciaux.

individus <- 
  individus %>% 
  tbl_df() %>% 
  mutate(nom = stri_replace_all_regex(nom, "\\((.*?)\\)", "") %>% 
           stri_trans_general(., "Latin-ASCII") %>% 
           stri_replace_all_regex(., "'|\\(|\\)|[[:blank:]]|\\?|_|\\*", "") %>% 
           stri_trim(.) %>% 
           str_to_lower()) %>% 
  mutate(prenom = 
           str_replace_all(prenoms, "\\((.*?)\\)", "") %>% 
           stri_replace_all_regex(., "'|\\(|\\)|\\?|_|\\*|&gt|\\^|\\+|[[:digit:]]|\\.|<|>", "") %>% 
           str_trim() %>% 
           gsub(" .*$", "", .) %>%
           stri_trans_general(., "Latin-ASCII") %>% 
           stri_replace_all_regex(., "-", "") %>% 
           stri_trim(.) %>% 
           str_to_lower()) %>% 
  data.table()

# Les individus uniquement
individus_simple <- copy(individus)

Ensuite, nous ajoutons les informations relatives aux parents de chaque individu, afin de pouvoir identifier les individus en fonction de leurs parents. Si Jeanne Michu dans l’arbre d’un utilisateur a pour parents Jean Michu et Irène Dupond, et que dans un autre arbre, on trouve une personne née au même moment avec des parents portant également les mêmes noms, cela nous permet de rapprocher non seulement Jeanne, mais aussi Jean et Irène dans les deux arbres pour les fusionner.

On prépare de fait des tableaux pour les informations des parents.

# On va creer une base large avec pour chaque ligne,
# les informations sur les parents egalement.
# De fait, chaque ligne representera un triplet (enfant ; mere ; pere)
ind_meres <- match(individus$id_mere, individus$id_personne)
ind_peres <- match(individus$id_pere, individus$id_personne)

# Les meres
individus_meres <- 
  individus[ind_meres,] %>% 
  magrittr::set_colnames(str_c(colnames(individus), "_mere")) %>% 
  select(-sourcename_mere, -id_mere_mere, -id_pere_mere, -id_personne_mere)

# Les peres
individus_peres <- 
  individus[ind_peres,] %>% 
  magrittr::set_colnames(str_c(colnames(individus), "_pere")) %>% 
  select(-sourcename_pere, -id_mere_pere, -id_pere_pere, -id_personne_pere)

Puis on rajoute ces informations aux lignes des individus.

# On ajoute les informations relatives aux parents
individus <- 
  individus %>% 
  cbind(individus_meres) %>% 
  cbind(individus_peres) %>% 
  select(-id_mere, -id_pere) %>% 
  tbl_df()
rm(individus_meres, individus_peres)

Nous supprimons les observations pour lesquelles il y a une anomalie concernant la date de naissance les concernant eux ou leurs parents. Si un des deux parents a son enfant avant 13 ans, nous considérons qu’il s’agit d’une anomalie.

age_min_nenf <- 13

individus <-
  individus %>% 
  mutate(date_N_annee = str_sub(date_N, 1, 4) %>% as.numeric(),
         date_N_annee_mere = str_sub(date_N_mere, 1, 4) %>% as.numeric(),
         date_N_annee_pere = str_sub(date_N_pere, 1, 4) %>% as.numeric()) %>% 
  mutate(anomalie_mere = date_N_annee <= (date_N_annee_mere+age_min_nenf),
         anomalie_mere = ifelse(is.na(anomalie_mere), yes = FALSE, no = anomalie_mere)) %>% 
  mutate(anomalie_pere = date_N_annee <= (date_N_annee_pere+age_min_nenf),
         anomalie_pere = ifelse(is.na(anomalie_pere), yes = FALSE, no = anomalie_pere)) %>% 
  mutate(anomalie = anomalie_mere + anomalie_pere) %>% 
  filter(anomalie == 0) %>% 
  select(-anomalie_mere, -anomalie_pere, -anomalie) %>% 
  data.table()

# On repercute sur la table des individus
individus_simple <- individus_simple[id_personne_num %in% individus$id_personne_num]

individus_simple <- individus_simple[, mget(c("id_personne_num", variables_names))]

Nous rajoutons une colonne indiquant le poids de chaque observation, en prévision des simplifications à venir.

individus$weight <- 1
individus_simple$weight <- 1

1.3.3.2 Création d’un registre

Nous créons un registre indiquant les liens de parentés entre les individus. Ce registre indique pour chaque personne identifiée par son identifiant numérique, les identifiants numériques respectifs de ses parents. Notons que ce registre sera mis à jour par la suite, et se révèlera utile afin d’effectuer les rapprochements entre individus.

registre <- 
  individus %>% 
  select(id_personne_num, id_personne_num_mere, id_personne_num_pere) %>% 
  tbl_df()

1.3.4 Première simplification

Pour commencer, nous allons proceder a un ecremage tres grossier de la base, en supprimant les doublons evidents, c’est-à-dire ayant :

• le même nom de famille ;
• le même prénom ;
• le même sexe ;
• la même date de naissance ;
• le même nom de la mère ;
• le même nom de père.

1.3.4.1 Repérage d’individus identiques

Avant de procéder à la fusion des individus dédoublonés en entités uniques, nous devons les repérer. Le registre devient utile ici. Dès que deux individus sont repérées comme désignant une seule et même personne, l’identifiant numérique le plus petit d’entre les deux devient le nouvel identifiant des deux. Il en va de même pour les parents.

corresp_meres <- 
  individus[, list(id_personne,id_personne_num, nom, prenom, sexe, date_N, nom_mere, prenom_mere, sexe_mere, date_N_mere)] %>% 
  group_by(nom, prenom, sexe, date_N, nom_mere, prenom_mere, sexe_mere, date_N_mere) %>% 
  mutate(id_personne_num_groupe_mere = min(id_personne_num)) %>% 
  ungroup() %>% 
  select(id_personne, id_personne_num, id_personne_num_groupe_mere) %>% 
  data.table()

corresp_peres <- 
  individus[, list(id_personne, id_personne_num, nom, prenom, sexe, date_N, nom_pere, prenom_pere, sexe_pere, date_N_pere)] %>% 
  group_by(nom, prenom, sexe, date_N, nom_pere, prenom_pere, sexe_pere, date_N_pere) %>% 
  mutate(id_personne_num_groupe_pere = min(id_personne_num)) %>% 
  ungroup() %>% 
  select(id_personne, id_personne_num, id_personne_num_groupe_pere) %>% 
  data.table()

corresp <- 
  corresp_meres %>% 
  left_join(corresp_peres, by = c("id_personne_num", "id_personne")) %>% 
  data.table()

corresp <- 
  corresp %>%
  # rowwise() %>% 
  mutate(id_personne_num_groupe = pmin(id_personne_num_groupe_mere, id_personne_num_groupe_pere, na.rm=T)) %>%
  select(-id_personne_num_groupe_mere, -id_personne_num_groupe_pere) %>% 
  data.table()

corresp <- unique(corresp)

Cette première simplification permet peut-être (et c’est le cas) de regrouper davantage d’individus entre eux, en fonction du registre. Nous mettons alors une procédure itérative visant à repérer les fusions qu’il est possible de faire. Dès qu’il n’y a plus de regroupement possibles, la boucle s’achève.

continuer <- TRUE
nb_iterations <- 1

while(continuer){
  registre_new <- 
    registre %>% 
    tbl_df() %>%
    left_join(corresp %>% select(-id_personne) %>% rename(id_personne_num_new = id_personne_num_groupe), by = "id_personne_num") %>% 
    left_join(
      corresp %>% select(-id_personne) %>% rename(id_personne_num_mere = id_personne_num,
                                                  id_personne_num_mere_new = id_personne_num_groupe),
      by = "id_personne_num_mere"
    ) %>%
    left_join(
      corresp %>% select(-id_personne) %>% rename(id_personne_num_pere = id_personne_num,
                                                  id_personne_num_pere_new = id_personne_num_groupe),
      by = "id_personne_num_pere"
    )
  
  
  registre_new <- 
    registre_new %>% 
    mutate(same_id = id_personne_num == id_personne_num_new,
           same_id_mother = id_personne_num_mere == id_personne_num_mere_new,
           same_id_father = id_personne_num_pere == id_personne_num_pere_new) %>% 
    mutate(same_id = ifelse(is.na(same_id), yes = TRUE, no = same_id),
           same_id_mother = ifelse(is.na(same_id_mother), yes = TRUE, no = same_id_mother),
           same_id_father = ifelse(is.na(same_id_father), yes = TRUE, no = same_id_father))
  
  continuer <- any(!registre_new$same_id | !registre_new$same_id_mother | !registre_new$same_id_father)
  
  if(!continuer) next
  
  # On a pu retrouver d'autres associations !
  # Celles des parents d'un enfant.
  
  # Pour les mères
  corresp_meres <-
    registre_new %>% 
    select(id_personne_num_new, id_personne_num_mere) %>% 
    unique() %>% 
    arrange(id_personne_num_new, id_personne_num_mere) %>% 
    filter(!is.na(id_personne_num_mere))
  
  corresp_meres <- data.table(corresp_meres)
  corresp_meres[, id_personne_num_mere_2 := my_min(id_personne_num_mere), by = id_personne_num_new]
  corresp_meres <- 
    corresp_meres %>% select(-id_personne_num_new) %>% 
    rename(id_personne_num = id_personne_num_mere, id_personne_num_groupe = id_personne_num_mere_2) %>% 
    unique()
  
  corresp_meres[, id_personne_num_groupe := my_min(id_personne_num_groupe), by = id_personne_num]
  
  # Idem pour les pères
  corresp_peres <-
    registre_new %>% 
    select(id_personne_num_new, id_personne_num_pere) %>% 
    unique() %>% 
    arrange(id_personne_num_new, id_personne_num_pere) %>% 
    filter(!is.na(id_personne_num_pere)) 
  corresp_peres <- data.table(corresp_peres)
  corresp_peres[, id_personne_num_pere_2 := my_min(id_personne_num_pere), by = id_personne_num_new]
  corresp_peres <- 
    corresp_peres %>% select(-id_personne_num_new) %>% 
    rename(id_personne_num = id_personne_num_pere, id_personne_num_groupe = id_personne_num_pere_2) %>% 
    unique()
  
  corresp_peres[, id_personne_num_groupe := my_min(id_personne_num_groupe), by = id_personne_num]
  
  corresp <-
    corresp %>% 
    tbl_df() %>% 
    left_join(
      corresp_meres %>% 
        rename(id_personne_num_groupe_new = id_personne_num_groupe),
      by = "id_personne_num"
    ) %>% 
    tbl_df() %>% 
    mutate(id_personne_num_groupe = pmin(id_personne_num_groupe, id_personne_num_groupe_new, na.rm=T)) %>% 
    select(-id_personne_num_groupe_new) %>% 
    left_join(
      corresp_peres %>% 
        rename(id_personne_num_groupe_new = id_personne_num_groupe),
      by = "id_personne_num"
    ) %>% 
    mutate(id_personne_num_groupe = pmin(id_personne_num_groupe, id_personne_num_groupe_new, na.rm=T)) %>% 
    select(-id_personne_num_groupe_new) %>% 
    unique()
  
  registre_new <- 
    registre_new %>% 
    mutate(id_personne_num = id_personne_num_new,
           id_personne_num_mere = id_personne_num_mere_new,
           id_personne_num_pere = id_personne_num_pere_new) %>% 
    select(-id_personne_num_new, -id_personne_num_mere_new, -id_personne_num_pere_new,
           -same_id, -same_id_mother, -same_id_father)
  
  registre_new <- data.table(registre_new)
  registre_new <- registre_new[, list(
    id_personne_num_mere = min(id_personne_num_mere, na.rm=T),
    id_personne_num_pere = min(id_personne_num_pere, na.rm=T)), by = id_personne_num] %>% 
    unique()
  
  registre_new$id_personne_num_mere[which(registre_new$id_personne_num_mere == Inf)] <- NA
  registre_new$id_personne_num_pere[which(registre_new$id_personne_num_pere == Inf)] <- NA
  
  # toto <- data.frame(id = c(1,1,2,3,3,4,3), val = c(1,0,2,NA, NA, 2,NA)) %>% data.table()
  # toto[, list(test = pmin(val)), by = id]
  # 
  # 
  # registre_new <- 
  #   registre_new %>% 
  #   tbl_df() %>% 
  #   group_by(id_personne_num) %>% 
  #   summarise(id_personne_num_mere = my_min(id_personne_num_mere),
  #             id_personne_num_pere = my_min(id_personne_num_pere)) %>% 
  #   ungroup()
  
  # registre_new <- unique(registre_new)
  
  registre <- copy(registre_new)
  nb_iterations <- nb_iterations+1
}# Fin de la boucle while

corresp <- data.table(corresp)
registre <- data.table(registre)

1.3.4.2 Fusion des informations des individus

Comme le repérage des dédoublonsages est effectué, il reste à fusionner les informations des doublons. Nous utilisons la fonction simplifier_personne() présentée précédemment pour le faire. Si trois Madame Michu sont repérées comme désignant la même personne, et que deux d’entre-elles ont la date de naissance suivante 1800-01-01 tandis que la troisième renseigne un jour différent, alors, du fait d’une fréquence plus élevée pour la valeur 1800-01-01, nous retiendrons cette dernière comme la date de naissance de Madame Michu.

Nous procédons à la simplifications en fonction des groupes de personnes identiques que l’on a identifiés.

corresp_nb <- 
  corresp %>% 
  tbl_df() %>% 
  group_by(id_personne_num_groupe) %>% 
  tally()

Il n’est pas nécessaire de s’attarder sur les singletons, qui ne nécessitent aucun traitement et qui ralentiraient l’exécution du code. Nous les mettons de côté ; ils seront réintégrés une fois la simplification effectuée.

groupes <- corresp_nb %>% filter(n > 1) %>% magrittr::extract2("id_personne_num_groupe")

Il reste à parcourir chacun de ces groupes afin d’effectuer la simplification. Comme ils sont nombreux, nous parallélisons l’exécution du code.

library(parallel)
ncl <- detectCores()-1
(cl <- makeCluster(ncl))

invisible(clusterEvalQ(cl, library(dplyr, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(stringr, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(tidyr, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(stringi, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(data.table, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(stringdist, warn.conflicts=FALSE, quietly=TRUE)))

# Envoi des fonctions/tableaux 
clusterExport(cl, c("most_probable_value", "variables_names", "individus_simple"), envir=environment())

individus_simple_enfants <- pblapply(groupes, simplifier_personne, df_corresp = corresp, cl = cl)
stopCluster(cl)
individus_simple_enfants <- rbindlist(individus_simple_enfants)

Nous rajoutons les singletons laissés de côté.

groupes_un <- corresp_nb %>% filter(n == 1) %>% magrittr::extract2("id_personne_num_groupe")
individus_simple_un <- individus_simple[id_personne_num %in% groupes_un]
individus_simple <- rbind(individus_simple_un, individus_simple_enfants)
rm(corresp_meres, corresp_peres, individus_simple_enfants, registre_new, groupes, groupes_un, individus_simple_un)

Nous disposons donc de trois tables :

1. La table de registres, qui indique pour chaque personne : son identifiant, ainsi que ceux de ses parents (si disponibles) ;
2. La table des individus, avec leur id et leurs renseignements ;
3. La table de correspondance (n’est plus utile, nous allons en construire une autre).

1.3.5 Seconde simplification

Il reste dans la base de données, des doublons moins évidents. Dans l’esprit, nous procédons de la même manière que pour la première simplification, en rajoutant les informations relatives au parent et en créan un registre, puis en simplifiant les doublons. Dans cette seconde simplification, nous prenons en compte les erreurs d’orthographe dans les noms et prénoms.

1.3.5.1 Repérage d’individus identiques

# On va creer une base large avec pour chaque ligne,
# les informations sur les parents egalement.
# De fait, chaque ligne representera un triplet (enfant ; mere ; pere)
individus_2 <- 
  copy(individus_simple)

# On rajoute les identifiants des parents
individus_2$id_personne_num_mere <- registre$id_personne_num_mere[match(individus_simple$id_personne_num, registre$id_personne_num)]
individus_2$id_personne_num_pere <- registre$id_personne_num_pere[match(individus_simple$id_personne_num, registre$id_personne_num)]


ind_meres <- match(individus_2$id_personne_num_mere, individus_2$id_personne_num)
ind_peres <- match(individus_2$id_personne_num_pere, individus_2$id_personne_num)

# Les meres
individus_meres <- 
  individus_2[ind_meres,] %>% 
  magrittr::set_colnames(str_c(colnames(individus_2), "_mere")) %>% 
  select(-id_personne_num_mere, -id_personne_num_mere_mere, -id_personne_num_pere_mere)

# Les peres
individus_peres <- 
  individus_2[ind_peres,] %>% 
  magrittr::set_colnames(str_c(colnames(individus_2), "_pere")) %>% 
  select(-id_personne_num_pere, -id_personne_num_mere_pere, -id_personne_num_pere_pere)


# On ajoute les informations relatives aux parents
individus_2 <- 
  individus_2 %>% 
  cbind(individus_meres) %>% 
  cbind(individus_peres) %>% 
  # select(-id_mere, -id_pere) %>% 
  tbl_df()

individus_2 <- data.table(individus_2)

individus_2$date_N_annee <- individus_2$date_N %>% str_sub(1,4) %>% as.numeric()

Nous repérons les doublons en tentant de prendre en compte les erreurs d’orthographe dans les noms et prénoms. Pour ce faire, nous nous appuyons sur la fonction stringdist() du package du même nom.

library(stringdisy)
#' calculer_distance
#' help function for ecremer()
# variable <- "nom"
calculer_distance <- function(personne_cour, candidats, variable, threshold){
  stringdist(personne_cour[,get(variable)], candidats[,get(variable)], method = "cosine") < threshold
}# End of calculer_distance()

Nous utilisons la fonction ecremer() définie ci-après pour repérer les doublons potentiels. Nous considérons à présent que deux observations avec à peu près les mêmes nom et prénom, à peu près le même nom de mère ou de père sont potentiellement fusionnables.

#' ecremer
#' Retourne un tableau de candidats possibles
#' @candidats: (data.frame) table dans laquelle chercher les candidats
ecremer <- function(candidats, threshold = 0.1, personne_cour, prenoms = FALSE){
  # Ecremer par le nom
  ind_nom <- which(calculer_distance(personne_cour = personne_cour, candidats = candidats, variable = "nom", threshold = threshold))
  candidats <- candidats[ind_nom,]
  
  if(nrow(candidats)>0){
    # Par le prenom
    if(prenoms == TRUE){
      ind_prenom <- which(calculer_distance(personne_cour = personne_cour, candidats = candidats, variable = "prenoms", threshold = threshold)) 
    }else{
      ind_prenom <- which(calculer_distance(personne_cour = personne_cour, candidats = candidats, variable = "prenom", threshold = threshold))
    }
    candidats <- candidats[ind_prenom,]
    
    if(nrow(candidats)>0 & !is.na(personne_cour$nom_mere)){
      ind_mere_nom <- calculer_distance(personne_cour = personne_cour, candidats = candidats, variable = "nom_mere", threshold = threshold)
      ind_mere_nom <- c(which(is.na(ind_mere_nom)), which(ind_mere_nom))
      candidats <- candidats[ind_mere_nom,]
      
    }# Fin if mere
    
    if(nrow(candidats)>0 & !is.na(personne_cour$nom_pere)){
      ind_pere_nom <- calculer_distance(personne_cour = personne_cour, candidats = candidats, variable = "nom_pere", threshold = threshold)
      ind_pere_nom <- c(which(is.na(ind_pere_nom)), which(ind_pere_nom)) %>% sort()
      candidats <- candidats[ind_pere_nom,]
    }# Fin if pere

  }# Fin du if
  candidats
}# End of ecremer()

Enfin, pour pouvoir effectuer le repérage, nous utilisons la fonction trouver_corresp_indice_naissance() qui s’appuie sur les résultats obtenus par la fonction ecremer() pour définir si un rapprochement peut être effectué entre plusieurs individus. Il faut toutefois que les personnes, outre le fait d’avoir un nom similaire, soient nées dans le même département la même année.

#' trouver_corresp_indice_naissance
#' Pour un numero de ligne du tableau individus
#' retourne sous forme de table les correspondances entre individus.
#' Chaque ligne du resultat correspond a l'identifiant d'une personne
#' et au plus petit identifiant d'une correspondance trouvee dans la base.
trouver_corresp_indice_naissance <- function(i){
  personne_cour <- individus_3[i, ]
  
  resul <- data.table(id_personne_num = personne_cour$id_personne_num,
                      id_personne_num_groupe = personne_cour$id_personne_num)
  
  # Les identifiants des parents
  personne_cour_id_num_mere <- personne_cour$id_personne_num_mere
  personne_cour_id_num_pere <- personne_cour$id_personne_num_pere
  
  # Si les identifiants des deux parents sont manquants, on ne va pas plus loin
  if(!(is.na(personne_cour_id_num_mere) & is.na(personne_cour_id_num_pere))){
    # On a un identifiant pour: mere U pere
    
    id_personne_num <- personne_cour$id_personne_num
    
    # Il faut que les candidats soient de la même année
    # Et du meme sexe
    
    # Il faut que les candidats soient nes la meme annee et soient du meme sexe
    ind <- which(individus_3$date_N_annee %in% personne_cour$date_N_annee & 
                   individus_3$sexe %in% personne_cour$sexe)
    
    candidats <- individus_3[ind]
    
    # individus %>% filter(id_personne == "yvono@le garrec|jean|1") %>% View()
    
    candidats <- candidats[!id_personne_num %in% personne_cour$id_personne_num]
    
    candidats <- ecremer(candidats, personne_cour = personne_cour)
    
    if(nrow(candidats)>0){
      id_retenir <- c(id_personne_num, candidats$id_personne_num)
      resul <- data.table(id_personne_num = id_retenir, id_personne_num_groupe = my_min(id_retenir))
    }# Fin du if()
  }# Fin du if sur les parents
  
  # }# Fin du if sur la date de naissance
  
  resul
}# End of trouver_corresp_indice_naissance()

Nous pouvons utiliser les fonctions précédemment introduite pour repérer les individus à fusionner. Afin de ne pas perdre de temps dans l’exécution du code pour des individus étant rapidement identifiés comme ne pouvant pas être fusionnés, nous repérons a priori ces individus et les mettons de côté.

ind_ne_pas_faire <- which(is.na(individus_2$date_N_annee) | 
                            is.na(individus_2$sexe) | 
                            (is.na(individus_2$id_personne_num_mere) & is.na(individus_2$id_personne_num_pere)))

corresp_naissance_ne_pas_faire <- 
  individus_2[ind_ne_pas_faire,] %>% 
  select(id_personne_num) %>% 
  mutate(id_personne_num_groupe = id_personne_num) %>% 
  data.table()

Reste à exécuter la fonction trouver_corresp_indice_naissance() pour chaque individu restant. À nouveau, cette exécution est parallélisée.

ind_faire <- seq(1, nrow(individus_2))
ind_faire <- ind_faire[-ind_ne_pas_faire]

library(parallel)
ncl <- detectCores()-1
(cl <- makeCluster(ncl))

invisible(clusterEvalQ(cl, library(dplyr, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(stringr, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(tidyr, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(stringi, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(data.table, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(stringdist, warn.conflicts=FALSE, quietly=TRUE)))

individus_3 <- individus_2[ind_faire]

clusterExport(cl, c("individus_3"), envir=environment())
clusterExport(cl, c("calculer_distance", "ecremer", "my_min"), envir=environment())
corresp_naissance <- pblapply(1:nrow(individus_3), trouver_corresp_indice_naissance, cl = cl)
stopCluster(cl)

corresp_naissance <- 
  corresp_naissance %>% 
  rbindlist()

corresp_naissance <- unique(corresp_naissance)

corresp_naissance <- 
  corresp_naissance %>% 
  bind_rows(corresp_naissance_ne_pas_faire) %>% 
  unique()

rm(ind_faire, individus_3)

Nous mettons ensuite la table de correspondance à jour.

corresp_naissance <- corresp_naissance[, list(id_personne_num_groupe = min(id_personne_num_groupe, na.rm=T)), by = id_personne_num] %>% unique()

corresp <- 
  corresp %>% 
  left_join(corresp_naissance %>% rename(id_personne_num_groupe_new = id_personne_num_groupe)) %>% 
  mutate(id_personne_num_groupe = pmin(id_personne_num_groupe, id_personne_num_groupe_new, na.rm=T)) %>% 
  select(-id_personne_num_groupe_new) %>% 
  data.table()

Puis, nous regardons dans cette table, de manière itérative comme lors de la première opération de simplification, s’il n’est pas possible d’effectuer d’autres rapprochements en fonction des liens de parentés.

continuer <- TRUE
nb_iterations <- 1

while(continuer){
  registre_new <- 
    registre %>% 
    tbl_df() %>%
    left_join(corresp %>% select(-id_personne) %>% rename(id_personne_num_new = id_personne_num_groupe), by = "id_personne_num") %>% 
    left_join(
      corresp %>% select(-id_personne) %>% rename(id_personne_num_mere = id_personne_num,
                                                  id_personne_num_mere_new = id_personne_num_groupe),
      by = "id_personne_num_mere"
    ) %>%
    left_join(
      corresp %>% select(-id_personne) %>% rename(id_personne_num_pere = id_personne_num,
                                                  id_personne_num_pere_new = id_personne_num_groupe),
      by = "id_personne_num_pere"
    )
  
  
  registre_new <- 
    registre_new %>% 
    mutate(same_id = id_personne_num == id_personne_num_new,
           same_id_mother = id_personne_num_mere == id_personne_num_mere_new,
           same_id_father = id_personne_num_pere == id_personne_num_pere_new) %>% 
    mutate(same_id = ifelse(is.na(same_id), yes = TRUE, no = same_id),
           same_id_mother = ifelse(is.na(same_id_mother), yes = TRUE, no = same_id_mother),
           same_id_father = ifelse(is.na(same_id_father), yes = TRUE, no = same_id_father))
  
  continuer <- any(!registre_new$same_id | !registre_new$same_id_mother | !registre_new$same_id_father)
  
  if(!continuer) next
  
  corresp_meres <-
    registre_new %>% 
    select(id_personne_num_new, id_personne_num_mere) %>% 
    unique() %>% 
    arrange(id_personne_num_new, id_personne_num_mere) %>% 
    filter(!is.na(id_personne_num_mere))
  
  corresp_meres <- data.table(corresp_meres)
  corresp_meres[, id_personne_num_mere_2 := my_min(id_personne_num_mere), by = id_personne_num_new]
  corresp_meres <- 
    corresp_meres %>% select(-id_personne_num_new) %>% 
    rename(id_personne_num = id_personne_num_mere, id_personne_num_groupe = id_personne_num_mere_2) %>% 
    unique()
  
  corresp_meres[, id_personne_num_groupe := my_min(id_personne_num_groupe), by = id_personne_num]
  
  corresp_peres <-
    registre_new %>% 
    select(id_personne_num_new, id_personne_num_pere) %>% 
    unique() %>% 
    arrange(id_personne_num_new, id_personne_num_pere) %>% 
    filter(!is.na(id_personne_num_pere)) 
  
  corresp_peres <- data.table(corresp_peres)
  corresp_peres[, id_personne_num_pere_2 := my_min(id_personne_num_pere), by = id_personne_num_new]
  corresp_peres <- 
    corresp_peres %>% select(-id_personne_num_new) %>% 
    rename(id_personne_num = id_personne_num_pere, id_personne_num_groupe = id_personne_num_pere_2) %>% 
    unique()
  
  corresp_peres[, id_personne_num_groupe := my_min(id_personne_num_groupe), by = id_personne_num]
  
  corresp <-
    corresp %>% 
    tbl_df() %>% 
    left_join(
      corresp_meres %>% 
        rename(id_personne_num_groupe_new = id_personne_num_groupe),
      by = "id_personne_num"
    ) %>% 
    tbl_df() %>% 
    mutate(id_personne_num_groupe = pmin(id_personne_num_groupe, id_personne_num_groupe_new, na.rm=T)) %>% 
    select(-id_personne_num_groupe_new) %>% 
    left_join(
      corresp_peres %>% 
        rename(id_personne_num_groupe_new = id_personne_num_groupe),
      by = "id_personne_num"
    ) %>% 
    mutate(id_personne_num_groupe = pmin(id_personne_num_groupe, id_personne_num_groupe_new, na.rm=T)) %>% 
    select(-id_personne_num_groupe_new) %>% 
    unique()
  
  registre_new <- 
    registre_new %>% 
    mutate(id_personne_num = id_personne_num_new,
           id_personne_num_mere = id_personne_num_mere_new,
           id_personne_num_pere = id_personne_num_pere_new) %>% 
    select(-id_personne_num_new, -id_personne_num_mere_new, -id_personne_num_pere_new,
           -same_id, -same_id_mother, -same_id_father)
  
  registre_new <- data.table(registre_new)
  registre_new <- registre_new[, list(
    id_personne_num_mere = min(id_personne_num_mere, na.rm=T),
    id_personne_num_pere = min(id_personne_num_pere, na.rm=T)), by = id_personne_num] %>% 
    unique()
  
  registre_new$id_personne_num_mere[which(registre_new$id_personne_num_mere == Inf)] <- NA
  registre_new$id_personne_num_pere[which(registre_new$id_personne_num_pere == Inf)] <- NA
  
  registre <- copy(registre_new)
  nb_iterations <- nb_iterations+1
}# Fin de la boucle while

corresp <- data.table(corresp)
registre <- data.table(registre)

1.3.5.2 Fusion des informations des individus

Le répérage des nouveaux doublons est effectué, il reste alors à fusionner les doublons, comme lors de la première fusion.

individus_simple <- individus_2[, mget(colnames(individus_simple))]

corresp_naissance_nb <- 
  corresp %>% 
  tbl_df() %>% 
  group_by(id_personne_num_groupe) %>% 
  tally()


groupes_naissance <- corresp_naissance_nb %>% filter(n > 1) %>% magrittr::extract2("id_personne_num_groupe")

library(parallel)
ncl <- detectCores()-1
(cl <- makeCluster(ncl))

invisible(clusterEvalQ(cl, library(dplyr, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(stringr, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(tidyr, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(stringi, warn.conflicts=FALSE, quietly=TRUE)))
invisible(clusterEvalQ(cl, library(data.table, warn.conflicts=FALSE, quietly=TRUE)))

clusterExport(cl, c("most_probable_value", "variables_names", "individus_simple"), envir=environment())

individus_simple_naissance <- pblapply(groupes_naissance, simplifier_personne, df_corresp = corresp, cl = cl, weights=TRUE)
individus_simple_naissance <- rbindlist(individus_simple_naissance)
stopCluster(cl)

Nous rajoutons ensuite les individus mis de côté précédemment.

groupes_un <- corresp_naissance_nb %>% filter(n == 1) %>% magrittr::extract2("id_personne_num_groupe")
individus_simple_naissance_un <- individus_simple[id_personne_num %in% groupes_un]
individus_simple_naissance <- rbind(individus_simple_naissance_un, individus_simple_naissance)
individus_simple <- individus_simple_naissance

1.3.6 Troisième simplification

La troisième simplification ressemble en tous points à la seconde, à l’exception de l’étape de repérage, dans lequel des conditions différentes sont utilisées pour indiquer si deux individus désignent une même personne. Nous ne présentons donc que les deux fonctions qui changent à cette étape, afin de faciliter la lecture.

Cette troisième simplification vise à corriger des erreurs de renseignements dans le sexe des individus.

Pour que deux individus soient désignés identiques, ils doivent avoir au moins un parent en commun, être nés, mariés ou décédés dans le même département. À nouveau, il faut que leur nom soit à peu près identique (nous faisons de nouveau appel à la fonction ecremer()).

#' trouver_corresp_indice_sexe
#' Pour un numero de ligne du tableau individus
#' retourne sous forme de table les correspondances entre individus.
#' Chaque ligne du resultat correspond a l'identifiant d'une personne
#' et au plus petit identifiant d'une correspondance trouvee dans la base.
# i <- 199 ; i <- 32
trouver_corresp_indice_sexe <- function(i){
  personne_cour <- individus_3[i, ]
  
  resul <- data.table(id_personne_num = personne_cour$id_personne_num,
                      id_personne_num_groupe = personne_cour$id_personne_num)
  
  # Les identifiants des parents
  personne_cour_id_num_mere <- personne_cour$id_personne_num_mere
  personne_cour_id_num_pere <- personne_cour$id_personne_num_pere
  
  # Si les identifiants des deux parents sont manquants, on ne va pas plus loin
  if(!(is.na(personne_cour_id_num_mere) & is.na(personne_cour_id_num_pere))){
    # On a un identifiant pour: mere U pere
    id_personne_num <- personne_cour$id_personne_num
    
    # Il faut que les candidats soient du meme departement
    if(!is.na(personne_cour$date_N)){
      ind_N <- which(individus_3$date_N %in% personne_cour$date_N)
    }else{
      ind_N <- NULL
    }
    
    if(!is.na(personne_cour$date_M)){
      ind_M <- which(individus_3$date_M %in% personne_cour$date_M)
    }else{
      ind_M <- NULL
    }
    
    if(!is.na(personne_cour$date_D)){
      ind_D <- which(individus_3$date_D %in% personne_cour$date_D)
    }else{
      ind_D <- NULL
    }
    
    ind <- c(ind_N, ind_M, ind_D) %>% unique()
    
    candidats <- individus_3[ind]
    
    if(!is.na(personne_cour$dept_N)){
      candidats <- 
        candidats[dept_N %in% personne_cour$dept_N]
    }
    
    candidats <- candidats[!id_personne_num %in% personne_cour$id_personne_num]
    
    candidats <- ecremer(candidats, personne_cour = personne_cour, threshold = 0.02)
    
    if(nrow(candidats)>0){
      id_retenir <- c(id_personne_num, candidats$id_personne_num)
      resul <- data.table(id_personne_num = id_retenir, id_personne_num_groupe = my_min(id_retenir))
    }# Fin du if()
  }# Fin du if sur les parents
  
  # }# Fin du if sur la date de naissance
  
  resul
}# End of trouver_corresp_indice_sexe()

Le reste du code est similaire à celui de la seconde simplification.

1.3.7 Quatrième simplification

La quatrième simplification vise à corriger les erreurs concernant les dates. Certaines dates indiquent une année, un mois et un jour. D’autres ne sont pas précises sur le mois ou le jour, et en lieu du nombre correspondant au mois ou au jour, la valeur 00 est indiquée. La simplification présente tente de prendre en compte cet aspect.

À nouveau, la fonction de repérage de doublons est l’élément changeant ici. Nous rajoutons comme contrainte pour le rapprochement, le fait que les individus doivent être nés relativement proches les uns des autres (<5km) tout en restant dans le même département (si ce dernier est indiqué).

#' trouver_corresp_indice_00
#' Pour un numero de ligne du tableau individus
#' retourne sous forme de table les correspondances entre individus.
#' Chaque ligne du resultat correspond a l'identifiant d'une personne
#' et au plus petit identifiant d'une correspondance trouvee dans la base.
trouver_corresp_indice_00 <- function(i){
  personne_cour <- individus_3[i, ]
  
  resul <- data.table(id_personne_num = personne_cour$id_personne_num,
                      id_personne_num_groupe = personne_cour$id_personne_num)
  
  # Si la date de naissance est entiere, on ne va pas chercher plus loin
  if(personne_cour$date_N_mois == 0 | personne_cour$date_N_jour == 0){
    # Les identifiants des parents
    personne_cour_id_num_mere <- personne_cour$id_personne_num_mere
    personne_cour_id_num_pere <- personne_cour$id_personne_num_pere
    
    # Si les identifiants des deux parents sont manquants, on ne va pas plus loin
    if(!(is.na(personne_cour_id_num_mere) & is.na(personne_cour_id_num_pere))){
      # On a un identifiant pour: mere U pere
      
      id_personne_num <- personne_cour$id_personne_num
      
      # Il faut que les candidats soient de la même année et si possible du même département
      candidats <- individus_3[date_N_annee %in% personne_cour$date_N_annee]
      
      if(!is.na(personne_cour$dept_N)){
        candidats <- 
          candidats[dept_N %in% personne_cour$dept_N]
      }
      
      candidats <- candidats[!id_personne_num %in% personne_cour$id_personne_num]
      
      if(nrow(candidats)>0){
        les_distances <- 
          geosphere::distGeo(c(as.numeric(personne_cour$long_N), as.numeric(personne_cour$lat_N)),
                             cbind(as.numeric(candidats$long_N), as.numeric(candidats$lat_N)))
        candidats <- candidats[which(les_distances < 5000),]
        
        candidats <- ecremer(candidats, personne_cour = personne_cour, threshold = 0.01, prenoms = TRUE)
      }
      
      
      if(nrow(candidats)>0){
        id_retenir <- c(id_personne_num, candidats$id_personne_num)
        resul <- data.table(id_personne_num = id_retenir, id_personne_num_groupe = my_min(id_retenir))
      }# Fin du if()
    }# Fin du if sur les parents
  }
  # }
  resul
}# End of trouver_corresp_indice_00()

1.3.8 Cinquième simplification

Dans cette cinquième simplification, nous effectuons le rapprochement par les prénoms, en utilisant le prénom uniquement, et plus le nom de famille. La fonction ecremer() est remplacée par la suivante :

#' ecremer
#' Retourne un tableau de candidats possibles
#' @candidats: (data.frame) table dans laquelle chercher les candidats
ecremer_2 <- function(candidats, threshold = 0.1, personne_cour){
  # Ecremer par le nom
  ind_nom <- which(calculer_distance(personne_cour = personne_cour, candidats = candidats, variable = "nom", threshold = threshold))
  candidats <- candidats[ind_nom,]
  
  if(nrow(candidats)>0){
    # Par le prenom
    pren_cour <- personne_cour$prenom
    ind_prenom <- which(str_detect(string = candidats$prenoms,
                                   pattern = fixed(pren_cour, ignore_case = TRUE)))
    candidats <- candidats[ind_prenom,]
    
    if(nrow(candidats)>0 & !is.na(personne_cour$nom_mere)){
      ind_mere_nom <- calculer_distance(personne_cour = personne_cour, candidats = candidats, variable = "nom_mere", threshold = threshold)
      ind_mere_nom <- c(which(is.na(ind_mere_nom)), which(ind_mere_nom))
      candidats <- candidats[ind_mere_nom,]
    }# Fin if mere
    
    if(nrow(candidats)>0 & !is.na(personne_cour$nom_pere)){
      ind_pere_nom <- calculer_distance(personne_cour = personne_cour, candidats = candidats, variable = "nom_pere", threshold = threshold)
      ind_pere_nom <- c(which(is.na(ind_pere_nom)), which(ind_pere_nom)) %>% sort()
      candidats <- candidats[ind_pere_nom,]
    }# Fin if pere
    
    
  }# Fin du if
  candidats
}# End of ecremer_2()

La fonction de repérage est quant à elle modifiée comme suit. Les personnes doivent être nées le même jour, pas trop loin l’une de l’autre (<5km), leurs parents doivent avoir des noms proches.

#' trouver_corresp_indice_fifth
#' Pour un numero de ligne du tableau individus
#' retourne sous forme de table les correspondances entre individus.
#' Chaque ligne du resultat correspond a l'identifiant d'une personne
#' et au plus petit identifiant d'une correspondance trouvee dans la base.
trouver_corresp_indice_fifth <- function(i){
  personne_cour <- individus_3[i, ]
  
  resul <- data.table(id_personne_num = personne_cour$id_personne_num,
                      id_personne_num_groupe = personne_cour$id_personne_num)
  
  id_personne_num <- personne_cour$id_personne_num
  
  # Il faut que les candidats soient nes le meme jour
  ind <- which(individus_3$date_N %in% personne_cour$date_N)
  candidats <- individus_3[ind]
  
  if(!is.na(personne_cour$dept_N)){
    candidats <- 
      candidats[dept_N %in% personne_cour$dept_N]
  }
  
  candidats <- candidats[!id_personne_num %in% personne_cour$id_personne_num]
  
  if(nrow(candidats)>0){
    les_distances <- 
      geosphere::distGeo(c(as.numeric(personne_cour$long_N), as.numeric(personne_cour$lat_N)),
                         cbind(as.numeric(candidats$long_N), as.numeric(candidats$lat_N)))
    candidats <- candidats[which(les_distances < 5000),]
    
    candidats <- ecremer_2(candidats, personne_cour = personne_cour, threshold = 0.01)
  }
  
  
  if(nrow(candidats)>0){
    id_retenir <- c(id_personne_num, candidats$id_personne_num)
    resul <- data.table(id_personne_num = id_retenir, id_personne_num_groupe = my_min(id_retenir))
  }# Fin du if()
  
  
  # }
  resul
}# End of trouver_corresp_indice_fifth()

1.3.9 Sixième simplification

La dernière simplification s’appuie sur les frères et soeurs.

Nous rajoutons d’abord les informations relatives aux parents des individus.

individus_2 <- individus_simple
individus_2$id_personne_num_mere <- registre$id_personne_num_mere[match(individus_simple$id_personne_num, registre$id_personne_num)]
  individus_2$id_personne_num_pere <- registre$id_personne_num_pere[match(individus_simple$id_personne_num, registre$id_personne_num)]

Nous définissons la fonction trouver_enfants() qui comme son nom l’indique, permet de trouver les identifiants des enfants pour un individu situé à la position i du tableau individus_2.

trouver_enfants <- function(i){
  id_cour <- individus_2$id_personne_num[i]
  ind_enfants <- c(which(individus_2$id_personne_num_mere == id_cour),
                   which(individus_2$id_personne_num_pere == id_cour)) %>%
    unique()
  # ids_enfants <- individus_2$id_personne_num[ind_enfants] %>% str_c(., collapse = "|") %>%
  #   str_c("|", ., "|")
  ids_enfants <- individus_2$id_personne_num[ind_enfants]
  ids_enfants
}# Fin de trouver_enfants()

Nous appliquons cette fonction à chaque observation de la base, en parallélisant l’exécution du code.

library(parallel)
  ncl <- detectCores()-1
  (cl <- makeCluster(ncl))

  invisible(clusterEvalQ(cl, library(dplyr, warn.conflicts=FALSE, quietly=TRUE)))
  invisible(clusterEvalQ(cl, library(stringr, warn.conflicts=FALSE, quietly=TRUE)))
  invisible(clusterEvalQ(cl, library(tidyr, warn.conflicts=FALSE, quietly=TRUE)))
  invisible(clusterEvalQ(cl, library(stringi, warn.conflicts=FALSE, quietly=TRUE)))
  invisible(clusterEvalQ(cl, library(data.table, warn.conflicts=FALSE, quietly=TRUE)))
  invisible(clusterEvalQ(cl, library(stringdist, warn.conflicts=FALSE, quietly=TRUE)))
  
  clusterExport(cl, c("individus_2"), envir=environment())
  
  ids_enfants <- pbsapply(seq(1, nrow(individus_2)), trouver_enfants, cl = cl)
  # stopCluster(cl)

Notons que nous ne fermons pas la connexion aux clusters, nous allons les réutiliser un peu plus loin.

Nous rajoutons à chaque individu de la base les identifiants trouvés pour leurs enfants. Puis, nous recréons une variable pour le prénom des individus, en ne gardant que le premier prénom. “Marie Chantal Jacqueline Michu”" aura comme prénom court “Marie”.

individus_2$ids_enfants <- ids_enfants
  
individus_2 <- 
  individus_2 %>% 
  mutate(date_N_annee = str_sub(date_N, 1, 4) %>% as.numeric(),
         date_D_annee = str_sub(date_D, 1, 4) %>% as.numeric(),
         prenom_court = str_replace_all(prenoms, "\\((.*?)\\)", ""),
         prenom_court = stri_replace_all_regex(prenom_court, "'|\\(|\\)|\\?|_|\\*|&gt|\\^|\\+|[[:digit:]]|\\.|<|>", "") %>% 
           str_trim(),
         prenom_court = gsub(" .*$", "", prenom_court),
         prenom_court = stri_trans_general(prenom_court, "Latin-ASCII"),
         prenom_court = str_extract(prenom_court, '\\w*') %>% 
           stri_trim() %>% 
           str_to_lower()) %>% 
  data.table()

Nous créons ensuite la fonction same_siblings() qui permet de repérer les personnes dédoublonnées, en fonction de leurs frères et soeurs.

#' same_siblings
#' Retourne les freres/soeurs en doublons pour un individu de la base
#' i <- ind_faire[55]
#' i <- 2407
#' rm(individu,ids_parents, parents,ids_siblings,res,siblings,prenom_individu,prenom_individu,candidats, i)
same_siblings <- function(i){
  individu <- individus_2[i, ]
  # Rechercher si cette personne a des freres et soeurs
  ids_parents <- c(individu$id_personne_num_mere, individu$id_personne_num_pere)
  ids_parents <- ids_parents[!is.na(ids_parents)]
  
  parents <- individus_2[id_personne_num %in% ids_parents]
  ids_siblings <- unlist(parents$ids_enfants) %>% unique()
  ids_siblings <- ids_siblings[-which(ids_siblings == individu$id_personne_num)]
  
  res <- NULL
  
  
  # Si la personne a des freres et soeurs
  if(length(ids_siblings)>0){
    siblings <- individus_2[id_personne_num %in% ids_siblings]
    nom_individu <- individu$nom
    prenom_court_individu <- individu$prenom_court
    # Reperer parmi eux ceux qui ont le meme prenom
    candidats <- siblings[prenom_court == prenom_court_individu & nom == nom_individu]
    
    # S'il y a une dates de deces indiquee, elle doit correspondre (meme annee)
    # Si la date est manquante, on emet l'hypothese qu'il peut toujours s'agir de la meme personne
    if(nrow(candidats)>0){
      annee_D_individu <- individu$date_D_annee
      if(is.na(annee_D_individu)){
        annee_D_potentielle <- NA
      }else{
        annee_D_potentielle <- c(NA, annee_D_individu)
      }
      
      candidats <- candidats[date_D_annee %in% annee_D_potentielle]
      
      if(nrow(candidats)>0){
        annee_N_individu <- individu$date_N_annee
        if(is.na(annee_N_individu)){
          annee_N_potentielle <- NA
        }else{
          annee_N_potentielle <- c(NA,annee_N_individu)
        }
        
        candidats <- candidats[date_N_annee %in% annee_N_potentielle]
        
        # A present, il reste des freres ou soeurs avec le meme prenom et la meme annee de naissance et de deces
        if(nrow(candidats)>0){
          res <- c(individu$id_personne_num, candidats$id_personne_num)
        }
        
      }# Fin du if pour la date de naissance
      
    }# Fin du if pour la date de deces
    
    
  }# Fin de s'il y a des siblings
  
  res
}# Fin de same_siblings()

Reste à appliquer cette fonction, en utilisant à nouveau les clusters que nous n’avons pas révoqués.

ind_ne_pas_faire <- c(which(is.na(individus_2$id_personne_num_mere) & is.na(individus_2$id_personne_num_pere)))

ind_faire <- seq(1, nrow(individus_2))
ind_faire <- ind_faire[-ind_ne_pas_faire]

# Calcul en parallele
clusterExport(cl, c("individus_2"), envir=environment())
fusion_siblings <- pblapply(ind_faire, same_siblings, cl = cl)
stopCluster(cl)

# Il faut fusionner ces personnes
fusion_siblings <- fusion_siblings[which(!sapply(fusion_siblings, is.null))]

corresp_siblings <- lapply(fusion_siblings, function(x) data.table(id_personne_num = x, id_personne_num_groupe = min(x))) %>% 
  rbindlist() %>% 
  unique()

corresp_siblings <- corresp_siblings[, list(id_personne_num_groupe = min(id_personne_num_groupe, na.rm=T)), by = id_personne_num] %>% unique()

# On repercute ces nouvelles informations sur corresp
corresp <- 
    corresp %>% 
    left_join(corresp_siblings %>% rename(id_personne_num_groupe_new = id_personne_num_groupe)) %>% 
    mutate(id_personne_num_groupe = pmin(id_personne_num_groupe, id_personne_num_groupe_new, na.rm=T)) %>%
    select(-id_personne_num_groupe_new) %>% 
    data.table()

Reste ensuite à boucler sur le registre et à appliquer la fonction de simplification, comme dans les étapes précédentes.

1.4 Mise en commun des données départementales

Dans la section précédente, nous avons rassemblé les arbres par département. Il est question dans cette nouvelle section de mettre en commun les bases départementales, afin de créer une base française.

1.4.1 Chargement des données

Chargeons la liste des départements que nous avions créée auparavant et rajoutons une colonne avec les numéros du département.

load("liste_depts.rda")
liste_depts <- 
  liste_depts %>% 
  mutate(dep_num = str_sub(dept, 2))

Nous pouvons ajouter le nom de la région, en s’appuyant sur des données contenues dans le package raster.

map_fr_sp <- raster::getData(name = "GADM", country = "FRA", level = 2)
liste_depts <- 
  liste_depts %>% 
  left_join(
    data.frame(dep_num = map_fr_sp$CCA_2, dep_name = map_fr_sp$NAME_2)
  ) %>% 
  tbl_df()

Nous chargeons ensuite les données de généalogie départementales obtenues à l’étape précédente.

#' charger_corresp
#' Charger la table de correspondance, pour un département donné
#' @i_departement: (int) ligne de liste_depts à traiter
charger_corresp <- function(i_departement){
  departement <- liste_depts$dept[i_departement]
  load(str_c("../data/dept/migration/", departement, "/resul.rda"))
  resul$corresp
}# Fin de charger_corresp()

# Charger les données
corresp <- pblapply(seq(1, nrow(liste_depts)), charger_corresp)
corresp <- bind_rows(corresp)

Les individus sont identifiés à l’aide d’une valeur numérique, propre à chaque département, ainsi qu’une chaîne de caractères composée du nom d’identifiant du généalogiste amateur, du nom de famille de l’individu et de ses prénoms. Tandis qu’une personne sera désignée de manière unique par son identifiant textuel, son identidiant numérique pourra être différent en fonction du département. Nous avons en effet travaillé par département, en définissant à ce moment un identifiant pour chaque individu présent dans l’extrait de la base. Maintenant que nous réunissons les départements, nous devons créer un nouvel identifiant numérique. Prenons un exemple : l’individu Pierre Dupond présent dans l’arbre de l’utilisateur batman36 a pour mère une femme née en 1801 dans le Finistère et pour père un homme né en 1800 dans le Morbihan. Aussi, cet individu a reçu un identifiant numérique lorsque nous avons traité le Finistère, et un autre identifiant numérique lorsque nous avons traité le Morbihan. Son identifiant textuel est en revanche identique dans les deux sous-bases. Nous allons ainsi attribuer un nouvel identifiant numérique à ce Pierre Dupond, en nous appuyant sur son identifiant textuel. Nous faisons cela puisque nous allons à nouveau tenter de repérer des associations qui auraient pu nous échaper sur les sous bases, et tenter de compléter certaines informations.

corresp_id <- data.frame(id_personne = unique(corresp$id_personne)) %>% tbl_df()
corresp_id <- corresp_id %>% mutate(id_personne_num_new = row_number())

Maintenant que nous disposons de nouveaux identifiants pour les individus, nous pouvons charger les bases départementales d’individus, mettre à jour leur identifiant numérique, et joindre les bases. Nous le faisons en faisant appel à une fonction que nous appliquons à chaque département.

#' charger_individus_simple
#' Retourne une liste pour un département contenant :
#' - registre : le registre mis à jour
#' - individus_simple : les informations sur les individus
#' - corresp : le tableau de correspondances mis à jour (à chaque personne présente initialement dans la base, l'identifiant numérique correspondant)
#' - departement : le numéro de département
#' @i_departement: (int) ligne de liste_depts à consulter pour récupérer le numéro de département
charger_individus_simple <- function(i_departement){
  departement <- liste_depts$dept[i_departement]
  load(str_c("../data/dept/migration/", departement, "/resul.rda"))
  
  # TABLE DE CORRESPONDANCE #
  
  corresp <- resul$corresp %>% tbl_df()
  # Rajout des identifiants sous forme de texte dans la table de correspondance
  corresp_2 <- 
    corresp %>% 
    left_join(
      corresp %>% 
        select(id_personne, id_personne_num) %>% 
        rename(id_personne_groupe = id_personne, id_personne_num_groupe = id_personne_num) %>% 
        unique(),
      by = "id_personne_num_groupe"
    )
  corresp_2 <- corresp_2 %>% filter(!is.na(id_personne_groupe))
  
  # Rajout des nouveaux identifiants numeriques
  corresp_2 <- 
    corresp_2 %>% left_join(corresp_id %>% mutate(id_personne = as.character(id_personne)), by = "id_personne")
  
  
  corresp_2 <- 
    corresp_2 %>% 
    left_join(
      corresp_id %>% rename(id_personne_groupe = id_personne, id_personne_num_groupe_new = id_personne_num_new) %>% 
        mutate(id_personne_groupe = as.character(id_personne_groupe)),
      by = "id_personne_groupe"
    )
  
  # INDIVIDUS #
  
  individus_simple <- resul$individus_simple
  individus_simple <- individus_simple[weight != 0,]
  individus_simple <- 
    individus_simple %>% tbl_df() %>% 
    left_join(
      corresp_2 %>% 
        select(id_personne_num, id_personne_num_new),
      by = "id_personne_num"
    ) %>% 
    tbl_df()
  
  # REGISTRE DES LIENS DE PARENTE #
  
  registre <- resul$registre %>% data.table()
  registre <- registre[id_personne_num %in% individus_simple$id_personne_num] %>% 
    tbl_df()
  
  registre <- 
    registre %>% 
    left_join(corresp_2 %>% select(id_personne_num, id_personne_num_new), by = "id_personne_num") %>% 
    # Ajout mere
    left_join(corresp_2 %>% select(id_personne_num, id_personne_num_new) %>% 
                rename(id_personne_num_mere = id_personne_num, id_personne_num_mere_new = id_personne_num_new),
              by = "id_personne_num_mere") %>% 
    # Ajout pere
    left_join(corresp_2 %>% select(id_personne_num, id_personne_num_new) %>% 
                rename(id_personne_num_pere = id_personne_num, id_personne_num_pere_new = id_personne_num_new),
              by = "id_personne_num_pere")
  
  registre <- 
    registre %>% 
    select(id_personne_num_new, id_personne_num_mere_new, id_personne_num_pere_new)
  
  corresp_2 <- corresp_2 %>% 
    select(id_personne, id_personne_num_new, id_personne_num_groupe_new) %>% 
    filter(!is.na(id_personne_num_new))
  
  registre <- 
    registre %>% 
    rename(id_personne_num = id_personne_num_new,
           id_personne_num_mere = id_personne_num_mere_new,
           id_personne_num_pere = id_personne_num_pere_new)
  
  individus_simple <- 
    individus_simple %>% 
    select(-id_personne_num) %>% 
    rename(id_personne_num = id_personne_num_new)
  
  corresp_2 <- 
    corresp_2 %>% 
    rename(id_personne_num = id_personne_num_new,
           id_personne_num_groupe = id_personne_num_groupe_new)
  
  list(registre = registre, individus_simple = individus_simple, corresp = corresp_2, departement = departement)
  
}# Fin de charger_individus_simple()

res <- pblapply(seq(1, nrow(liste_depts)), charger_individus_simple)

Nous pouvons extraire et joindre entre-eux les tableaux de correspondance.

# Table des correspondances
corresp <- pblapply(res, function(x) x$corresp)
corresp <- corresp %>% rbindlist() %>% unique()

corresp <- corresp[, list(id_personne_num_groupe = min(id_personne_num_groupe, na.rm=T)), by = id_personne_num] %>% unique()
corresp_id <- 
  corresp_id %>% 
  rename(id_personne_num = id_personne_num_new) %>% 
  data.table()

corresp <- corresp_id[corresp, on = "id_personne_num"]

Idem pour les tableaux contenant les informations relatives à chaque individu.

individus_simple <- pblapply(res, function(x) x$individus_simple)
individus_simple <- individus_simple %>% rbindlist() %>% unique()
individus_simple <- individus_simple[!is.na(id_personne_num)]

1.4.2 Nouvelle simplification

Nous recherchons les associations évidentes qu’il est possible de faire avec ces nouvelles données. Pour cela, nous regardons l’identifiant textuel des individus et voyons si nous pouvons dénombrer des valeurs dupliquées.

corresp_2 <- corresp[, list(id_personne_num, id_personne_num_groupe)] %>% 
  rename(id_personne_num_groupe_new = id_personne_num_groupe,
    id_personne_num_groupe = id_personne_num)

setkey(corresp, id_personne_num_groupe)
setkey(corresp_2, id_personne_num_groupe)
corresp <- corresp_2[corresp]
corresp[, id_personne_num_groupe := pmin(id_personne_num_groupe, id_personne_num_groupe_new, na.rm=T)]
corresp$id_personne_num_groupe_new <- NULL
rm(corresp_2)

Nous faisons de nouveau appel à une fonction pour simplifier les valeurs, basée sur la fréquence observée des valeurs.

#' most_probable_value
#' Find the most probable value for a variable
#' using the values found within all candidates
#' (excluding NAs)
#' @variable_name: (string) name of the variable
#' @df: (data.table)
#' @weights: (logic) should the simplification consider the weights? (default: FALSE)
#' variable_name <- "date_N"
most_probable_value <- function(df, variable_name, weights = FALSE){
  valeurs <- df[[variable_name]]
  # valeurs <- individu_cour_a_fusionner %>% magrittr::extract2(variable_name)
  if(!all(is.na(valeurs))){
    if(weights){
      poids <- df[["weight"]]
      res <- data.frame(valeurs, poids, stringsAsFactors = FALSE) %>% 
        group_by(valeurs) %>% 
        summarise(poids = sum(poids)) %>% 
        ungroup()
      
      # Si la variable est une date
      # on va privilegier les informatins relatives aux dates completes
      if(variable_name %in% c("date_N", "date_M", "date_D")){
        tmp <- res %>% 
          mutate(mois = str_sub(valeurs, 5, 6),
                 jour = str_sub(valeurs, 7, 8)) %>% 
          filter(!(mois == "00" | jour == "00"))
        # S'il reste des informations, on se base sur elles
        if(nrow(tmp)>0) res <- tmp
      }
      
      res <- 
        res %>%
        arrange(desc(poids)) %>% 
        slice(1) %>% 
        magrittr::extract2("valeurs")
    }else{
      table_freq <- sort(table(valeurs))
      res <- names(table_freq)[1]
    }
  }else{
    res <- NA
  }
  res
}# End of most_probable_value()

#' simplifier_personne
#' À partir d'un identifiant numérique, retourne un tableau à une seule ligne contenant les valeurs de chaque variable relative à cette personne.
#' @weights: (logic) should the simplification consider the weights? (default: FALSE)
simplifier_personne <- function(un_groupe_num, weights = TRUE){
  groupe_cour <- df_corresp[id_personne_num_groupe %in% un_groupe_num]
  
  individus_cour <- df_individus_simple[id_personne_num %in% groupe_cour$id_personne_num]
  
  weight <- sum(individus_cour$weight)
  
  if(nrow(individus_cour) == 1){
    nouvelles_valeurs <- individus_cour[, mget(c(variables_names))]
  }else{
    # Les valeurs probables pour les variables d'interet
    nouvelles_valeurs <- 
      variables_names %>% 
      sapply(., most_probable_value, df = individus_cour, weights = weights) %>% 
      t() %>% 
      data.table(stringsAsFactors = FALSE)
  }
  
  
  cbind(id_personne_num = un_groupe_num, nouvelles_valeurs, weight = weight)
}# End of simplifier_personne()

Nous utilisons à nouveau la fonction my_min(), qui permet de retourner la valeur NA plutôt que Inf quand toutes les valeurs comparées sont manquantes.

my_min <- function(x) ifelse( !all(is.na(x)), min(x, na.rm=T), NA)

Les variables que nous souhaitons extraire sont les suivantes :

# Les variables qui nous interessent
variables_names <- c("nom", "prenoms",
                     "sexe",
                     "lieu_N", "dept_N", "region_N", "pays_N", "lat_N", "long_N", "date_N",
                     "lieu_M", "dept_M", "region_M", "pays_M", "lat_M", "long_M", "date_M",
                     "lieu_D", "dept_D", "region_D", "pays_D", "lat_D", "long_D", "date_D",
                     "prenom")

Puis, nous effectuons la simplification des individus, de la même manière que ce que nous avons effectué dans la section précédente.

Reste alors à sauvegarder le résultat :

res <- list(nb_obs = nb_obs,
            corresp = corresp,
            registre = registre,
            individus_simple = individus_simple)

save(res, file = "migration_1800_04_dep.rda")

1.5 Création du tableau de générations

Les étapes précédents ont permis d’obtenir une base de données de généalogie contenant des informations relatives à des individus nés en France entre 1800 et 1804, leur parents et leur descendants. Nous pouvons, à partir de ces données, créer une nouvelle base permettant d’étudier les générations. Cette section donne les étapes utilisées à cet effet et présente les codes utilisés pour générer les graphiques.

Nous commençons par charger certains packages nécessaires.

library(tidyverse)
library(lubridate)
library(stringr)
library(stringi)
library(data.table)
library(pbapply)
library(dplyrExtras)
library(stringdist)
library(sp)

Nous chargeons également la liste des départements en y ajoutant leur nom.

load("liste_depts.rda")
liste_depts <- 
  liste_depts %>% 
  mutate(dep_num = str_sub(dept, 2))

map_fr_sp <- raster::getData(name = "GADM", country = "FRA", level = 2)

liste_depts <- 
  liste_depts %>% 
  left_join(
    data.frame(dep_num = map_fr_sp$CCA_2, dep_name = map_fr_sp$NAME_2)
  ) %>% 
  tbl_df()

Nous récupérons les données de généalogie obtenues à l’issue de l’étape précédente.

load("migration_1800_04_dep.rda")
registre <- res$registre
individus_simple <- res$individus_simple
individus_2 <- individus_simple

À partir des identifiants des parents présents dans chaque ligne du tableau des données individuelles, nous pouvons rajouter à chaque individu les informations relatives à leurs parents, si toutefois elles sont existantes.

# On rajoute les identifiants des parents
individus_2$id_personne_num_mere <- registre$id_personne_num_mere[match(individus_simple$id_personne_num, registre$id_personne_num)]
individus_2$id_personne_num_pere <- registre$id_personne_num_pere[match(individus_simple$id_personne_num, registre$id_personne_num)]

ind_meres <- match(individus_2$id_personne_num_mere, individus_2$id_personne_num)
ind_peres <- match(individus_2$id_personne_num_pere, individus_2$id_personne_num)

# Les meres
individus_meres <- 
  individus_2[ind_meres,] %>% 
  magrittr::set_colnames(str_c(colnames(individus_2), "_mere")) %>% 
  select(-id_personne_num_mere, -id_personne_num_mere_mere, -id_personne_num_pere_mere)

# Les peres
individus_peres <- 
  individus_2[ind_peres,] %>% 
  magrittr::set_colnames(str_c(colnames(individus_2), "_pere")) %>% 
  select(-id_personne_num_pere, -id_personne_num_mere_pere, -id_personne_num_pere_pere)

# On ajoute les informations relatives aux parents
individus_2 <- 
  individus_2 %>% 
  cbind(individus_meres) %>% 
  cbind(individus_peres) %>% 
  # select(-id_mere, -id_pere) %>% 
  data.table()

Concentrons-nous sur des individus nés entre 1680 et 2017.

# On retire les individus nes apres 2017 et avant 1680
individus_2 <- individus_2[date_N_annee > 1680 & date_N_annee <= 2017]

1.5.1 Retrait des anomalies

Il existe des observations avec des valeurs aberrantes. Nous alons les identifier et les retirer.

Une source d’anomalie provient de l’âge auquel une personne peut avoir un enfant. Si nous constations qu’un enfant est né avant qu’un de ses parents a atteint l’âge de 13 ans, nous considérons qu’il s’agit d’une erreur de saisie (ou de fusion, ce qui est plus problématique). Nous considérons qu’un individu est en capacité de se reproduire entre ses 13 et 70 ans.

# On retire les anomalies concernant les dates de naissances
age_min_nenf <- 13
age_max_enf <- 70

individus_2 <-
  individus_2 %>% 
  mutate(date_N_annee = str_sub(date_N, 1, 4) %>% as.numeric(),
         date_D_annee = str_sub(date_D, 1, 4) %>% as.numeric(),
         date_N_annee_mere = str_sub(date_N_mere, 1, 4) %>% as.numeric(),
         date_N_annee_pere = str_sub(date_N_pere, 1, 4) %>% as.numeric()) %>% 
  # Anomalie si la mere a eu son enfant a moins de 13 ans ou a plus de 70 ans
  mutate(anomalie_mere = (date_N_annee <= (date_N_annee_mere+age_min_nenf)) | (date_N_annee >= (date_N_annee_mere + age_max_enf)),
         anomalie_mere = ifelse(is.na(anomalie_mere), yes = FALSE, no = anomalie_mere)) %>% 
  # Idem pour le pere
  mutate(anomalie_pere = (date_N_annee <= (date_N_annee_pere+age_min_nenf)) | (date_N_annee >= (date_N_annee_pere + age_max_enf)),
         anomalie_pere = ifelse(is.na(anomalie_pere), yes = FALSE, no = anomalie_pere)) %>% 
  mutate(anomalie = anomalie_mere + anomalie_pere) %>% 
  filter(anomalie == 0) %>% 
  select(-anomalie_mere, -anomalie_pere, -anomalie) %>% 
  data.table()

Nous retirons les observations pour lesquelles des âges négatifs ou supérieurs à 122 ans sont renseignés.

individus_2 <-
  individus_2 %>% 
  # mutate(age = interval(ymd(date_N), ymd(date_D)) / years(1)) %>% 
  mutate(age = date_D_annee - date_N_annee) %>% 
  mutate(keep = age >= 0 & age <= 122,
         keep = ifelse(is.na(keep), yes = TRUE, no = keep)) %>% 
  filter(keep == TRUE) %>% 
  select(-keep)

Nous retirons les observations pour lesquelles l’année de décès survient avant cette de la naissance

individus_2 <- 
  individus_2 %>% 
  mutate(keep = (str_sub(date_D, 1, 4) %>% as.numeric) >= (str_sub(date_N, 1, 4) %>% as.numeric),
         keep = ifelse(is.na(keep), yes = TRUE, no = keep)) %>% 
  filter(keep == TRUE) %>% 
  select(-keep)

Filtrons les individus de la base individus en ne retenant que ceux qui ont été obtenus dans individus_2.

individus_simple <- 
  individus_simple[id_personne_num %in% individus_2$id_personne_num]

1.5.2 Les générations

1.5.2.1 Aïeux

Nous pouvons à présent construire le tableau de générations, en commençant par les individus à l’origine de l’étude, c’est-à-dire ceux qui sont nés entre 1800 et 1804, en France.

gen_0 <- individus_2[date_N_annee %in% seq(1800, 1804)]
gen_0 <- registre[id_personne_num %in% gen_0$id_personne_num]

Nous allons construire cet tableau génération par génération. Nous définissons à cet effet un tableau contenant la génération courante, qui sera modifé quand on s’intéressera aux enfants, petits-enfants, etc. Il s’agit de l’initialiser.

gen_n <- copy(gen_0)
gen_n$generation <- 0

Nous retirons de notre tableau d’individus les individus de la génération courante.

gen_restants <- registre[id_personne_num %in% individus_2$id_personne_num]
gen_restants <- gen_restants[!id_personne_num %in% gen_n$id_personne_num]

Nous initialisons notre tableau de génération, composé pour l’instant des individus de la génération courante, i.e., la génération 0.

generations <- 
  data.frame(id_personne_num = gen_0$id_personne_num) %>% 
  rename(id_aieul = id_personne_num) %>% 
  data.table()

1.5.2.2 Enfants

Identifionts à présent les enfants de la génération courante.

peres_potentiels <- 
  gen_restants %>% 
  select(id_personne_num, id_personne_num_pere) %>% 
  rename(id_aieul = id_personne_num_pere, id_enfant = id_personne_num) %>% 
  filter(!is.na(id_aieul))

meres_potentiels <- 
  gen_restants %>% 
  select(id_personne_num, id_personne_num_mere) %>% 
  rename(id_aieul = id_personne_num_mere, id_enfant = id_personne_num) %>% 
  filter(!is.na(id_aieul))

parents_potentiels <- 
  peres_potentiels %>% 
  rbind(meres_potentiels)

setkey(parents_potentiels, id_aieul)
setkey(generations, id_aieul)

generations <- parents_potentiels[generations]

1.5.2.3 Petits-enfants

Puis les petits-enfants.

peres_potentiels <- 
  gen_restants %>% 
  select(id_personne_num, id_personne_num_pere) %>% 
  rename(id_enfant = id_personne_num_pere, id_p_enfant = id_personne_num) %>% 
  filter(!is.na(id_enfant))

meres_potentiels <- 
  gen_restants %>% 
  select(id_personne_num, id_personne_num_mere) %>% 
  rename(id_enfant = id_personne_num_mere, id_p_enfant = id_personne_num) %>% 
  filter(!is.na(id_enfant))

parents_potentiels <- 
  peres_potentiels %>% 
  rbind(meres_potentiels)


setkey(parents_potentiels, id_enfant)
setkey(generations, id_enfant)

generations <- parents_potentiels[generations]

1.5.2.4 Arrière-petits-enfants

Et enfin, les arrière-petits-enfants.

peres_potentiels <- 
  gen_restants %>% 
  select(id_personne_num, id_personne_num_pere) %>% 
  rename(id_p_enfant = id_personne_num_pere, id_ap_enfant = id_personne_num) %>% 
  filter(!is.na(id_p_enfant))

meres_potentiels <- 
  gen_restants %>% 
  select(id_personne_num, id_personne_num_mere) %>% 
  rename(id_p_enfant = id_personne_num_mere, id_ap_enfant = id_personne_num) %>% 
  filter(!is.na(id_p_enfant))

parents_potentiels <- 
  peres_potentiels %>% 
  rbind(meres_potentiels)

setkey(parents_potentiels, id_p_enfant)
setkey(generations, id_p_enfant)

generations <- parents_potentiels[generations]

1.5.3 Conversion des départements en polygones

Nous avons besoin d’identifier dans quel département tombent certains points identifiés par des coordonnées (longitude,latitude). Pour cela, une manière simple consiste à utiliser la fonction inside.gpc.poly() du package surveillance. Il faut toutefois fournir un polygone de la classe gpc.poly(). Pour l’heure, nous disposons de données pour tracer les départements sous forme de tableaux de données. Nous convertissons de fait ces tableaux en objets gpc.poly().

library(rgeos)
#' group_to_gpc_poly
#' Returns the polygon from a group within the data frame containing
#' the coordinates of regions
#' @df: data frame with the coordinates of the polygons
#' @group: (string) group ID
group_to_gpc_poly <- function(df, group){
  group_coords <- df[df$group == group, c("long", "lat")]
  as(group_coords, "gpc.poly")
}

#' departement_polygone
#' Returns a French region as a gpc.poly object
#' un_dep: (string) code of the region
departement_polygone <- function(un_dep){
  dep_tmp <-
    france_1800 %>%
    filter(departement == un_dep)
  df_tmp_polys <- lapply(unique(dep_tmp$group), group_to_gpc_poly, df = dep_tmp)
  df_tmp_polys_tot <- df_tmp_polys[[1]]
  if(length(df_tmp_polys)>1){
    for(i in 2:length(df_tmp_polys)){
      df_tmp_polys_tot <- gpclib::union(df_tmp_polys_tot,  df_tmp_polys[[i]])
    }
  }
  df_tmp_polys_tot
}# End of departement_polygone()

france_1800_gpc <- lapply(unique(france_1800$departement), departement_polygone)
names(france_1800_gpc) <- unique(france_1800$departement)
save(france_1800_gpc, file = "france_1800_gpc.rda")

Maintenant que nous disposons des départements sous forme de gpc.poly, nous pouvons créer une fonction permettant d’indiquer si un point (lon,lat) se situe à l’intérieur d’un département donné.

est_dans_poly_naissance <- function(longitude, latitude, departement){
  surveillance::inside.gpc.poly(x = longitude,
                                y = latitude,
                                departement,
                                mode.checked = FALSE)
}

1.5.4 Repérage des départements de naissance des individus

Pour repérer le département de naissance des individus, nous convertissons les coordonnées en valeurs numériques. Elles sont pour l’heure en chaînes de caractères.

individus_simple_reste <- individus_simple
individus_simple_reste <- individus_simple_reste[!is.na(lat_N)]
individus_simple_reste <- 
  individus_simple_reste %>% 
  mutate(long_N = as.numeric(long_N),
         lat_N = as.numeric(lat_N))

individus_simple_reste <- data.table(individus_simple_reste)

Nous allons effectuer une boucle sur les départements pour déterminer, à chaque itération, quels sont les individus nés dans ce département. Comme nous avons un nombre conséquent d’invididus et que le territoire français est relativement grand à l’échelle des départements, nous procédons par un filtrage grossier de nos données à chaque tour. Pour un département donné, nous récupérons les coordonnées d’un rectangle le contenant, à l’aide de la fonction get.bbox(). Nous filtrons ensuite les individus en fonction de ce rectangle, pour écarter ceux dont les coordonnées sont à l’extérieur. Reste alors à appliquer la fonction est_dans_poly_naissance() sur les points filtrés, qui sont bien moins nombreux qu’à l’origine.

dept_n_indiv <- rep(list(NA), length(france_1800_gpc))

pb <- txtProgressBar(min = 1, max = length(france_1800_gpc), style = 3)
for(i in 1:length(france_1800_gpc)){
  poly <- france_1800_gpc[[i]]
  nom_dep <- names(france_1800_gpc)[i]
  
  individus_simple_reste_2 <- 
    individus_simple_reste %>% 
    select(long_N, lat_N, id_personne_num) %>% 
    filter(long_N >= get.bbox(poly)$x[1],
           long_N <= get.bbox(poly)$x[2],
           lat_N >= get.bbox(poly)$y[1],
           lat_N <= get.bbox(poly)$y[2]) %>% 
    rowwise() %>% 
    mutate(in_poly = est_dans_poly_naissance(long_N, lat_N, poly)) %>% 
    select(id_personne_num, in_poly) %>% 
    mutate(dept_N = nom_dep) %>% 
    filter(in_poly == TRUE) %>% 
    select(-in_poly)
  
  
  individus_simple_reste <- 
    individus_simple_reste[! id_personne_num %in% individus_simple_reste_2$id_personne_num]
  
  dept_n_indiv[[i]] <- individus_simple_reste_2
  setTxtProgressBar(pb, i)
  rm(individus_simple_reste_2)
  # On appelle le garbage collector toutes les 10 itérations
  if(i %% 10 == 0) gc()
}

dept_n_indiv <- 
  dept_n_indiv %>% 
  rbindlist()

Reste alors à mettre à jour individus_simple pour répercuter les informations obtenues sur le département de naissance des individus.

dept_n_indiv <- dept_n_indiv %>% rename(dept_N_2 = dept_N)

setkey(dept_n_indiv, id_personne_num)
setkey(individus_simple, id_personne_num)

individus_simple <- dept_n_indiv[individus_simple]
head(individus_simple)

individus_simple <- 
  individus_simple %>% 
  mutate(dept_N = ifelse(!is.na(dept_N_2), yes = dept_N_2, no = dept_N))

individus_simple <- data.table(individus_simple)

1.5.5 Rajout de l’âge des individus

Nous pouvons recalculer l’âge des individus.

individus_simple[, date_N_annee := as.numeric(str_sub(date_N, 1, 4))]
individus_simple[, date_D_annee := as.numeric(str_sub(date_D, 1, 4))]
individus_simple[, age := as.numeric(date_D_annee) - as.numeric(date_N_annee)]

1.5.6 Toutes les informations sur une même ligne

Pour notre tableau de génération, chaque ligne doit permettre de suivre un descendant. Chaque ligne doit contenir les informations d’une personne et de ses aïeux jusqu’à celui né entre 1800 et 1804.

var_recup <- c("id_personne_num", "long_N", "lat_N", "sexe", "date_N", "date_D", "date_M", "dept_N", "date_N_annee", "date_D_annee", "age")

infos_aieul <- 
  individus_simple %>% 
  s_select(var_recup) %>% 
  magrittr::set_colnames(str_c(var_recup, "_aieul")) %>% 
  rename(id_aieul = id_personne_num_aieul)

infos_enfants <- 
  individus_simple %>% 
  s_select(var_recup) %>% 
  magrittr::set_colnames(str_c(var_recup, "_enfants")) %>% 
  rename(id_enfant = id_personne_num_enfants)

infos_p_enfant <- 
  individus_simple %>% 
  s_select(var_recup) %>% 
  magrittr::set_colnames(str_c(var_recup, "_p_enfants")) %>% 
  rename(id_p_enfant = id_personne_num_p_enfants)

infos_ap_enfant <- 
  individus_simple %>% 
  s_select(var_recup) %>% 
  magrittr::set_colnames(str_c(var_recup, "_ap_enfants")) %>% 
  rename(id_ap_enfant = id_personne_num_ap_enfants)

setkey(infos_aieul, id_aieul)
setkey(infos_enfants, id_enfant)
setkey(infos_p_enfant, id_p_enfant)
setkey(infos_ap_enfant, id_ap_enfant)


setkey(generations, id_aieul)
generations <- infos_aieul[generations]

setkey(generations, id_enfant)
generations <- infos_enfants[generations]

setkey(generations, id_p_enfant)
generations <- infos_p_enfant[generations]

setkey(generations, id_ap_enfant)
generations <- infos_ap_enfant[generations]

rm(infos_aieul, infos_enfants, infos_p_enfant, infos_ap_enfant)

generations

1.5.7 Distances avec l’aieul

En fonction des coordonnées GPS fournies, nous sommes à même de calculer la distance séparant un enfant du lieu de naissance de ses aïeux. Passons d’abord toutes les variables de coordonnées en valeurs numériquues.

library(geosphere)
generations <- 
  generations %>% 
  mutate(long_N_aieul = as.numeric(long_N_aieul),
         lat_N_aieul = as.numeric(lat_N_aieul),
         long_N_enfants = as.numeric(long_N_enfants),
         lat_N_enfants = as.numeric(lat_N_enfants),
         long_N_p_enfants = as.numeric(long_N_p_enfants),
         lat_N_p_enfants = as.numeric(lat_N_p_enfants),
         long_N_ap_enfants = as.numeric(long_N_ap_enfants),
         lat_N_ap_enfants = as.numeric(lat_N_ap_enfants))

generations <- data.table(generations)

Nous allons utiliser la fonction distm() du package geosphere pour calculer les distances (en mètres) entre deux points.

calcul_distance <- function(long_1, lat_1, long_2, lat_2){
  distm(c(long_1, lat_1), c(long_2, lat_2), fun = distHaversine) %>% "["(1)
}

Pour chaque ligne du tableau de génération, nous calculons les trois distances suivantes :

1. dist_enfants : distance entre un enfant et l’aïeul de 1800 ;
2. dist_p_enfants : distance entre un petit-enfant et l’aïeul de 1800 ;
3. dist_ap_enfants : distance entre un arrière-petit-enfant et l’aïeul de 1800.

distances_l <- 
  pblapply(1:nrow(generations), function(y){
      x <- generations[y,]
      dist_enfants = calcul_distance(x$long_N_aieul, x$lat_N_aieul, x$long_N_enfants, x$lat_N_enfants)
      dist_p_enfants = calcul_distance(x$long_N_aieul, x$lat_N_aieul, x$long_N_p_enfants, x$lat_N_p_enfants)
      dist_ap_enfants = calcul_distance(x$long_N_aieul, x$lat_N_aieul, x$long_N_ap_enfants, x$lat_N_ap_enfants)
      data.frame(dist_enfants = dist_enfants,
                 dist_p_enfants = dist_p_enfants,
                 dist_ap_enfants = dist_ap_enfants)
    }) %>% 
    rbindlist() %>% 
    data.table()

Il reste à rajouter ces informations au tableau generation.

generations <- 
  generations %>% 
  cbind(distances_l)

Nous pouvons convertir ces distances en kilomètres.

generations <- 
  generations %>% 
  mutate(dist_enfants_km =  dist_enfants / 1000,
         dist_p_enfants_km = dist_p_enfants / 1000,
         dist_ap_enfants_km = dist_ap_enfants / 1000)

Effectuons quelques regroupements de départements, pour refléter l’état approximatif des départements au XIXe siècle.

liste_depts <- 
  liste_depts %>% 
  mutate(dept2 = dept,
         dept2 = ifelse(dept2 %in% c("F78", "F91", "F92", "F95"), yes = "F78", no = dept2),
         dept2 = ifelse(dept2 %in% c("F93", "F94", "F77"), yes = "F77", no = dept2)
  )

Faisons de même sur la base des générations.

generations <- 
  generations %>% 
  mutate(dept_N_aieul = ifelse(dept_N_aieul %in% c("F78", "F91", "F92", "F95"), yes = "F78", no = dept_N_aieul),
         dept_N_enfants = ifelse(dept_N_enfants %in% c("F78", "F91", "F92", "F95"), yes = "F78", no = dept_N_enfants),
         dept_N_p_enfants = ifelse(dept_N_p_enfants %in% c("F78", "F91", "F92", "F95"), yes = "F78", no = dept_N_p_enfants),
         dept_N_ap_enfants = ifelse(dept_N_ap_enfants %in% c("F78", "F91", "F92", "F95"), yes = "F78", no = dept_N_ap_enfants),
         #
         dept_N_aieul = ifelse(dept_N_aieul %in% c("F93", "F94", "F77"), yes = "F77", no = dept_N_aieul),
         dept_N_enfants = ifelse(dept_N_enfants %in% c("F93", "F94", "F77"), yes = "F77", no = dept_N_enfants),
         dept_N_p_enfants = ifelse(dept_N_p_enfants %in% c("F93", "F94", "F77"), yes = "F77", no = dept_N_p_enfants),
         dept_N_ap_enfants = ifelse(dept_N_ap_enfants %in% c("F93", "F94", "F77"), yes = "F77", no = dept_N_ap_enfants)
  )

generations <- data.table(generations)

1.6 Adelphies

Nous pouvons extraire des informations relatives aux adelphies (fratries et sororités).

Le nombre d’enfants par femme de la génération 1 (enfants).

nb_enf_aieux <- 
  generations %>% 
  filter(sexe_aieul == 2) %>% 
  mutate(keep = TRUE) %>% 
  mutate(keep = ifelse(age_aieul < 13, yes = FALSE, no = keep)) %>% 
  mutate(keep = ifelse(is.na(dept_M_aieul), yes = FALSE, no = keep)) %>% 
  filter(keep) %>% 
  select(id_aieul, id_enfant, dept_N_aieul) %>% 
  unique() %>% 
  mutate(nb_enf = ifelse(is.na(id_enfant), yes = 0, no = 1)) %>% 
  group_by(id_aieul, dept_N_aieul) %>% 
  summarise(nb_enf = sum(nb_enf))

summary(nb_enf_aieux$nb_enf)

Le nombre d’enfants par femme de la generation 2 (petits enfants).

nb_enf_enfants <- 
  generations %>% 
  filter(sexe_enfants == 2) %>% 
  mutate(keep = TRUE) %>% 
  mutate(keep = ifelse(age_enfants < 13, yes = FALSE, no = keep)) %>% 
  mutate(keep = ifelse(is.na(dept_M_enfant), yes = FALSE, no = keep)) %>% 
  filter(keep) %>% 
  select(id_enfant, id_p_enfant, dept_N_p_enfants) %>% 
  unique() %>% 
  mutate(nb_enf = ifelse(is.na(id_p_enfant), yes = 0, no = 1)) %>% 
  group_by(id_enfant, dept_N_p_enfants) %>% 
  summarise(nb_enf = sum(nb_enf))

summary(nb_enf_enfants$nb_enf)

Le nombre d’enfants par femme de la generation 3 (arrieres petits enfants).

nb_enf_p_enfants <- 
  generations %>% 
  filter(sexe_p_enfants == 2) %>% 
  mutate(keep = TRUE) %>% 
  mutate(keep = ifelse(age_p_enfants < 13, yes = FALSE, no = keep)) %>% 
  mutate(keep = ifelse(is.na(dept_M_p_enfant), yes = FALSE, no = keep)) %>% 
  filter(keep) %>% 
  select(id_p_enfant, id_ap_enfant, dept_N_ap_enfants) %>% 
  unique() %>% 
  mutate(nb_enf = ifelse(is.na(id_ap_enfant), yes = 0, no = 1)) %>% 
  group_by(id_p_enfant, dept_N_ap_enfants) %>% 
  summarise(nb_enf = sum(nb_enf))

summary(nb_enf_p_enfants$nb_enf)

1.7 Base des individus mobiles

Nous nous intéressons dans l’étude à la mobilité des individus, entre les générations. Il faut de fait se concentrer sur des individus en capacité de se déplacer. Nous plaçons une limite inférieure arbitraire sur les âges pour nos individus, de 16 ans. Nous émettons l’hypothèse qu’avant 16 ans, les individus ne changeront pas de village.

generations_mobilite <- 
  generations %>% 
  filter(age_aieul >= 16,
         age_enfants >= 16,
         age_p_enfants >= 16,
         age_ap_enfants >= 16) %>% 
  data.table()

Enfin, nous sauvegardons les deux bases de données generations et generations_mobilite.

save(generations, file = "generations_migration.rda")
save(generations_mobilite, file = "generations_mobilite.rda")

2 Représentativité des données

Maintenant, nous disposons d’une base nettoyée, bien qu’elle contienne encore sûrement des doublons que les exécutions de nos algorithmes n’ont pas su déceler. Nous pouvons alors l’explorer. Nous nous intéressons dans cette section à quelques comparaisons de statistiques descriptives des données de généalogie avec certains résultats de la littérature.

2.1 Naissances

Un premier élément permettant de caractériser notre échantillon est sa taille. Nous regardons le nombre de naissances par année en différentiant les générations. Puis, nous comparons le nombre de naissances de l’échantillon initial avec des données officielles, pour avoir une idée de la représentativité de la génération à l’origine de notre étude.

2.1.1 Nombre de naissances par génération

Nous chargeons quelques packages, ainsi que les données de généalogie.

library(tidyverse)
library(stringr)
library(data.table)
library(pbapply)
library(dplyrExtras)
library(sp)

load("generations_migration.rda")

Pour chacune des générations, nous calculons le nombre de naissance par année.

nb_nais_aieux <- 
  generations %>% 
  tbl_df() %>% 
  select(id_aieul, date_N_aieul) %>% 
  unique() %>% 
  mutate(date_N_aieul_annee = str_sub(date_N_aieul, 1, 4) %>% as.numeric()) %>% 
  group_by(date_N_aieul_annee) %>% 
  tally() %>% 
  arrange(date_N_aieul_annee)

nb_nais_enfants <- 
  generations %>% 
  tbl_df() %>% 
  select(id_enfant, date_N_enfants) %>% 
  unique() %>% 
  mutate(date_N_enfants_annee = str_sub(date_N_enfants, 1, 4) %>% as.numeric()) %>% 
  group_by(date_N_enfants_annee) %>% 
  tally() %>% 
  arrange(date_N_enfants_annee)

nb_nais_p_enfants <- 
  generations %>% 
  tbl_df() %>% 
  select(id_p_enfant, date_N_p_enfants) %>% 
  unique() %>% 
  mutate(date_N_p_enfants_annee = str_sub(date_N_p_enfants, 1, 4) %>% as.numeric()) %>% 
  group_by(date_N_p_enfants_annee) %>% 
  tally() %>% 
  arrange(date_N_p_enfants_annee)

nb_nais_ap_enfants <- 
  generations %>% 
  tbl_df() %>% 
  select(id_ap_enfant, date_N_ap_enfants) %>% 
  unique() %>% 
  mutate(date_N_ap_enfants_annee = str_sub(date_N_ap_enfants, 1, 4) %>% as.numeric()) %>% 
  group_by(date_N_ap_enfants_annee) %>% 
  tally() %>% 
  arrange(date_N_ap_enfants_annee)

Nous mettons en commun ces observations.

nb_nais <- 
  nb_nais_aieux %>% 
  rename(annee_N = date_N_aieul_annee) %>% 
  mutate(generation = 0) %>% 
  bind_rows(nb_nais_enfants %>% 
              rename(annee_N = date_N_enfants_annee) %>% 
              mutate(generation = 1)
  ) %>% 
  bind_rows(
    nb_nais_p_enfants %>% 
      rename(annee_N = date_N_p_enfants_annee) %>% 
      mutate(generation = 2)
  ) %>% 
  bind_rows(
    
    nb_nais_ap_enfants %>% 
      rename(annee_N = date_N_ap_enfants_annee) %>% 
      mutate(generation = 3)
  ) %>% 
  mutate(generation = factor(generation, levels = seq(0,3), labels = c('Aïeux', "Enfants",
                                                                       "Petits-enfants", "Arrière-petits-enfants")))

Puis, nous pouvons tracer un histogramme, avec une échelle log pour les ordonnées, le nombre de naissance d’aïeux étant relativement explosif vis-à-vis du nombre de naissances annuelles pour les descendants.

p_distrib_naissances_gen <- 
  ggplot(data = nb_nais, aes(x = annee_N, y = log(n))) +
  geom_bar(stat="identity", position = "identity", aes(fill =generation , group = generation), alpha = .5) +
  scale_x_continuous(breaks = scales::pretty_breaks(n = 10), limits = c(1800, 1940)) +
  xlab("Nombre d'individus") + ylab("Année de naissance") +
  scale_fill_discrete("")

p_distrib_naissances_gen

2.1.2 Comparaison des nombres de naissances avec les données de l’INSEE

Pour avoir une idée de la représentativité des données, nous regardons, par département, combien de naissances sont renseignées dans nos données et comparons le résultat avec les enregistrements de la Statistique Générale de la France (SGF) proposés sur le site de l’INSEE.

Nous chargeons nos données ainsi que celles de la SGF Nous allons nous intéresser seulement à quelques variables parmi les 228 que propose le fichier de la SGF.

load("generations_migration.rda")
pop_insee <- readxl::read_excel("../data/population_insee/TERR_T86.xls", skip = 8, col_names = FALSE)
pop_insee_col <- readxl::read_excel("../data/population_insee/TERR_T86.xls", skip = 6) %>% 
colnames()
colnames(pop_insee) <- pop_insee_col
# V3 : Numéro de département
# V6 : Nom de l'unité d'analyse
# V221 : Naissances totales; légitimes et naturels; garçons, 1800 à 1801
# V222 : Naissances totales; légitimes et naturels: filles 1800 à 1801

Une petite correction de numéro de département, pour pouvoir apparier nos données à celles de la SGF.

pop_insee <- 
pop_insee %>% 
mutate(V3 = ifelse(V3 == "99", yes = "54", no = V3))

Nous retenons uniquement les variables qui nous intéressent.

nb_nais_region <- 
pop_insee %>% 
select(V3, V221, V222) %>% 
rename(num_dep = V3, nb_nais_f = V222, nb_nais_h = V221) %>% 
filter(num_dep != "00")

head(nb_nais_region)

# A tibble: 6 x 3
  num_dep nb_nais_h nb_nais_f
  <chr>       <dbl>     <dbl>
1 01           4892      4811
2 02           6801      6534
3 03           4108      3883
4 04           2813      2739
5 05           2070      1945
6 07           4496      4288

Nous sommons les naissances pour les femmes et les hommes.

nb_nais_1801 <- 
nb_nais_region %>% 
summarise(nb_nais_h = sum(nb_nais_h, na.rm=T), nb_nais_f = sum(nb_nais_f, na.rm=T))

head(nb_nais_1801)

# A tibble: 1 x 2
  nb_nais_h nb_nais_f
      <dbl>     <dbl>
1    464562    439126

Nous calculons à présent le nombre de naissances de filles et de garçons pour chaque génération, pour chaque année, dans les données de généalogie.

nb_naiss_par_annee <- 
generations %>% 
tbl_df() %>% 
select(id_aieul, sexe_aieul, date_N_annee_aieul) %>% 
unique() %>% 
rename(id_personne = id_aieul, sexe = sexe_aieul, date_N_annee = date_N_annee_aieul) %>% 
bind_rows(
generations %>% 
tbl_df() %>% 
select(id_enfant, sexe_enfants, date_N_annee_enfants) %>% 
unique() %>% 
rename(id_personne = id_enfant, sexe = sexe_enfants, date_N_annee = date_N_annee_enfants)
) %>% 
bind_rows(
generations %>% 
tbl_df() %>% 
select(id_p_enfant, sexe_p_enfants, date_N_annee_p_enfants) %>% 
unique() %>% 
rename(id_personne = id_p_enfant, sexe = sexe_p_enfants, date_N_annee = date_N_annee_p_enfants)
) %>% 
bind_rows(
generations %>% 
tbl_df() %>% 
select(id_ap_enfant, sexe_ap_enfants, date_N_annee_ap_enfants) %>% 
unique() %>% 
rename(id_personne = id_ap_enfant, sexe = sexe_ap_enfants, date_N_annee = date_N_annee_ap_enfants)
) %>% 
unique() %>% 
group_by(sexe, date_N_annee) %>% 
tally()

L’idée est ensuite de construire un tableau (representativite) uniquement pour les aïeux nés en 1801, pour comparer avec les valeurs de la SGF. Le tableau que nous construisons indique, pour chaque département, le nombre de naissances enregistrées dans les données de généalogie, pour les femmes et pour les hommes. Nous rajoutons ensuite les données de la SGF contenues dans le tableau nb_nais_region. Il s’agit ensuite de comparer les valeurs, en faisant un ratio du nombre de naissances dans les données de généalogies sur le nombre de naissances reporté par la SGF.

representativite <- 
generations %>% 
tbl_df() %>% 
filter(date_N_annee_aieul == 1801) %>% 
select(id_aieul, sexe_aieul, date_N_annee_aieul, dept_N_aieul) %>% 
unique() %>% 
rename(dept_N = dept_N_aieul,
date_N_annee = date_N_annee_aieul,
sexe = sexe_aieul) %>% 
group_by(sexe, date_N_annee, dept_N) %>% 
summarise(nb_nais = n()) %>% 
left_join(
nb_nais_region %>% gather(sexe, nb_nais_1801, -num_dep) %>% 
mutate(sexe = ifelse(sexe == "nb_nais_h", yes = "1", no = "2")) %>% 
mutate(dept_N = str_c("F", num_dep))
) %>% 
filter(!is.na(sexe)) %>% 
mutate(pct_nais = nb_nais / nb_nais_1801 * 100,
pct_nais = round(pct_nais, 2)) %>% 
ungroup() %>% 
mutate(sexe = factor(sexe, levels = c(2,1), labels = c("Femmes", "Hommes"))) %>% 
arrange(date_N_annee, sexe)

head(representativite)

# A tibble: 6 x 7
  sexe   date_N_annee dept_N nb_nais num_dep nb_nais_1801 pct_nais
  <fct>         <dbl> <chr>    <int> <chr>          <dbl>    <dbl>
1 Femmes         1801 F01       1332 01              4811     27.7
2 Femmes         1801 F02       2341 02              6534     35.8
3 Femmes         1801 F03       1659 03              3883     42.7
4 Femmes         1801 F04        668 04              2739     24.4
5 Femmes         1801 F05        842 05              1945     43.3
6 Femmes         1801 F06        779 <NA>              NA     NA  

Nous retirons les départements pour lesquels la SGF ne communique pas de données.

representativite <- 
representativite %>% 
filter(!is.na(dept_N))

Nous extrayons les numéros du département.

representativite <- 
representativite %>% 
select(dept_N, sexe, nb_nais, nb_nais_1801, pct_nais) %>% 
mutate(num_dep = str_sub(dept_N, 2))

head(representativite)

# A tibble: 6 x 6
  dept_N sexe   nb_nais nb_nais_1801 pct_nais num_dep
  <chr>  <fct>    <dbl>        <dbl>    <dbl> <chr>  
1 F01    Femmes    1332         4811     27.7 01     
2 F01    Hommes    1474         4892     30.1 01     
3 F02    Femmes    2341         6534     35.8 02     
4 F02    Hommes    2836         6801     41.7 02     
5 F03    Femmes    1659         3883     42.7 03     
6 F03    Hommes    1866         4108     45.4 03     

Nous pouvons ensuite proposer une représentation graphique sous forme de carte choroplèthe la représentativité régionale de notre échantillon d’aïeux. Pour ce faire, nous chargons les données permettant de tracer la France. Plus de détails sur l’obtention de ce fichier sont donnés en annexe.

load("france_1800.rda")

representativite_map <- 
france_1800 %>% 
left_join(
representativite %>% 
filter(!is.na(sexe)) %>% 
select(num_dep, sexe, pct_nais) %>% 
spread(sexe, pct_nais)
)

p_repres <- ggplot() +
geom_polygon(data = representativite_map %>% gather(sexe, pct_nais, Femmes, Hommes),
aes(x = long, y = lat, group = group, fill = pct_nais), size = .4) +
coord_quickmap() +
scale_fill_gradientn("%", colours = rev(heat.colors(10)), na.value = "grey") +
facet_wrap(~sexe)

p_repres

2.2 Ratio Femmes/Hommes

On peut également calculer le pourcentages d’hommes à partir des données de généalogie de la manière suivante : on créé un tableau avec tous les individus uniques de l’échantillon, puis ont somme le nombre d’individus par sexe.

sexes <- generations %>% 
  select(id_aieul, sexe_aieul, dept_N_aieul) %>% 
  unique() %>% 
  rename(id_personne = id_aieul, sexe = sexe_aieul, dept_N = dept_N_aieul) %>% 
  mutate(sexe = as.character(sexe)) %>% 
  bind_rows(
    generations %>% 
      select(id_enfant, sexe_enfants, dept_N_enfants) %>% 
      unique() %>% 
      rename(id_personne = id_enfant, sexe = sexe_enfants, dept_N = dept_N_enfants) %>% 
      mutate(sexe = as.character(sexe))
  ) %>% 
  bind_rows(
    generations %>% 
      select(id_p_enfant, sexe_p_enfants, dept_N_p_enfants) %>% 
      unique() %>% 
      rename(id_personne = id_p_enfant, sexe = sexe_p_enfants, dept_N = dept_N_p_enfants) %>% 
      mutate(sexe = as.character(sexe))
  ) %>% 
  bind_rows(
    generations %>% 
      select(id_ap_enfant, sexe_ap_enfants,dept_N_ap_enfants ) %>% 
      unique() %>% 
      rename(id_personne = id_ap_enfant, sexe = sexe_ap_enfants, dept_N = dept_N_ap_enfants) %>% 
      mutate(sexe = as.character(sexe))
  ) %>% 
  unique()

sexes %>% 
  tbl_df() %>% 
  group_by(sexe) %>% 
  tally() %>% 
  ungroup() %>% 
  filter(!is.na(sexe)) %>% 
  mutate(pct = n/sum(n)) %>% 
  head()

# A tibble: 2 x 3
  sexe        n   pct
  <chr>   <int> <dbl>
1 1     1316749 0.538
2 2     1132795 0.462

On peut également calculer le taux de masculinité :

r sexes %>% tbl_df() %>% group_by(sexe) %>% tally() %>% ungroup() %>% filter(!is.na(sexe)) %>% spread(sexe, n) %>% mutate(sex_ratio = `1` / `2` * 100)

# A tibble: 1 x 3 `1` `2` sex_ratio <int> <int> <dbl> 1 1316749 1132795 116

Il n’est pas compliqué de calculer le pourcentage d’hommes par départements :

```r representativite_fh <- sexes %>% filter(str_detect(dept_N, “^F[[:digit:]]{2}”)) %>% filter(!is.na(sexe)) %>% group_by(dept_N, sexe) %>% summarise(n = n()) %>% group_by(dept_N) %>% mutate(pct = n / sum(n)) %>% rename(departement = dept_N)

representativite_fh_sex_ratio <- sexes %>% filter(str_detect(dept_N, “^F[[:digit:]]{2}”)) %>% filter(!is.na(sexe)) %>% group_by(dept_N, sexe) %>% summarise(n = n()) %>% spread(sexe, n) %>% mutate(sex_ratio = 1 / 2 * 100) ```

On peut alors préparer un tableau de données pour effectuer une représentation graphique :

r representativite_map_fh <- france_1800 %>% left_join( representativite_fh %>% filter(sexe == 1) )

p_repres_fh <- ggplot() +
  geom_polygon(data = representativite_map_fh,
               aes(x = long, y = lat, group = group, fill = pct), col = "black", size = .4) +
  coord_quickmap() +
  scale_fill_gradientn("%", colours = rev(heat.colors(10)), na.value = "grey")
p_repres_fh

2.3 Taux de fécondité

Nous allons nous intéresser au taux de fécondité. Aussi, nous devons calculer le nombre d’enfants par femme. Nous le faisons pour chacune des générations.

# Nombre d'enfants par femme : generation 1 (enfants)
nb_enf_aieux <- 
  generations %>% 
  tbl_df() %>% 
  filter(sexe_aieul == 2) %>% 
  mutate(keep = TRUE) %>% 
  mutate(keep = ifelse(age_aieul < 15, yes = FALSE, no = keep)) %>% 
  mutate(keep = ifelse(is.na(dept_M_aieul), yes = FALSE, no = keep)) %>% 
  filter(keep) %>% 
  select(id_aieul, id_enfant, dept_N_aieul) %>% 
  unique() %>% 
  mutate(nb_enf = ifelse(is.na(id_enfant), yes = 0, no = 1)) %>% 
  group_by(id_aieul, dept_N_aieul) %>% 
  summarise(nb_enf = sum(nb_enf))

# Nombre d'enfants par femme : generation 2 (petits enfants)
nb_enf_enfants <- 
  generations %>% 
  filter(sexe_enfants == 2) %>% 
  mutate(keep = TRUE) %>% 
  mutate(keep = ifelse(age_enfants < 15, yes = FALSE, no = keep)) %>% 
  mutate(keep = ifelse(is.na(dept_M_enfants), yes = FALSE, no = keep)) %>% 
  filter(keep) %>% 
  select(id_enfant, id_p_enfant, dept_N_p_enfants) %>% 
  unique() %>% 
  mutate(nb_enf = ifelse(is.na(id_p_enfant), yes = 0, no = 1)) %>% 
  group_by(id_enfant, dept_N_p_enfants) %>% 
  summarise(nb_enf = sum(nb_enf))

# Nombre d'enfants par femme : generation 3 (arrieres petits enfants)
nb_enf_p_enfants <- 
  generations %>% 
  filter(sexe_p_enfants == 2) %>% 
  mutate(keep = TRUE) %>% 
  mutate(keep = ifelse(age_p_enfants < 15, yes = FALSE, no = keep)) %>% 
  mutate(keep = ifelse(is.na(dept_M_p_enfants), yes = FALSE, no = keep)) %>% 
  filter(keep) %>% 
  select(id_p_enfant, id_ap_enfant, dept_N_ap_enfants) %>% 
  unique() %>% 
  mutate(nb_enf = ifelse(is.na(id_ap_enfant), yes = 0, no = 1)) %>% 
  group_by(id_p_enfant, dept_N_ap_enfants) %>% 
  summarise(nb_enf = sum(nb_enf))

summary(nb_enf_aieux$nb_enf)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  0.000   0.000   1.000   1.357   2.000  33.000 

summary(nb_enf_enfants$nb_enf)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  0.000   0.000   1.000   1.605   2.000  19.000 

summary(nb_enf_p_enfants$nb_enf)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  0.000   0.000   1.000   1.615   2.000  21.000 

Ensuite, nous rejoignons toutes ces valeurs dans un seul tableau de données.

meres <- 
  nb_enf_aieux %>% 
  ungroup() %>% 
  rename(id_personne_num = id_aieul, dept_N = dept_N_aieul) %>% 
  mutate(generation = 0) %>% 
  bind_rows(
    nb_enf_enfants %>% 
      ungroup() %>% 
      rename(id_personne_num = id_enfant, dept_N = dept_N_p_enfants) %>% 
      mutate(generation = 1)
  ) %>% 
  bind_rows(
    nb_enf_p_enfants %>% 
      ungroup() %>% 
      rename(id_personne_num = id_p_enfant, dept_N = dept_N_ap_enfants) %>% 
      mutate(generation = 2)
  ) %>% 
  tbl_df()

Le nombre d’enfants moyen, sans se soucier des générations :

mean(meres$nb_enf)

[1] 1.467611

La répartition du nombre d’enfants, par génération (en lignes) est donnée dans le tableau ci-après. Les colonnes indiquent le nombre d’enfants. La colonne intitulée 11 concerne le pourcentage de femmes de l’échantillon a avoir donné naissance à 11 enfants et plus.

meres %>% 
  mutate(nb_enf = ifelse(nb_enf > 10, yes = 11, no = nb_enf)) %>% 
  select(generation, nb_enf) %>% 
  group_by(generation, nb_enf) %>% 
  summarise(nn = n()) %>% 
  ungroup () %>% 
  group_by(generation) %>% 
  mutate(pct_nn = nn / sum(nn)*100) %>% 
  select(generation, nb_enf, pct_nn) %>% 
  spread(nb_enf, pct_nn) %>% 
  knitr::kable(digits=2)

generation	0	1	2	3	4	5	6	7	8	9	10	11
0	46.98	24.54	10.26	6.14	4.02	2.62	1.89	1.20	0.85	0.53	0.36	0.63
1	27.36	39.18	13.18	7.48	4.80	2.86	1.87	1.22	0.84	0.51	0.27	0.43
2	25.81	40.33	13.72	7.88	4.50	2.68	1.89	1.17	0.78	0.47	0.30	0.46

3 Dimension temporelle des données : mortalité

Dans cette section, nous présentons les codes développés pour étudier les aspects temporels relatifs aux données. Nous nous concentrons principalement sur les aïeux de la base, à savoir les individus nés entre 1800 et 1804.

Nous pouvons explorer les questions de mortalité grâce aux dates de naissance et de décès des individus. Nous allons principalement nous concentrer sur les aïeux, mais les méthodes peuvent être appliquées à l’ensemble de la base.

aieux <-
  generations %>%
  select(id_aieul, date_N_annee_aieul, date_D_annee_aieul, age_aieul, sexe_aieul, dept_N_aieul) %>%
  unique() %>%
  filter(!is.na(date_N_annee_aieul), !is.na(date_D_annee_aieul)) %>%
  mutate(end_age = date_D_annee_aieul - date_N_annee_aieul) %>%
  tbl_df()

Nous avons besoin de charger quelques packages.

library(tidyverse)
library(demography)
library(survival)
library(data.table)

3.1 Estimation de la survie des individus

3.1.1 Estimations à partir des données issues de tables de mortalité

Nous avons également besoin de récupérer des données de mortalité. Nous téléchargeons les tables de mortalités de Vallin et al. (2001).

lifetable_fr <- read_csv("lifetable_FRA.csv")

Nous définissons une fonction pour extraire, pour une cohorte donnée, la probabilité de décès instantanné \(\mu_x\) (qui est équivalente à \(q_x\), c’est à dire la probabilité qu’un individu âgé de \(x\) ne survive pas en \(x+1\)).

#' qx_annee_age
#' Obtenir la proba qu'un individu âgé de @annee - @annee_naiss_depart
#' ne survive pas à l'âge de x
#' @annee (int) : Année pour laquelle extraire la table de mortalité
#' @annee_naiss_depart (int) : Annee de naissance des individus que l'on suit
# annee = 1806 ; annee_naiss_depart <- 1804
qx_annee_age <- function(annee, annee_naiss_depart, sexe){
  # Âge des individus que l'on suit
  age_courant <- annee - annee_naiss_depart
  # La table de mortalite en France pour l'année annee et le sexe sexe
  table_mortalite_annee <- lifetable_fr %>% filter(Year1 == annee, Sex == sexe, TypeLT == 1)
  # Extraire la ligne de cette table pour les individus que l'on suit
  res <- 
    table_mortalite_annee %>% 
    filter(Age == age_courant) %>% 
    slice(1) %>% 
    select(`q(x)`, `m(x)`) %>% 
    rename(qx = `q(x)`, mx = `m(x)`)
    # magrittr::extract2("q(x)")
  if(nrow(res) == 0) res <- data.frame(qx = NA, mx = NA)
  res
}# Fin de qx_annee_age()

Nous définissons une autre fonction pour calculer la probabilité qu’un individu âgé de \(x\) atteigne \(x+1\), soit \(p_x\), et la probabilité de survie conditionnelle à l’âge \(S_x(t)\). Ces valeurs sont estimées par cohorte.

#' survie_diag
#' On parcourt chaque année jusqu'à annee_naiss_depart+105
#' Mais on commence à 1806 (pas de données avant)
#' annee_naiss_depart <- 1804 ; sexe <- 1 ; age_min = 6
survie_diag <- function(annee_naiss_depart, sexe = c(1,2), age_min = 6){
  nqx <- lapply(seq(1806, annee_naiss_depart+105), qx_annee_age, annee_naiss_depart = annee_naiss_depart, sexe = sexe) %>% 
    bind_rows()

  ages <- seq(1806, annee_naiss_depart+105) - annee_naiss_depart
  res <- nqx
  res$age <- ages
  res$sexe <- as.character(sexe)
  
  res <- res %>% filter(age >= age_min)
  # Calcul de p_x (proba que (x) atteigne x+1)
  res$px = 1-res$qx
  # Rajout de la survie
  res$sx <- cumprod(res$px)
  res$cohorte <- annee_naiss_depart

  res
}# Fin de survie_diag

On applique la fonction survie_diag() aux années de naissances des individus des cohortes de 1800 à 1804, pour les femmes, puis pour les hommes.

survies_diagonale_1800_1804_f <- 
  pblapply(seq(1800, 1804), survie_diag, sexe = "2") %>% 
  bind_rows()


survies_diagonale_1800_1804_h <- 
  pblapply(seq(1800, 1804), survie_diag, sexe = "1") %>% 
  bind_rows()


survies_diagonale_1800_1804 <- 
  survies_diagonale_1800_1804_f %>% 
  bind_rows(survies_diagonale_1800_1804_h)

head(survies_diagonale_1800_1804)

# A tibble: 6 x 7
       qx      mx   age sexe     px    sx cohorte
    <dbl>   <dbl> <int> <chr> <dbl> <dbl>   <int>
1 0.0138  0.0139      6 2     0.986 0.986    1800
2 0.0118  0.0119      7 2     0.988 0.975    1800
3 0.0101  0.0102      8 2     0.990 0.965    1800
4 0.00808 0.00811     9 2     0.992 0.957    1800
5 0.00666 0.00668    10 2     0.993 0.950    1800
6 0.00591 0.00593    11 2     0.994 0.945    1800

survies_diagonale_1800_1804 %>% 
  mutate(sexe = factor(sexe, levels = c(2,1), labels = c("Femmes", "Hommes"))) %>% 
  ggplot(data = .,
       aes(x = age, y = sx, colour = factor(cohorte), linetype = sexe)) +
  geom_line() +
  ggtitle("Fonctions de survie des femmes et des hommes par cohortes")

Les estimations sont données par cohorte, nous les aggrégeons pour le groupe 1800-1804, en prenant la moyenne des probabilités de survie et de décès instantanné, par âge et sexe.

surv_hist <- 
  survies_diagonale_1800_1804 %>%
  group_by(age, sexe) %>%
  summarise(sx = mean(sx, na.rm=T),
            mx = mean(mx))
head(surv_hist)

# A tibble: 6 x 4
# Groups: age [3]
    age sexe     sx      mx
  <int> <chr> <dbl>   <dbl>
1     6 1     0.986 0.0137 
2     6 2     0.987 0.0135 
3     7 1     0.975 0.0112 
4     7 2     0.976 0.0112 
5     8 1     0.966 0.0102 
6     8 2     0.966 0.00983

3.1.2 Estimations à partir des données de généalogie

Maintenant que nous disposons d’estimations de la probabilité de survie à chaque âge à partir des données issues de tables de mortalité, nous proposons nos estimations réalisées à partir des données de généalogie de Geneanet. Nous définissons une fonction, calcul_survie_geneanet(), qui à partir d’une année de naissance, estime les probabilités de survie, de décès instantanné pour les hommes et les femmes, pour une la cohorte des individus nés cette année là.

#' calcul_survie_geneanet
#' Estime les probabilités de survie, de décès instantanné (en log) pour les hommes et les femmes
#' pour des individus d'une cohorte donnée
#' @une_annee (int) : année de naissance des individus à suivre
calcul_survie_geneanet <- function(une_annee){
  aieux_annee <- aieux %>% 
    tbl_df() %>% 
    filter(date_N_annee_aieul == une_annee) %>% 
    filter(age_aieul >= 6)
  # On compte le nombre de personne encore en vie à chaque année
  calcul_exposition <- function(une_annee, sexe){
    aieux_annee %>% filter(date_D_annee_aieul >= une_annee, sexe_aieul == sexe) %>% nrow()
  }
  
  annees <- seq(une_annee+6, une_annee+110)
  
  exposition_annee_h <- data.frame(annee = annees, exposition = sapply(annees, calcul_exposition,  sexe = "1")) %>% tbl_df()
  exposition_annee_f <- data.frame(annee = annees, exposition = sapply(annees, calcul_exposition, sexe = "2")) %>% tbl_df()
  exposition_annee <- 
    exposition_annee_h %>% rename(exposition_h = exposition) %>% 
    left_join(exposition_annee_f %>% rename(exposition_f = exposition), by = "annee") %>% 
    mutate(age = annees-une_annee)
  
  
  exposition_annee <- 
    exposition_annee %>% 
    mutate(exposition_h = ifelse(exposition_h <= 0, yes = 1, no = exposition_h),
           exposition_f = ifelse(exposition_h <= 0, yes = 1, no = exposition_f)) %>% 
    # Calcul du nombre de décès
    # mutate(expo_moins_1_h = lag(exposition_h),
    #        expo_moins_1_f = lag(exposition_f)) %>% 
    mutate(expo_plus_1_h = lead(exposition_h),
           expo_plus_1_f = lead(exposition_f)) %>%
    mutate(deces_h = exposition_h - expo_plus_1_h,
           deces_f = exposition_f - expo_plus_1_f) %>% 
    mutate(sx_h = (exposition_h - deces_h) / exposition_h,
           sx_f = (exposition_f - deces_f) / exposition_f) %>% 
    select(-expo_plus_1_h, -expo_plus_1_f) %>% 
    # Proba de déceder dans l'année
    mutate(qxt_h = deces_h / exposition_h,
           qxt_f = deces_f / exposition_f) %>% 
    mutate(qxt_h = ifelse(qxt_h <= 0, yes = NA, no = qxt_h),
           qxt_f = ifelse(qxt_f <= 0, yes = NA, no = qxt_f)) %>% 
    # Calcul de la force de mortalité
    mutate(fm_h = log(qxt_h),
           fm_f = log(qxt_f)) %>% 
    mutate(annee_naissance = une_annee)
  
  # Survie
  # exposition_annee$lx_f[1] <- exposition_annee$lx_h[1] <- 1
  exposition_annee$sx_h <- cumprod(exposition_annee$sx_h)
  exposition_annee$sx_f <- cumprod(exposition_annee$sx_f)
  exposition_annee
}# Fin de calcul_survie_geneanet()

On calcule donc les probabilités de survie et de décès instantanné pour les individus des cohortes 1800 à 1804, à partir des données de généalogie :

library(pbapply)
surv_geneanet_tot <- 
  pblapply(1800:1804, calcul_survie_geneanet) %>% 
  bind_rows()

head(surv_geneanet_tot)

# A tibble: 6 x 13
  annee exposition_h exposition_f   age deces_h deces_f  sx_h  sx_f
  <int>        <int>        <int> <int>   <int>   <int> <dbl> <dbl>
1  1806        72587        61392     6     550     570 0.992 0.991
2  1807        72037        60822     7     492     427 0.986 0.984
3  1808        71545        60395     8     376     397 0.980 0.977
4  1809        71169        59998     9     307     277 0.976 0.973
5  1810        70862        59721    10     245     233 0.973 0.969
6  1811        70617        59488    11     272     263 0.969 0.965
# ... with 5 more variables: qxt_h <dbl>, qxt_f <dbl>, fm_h <dbl>, fm_f
#   <dbl>, annee_naissance <int>

Comme pour les données issues des tables de mortalité, nous calculons la moyenne par age et par sexe :

surv_geneanet <- 
  surv_geneanet_tot %>% 
  group_by(age) %>% 
  summarise(fm_h = mean(fm_h, na.rm=T), fm_f = mean(fm_f, na.rm=T),
            sx_h = mean(sx_h, na.rm=T), sx_f = mean(sx_f, na.rm=T))

head(surv_geneanet)

# A tibble: 6 x 5
    age  fm_h  fm_f  sx_h  sx_f
  <int> <dbl> <dbl> <dbl> <dbl>
1     6 -4.77 -4.63 0.992 0.990
2     7 -4.96 -4.82 0.985 0.982
3     8 -5.18 -5.07 0.979 0.976
4     9 -5.32 -5.22 0.974 0.971
5    10 -5.44 -5.32 0.970 0.966
6    11 -5.51 -5.46 0.966 0.962

Nous modifions la forme du dernier tableau pour obtenir une disposition propice à l’utilisation par les fonctions de graphiques de ggplot2.

surv_geneanet_sx <- 
  surv_geneanet %>% 
  select(age, sx_h, sx_f) %>% 
  gather(sexe, sx, sx_h, sx_f) %>% 
  mutate(sexe = factor(sexe, levels = c("sx_h", "sx_f"), labels = c(1,2))) %>% 
  mutate(sexe = as.character(sexe))

Nous mettons en commun les données historiques issues des tables avec celles de notre échantillon de données généalogiques.

df_plot_survie <- 
  surv_hist %>% mutate(type = "Historique") %>%
  bind_rows(
    surv_geneanet_sx %>% 
      mutate(type = "Geneanet")
  ) %>% 
  mutate(sexe = factor(sexe, levels = c("2","1"), labels = c("Femmes", "Hommes")))

Il est alors possible de représenter les fonctions de survie par sexe pour chacune des deux bases.

p <- 
  ggplot(data = df_plot_survie,
         aes(x = age, y = sx, colour = sexe, linetype = type)) +
  geom_line() +
  scale_colour_manual("", values = c("Femmes" = "red", "Hommes" = "blue")) +
  scale_linetype_manual("", values = c("Historique" = "dashed", "Geneanet" = "solid")) +
  xlab("âge") + ylab("Probabilité") +
  scale_x_continuous(breaks = scales::pretty_breaks(n = 10))
p

3.1.3 Estimations à partir des données de généalogie, pour l’ensemble des individus

Cette fois, on peut effectuer les calculs pour l’ensemble des individus, non plus en se restreignant à ceux ayant au moins 6 ans.

#' calcul_survie_geneanet_tot
#' Estime les probabilités de survie, de décès instantanné (en log) pour les hommes et les femmes
#' pour des individus d'une cohorte donnée, dans un département donné
#' @une_annee (int) : année de naissance des individus à suivre
# une_annee <- 1800
calcul_survie_geneanet_tot <- function(une_annee){
  aieux_annee <- 
    aieux %>% 
    tbl_df() %>% 
    filter(date_N_annee_aieul == une_annee)
  
  # On compte le nombre de personne encore en vie à chaque année
  calcul_exposition <- function(une_annee, sexe){
    aieux_annee %>% filter(date_D_annee_aieul >= une_annee, sexe_aieul == sexe) %>% nrow()
  }
  
  annees <- seq(une_annee, une_annee+110)
  
  exposition_annee_h <- data.frame(annee = annees, exposition = sapply(annees, calcul_exposition,  sexe = "1")) %>% tbl_df()
  exposition_annee_f <- data.frame(annee = annees, exposition = sapply(annees, calcul_exposition, sexe = "2")) %>% tbl_df()
  exposition_annee <- 
    exposition_annee_h %>% rename(exposition_h = exposition) %>% 
    left_join(exposition_annee_f %>% rename(exposition_f = exposition), by = "annee") %>% 
    mutate(age = annees-une_annee)
  
  
  exposition_annee <- 
    exposition_annee %>% 
    mutate(exposition_h = ifelse(exposition_h <= 0, yes = 1, no = exposition_h),
           exposition_f = ifelse(exposition_h <= 0, yes = 1, no = exposition_f)) %>% 
    # Calcul du nombre de décès
    # mutate(expo_moins_1_h = lag(exposition_h),
    #        expo_moins_1_f = lag(exposition_f)) %>% 
    mutate(expo_plus_1_h = lead(exposition_h),
           expo_plus_1_f = lead(exposition_f)) %>%
    mutate(deces_h = exposition_h - expo_plus_1_h,
           deces_f = exposition_f - expo_plus_1_f) %>% 
    mutate(sx_h = (exposition_h - deces_h) / exposition_h,
           sx_f = (exposition_f - deces_f) / exposition_f) %>% 
    select(-expo_plus_1_h, -expo_plus_1_f) %>% 
    # Proba de déceder dans l'année
    mutate(qxt_h = deces_h / exposition_h,
           qxt_f = deces_f / exposition_f) %>% 
    mutate(qxt_h = ifelse(qxt_h <= 0, yes = NA, no = qxt_h),
           qxt_f = ifelse(qxt_f <= 0, yes = NA, no = qxt_f)) %>% 
    # Calcul de la force de mortalité
    mutate(fm_h = log(qxt_h),
           fm_f = log(qxt_f)) %>% 
    mutate(annee_naissance = une_annee)
  
  # Survie
  # exposition_annee$lx_f[1] <- exposition_annee$lx_h[1] <- 1
  exposition_annee$sx_h <- cumprod(exposition_annee$sx_h)
  exposition_annee$sx_f <- cumprod(exposition_annee$sx_f)
  exposition_annee
}# Fin de calcul_survie_geneanet_tot()

surv_geneanet_2_tot_c <- pblapply(1800:1804, calcul_survie_geneanet_tot)
surv_geneanet_2_tot_c <- bind_rows(surv_geneanet_2_tot_c)

surv_geneanet_2_tot <- 
  surv_geneanet_2_tot_c %>% 
  group_by(age) %>% 
  summarise(fm_h = mean(fm_h, na.rm=T), fm_f = mean(fm_f, na.rm=T),
            sx_h = mean(sx_h, na.rm=T), sx_f = mean(sx_f, na.rm=T))

3.1.4 Estimations à partir des données de généalogie, par département

Nous voulons également explorer la mortalité au travers de ses différences régionales. Nous adaptons le code précédent pour que les estimations s’effectuent à l’échelle des régions.

#' calcul_survie_geneanet_2
#' Estime les probabilités de survie, de décès instantanné (en log) pour les hommes et les femmes
#' pour des individus d'une cohorte donnée, dans un département donné
#' @une_annee (int) : année de naissance des individus à suivre
#' @un_dep (string) : identifiant du département (lu dans le tableau `aieux`)
calcul_survie_geneanet_2 <- function(une_annee, un_dep){
  aieux_annee <- 
    aieux %>% 
    tbl_df() %>% 
    filter(date_N_annee_aieul == une_annee) %>% 
    filter(dept_N_aieul %in% un_dep)
  
  # On compte le nombre de personne encore en vie à chaque année
  calcul_exposition <- function(une_annee, sexe){
    aieux_annee %>% filter(date_D_annee_aieul >= une_annee, sexe_aieul == sexe) %>% nrow()
  }
  
  annees <- seq(une_annee, une_annee+110)
  
  exposition_annee_h <- data.frame(annee = annees, exposition = sapply(annees, calcul_exposition,  sexe = "1")) %>% tbl_df()
  exposition_annee_f <- data.frame(annee = annees, exposition = sapply(annees, calcul_exposition, sexe = "2")) %>% tbl_df()
  exposition_annee <- 
    exposition_annee_h %>% rename(exposition_h = exposition) %>% 
    left_join(exposition_annee_f %>% rename(exposition_f = exposition), by = "annee") %>% 
    mutate(age = annees-une_annee)
  
  
  exposition_annee <- 
    exposition_annee %>% 
    mutate(exposition_h = ifelse(exposition_h <= 0, yes = 1, no = exposition_h),
           exposition_f = ifelse(exposition_h <= 0, yes = 1, no = exposition_f)) %>% 
    # Calcul du nombre de décès
    mutate(expo_plus_1_h = lead(exposition_h),
           expo_plus_1_f = lead(exposition_f)) %>%
    mutate(deces_h = exposition_h - expo_plus_1_h,
           deces_f = exposition_f - expo_plus_1_f) %>% 
    mutate(sx_h = (exposition_h - deces_h) / exposition_h,
           sx_f = (exposition_f - deces_f) / exposition_f) %>% 
    select(-expo_plus_1_h, -expo_plus_1_f) %>% 
    # Proba de déceder dans l'année
    mutate(qxt_h = deces_h / exposition_h,
           qxt_f = deces_f / exposition_f) %>% 
    mutate(qxt_h = ifelse(qxt_h <= 0, yes = NA, no = qxt_h),
           qxt_f = ifelse(qxt_f <= 0, yes = NA, no = qxt_f)) %>% 
    # Calcul de la force de mortalité
    mutate(fm_h = log(qxt_h),
           fm_f = log(qxt_f)) %>% 
    mutate(annee_naissance = une_annee)
  
  # Survie
  exposition_annee$sx_h <- cumprod(exposition_annee$sx_h)
  exposition_annee$sx_f <- cumprod(exposition_annee$sx_f)
  exposition_annee <- 
    exposition_annee %>% 
    mutate(dept = un_dep)
  exposition_annee
}# Fin de calcul_survie_geneanet_2()

Nous allons parcourir chacun des départements disponibles. Nous dressons la liste de ces derniers.

les_depts <- unique(aieux$dept_N_aieul)
les_depts <- les_depts[!is.na(les_depts)]

Nous parcourons chacun de ces départements, et pour chaque département, nous parcourons chaque cohorte pour lui appliquer la fonction calcul_survie_geneanet_2().

surv_geneanet_2_depts_c <- 
  pblapply(les_depts, function(un_dep){
    lapply(1800:1804, calcul_survie_geneanet_2, un_dep = un_dep) %>% 
      bind_rows()
  })

surv_geneanet_2_depts_c <- surv_geneanet_2_depts_c %>% bind_rows()

head(surv_geneanet_2_depts_c)

# A tibble: 6 x 14
  annee exposition_h exposition_f   age deces_h deces_f  sx_h  sx_f
  <int>        <dbl>        <int> <int>   <dbl>   <int> <dbl> <dbl>
1  1800         2231         2268     0   216       161 0.903 0.929
2  1801         2015         2107     1    67.0      90 0.873 0.889
3  1802         1948         2017     2    31.0      34 0.859 0.874
4  1803         1917         1983     3    29.0      42 0.846 0.856
5  1804         1888         1941     4    22.0      29 0.836 0.843
6  1805         1866         1912     5    13.0      15 0.831 0.836
# ... with 6 more variables: qxt_h <dbl>, qxt_f <dbl>, fm_h <dbl>, fm_f
#   <dbl>, annee_naissance <int>, dept <chr>

Puis, on agrège par âge et département, en calculant les moyennes pour la force de mortalité et pour la fonction de survie.

surv_geneanet_2_depts <- 
  surv_geneanet_2_depts_c %>% 
  group_by(age, dept) %>% 
  summarise(fm_h = mean(fm_h, na.rm=T), fm_f = mean(fm_f, na.rm=T),
            sx_h = mean(sx_h, na.rm=T), sx_f = mean(sx_f, na.rm=T))

On prépare les données pour que le format corresponde à ce qui est attendu par les fonctions de graphiques du package ggplot2.

surv_geneanet_2_dept_sx <- 
  surv_geneanet_2_depts %>% 
  select(dept, age, sx_h, sx_f) %>% 
  gather(sexe, sx, sx_h, sx_f) %>% 
  mutate(sexe = factor(sexe, levels = c("sx_h", "sx_f"), labels = c(1,2))) %>% 
  mutate(sexe = as.character(sexe))

df_plot_survie_2 <- 
  surv_geneanet_2_dept_sx %>% 
  mutate(sexe = factor(sexe, levels = c("2","1"), labels = c("Femmes", "Hommes")))

p <- 
  ggplot(data = df_plot_survie_2,
         aes(x = age, y = sx, colour = sexe, group = interaction(dept, sexe))) +
  geom_line(alpha = .2) +
  xlab("âge") + ylab("Probabilité") +
  scale_x_continuous(breaks = scales::pretty_breaks(n = 10))
p

3.2 Estimation de la force de mortalité

À partir des estimations réalisées dans la section précédente, on peut se pencher sur la force de mortalité, c’est-à-dire la probabilité, pour un individu ayant atteint un âge donné, de décéder dans le courant de l’année.

3.2.1 Pour les données de Vallin & Meslé comparées à celles de Geneanet

df_plot_mx <- 
  surv_hist %>% 
  select(age, sexe, mx) %>% 
  mutate(mx = log(mx)) %>% 
  mutate(type = "Historique") %>% 
  bind_rows(
    surv_geneanet %>% 
      select(age, fm_h, fm_f) %>% 
      gather(sexe, mx, fm_h, fm_f) %>% 
      mutate(type = "Geneanet") %>% 
      mutate(sexe = ifelse(sexe == "fm_h", yes = "1", no = sexe),
             sexe = ifelse(sexe == "fm_f", yes = "2", no = sexe))
  ) %>% 
  mutate(sexe = factor(sexe, levels = c("2","1"), labels = c("Femmes", "Hommes")))

head(df_plot_mx)

# A tibble: 6 x 4
# Groups: age [3]
    age sexe      mx type      
  <int> <fct>  <dbl> <chr>     
1     6 Hommes -4.29 Historique
2     6 Femmes -4.31 Historique
3     7 Hommes -4.49 Historique
4     7 Femmes -4.49 Historique
5     8 Hommes -4.58 Historique
6     8 Femmes -4.62 Historique

p_mx <- 
  ggplot(data = df_plot_mx, aes(x = age, y = mx, colour = sexe, linetype = type)) +
  geom_line() +
  scale_colour_manual("", values = c("Femmes" = "red", "Hommes" = "blue")) +
  scale_linetype_manual("", values = c("Historique" = "dashed", "Geneanet" = "solid")) +
  xlab("âge") + ylab("log(mu_x)") +
  scale_x_continuous(breaks = scales::pretty_breaks(n = 10))
p_mx

Et enfin, on peut tracer ces deux graphiques sur une seule figure :

df_plot <- 
  df_plot_survie %>% select(-mx) %>% rename(value = sx) %>% mutate(variable = "sx") %>% 
  bind_rows(df_plot_mx %>% rename(value = mx) %>% mutate(variable = "mx"))

p_survie <- 
  ggplot(data = df_plot %>% mutate(variable = factor(variable, levels = c("sx", "mx"),
labels = c("Probabilité de survie", "Force de mortalité (en log)"))),
         aes(x = age, y = value, colour = sexe, linetype = type)) +
  geom_line() +
  scale_colour_manual("", values = c("Femmes" = "red", "Hommes" = "blue")) +
  scale_linetype_manual("", values = c("Historique" = "dashed", "Geneanet" = "solid")) +
  xlab("age") + ylab(NULL) +
  scale_x_continuous(breaks = scales::pretty_breaks(n = 10)) +
  facet_wrap(~variable, scales = "free_y", nrow = 1)

p_survie

3.2.2 Pour l’ensemble des individus de Geneanet, par département

df_plot_mx_2 <- 
  surv_geneanet_2_depts %>% 
  select(dept, age, fm_h, fm_f) %>% 
  gather(sexe, mx, fm_h, fm_f) %>% 
  mutate(sexe = ifelse(sexe == "fm_h", yes = "1", no = sexe),
         sexe = ifelse(sexe == "fm_f", yes = "2", no = sexe)) %>% 
  mutate(sexe = factor(sexe, levels = c("2","1"), labels = c("Femmes", "Hommes")))



p_mx_2 <- 
  ggplot(data = df_plot_mx_2, aes(x = age, y = mx, colour = sexe, group = interaction(dept, sexe))) +
  geom_line(alpha = .2) +
  scale_colour_manual("", values = c("Femmes" = "red", "Hommes" = "blue")) +
  xlab("\\^age") + ylab("$log(\\mu\\_x$)") +
  scale_x_continuous(breaks = scales::pretty_breaks(n = 10))

p_mx_2+theme_paper()


# Les deux graphiques combinés

df_plot_2 <- 
  df_plot_survie_2 %>% rename(value = sx) %>% mutate(variable = "sx") %>% 
  bind_rows(df_plot_mx_2 %>% rename(value = mx) %>% mutate(variable = "mx"))

p_survie_2 <- 
  ggplot(data = df_plot_2 %>% 
           mutate(variable = factor(variable, levels = c("sx", "mx"),
                                    labels = c("Probabilit\\'e de survie", "Force de mortalit\\'e (en log)"))),
         aes(x = age, y = value, colour = sexe, group = interaction(dept, sexe))) +
  geom_line(alpha=.1) +
  scale_colour_manual("", values = c("Femmes" = "red", "Hommes" = "blue")) +
  xlab("\\^age") + ylab(NULL) +
  scale_x_continuous(breaks = scales::pretty_breaks(n = 10)) +
  facet_wrap(~variable, scales = "free_y", nrow = 1)


p_survie_2+theme_paper()

3.3 Mortalité régionale par cohorte

df_plot_mx_2 <- 
  surv_geneanet_2_depts %>% 
  select(dept, age, fm_h, fm_f) %>% 
  gather(sexe, mx, fm_h, fm_f) %>% 
  mutate(sexe = ifelse(sexe == "fm_h", yes = "1", no = sexe),
         sexe = ifelse(sexe == "fm_f", yes = "2", no = sexe)) %>% 
  mutate(sexe = factor(sexe, levels = c("2","1"), labels = c("Femmes", "Hommes")))

p_mx_2 <- 
  ggplot(data = df_plot_mx_2, aes(x = age, y = mx, colour = sexe, group = interaction(dept, sexe))) +
  geom_line(alpha = .2) +
  scale_colour_manual("", values = c("Femmes" = "red", "Hommes" = "blue")) +
  xlab("âge") + ylab("log(mu_x)") +
  scale_x_continuous(breaks = scales::pretty_breaks(n = 10))

p_mx_2

En combinant les deux graphiques :

df_plot_2 <- 
  df_plot_survie_2 %>% rename(value = sx) %>% mutate(variable = "sx") %>% 
  bind_rows(df_plot_mx_2 %>% rename(value = mx) %>% mutate(variable = "mx"))

p_survie_2 <- 
  ggplot(data = df_plot_2 %>% 
           mutate(variable = factor(variable, levels = c("sx", "mx"),
                                    labels = c("Probabilité de survie", "Force de mortalité (en log)"))),
         aes(x = age, y = value, colour = sexe, group = interaction(dept, sexe))) +
  geom_line(alpha=.1) +
  scale_colour_manual("", values = c("Femmes" = "red", "Hommes" = "blue")) +
  xlab("âge") + ylab(NULL) +
  scale_x_continuous(breaks = scales::pretty_breaks(n = 10)) +
  facet_wrap(~variable, scales = "free_y", nrow = 1)


p_survie_2

3.4 Espérance de vie à la naissance

À l’aide des estimations des fonctions de survie pour les femmes et les hommes, il est aisé de calculer les espérances de vie à la naissance des individus.

3.4.1 Pour Vallin et Meslé (2010)

Calculons l’espérance de vie à la naissance avec les données de Vallin et Meslé.

# Espérance de vie à la naissance pour la cohorte de 1806, données de Vallin et Meslé.
vallin_1806 <- 
  survie_diag(annee_naiss_depart = 1806, sexe = "1", age_min = 0) %>% 
  bind_rows(survie_diag(annee_naiss_depart = 1806, sexe = "2", age_min = 0))
vallin_1806 %>% 
  group_by(sexe) %>% 
  summarise(e = sum(cumprod(exp(-mx)), na.rm=T))

3.4.2 Pour l’ensemble des aïeux

Nous construisons une fonction effectuant le calcul.

calcul_esperance_vie_naissance_aieux <- function(df, summ = FALSE, age_min = 0){
  res <- 
    df %>% 
    filter(age_min >= age_min) %>% 
    left_join(
      df %>% 
        # group_by(annee_naissance, dept) %>% 
        filter(age == 0) %>% 
        select(annee_naissance, exposition_h, exposition_f) %>% 
        rename(exposition_f1 = exposition_f, exposition_h1 = exposition_h),
      by = c("annee_naissance")
    ) %>% 
    mutate(e_h = exposition_h / exposition_h1,
           e_f = exposition_f / exposition_f1) %>% 
    group_by(annee_naissance) %>% 
    summarise(e_h = sum(e_h)-0.5,
              e_f = sum(e_f)-0.5) %>% 
    ungroup()
  if(summ){
    res <- res %>% summarise(e_h = mean(e_h), e_f = mean(e_f))
  }
  res
}

Nous appliquons cette fonction à notre table de mortalité concernant la totalité des aïeux.

esperance_vie_aieux <- 
  calcul_esperance_vie_naissance_aieux(surv_geneanet_2_tot_c)

esperance_vie_aieux %>% 
  knitr::kable(digit=2)

annee_naissance	e_h	e_f
1800	45.27	43.31
1801	43.66	41.88
1802	43.06	41.32
1803	42.60	40.93
1804	43.12	41.62

esperance_vie_aieux_mean <- 
  esperance_vie_aieux %>% 
  summarise(e_h = mean(e_h), e_f = mean(e_f))
esperance_vie_aieux_mean

# A tibble: 1 x 2
    e_h   e_f
  <dbl> <dbl>
1  43.5  41.8

3.4.3 Par département

calcul_esperance_vie_naissance <- function(df, summ = FALSE, age_min = 0){
  res <- 
    df %>% 
    filter(age_min >= age_min) %>% 
    left_join(
      df %>% 
        # group_by(annee_naissance, dept) %>% 
        filter(age == 0) %>% 
        select(annee_naissance, dept, exposition_h, exposition_f) %>% 
        rename(exposition_f1 = exposition_f, exposition_h1 = exposition_h),
      by = c("annee_naissance", "dept")
    ) %>% 
    mutate(e_h = exposition_h / exposition_h1,
           e_f = exposition_f / exposition_f1) %>% 
    group_by(annee_naissance, dept) %>% 
    summarise(e_h = sum(e_h)-0.5,
              e_f = sum(e_f)-0.5) %>% 
    ungroup()
  if(summ){
    res <- res %>% summarise(e_h = mean(e_h), e_f = mean(e_f))
  }
  res
}

Les espérances de vie à la naissance pour chaque cohorte s’obtient de la manière suivante :

esperance_vie <- 
  calcul_esperance_vie_naissance(surv_geneanet_2_depts_c)
head(esperance_vie)

# A tibble: 6 x 4
  annee_naissance dept    e_h   e_f
            <int> <chr> <dbl> <dbl>
1            1800 F01    42.9  40.3
2            1800 F02    43.8  42.5
3            1800 F03    44.7  38.2
4            1800 F04    44.6  38.3
5            1800 F05    39.2  32.3
6            1800 F06    44.0  48.0

Il suffit alors d’effectuer la moyenne par département pour obtenir une estimation pour l’ensemble des cohortes de 1800 à 1804 :

esperance_vie_mean <- 
  esperance_vie %>% 
  group_by(dept) %>% 
  summarise(e_h = mean(e_h), e_f = mean(e_f))

head(esperance_vie_mean)

# A tibble: 6 x 3
  dept    e_h   e_f
  <chr> <dbl> <dbl>
1 F01    41.9  40.5
2 F02    41.1  39.9
3 F03    41.3  36.7
4 F04    42.7  37.3
5 F05    36.5  31.6
6 F06    42.7  42.5

Une représentation par le biais d’une carte choroplèthe est réalisable, en utilisant à nouveau la carte de france france_1800 dont la création est détaillée en annexe.

esperance_vie_map <- 
  france_1800 %>% 
  left_join(esperance_vie_mean, by = c("departement" = "dept"))

p_map <- 
  ggplot(esperance_vie_map %>% 
           gather(sexe, esperance, e_h, e_f) %>% 
           mutate(sexe = factor(sexe, levels = c("e_f", "e_h"), labels = c("Femmes", "Hommes"))),
         aes(x = long, y = lat, group = group, fill = esperance)) +
  geom_polygon() +
  facet_grid(~sexe) +
  scale_x_continuous(limits = c(-6, 10)) +
  scale_y_continuous(limits = c(41, 51.5)) +
  scale_colour_manual("Sexe", values = c("Femmes" = "blue", "Hommes" = "red"), na.value = "#377eb8") +
  scale_fill_gradient(low = "yellow", high = "red") +
  coord_quickmap() +
  theme(legend.position = "bottom", text = element_text(size = 6),
        legend.key.height   = unit(1, "line"),
        legend.key.width    = unit(1.5, "line"))


p_map

Il est intéressant de comparer les valeurs avec celles qui existent dans la littérature. Nous reportons les valeurs trouvées dans Walle (1973).

esperances_vandewalle <- 
  c('F01', 'Ain', '32.6',
    'F02', 'Aisne', '35.9',
    'F03', 'Allier', '31.7',
    'F04', 'Alpes-de-Haute-Provence', '34.2',
    'F05', 'Hautes-Alpes', '38.3',
    'F06', 'Alpes-Maritimes', NA,
    'F07', 'Ardèche', '41.9',
    'F08', 'Ardennes', '40.7',
    'F09', 'Ariège', '45.2',
    'F10', 'Aube', '34.2',
    'F11', 'Aude', '39.4',
    'F12', 'Aveyron', '44.6',
    'F13', 'Bouches-du-Rhône', NA,
    'F14', 'Calvados', '47.7',
    'F15', 'Cantal', '43.0',
    'F16', 'Charente', '35.8',
    'F17', 'Charente-Maritime', '32.8',
    'F18', 'Cher', '32.3',
    'F19', 'Corrèze', '34.8',
    'F20', 'Corse-du-Sud', '39',
    'F21', "Côte-d'Or", '35.4',
    'F22', "Côtes-d'Armor", '37.9',
    'F23', 'Creuse', '35.2',
    'F24', 'Dordogne', '34.3',
    'F25', 'Doubs', '39.6',
    'F26', 'Drôme', '37.3',
    'F27', 'Eure', '40.6',
    'F28', 'Eure-et-Loir', '30.1',
    'F29', 'Finistère', '29.2',
    'F30', 'Gard', '36.9',
    'F31', 'Haute-Garonne', '38.0',
    'F32', 'Gers', '40.4',
    'F33', 'Gironde', '37.7',
    'F34', 'Hérault', '37.8',
    'F35', 'Ille-et-Vilaine', '32.4',
    'F36', 'Indre', '29.1',
    'F37', 'Indre-et-Loire', '35.2',
    'F38', 'Isère', '37.9',
    'F39', 'Jura', '37.5',
    'F40', 'Landes', '30.7',
    'F41', 'Loir-et-Cher', '24.9',
    'F42', 'Loire', '37.0',
    'F43', 'Haute-Loire', '38.8',
    'F44', 'Loire-Atlantique', '39.4',
    'F45', 'Loiret', '24.7',
    'F46', 'Lot', '35.8',
    'F47', 'Lot-et-Garonne', '41.6',
    'F48', 'Lozère', '45.3',
    'F49', 'Maine-et-Loire', '38.4',
    'F50', 'Manche', '45.9',
    'F51', 'Marne', '38.1',
    'F52', 'Haute-Marne', '33.2',
    'F53', 'Mayenne', '34.7',
    'F54', 'Meurthe-et-Moselle', '36.0',
    'F55', 'Meuse', '35.0',
    'F56', 'Morbihan', '32.1',
    'F57', 'Moselle', '39.2',
    'F58', 'Nièvre', '26.9',
    'F59', 'Nord', '36.3',
    'F60', 'Oise', '35.8',
    'F61', 'Orne', '39.1',
    'F62', 'Pas-de-Calais', '38.0',
    'F63', 'Puy-de-Dôme', '35.9',
    'F64', 'Pyrénées-Atlantiques', '44.2',
    'F65', 'Hautes-Pyrénées', '46.5',
    'F66', 'Pyrénées-Orientales', '32.1',
    'F67', 'Bas-Rhin', '37.1',
    'F68', 'Haut-Rhin', '39.9',
    'F69', 'Rhône', NA,
    'F70', 'Haute-Saône', '34.8',
    'F71', 'Saône-et-Loire', '32.8',
    'F72', 'Sarthe', '32.4',
    'F73', 'Savoie', NA,
    'F74', 'Haute-Savoie', NA,
    'F75', 'Paris', NA,
    'F76', 'Seine-Maritime', '44.0',
    'F77', 'Seine-et-Marne', '28.6',
    'F78', 'Essonne', NA,
    'F79', 'Deux-Sèvres', '41.0',
    'F80', 'Somme', '41.0',
    'F81', 'Tarn', '39.6',
    'F82', 'Tarn-et-Garonne', '37.3',
    'F83', 'Var', '33.9',
    'F84', 'Vaucluse', '33.4',
    'F85', 'Vendée', '36.6',
    'F86', 'Vienne', '46.0',
    'F87', 'Haute-Vienne', '29.5',
    'F88', 'Vosges', '38.3',
    'F89', 'Yonne', '29.2') %>% 
  matrix(ncol = 3, byrow = TRUE) %>% 
  data.frame(stringsAsFactors = FALSE) %>% 
  tbl_df()

colnames(esperances_vandewalle) <- c("dept", "nom_sp", "e_f_walle")

esperances_vandewalle <- 
  esperances_vandewalle %>% 
  mutate(e_f_walle = as.numeric(e_f_walle))

Nous joignons nos estimations régionalles avec celles d’Etienne van de Walle.

esperance_vie_mean_2 <- 
  esperance_vie_mean %>% 
  left_join(esperances_vandewalle)

head(esperance_vie_mean_2)

# A tibble: 6 x 5
  dept    e_h   e_f nom_sp                  e_f_walle
  <chr> <dbl> <dbl> <chr>                       <dbl>
1 F01    41.9  40.5 Ain                          32.6
2 F02    41.1  39.9 Aisne                        35.9
3 F03    41.3  36.7 Allier                       31.7
4 F04    42.7  37.3 Alpes-de-Haute-Provence      34.2
5 F05    36.5  31.6 Hautes-Alpes                 38.3
6 F06    42.7  42.5 Alpes-Maritimes              NA  

L’espérance de vie moyenne à la naissance des femmes vaut, pour notre échantillon et pour celui de van de Walle :

mean(esperance_vie_mean_2$e_f, na.rm=T)

[1] 42.20441

mean(esperance_vie_mean_2$e_f_walle, na.rm=T)

[1] 36.67805

Notre valeur moyenne estimée est supérieure. Toutefois, nous souhaitons voir si les mêmes effets géographiques apparaissent avec nos données. Nous calculons alors les écarts par départements par rapport à la valeur moyenne des départements, pour nos données, ainsi que pour celles de la littérature.

esperance_vie_mean_2 <- 
  esperance_vie_mean_2 %>% 
  mutate(e_f_3 = e_f - mean(e_f, na.rm=T),
         e_f_walle_3 = e_f_walle - mean(e_f_walle, na.rm=T))

head(esperance_vie_mean_2)

# A tibble: 6 x 7
  dept    e_h   e_f nom_sp                  e_f_walle   e_f_3 e_f_walle_3
  <chr> <dbl> <dbl> <chr>                       <dbl>   <dbl>       <dbl>
1 F01    41.9  40.5 Ain                          32.6 - 1.70      - 4.08 
2 F02    41.1  39.9 Aisne                        35.9 - 2.30      - 0.778
3 F03    41.3  36.7 Allier                       31.7 - 5.49      - 4.98 
4 F04    42.7  37.3 Alpes-de-Haute-Provence      34.2 - 4.92      - 2.48 
5 F05    36.5  31.6 Hautes-Alpes                 38.3 -10.6         1.62 
6 F06    42.7  42.5 Alpes-Maritimes              NA     0.325      NA    

Puis, nous représentons le résultats sur des cartes choroplèthes.

esperance_vie_map <- 
  france_1800 %>% 
  left_join(esperance_vie_mean_2 %>% 
              select(dept, e_f_3, e_f_walle_3), by = c("departement" = "dept"))

p_map <- 
  ggplot(esperance_vie_map %>% 
           gather(type, esperance, e_f_3, e_f_walle_3) %>% 
           mutate(type = factor(type, levels = c("e_f_3", "e_f_walle_3"), labels = c("Geneanet", "van de Walle"))),
         aes(x = long, y = lat, group = group, fill = esperance)) +
  geom_polygon(col = "grey") +
  facet_wrap(~type) +
  scale_x_continuous(limits = c(-6, 10)) +
  scale_y_continuous(limits = c(41, 51.5)) +
  scale_fill_gradient2("Écart à la moyenne départementale (années)", low = "red", mid = "white", high = "#71BC78", midpoint = 0, na.value = "grey",
                       guide = guide_colourbar(direction = "horizontal",
                                               barwidth = 30,
                                               title.position = "bottom",
                                               title.hjust = .5)) +
  coord_quickmap() +
  theme(legend.position = "bottom", text = element_text(size = 6),
        legend.key.height   = unit(1, "line"),
        legend.key.width    = unit(1.5, "line")) +
  theme(text = element_text(size = 16), legend.position = "bottom")


p_map

3.5 Intervalles de confiance pour les espérances de vie à la naissance

Nous allons estimer des intervalles de confiance pour les espérances de vie à la naissance. Deux méthodes sont envisagées. Une première est d’utiliser ce que propose Carlo Giovanni Camarda sur son site (https://sites.google.com/site/carlogiovannicamarda/r-stuff/life-expectancy-confidence-interval). Une seconde est fortement inspirée de cette méthode, mais s’avère plus simple à mettre en oeuvre. Les deux méthodes sont fondées sur du bootstrap.

Nous récupérons la fonction lifetable.qx() servant à constuire une table de mortalité à partir des probabilités de décès.

## function for constructing a lifetable starting from probabilities
lifetable.qx <- function(x, qx, sex="M", ax=NULL, last.ax=5.5){
  m <- length(x)
  n <- c(diff(x), NA)
  qx[is.na(qx)] <- 0
  if(is.null(ax)){
    ax <- rep(0,m)
    if(x[1]!=0 | x[2]!=1){
      ax <- n/2
      ax[m] <- last.ax
    }else{    
      if(sex=="F"){
        if(qx[1]>=0.1){
          ax[1] <- 0.350
        }else{
          ax[1] <- 0.05 + 3*qx[1]
        }
      }
      if(sex=="M"){
        if(qx[1]>=0.1){
          ax[1] <- 0.33
        }else{
          ax[1] <- 0.0425 + 2.875*qx[1]
        }
      }
      ax[-1] <- n[-1]/2
      ax[m] <- last.ax
    }
  }
  px  <- 1-qx
  lx  <- cumprod(c(1,px))*100000
  dx  <- -diff(lx)
  Lx  <- n*lx[-1] + ax*dx
  lx <- lx[-(m+1)]
  Lx[m] <- lx[m]*last.ax
  Lx[is.na(Lx)] <- 0 ## in case of NA values
  Lx[is.infinite(Lx)] <- 0 ## in case of Inf values
  Tx  <- rev(cumsum(rev(Lx)))
  ex  <- Tx/lx
  data.frame(x, n, ax, qx, px, lx, dx, Lx, Tx, ex)
}

Cette fonction est appelée par celle qui suit, cette dernière permettant quant à elle de calculer les intervalles de confiance pour l’espérance de vie à un âge donné. Elle est légèrement adaptée des codes proposés par Camarda.

#' ex_ci
#' @description Compute a confidence interval for life expectancy
#' https://sites.google.com/site/carlogiovannicamarda/r-stuff/life-expectancy-confidence-interval
#' @param LT Life table
#' @param which.x Which age for the life expectancy
#' @param level level of the confidence interval
#' @param ns number of samples
#' @param sex 'F' (female) of 'M' (male)
#' @param year_birth year of birth of the cohort
#' @param dept_birth departement of birth
#' which.x <- 0; level <- 0.95; ns <- 1000; sex="M"; year_birth <- 1800 ; dept_birth <- "F29"
ex_ci <- function(LT, which.x=0, level=0.95, ns=1000, sex="M", year_birth = 1800, dept_birth = "F29"){
  ## ages
  x <- LT$age
  if(sex == "M"){
    ## number of deaths
    Dx <- LT$deces_h
    ## estimated probs
    qx <- LT$qxt_h
  }else{
    ## number of deaths
    Dx <- LT$deces_f
    ## estimated probs
    qx <- LT$qxt_f
  }
  
  ## number of ages
  m <- nrow(LT)
  
  ## trials for binomial, rounded
  Ntil <- round(Dx/qx)
  ## ax for last age
  # last.ax <- LT$ax[m]
  ## simulated death counts
  ## from Binomial distribution
  Y <- suppressWarnings(matrix(rbinom(m*ns,
                                      Ntil,
                                      qx),
                               m,ns))
  ## simulated probabilities
  QX <- Y/Ntil
  ## which age?
  wh <- which(x==which.x)
  
  ## for all replicates, compute life expectancy
  fun.ex <- function(qx){
    return(lifetable.qx(x=x, qx, sex=sex)$ex[wh])
  }
  exsim <- apply(QX, 2, fun.ex)
  
  ## confidence interval
  CI <- quantile(exsim,
                 probs = c((1-level)/2,
                           1 - (1-level)/2))
  
  data.frame(ex_mean = mean(exsim),
             ex_inf = CI[[1]],
             ex_sup = CI[[2]],
             age = which.x,
             sexe = sex,
             annee_naissance = year_birth,
             dept = dept_birth,
             stringsAsFactors = FALSE)
  
} # Fin de ex_ci()

Nous introduisons également une fonction permettant d’utiliser la seconde méthode. Cette méthode consiste à considérer une population \(E_0\) née une année donnée (e.g., 1800) ayant une probabilité de décéder dans l’année (\(\widehat{q}_x\), estimée sur nos données pour cette cohorte). Il s’agit de procéder de manière itérative. À l’âge \(x\), s’il reste \(E_x\) personnes en vie, le nombre de décès \(D_x\) est tiré suivant une loi binomiale \(\mathcal{B}(E_x,\widehat{q}_x)\). Ensuite, nous posons \(E_{x+1}=E_x - D_x\), c’est-à-dire que nous considérons le nombre de personnes en vie à l’âge suivant égal au nombre de personnes en vie à l’âge \(x\) auquel nous retranchons le nombre de décès qui vient d’être tiré. Nous tirons alors un nouveau nombre de décès pour l’âge suivant, et ainsi de suite. Ces étapes itératives sont reproduites de manière indépendante dans \(1~000\) échantillons. Dans chacun de ces échantillons, nous calculons l’espérance de vie à la naissance selon les tirages effectués. Les quantiles empiriques d’ordres \(0,025\) et \(0,975\) deviennent alors les bornes de l’intervalle de confiance à \(95\%\) pour l’espérance de vie à la naissance des individus de la cohorte étudiée. La fonction simu() ci-après permet d’effectuer un tirage. Il s’agira de faire appel à simu() de manière itérative pour obtenir différents échantillons.

#' simu
#' @description Estime des IC par bootstrap pour les esperances de vie a la naissance
#' a partir d'une table de vie
#' @param i valeur seed
#' @param LT Table de mortalite
simu <- function(i, LT){
  set.seed(i)
  d_h=d_f=rep(0,111)
  e_h=LT$exposition_h
  e_f=LT$exposition_f
  q_h=LT$qxt_h
  q_f=LT$qxt_f
  
  # Last not NA value :
  ind_max_not_na_h <- max(which(!is.na(q_h)))
  ind_na_h <- which(is.na(q_h))
  q_h[ind_na_h[ind_na_h> ind_max_not_na_h]]=1
  
  ind_max_not_na_f <- max(which(!is.na(q_f)))
  ind_na_f <- which(is.na(q_f))
  q_f[ind_na_f[ind_na_f> ind_max_not_na_f]]=1
  
  if(any(is.na(q_h))){
    loess_h <- loess(q_h~age, data = data.frame(q_h = q_h, age = 1:length(q_h)))
    ind_missing_h <- which(is.na(q_h))
    predicted_h <- predict(loess_h, newdata = data.frame(age = ind_missing_h))
    predicted_h[predicted_h>1] <- 1
    predicted_h[predicted_h<0] <- 0
    q_h[ind_missing_h] <- predicted_h
  }
  
  if(any(is.na(q_f))){
    loess_f <- loess(q_f~age, data = data.frame(q_f = q_f, age = 1:length(q_f)))
    ind_missing_f <- which(is.na(q_f))
    predicted_f <- predict(loess_f, newdata = data.frame(age = ind_missing_f))
    predicted_f[predicted_f>1] <- 1
    predicted_f[predicted_f<0] <- 0
    q_f[ind_missing_f] <- predicted_f
  }
  
  for(t in 1:110){
    D_h = rbinom(1,size=e_h[t],prob=q_h[t])
    d_h[t] = D_h
    e_h[t+1] = e_h[t]-d_h[t]
    D_f = rbinom(1,size=e_f[t],prob=q_f[t])
    d_f[t] = D_f
    e_f[t+1] = e_f[t]-d_f[t]  
  }
  s_h=e_h/e_h[1]
  s_f=e_f/e_f[1]
  
  list(E_h=sum(s_h)-.5,E_f=sum(s_f)-.5)
}# Fin de simu()

3.5.1 Pour Vallin et Meslé

LT_h <- lifetable_fr %>% filter(Year1 == 1806, Sex == 1, TypeLT == 1)
LT_f <- lifetable_fr %>% filter(Year1 == 1806, Sex == 2, TypeLT == 1)

#'
#'@param LT Life table
#'@param ns number of samples to draw
#'@param which.x age for which to compute life expectancy
#'@param level level for the confidence interval
ex_ci_vallin <- function(LT, ns=1000, which.x=0, level=0.95){
  ## sex
  sex <- LT$Sex %>% unique()
  ## ages
  x <- LT$Age
  ## number of deaths
  Dx <- LT$`d(x)`
  ## estimated probs
  qx <- LT$`q(x)`
  
  ## number of ages
  m <- nrow(LT)
  
  ## trials for binomial, rounded
  Ntil <- round(Dx/qx)
  ## ax for last age
  # last.ax <- LT$ax[m]
  ## simulated death counts
  ## from Binomial distribution
  Y <- suppressWarnings(matrix(rbinom(m*ns,
                                      Ntil,
                                      qx),
                               m,ns))
  ## simulated probabilities
  QX <- Y/Ntil
  ## which age?
  wh <- which(x==which.x)
  
  ## for all replicates, compute life expectancy
  fun.ex <- function(qx){
    return(lifetable.qx(x=x, qx, sex=sex)$ex[wh])
  }
  exsim <- apply(QX, 2, fun.ex)
  
  ## confidence interval
  CI <- quantile(exsim,
                 probs = c((1-level)/2,
                           1 - (1-level)/2))
  
  data.frame(ex_mean = mean(exsim),
             ex_inf = CI[[1]],
             ex_sup = CI[[2]],
             age = which.x,
             sexe = sex,
             annee_naissance = which.x,
             stringsAsFactors = FALSE)
  
}# End of ex_ci_vallin()

ex_vallin <- ex_ci_vallin(LT = LT_h) %>% 
  bind_rows(
    ex_ci_vallin(LT = LT_f)
  )

ex_vallin

3.5.2 Pour les aïeux

Nous proposons une fonction, e_ic_cohorte(), qui à partir d’une année de naissance, estime par bootstrap les intervalles de confiance pour l’espérance de vie à la naissance de la cohorte d’individus nés cette année.

#' e_ic_cohorte
#' @param year_birth annee de naissance de la cohorte
#' @param ns nombre d'échantillons
#' #' @param level niveau de l'intervalle de confiance
e_ic_cohorte <- function(year_birth, ns=1000, level=0.95){
  LT <- surv_geneanet_2_tot_c %>% filter(annee_naissance == year_birth)
  E=lapply(1:ns,simu, LT=LT)
  BE=matrix(unlist(E),nrow=2)
  
  ex_mean <- apply(BE,1,mean)
  ex_inf <- apply(BE,1,function(x) quantile(x,(1-level)/2))
  ex_sup <- apply(BE,1,function(x) quantile(x,1 - (1-level)/2))
  
  data.frame(ex_mean = ex_mean,
             ex_inf = ex_inf,
             ex_sup = ex_sup,
             sexe = c("M", "F"),
             annee_naissance = year_birth,
             stringsAsFactors = FALSE)
}# Fin de e_ic_cohorte()

ex_aieux <- pblapply(1800:1804, e_ic_cohorte)
ex_aieux <- bind_rows(ex_aieux)

ex_aieux_avg <- 
  ex_aieux %>% 
  group_by(sexe) %>% 
  summarise(ex_mean = mean(ex_mean),
            ex_inf = mean(ex_inf),
            ex_sup = mean(ex_sup))

ex_aieux_avg

# A tibble: 2 x 4
  sexe  ex_mean ex_inf ex_sup
  <chr>   <dbl>  <dbl>  <dbl>
1 F        41.8   41.6   42.0
2 M        43.5   43.3   43.7

3.5.3 Pour les aïeux, par département

Ici, nous proposons d’utiliser les deux méthodes d’estimation des intervalles de confiance.

3.5.3.1 Première méthode

Nous appliquons la méthode de Carlo Giovanni Camarda (https://sites.google.com/site/carlogiovannicamarda/r-stuff/life-expectancy-confidence-interval).

Nous bouclons sur les années de naissance des aïeux, pour chaque sexe, et pour chaque département.

annees <- 1800:1804
sexes <- c("F", "M")
depts <- sort(unique(surv_geneanet_2_depts_c$dept))

n_tot <- length(annees)*length(sexes)*length(depts)
pb <- txtProgressBar(min = 1, max = n_tot, style = 3)
i <- 1
ex_depts <- NULL
for(annee in annees){
 for(sexe in sexes) {
   for(dept in depts){
     LT <- surv_geneanet_2_depts_c %>% filter(annee_naissance == annee, dept == dept)
     tmp <- ex_ci(LT=LT, year_birth = annee, sex = sexe, dept_birth = dept)
     ex_depts <- bind_rows(ex_depts, tmp)
     rm(tmp, LT)
     setTxtProgressBar(pb, i)
     i <- i+1
   }
 }
}

head(ex_depts)

   ex_mean   ex_inf   ex_sup age sexe annee_naissance dept
1 40.27501 38.32497 42.34023   0    F            1800  F01
2 42.42850 40.84798 43.98664   0    F            1800  F02
3 38.24195 36.38821 40.12682   0    F            1800  F03
4 38.34529 34.82395 41.78034   0    F            1800  F04
5 32.26114 29.72819 35.05179   0    F            1800  F05
6 47.84448 44.62942 51.02613   0    F            1800  F06

3.5.3.2 Seconde méthode

En utilisant notre méthode à présent.

#' ex_ci
#' @description Compute a confidence interval for life expectancy
#' From an initial population $E_0$ born in `year_birth`, et a probability of dying within the current year
#' $\widehat{q}_x$ previously estimated, proceeds as follows:
#' For age $x$, if $E_x$ people are still alive, the number of deaths $D_x$ is drawn from a Binomial $B(E_x,q_x)$ distribution
#' and we then posit $E_{x+1}=E_x-D_x$.
#' @param level level of the confidence interval
#' @param ns number of samples
#' @param year_birth year of birth of the cohort
#' @param dept_birth departement of birth
#' which.x <- 0; level <- 0.95; ns <- 1000; sex="M"; year_birth <- 1800 ; dept_birth <- "F29"
ex_ci_2 <- function(level=0.95, ns=1000, year_birth = 1800, dept_birth = "F29"){
  
  LT <- surv_geneanet_2_depts_c %>% filter(annee_naissance == year_birth, dept == dept_birth)
  
  E=lapply(1:ns,simu, LT=LT)
  BE=matrix(unlist(E),nrow=2)
  
  ex_mean <- apply(BE,1,mean)
  ex_inf <- apply(BE,1,function(x) quantile(x,(1-level)/2))
  ex_sup <- apply(BE,1,function(x) quantile(x,1 - (1-level)/2))
  
  data.frame(ex_mean = ex_mean,
             ex_inf = ex_inf,
             ex_sup = ex_sup,
             sexe = c("M", "F"),
             annee_naissance = year_birth,
             dept = dept_birth,
             stringsAsFactors = FALSE)
  
} # Fin de ex_ci_2()

Comme cette méthode effectue les estimations pour chaque sexe dans la fonction simu(), nous itérons uniquement sur les années de naissance et sur les départements.

n_tot <- length(annees)*length(depts)
pb <- txtProgressBar(min = 1, max = n_tot, style = 3)
i <- 1
ex_depts_2 <- NULL
for(annee in annees){
  for(dept in depts){
    tmp <- ex_ci_2(year_birth = annee, dept_birth = dept)
    ex_depts_2 <- bind_rows(ex_depts_2, tmp)
    rm(tmp)
    setTxtProgressBar(pb, i)
    i <- i+1
  }
}

head(ex_depts_2)

   ex_mean   ex_inf   ex_sup sexe annee_naissance dept
1 41.02150 39.07937 42.96377    M            1800  F01
2 39.66530 37.69183 41.63957    F            1800  F01
3 43.74457 42.35724 45.04229    M            1800  F02
4 42.40019 40.97947 43.99872    F            1800  F02
5 44.09319 42.46796 45.83167    M            1800  F03
6 37.91267 36.07012 39.75813    F            1800  F03

ex_depts_avg <- 
  ex_depts_2 %>% 
  group_by(dept, sexe) %>% 
  summarise(ex_mean = mean(ex_mean),
            ex_inf = mean(ex_inf),
            ex_sup = mean(ex_sup))

head(ex_depts_avg)

# A tibble: 6 x 5
# Groups: dept [3]
  dept  sexe  ex_mean ex_inf ex_sup
  <chr> <chr>   <dbl>  <dbl>  <dbl>
1 F01   F        39.7   37.6   41.8
2 F01   M        40.9   38.9   42.9
3 F02   F        39.9   38.3   41.5
4 F02   M        41.0   39.5   42.4
5 F03   F        35.0   33.0   37.1
6 F03   M        39.4   37.4   41.4

Nous pouvons alors grapher les estimations des espérances de vie à la naissance accompagnées des intervalles de confiance.

esperance_vie_mean_ic <- 
esperance_vie_mean_2 %>% 
  left_join(
    ex_depts_avg %>% 
      filter(sexe == "F") %>% 
      rename(ex_mean_F = ex_mean, ex_inf_F = ex_inf, ex_sup_F = ex_sup) %>% 
      select(-sexe) %>% 
      left_join(
        ex_depts_avg %>% 
          filter(sexe == "M") %>% 
          rename(ex_mean_M = ex_mean, ex_inf_M = ex_inf, ex_sup_M = ex_sup) %>% 
          select(-sexe)
      )
  )

esperance_vie_mean_ic %>% 
  mutate(outside = e_f_walle <= ex_inf_F,
         inside = e_f_walle > ex_inf_F & e_f_walle < ex_sup_F) %>% 
  summarise(outside = sum(outside, na.rm=T),
            inside = sum(inside, na.rm=T))

p_error_bar_e <- 
  esperance_vie_mean_ic %>% 
  # mutate(departement = str_c(nom_sp, " (", str_sub(dept, 2), ")")) %>% 
  mutate(departement = str_sub(dept, 2)) %>% 
  select(departement, e_h, e_f, ex_inf_M, ex_sup_M, ex_inf_F, ex_sup_F) %>% 
  unite(values_f, e_f, ex_inf_F, ex_sup_F, sep = "@") %>% 
  unite(values_h, e_h, ex_inf_M, ex_sup_M, sep = "@") %>% 
  gather(sexe, value, values_h, values_f) %>% 
  separate(value, c("e", "e_inf", "e_sup"), sep = "@", convert = TRUE) %>% 
  mutate(echantillon = "Geneanet") %>% 
  filter(!is.na(departement)) %>% 
  ggplot(data = ., aes(x = departement, y = e, colour = sexe)) +
  geom_point() +
  geom_errorbar(aes(ymin = e_inf, ymax = e_sup), width=1, size = .8) +
  scale_colour_manual("Sexe", values = c("values_f" = "red", "values_h" = "blue"),
                      labels = c("values_f" = "Femmes", "values_h" = "Hommes")) +
  geom_point(data = 
               esperance_vie_mean_ic %>% 
                   # mutate(departement = str_c(nom_sp, " (", str_sub(dept, 2), ")")) %>% 
                   mutate(departement = str_sub(dept, 2)) %>% 
                   select(departement, e_f_walle) %>% 
                   mutate(sexe = "values_f") %>% 
                   rename(e = e_f_walle) %>% 
               mutate(echantillon = "van de Walle"),
             aes(x = departement, y = e, shape = echantillon), colour = "black") +
  scale_shape_manual("", values = c("van de Walle" = 17)) +
  theme(axis.text.x = element_text(angle = 90, size = .6*20)) +
  xlab("Département") + ylab("Espérance de vie à la naissance")

p_error_bar_e