Licence

Du post : CC-BY-NC 2.0 FR : https://creativecommons.org/licenses/by-nc/2.0/fr/
Dépôt Github et Licence du code : http://github.com/JGravier/sitRada/tree/master/Effets2Sources

Objectifs

Visualisation des effets de sources et gestion de ces derniers sur des tableaux typo-chronologiques grâce à l’utilisation d’une visualisation heatmap des données. Plusieurs jeux de données sont ici testés. Ce sont des tableaux typo-chronologiques, au sens où les individus (lignes) sont relatifs à des périodes chronologiques (phases de sites, phases d’occupation d’un territoire…).

Pour plus d’information, voir le support de présentation lors de l’atelier Sitrada

Données

  1. “Compiegne.csv” : La céramique du site des “Halettes” à Compiègne (Oise). Le tableau est extrait de Desachy 2004, p. 44
    • lignes : périodes d’occupation du site archéologique, construites selon les observations de terrain,
    • colonnes : types de céramique (16 types)
  2. “Data_PCR_Troyes.csv” : Données du PCR de La Plaine de Troyes. Le tableau est extrait du compte rendu provisoire de l’atelier SITRADA : Sanson, Riquier 2016 (actuellement p. 20)
    • lignes : 7 phases d’occupation du territoire de la Plaine de Troyes, construites selon les observations de terrain,
    • colonnes : 3 variables quantitatives (surface, nombre de structures, NR récipients céramiques)
  3. “Noyon_50ans_Fonctions.csv” : Données du travail de thèse de J. Gravier sur la ville de Noyon. Tableau de contingence du nombre d’entités urbaines qui existent par fonction urbaine dans le temps long de 1800 ans. Voir notamment Gravier 2015
    • lignes : pas de temps chronologique arbitraire de 50 ans ([1;50], [51;100], [101;150]… jusqu’à la fin du 18e s.),
    • colonnes : fonctions urbaines définies selon l’adapatation du thésaurus du CNAU

Références

Outils

Réalisation avec R version 3.3.2 – “Sincere Pumpkin Patch” sur plateforme Windows 7 (64-bit)

Traitements

Packages nécessaires

library(tidyverse)
## Loading tidyverse: ggplot2
## Loading tidyverse: tibble
## Loading tidyverse: tidyr
## Loading tidyverse: readr
## Loading tidyverse: purrr
## Loading tidyverse: dplyr
## Conflicts with tidy packages ----------------------------------------------
## filter(): dplyr, stats
## lag():    dplyr, stats
library(gplots)
## 
## Attaching package: 'gplots'
## The following object is masked from 'package:stats':
## 
##     lowess

Importer les données

La céramique des Halettes à Compiègne

Compiegne <- read.csv(url("https://raw.githubusercontent.com/JGravier/sitRada/master/Heatmap/Compiegne.csv"),
                      header = TRUE,
                      sep = ";",
                      stringsAsFactors = F)

PCR Plaine de Troyes

PCR <- read.csv(url("https://raw.githubusercontent.com/JGravier/sitRada/master/Effets2Sources/Data_PCR_Troyes.csv"),
                header = TRUE,
                sep = ";",
                dec =",",
                stringsAsFactors = F)

Les entités urbaines de Noyon

Noyon <- read.csv(url("https://raw.githubusercontent.com/JGravier/These/master/Comparaison3Villes/Noyon_50ans_Fonctions.csv"),
                  header = T,
                  sep = ";",
                  stringsAsFactors = F)

On renomme les lignes, pour les avoir sur les plots finaux

row.names(Compiegne) <- Compiegne$Periode
row.names(PCR) <- PCR$ID
row.names(Noyon) <- Noyon$X

Les tableaux ont tous des variables quantitatives, mais ils sont de nature différente :
- tableau hétérogène (mesure et comptage) pour le cas du PCR
- tableau de contingence (dit également de comptage) pour Compiègne et Noyon

On a ainsi identifié deux chaînes de traitements pour visualiser et gérer les effets de sources sur des tableaux typo-chronologiques.

Chaîne de traitement 1

Application au cas du PCR de la Plaine de Troyes

Afin d’obtenir des heatmap grâce au package gplots, il est nécessaire de transformer les data.frames en matrices

PCRM <- PCR %>%
  select(surface:Nrceram) %>%
  as.matrix()
head(PCRM)
##       surface NbSutructures Nrceram
## CP001  133.86             1       0
## CP002  272.63             1       0
## CP003  601.38            17     169
## CP004   97.91             1       2
## CP005 2992.34             1       0
## CP006  305.19             2       0

Pour visualiser les effets de sources, on travaille directement sur le tableau de données standardisé. La standardisation permet en effet de rendre comparable des variables de nature ou d’amplitude différentes.

Création d’une fonction de standardisation

Standar <- function(x){
  scale(x, center = TRUE,  # center T = centrée (moins moyenne)
        scale = TRUE)  # scale T = réduite (div. par écart-type)
}
# Note : dans la fonction scale() la variance est exprimée en 1/n-1 (et non en 1/n)

Choix de palette de couleur

my_palette <- colorRampPalette(c("#018571", "#f5f5f5", "#a6611a"))(n = 50)

Visualisation sous forme de heatmap du tableau standardisé : les effets de sources

PCRM %>%
  Standar() %>%
  heatmap.2(Rowv = FALSE,
            dendrogram = "column", # dendrogram sur colonne + reoder
            col = my_palette,
            denscol = "black",
            trace = "none")

Il apparaît ici que le nombre de céramiques et de structures sont particulièrement associés à la phase d’occupation du territoire CP003 puisque les valeurs standardisées sont supérieures à 2. Compte tenu du fait que la valeur standardisée de la surface cumulée des sites archéologiques est également > à 2 pour CP007, cette variable est très fortement supérieure à la moyenne pour cette phase chronologique.
Les autres valeurs standardisées sont majoritairement négatives et comprises entre [-1 ; 0]. Ainsi, les valeurs initiales des variables sont relativement proches de la moyenne.

Avec cette méthode de standardisation, on a donc rendu comparable les variables mais pas les lignes (à savoir les périodes chronologiques). Or, il existe de fortes amplitudes des variables, et donc de fortes différences entre les périodes. C’est notamment ce qu’on remarque quand on observe les totaux des lignes :

PCRM %>% 
  margin.table(margin = 1)
##    CP001    CP002    CP003    CP004    CP005    CP006    CP007 
##   134.86   273.63   787.38   100.91  2993.34   307.19 19877.63

Cette somme n’a pas d’unité (puisqu’on à la fois des mesures et des effectifs). Cependant on note qu’entre la valeur minimale (CP004 ~ 100) et maximale (CP007 ~ 20 000), on a une différence énorme (CP007 est 200 fois plus conséquent que CP004). Dès lors, il s’avère nécessaire de rendre comparable les périodes chronologiques en réalisant, préalablement à l’opération de standardisation, une opération de mise en fréquence des variables. Ainsi, chaque ligne i du tableau vaudra 100% et chaque variable j sera connue selon sa proportion dans chaque période chronologique.

Création d’une fonction de mise en fréquence des variables

FrceVar <- function(x){
  prop.table(x, margin = 1)*100  # *100 pour obtenir la fréquence en pourcentage
}
# Note : margin = 1 car on s'intéresse aux fréquences des variables par ligne (période chronologique)

Transformation du tableau initial en fréquence des variables

PCRM %>%
  FrceVar()
##        surface NbSutructures    Nrceram
## CP001 99.25849    0.74150971  0.0000000
## CP002 99.63454    0.36545700  0.0000000
## CP003 76.37735    2.15905916 21.4635881
## CP004 97.02705    0.99098206  1.9819641
## CP005 99.96659    0.03340750  0.0000000
## CP006 99.34894    0.65106286  0.0000000
## CP007 99.61263    0.01509234  0.3722778

Puis applicarion de la standardisation sur ce nouveau tableau

PCRM %>%
  FrceVar() %>%
  Standar()
##          surface NbSutructures    Nrceram
## CP001  0.3890992    0.04547213 -0.4254785
## CP002  0.4325295   -0.46606968 -0.4254785
## CP003 -2.2534408    1.97375452  2.2584799
## CP004  0.1313908    0.38482759 -0.1776397
## CP005  0.4708778   -0.91775426 -0.4254785
## CP006  0.3995449   -0.07756208 -0.4254785
## CP007  0.4299987   -0.94266823 -0.3789262
## attr(,"scaled:center")
##       surface NbSutructures       Nrceram 
##    95.8893713     0.7080815     3.4025471 
## attr(,"scaled:scale")
##       surface NbSutructures       Nrceram 
##     8.6587667     0.7351358     7.9969900

Visualisation sous forme de heatmap du tableau mis en fréquence, puis standardisé : gérer les effets de sources

PCRM %>%
  FrceVar() %>%
  Standar() %>%
  heatmap.2(Rowv = F,
            dendrogram = "column",
            col = my_palette,
            denscol = "black",
            trace = "none")

Il est maintenant clair que la période CP003 est très différente des autres : sur-représentation du nombre de céramiques & de structures, et sous-représentation de la surface cumulée des sites. Au contraire, l’impression précédente selon laquelle CP007 se distinguait par une sur-représentation de la surface s’est avérée fausse. D’un point de vue proportionnel, la surface cumulée des sites archéologiques est semblable entre CP007 et les autres périodes.

Application au cas de Compiègne

CompM <- Compiegne %>%
  select(TypeA:TypeP) %>%
  as.matrix()
head(CompM)
##    TypeA TypeB TypeC TypeD TypeE TypeF TypeG TypeH TypeI TypeJ TypeK TypeL
## P5 13740   375  8270   250    20    10    40    80     0     5  1740   460
## P4 13540  1520 10110   740  1230   635   265   400     0   105  7210  1785
## P3 12490  5255  4220   980  1395  1415   440   680    30   580  6750  5930
## P2  6940  2880  5800  3400  1510  2280  1080  2840   620  2075  2130  2410
## P1  6490  2350  6900  2745   670  3020   950  6700   340  2660  1080   570
##    TypeM TypeN TypeO TypeP
## P5     0  1510   350     0
## P4   410   565   310   450
## P3   350   160    10   275
## P2   910   410   310   410
## P1   740   190   985    50

Visualisation des effets de sources

CompM %>%
  Standar() %>%
  heatmap.2(Rowv = FALSE, Colv = F,
            dendrogram = "none",
            col = my_palette,
            denscol = "black",
            trace = "none")

Gérer les effets de sources

CompM %>%
  FrceVar() %>%
  Standar() %>%
  heatmap.2(Rowv = FALSE, Colv = F,
            dendrogram = "none",
            col = my_palette,
            denscol = "black",
            trace = "none")

Dans le cas de Compiègne, il est notable que les deux visualisations sont semblables. Ceci se comprend au regard du fait que les différences entre les périodes ne sont pas très importantes. En effet, on a en effectif total des lignes :

margin.table(CompM, margin = 1)
##    P5    P4    P3    P2    P1 
## 26850 39275 40960 36005 36440

Application au cas de Noyon

NoyonM <- Noyon %>%
  select(F1:F10) %>%  # on a retiré la fonction F11 relative aux "formations naturelles", fonction donc non spécifiquement urbaine
  as.matrix()
head(NoyonM)
##    F1 F2 F3 F4 F5 F6 F7 F8 F9 F10
## P1  6  2  1  0  2  0  6  0  2   2
## P2  6  2  1  0  2  0  6  0  2   3
## P3  6  2  1  0  2  0  7  0  2   3
## P4  6  2  1  0  2  0  7  0  2   2
## P5  6  2  0  0  2  0  6  0  2   4
## P6  6  2  0  7  2  0  6  0  2   4

Visualisation des effets de sources

NoyonM %>%
  Standar() %>%
  heatmap.2(Rowv = FALSE,
            dendrogram = "column",
            col = my_palette,
            denscol = "black",
            trace = "none")

Il est très clair qu’il existe un changement entre P22 et P24 environ. Presque toutes les valeurs standardisées de P1 à P21 sont inférieures à 0, tandis que celles après P24 sont supérieures à 0. Ainsi, l’essentiel des valeurs des fonctions urbaines (variables) sont inférieures à la moyenne pour les périodes anciennes, tandis que c’est l’inverse pour les périodes plus récentes. Cette visualisation fait très clairement ressortir un effet de source : il existe beaucoup plus d’entités urbaines dans le tableau pour les périodes récentes puisque la documentation textuelle est nettement plus conséquente à partir du 11e-12e s. (équivalent aux lignes P22-P24). Ainsi, on connaît mieux la ville à partir de cette époque et l’on est en mesure de créer de nombreuses “entités urbaines” (= unité/brique d’analyse de l’espace intra-urbain).

Gérer les effets de sources

NoyonM %>%
  FrceVar() %>%
  Standar() %>%
  heatmap.2(Rowv = FALSE,
            dendrogram = "column",
            col = my_palette,
            denscol = "black",
            trace = "none")

Une fois la mise en fréquence réalisée, préalablement à la standardisation, les résultats sont très différents. Certaines fonctions sont plus représentées pour les périodes récentes (F8, F6 et F1) et d’autres pour les périodes plus anciennes (F7 à F5). Surtout, chaque variable peut être observée pour elle-même ou mise en comparaison avec les autres – en repérant visuellement des classes ou en opérant de nouvelles fonctions de classification.

Chaîne de traitement 2 : gérer les effets de sources pour des tableaux de contingence

Un autre moyen est possible pour gérer les effets de sources, mais cette chaîne de traitement ne s’applique qu’aux tableaux de contingence puisqu’elle se réfère à l’écart à l’indépendance.
L’objectif de ce traitement est de visualiser les écarts à l’indépendance en pourcentage et s’inspire directement des travaux réalisés par Bruno Desachy sur les écarts positifs aux pourcentage moyen, et son outil développé en VBA le sériographe EPPM. Contrairement au sériographe EPPM, la chaîne de traitement proposée ici prend en considération les écarts positifs et négatifs à l’indépendance.

Application au cas de Noyon

Création d’une fonction permettant d’obtenir les écarts à l’indépendance

TabEcart <- function(x){
  nindividus <- margin.table(x, margin = 1)  # totaux des lignes i
  njvariables <- margin.table(x, margin = 2)  # totaux des colonnes j
  ntotal <- sum(nindividus)  # effectif total
  TabInde <- (round(njvariables%*%t(nindividus)/ntotal, 2))  # tableau d'indépendance
  TabInde <- t(TabInde)
  TabEcart <- x - TabInde  # tableau des écarts à l'indépendance
}

Visualisation sous forme de heatmap du tableau mis en fréquence, puis des écarts à l’indépendance : gérer les effets de sources

NoyonM %>%
  FrceVar() %>%  # mise en fréquence
  TabEcart() %>%  # écarts à l'indépendance
  heatmap.2(Rowv = FALSE,
            dendrogram = "column",
            col = my_palette,
            denscol = "black",
            trace = "none")

Comme pour les variables standardisées, les écarts positifs sont en marron et les écarts négatifs en vert d’eau foncé. La distribution des écarts est globalement symétrique autour de 0 (cf. histogramme). Du fait qu’on étudie les écarts à l’indépendance en pourcentage, les valeurs sont plus facilement interprétables que dans le cas de la standardisation.
Cette deuxième chaîne de traitement permet de voir que les fonctions F7, F4, F1 et F8 se démarquent des autres fonctions. Les fonctions F1 et F8 sont particulièrement sous-représentées pour les périodes anciennes, tandis que c’est l’inverse pour F7 et F4.

Quelles différences entre ces deux chaînes de traitement ?

Ces deux chaînes de traitement ne doivent pas être considérées comme interchangeables, tout d’abord puisque la standardisation et les écarts à l’indépendance ne se fondent pas sur les mêmes hypothèses (bien que l’on se situe globalement dans le monde des statistiques probabilistes). Compte tenu du fait que les écarts à l’indépendance en pourcentage sont plus habituellement utilisés par les archéologues, la deuxième chaîne de traitement est sans doute plus accessible de prime abord. En revanche, la première chaîne de traitement permet deux choses supplémentaires : 1) travailler sur des tableaux hétérogènes, et 2) faciliter les comparaisons entre différents jeux de données (par exemple sur deux jeux données ayant les mêmes variables pour deux sites archéologiques synchrones ou deux jeux de données relatifs à un même site mais composés de variables différentes).
On peut en outre se demander si les différences de traitements jouent un rôle important quant aux classifications que l’on peut faire émerger, en particulier en matière de ressemblances et de différences entre périodes chronologiques.

Comparaison des classifications des périodes chronologiques de Noyon à partir de la standardisation & des écarts à l’indépendance

CAH sur chaîne de traitement 1

CAH1 <- NoyonM %>%
  FrceVar() %>%  # mise en fréquence
  Standar() %>%  # standardisation
  dist() %>%
  hclust(method = "ward.D2")  # méthode qui équivaut à la classification selon la méthode de Ward

CAH sur chaîne de traitement 2

CAH2 <- NoyonM %>%
  FrceVar() %>%  # mise en fréquence
  TabEcart() %>%  # écarts à l'indépendance
  dist() %>%
  hclust(method = "ward.D2")

Visualisation des dendrogrammes

par(mfrow = c(1,2))
plot(CAH1, hang = -1, cex = 0.6, main = "CAH 1", xlab = "Périodes")
plot(CAH2, hang = -1, cex = 0.6, main = "CAH 2", xlab = "Périodes")

De manière globale, les classifications sont relativement semblables. Il existe cependant des différences notables pour les périodes anciennes (~ P1 à P14) qui sont plus semblables quand on se fonde sur les écarts à l’indépendance (distances visiblement moins conséquentes).

On peut en outre se demander si l’on retrouve des différences notables lorsque le nombre de périodes étudiées est moins important, comme c’est souvent le cas pour des études intra-sites.

Comparaison des classifications des périodes d’occupation du site des Halettes de Compiègne

CAH

CAH1C <- CompM %>%
  FrceVar() %>%
  Standar() %>%
  dist() %>%
  hclust(method = "ward.D2")
CAH2C <- CompM %>%
  FrceVar() %>%
  TabEcart() %>%
  dist() %>%
  hclust(method = "ward.D2")

Visualisation des dendrogrammes

par(mfrow = c(1,2))
plot(CAH1C, hang = -1, cex = 0.6, main = "CAH 1", xlab = "Périodes")
plot(CAH2C, hang = -1, cex = 0.6, main = "CAH 2", xlab = "Périodes")

Pour le cas de Compiègne, on note qu’il n’existe pas de différences entre les classifications que l’on a pu faire émerger du tableau de données initial. On remarque en effets que P1 et P2 se distinguent des autres périodes pour les deux CAH. Par ailleurs, la période 5 se distingue à l’intérieur du groupe compopsé de P3 à P5.
Il nous semble que les différences entres les deux chaînes de traitements sont relativement faibles lorsqu’on réalise des CAH sur des jeux de données où les lignes sont peu nombreuses (ici des périodes).

Pour aller plus loin