6.2 Distance entre individus

Pour partir d’un exemple concret, imaginez que vous êtes en train d’analyser des données concernant les échantillons de plancton que vous avez prélevé sur votre lieu de recherche. Ce plancton a été numérisé (photo de chaque organisme) et les images ont été traitées avec un logiciel qui mesure automatiquement une vingtaine de variables telles que la surface de l’objet sur l’image, son périmètre, sa longueur … Vous vous trouvez donc face à un jeu de données qui a une taille non négligeable : 20 colonnes par 1262 lignes, soit le nombre d’individus mesurés dans vos échantillons.

SciViews vous propose, comme d’habitude, une série de fonctions qui simplifient et homogénéisent l’utilisation de R (avec la part belle à l’interface “formule”) pour votre exploration de données multivariées, notamment via le package {exploreit}. Tout ce que vous aurez besoin est chargé automatiquement si vous spécifiez "explore" comme section dans SciViews::R. Vous commencerez donc vos scripts R et vos documents R Markdowns avec l’instruction suivante :

SciViews::R("explore")
# Registered S3 methods overwritten by 'svBase':
#   method           from      
#   print.data.frame base      
#   print.data.table data.table
# ✓ Default data frame object (dtx): data.table
# ── Attaching packages ───────────────────────────────────── SciViews::R 1.5.0 ──
# ✓ exploreit  1.0.2       ✓ broom      0.8.0  
# ✓ SciViews   1.5.0       ✓ dtplyr     1.2.1  
# ✓ chart      1.4.0       ✓ dplyr      1.0.8  
# ✓ data.io    1.4.1       ✓ tidyr      1.2.0  
# ✓ svFlow     1.2.1       ✓ tibble     3.1.6  
# ✓ svBase     1.3.0       ✓ ggplot2    3.3.5  
# ✓ svMisc     1.3.1       ✓ data.table 1.14.2 
# ✓ fs         1.5.2       ✓ lattice    0.20.45
# ✓ collapse   1.7.6       ✓ MASS       7.3.56 
# ✓ forcats    0.5.1       ✓ rlang      1.0.2
# ── Conflicts ─────────────────────────────────────────── SciViews_conflicts() ──
# x data.table:::=() masks rlang:::=()
# x svFlow::!!()     masks rlang::!!()
# x svMisc::?()      masks utils::?()
# x dplyr::between() masks data.table::between()
# x collapse::D()    masks stats::D()
# x dplyr::filter()  masks stats::filter()
# x dplyr::first()   masks data.table::first()
# x dplyr::lag()     masks stats::lag()
# x dplyr::last()    masks data.table::last()
# x dplyr::select()  masks MASS::select()
# x data.io::write() masks base::write()

Ensuite, comme à notre habitude, nous lisons les données avec read() :

zoo <- read("zooplankton", package = "data.io")
zoo
# # A data.table: 1262 x 20
# # Language:     en
#      ecd  area perimeter feret major minor  mean  mode   min   max std_dev range
#    <dbl> <dbl>     <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>   <dbl> <dbl>
#  1 0.770 0.465      4.45 1.32  1.16  0.509 0.363 0.036 0.004 0.908   0.231 0.904
#  2 0.700 0.385      2.32 0.728 0.713 0.688 0.361 0.492 0.024 0.676   0.183 0.652
#  3 0.815 0.521      4.15 1.33  1.11  0.598 0.308 0.032 0.008 0.696   0.204 0.688
#  4 0.785 0.484      4.44 1.78  1.56  0.394 0.332 0.036 0.004 0.728   0.218 0.724
#  5 0.361 0.103      1.71 0.739 0.694 0.188 0.153 0.016 0.008 0.452   0.110 0.444
#  6 0.832 0.544      5.27 1.66  1.36  0.511 0.371 0.02  0.004 0.844   0.268 0.84 
#  7 1.23  1.20      15.7  3.92  1.37  1.11  0.217 0.012 0.004 0.784   0.214 0.78 
#  8 0.620 0.302      3.98 1.19  1.04  0.370 0.316 0.012 0.004 0.756   0.246 0.752
#  9 1.19  1.12      15.3  3.85  1.34  1.06  0.176 0.012 0.004 0.728   0.172 0.724
# 10 1.04  0.856      7.60 1.89  1.66  0.656 0.404 0.044 0.004 0.88    0.264 0.876
# # … with 1,252 more rows, and 8 more variables: size <dbl>, aspect <dbl>,
# #   elongation <dbl>, compactness <dbl>, transparency <dbl>, circularity <dbl>,
# #   density <dbl>, class <fct>

Vous voulez regrouper votre plancton en fonction de la ressemblance entre les organismes, c’est-à-dire, en fonction des écarts entre les mesures effectuées pour les 19 variables quantitatives, à l’exclusion de la vingtième colonne class qui est une variable factor). En raison de la taille du tableau, il est évident que cela ne pourra pas se faire de manière manuelle. Nous pouvons raisonnablement considérer que plus les mesures sont similaires entre deux individus, plus ils ont des chances d’être semblables, c’est-à-dire, d’appartenir au même groupe taxonomique. Mais comment faire pour synthétiser l’information de similarité ou différence contenue dans 19 paires de valeurs (une paire par variable) ? Nous avons besoin d’une mesure de distance qui quantifie la similarité (ou à l’inverse la dissimilarité) en un seul nombre. Celle qui vient naturellement à l’esprit est la distance euclidienne. Prenons un cas simplifié. Quelle est la distance qui sépare deux individus A et B par rapport à trois variables x, y, z ? Ici, nous pouvons représenter l’information graphiquement dans un espace à trois dimensions. La distance qui nous intéresse est la distance linéaire entre les deux points dans l’espace. Autrement dit, c’est la longueur du segment de droite qui relie les deux points dans l’espace. Cette distance, nous pouvons la calculer à l’aide de la formule suivante (voir par exemple ici pour une résolution dans le plan) :

\[\mathrm{D_{Euclidean}}_{A, B} = \sqrt{(x_A - x_B)^2 + (y_A - y_B)^2 + (z_A - z_B)^2}\]

Notez que cette formule se généralise à n dimensions et s’écrit alors, pour n’importe quelle paire d’individus indicés j et k dans notre tableau et pour les différentes mesures de 1 à i notés yi :

\[\mathrm{D_{Euclidean}}_{j, k} = \sqrt{\sum_{i=1}^{n}(y_{ij}-y_{ik})^2}\]

C’est la racine carrée de la somme dans les n dimensions des écarts entre les valeurs au carré pour toutes les variables yi. Plus sa valeur est grande, plus les individus sont éloignés (différents). Pour cette raison, nous appellerons cette distance, une mesure de dissimilarité.

6.2.1 Matrice de distances

Nous avons maintenant la possibilité de quantifier la dissimilarité entre nos organismes planctoniques… mais nous en avons un grand nombre. Cela va être impossible à gérer autant de mesures qu’il y a de paires possibles parmi 1262 individus12. La matrice de distance est une matrice ici 1262 par 1262 qui rassemble toutes les valeurs possibles. Notez que sur la diagonale, nous comparons chaque individu avec lui-même. La distance euclidienne vaut donc systématiquement zéro sur la diagonale.

\[\mathrm{D_{Euclidean}}_{j, j} = 0\]

De plus, de part et d’autre de cette diagonale, nous trouvons les paires complémentaires (j versus k d’un côté et k versus j de l’autre). Or qu’elle soit mesurée dans un sens ou dans l’autre, la distance du segment de droite qui relie deux points dans l’espace est toujours la même.

\[\mathrm{D_{Euclidean}}_{j, k} = \mathrm{D_{Euclidean}}_{k, j}\]

Par conséquent, seulement une portion (soit le triangle inférieur, soit le triangle supérieur hors diagonale) est informative. La diagonale ne porte aucune information utile, et l’autre triangle est redondant. Nous avons donc pour habitude de ne calculer et représenter que le triangle inférieur de cette matrice.

À vous de jouer !
h5p

Le package {exploreit} met à votre disposition les fonctions suivantes pour créer une matrice de dissimilarité et pour effectuer des regroupements (cluster en anglais) sur base de l’information qu’elle contient (nous verrons leur usage progressivement dans ce module) :

  • dissimilarity(data, method = "euclidean", scale = FALSE, transpose = FALSE) : matrice de dissimilarité, voir ?vegan::vegdist
    • print(), chart(), labels() et nobs() sont des méthodes disponibles pour les objets Dissimilarity/dist
  • cluster(x, method = "complete") : CAH (nous verron,s plus loin ce que c’est) à partir d’un objet Dissimilarity/*dist**, voir ?stats:hclust
    • print(), str(), plot(), chart(), labels(), nobs(), predict() et augment()
  • chart(cluster) : visualise un cluster (dendrogramme). Utiliser + geom_dendroline(h = XX, color = "red") pour y visualiser la coupure
  • predict(cluster, h = 5) ou predict(cluster, k = 5) extrait les groupes
  • augment(data = df, cluster, h =|k = ) ajoute les groupes dans le tableau df

Voici comment nous calculons une matrice de distances zoo6_dist, ici sur un petit sous-ensemble de six lignes de notre jeu de données (nous devons aussi éliminer la colonne class qui ne contient pas de données numériques et qui ne nous intéresse pas pour le moment) :

zoo %>.%
  sselect(., -class) %>.% # Élimination de la colonne class
  head(., n = 6) %->%
  zoo6 # Récupération des 6 premiers individus
zoo6_dist <- dissimilarity(zoo6, method = "euclidean")
zoo6_dist
# Dissimilarity matrix with metric: euclidean
# # A data.frame: 6 x 6
#   labels `1`    `2`    `3`    `4`    `5`   
#   <chr>  <dst>  <dst>  <dst>  <dst>  <dst> 
# 1 1                                        
# 2 2       8.219                            
# 3 3       2.565  5.732                     
# 4 4       0.858  7.881  2.222              
# 5 5       4.848  4.295  3.012  4.615       
# 6 6       2.427 10.620  4.923  2.848  7.177

Nous voyons bien ici que R n’imprime que le triangle inférieur de notre matrice 6 par 6. Notez aussi que les objets Dissimilarity de tailles plus réalistes que vous utiliserez dans vos analyses ne sont pas prévus pour être imprimés et visualisés telles quelles. Il s’agit seulement de la première étape vers une représentation utile qui sera réalisée à la page suivante, à l’aide de la classification hiérarchisée. Ne les imprimez donc jamais telles quelles dans vos scripts ou rapports !

Nous verrons plus loin comment nous pouvons utiliser l’information que cette matrice de distances contient pour regrouper les individus de manière pertinente à l’aide de la technique CAH. Mais avant cela, nous avons besoin d’un peu de théorie pour bien comprendre quelle métrique choisir pour calculer nos distances et pourquoi. On parle aussi d’indices de similarité ou dissimilarité.

Attention : nous n’avons pas considéré ici les unités respectives de nos variables. Une surface (mm2) ou une longeur (mm) ne sont pas mesurées dans les mêmes unités. Nous risquons alors de donner plus de poids dans nos calculs aux variables qui présentent des valeurs élevées. Nous aurions le même effet si nous décidions par exemple d’exprimer une mesure longitudinale en µm au lieu de l’exprimer en mm. Dans ce cas, il vaut mieux standardiser d’abord le tableau (moyenne de zéro et écart type de un) selon les colonnes avant d’effectuer le calcul. Une illustration de cette approche sera discutée plus loin.

À vous de jouer !
h5p

6.2.2 Indices de (dis)similarité

Un indice de similarité (similarity index en anglais) est un descripteur statistique (nombre unique) de la similitude de deux échantillons ou individus représentés par plusieurs variables dans un échantillon multivarié. Un indice de similarité prend une valeur comprise entre 0 (différence totale) et 1 ou 100% (similitude totale). Un indice de dissimilarité} est le complément d’un indice de similarité (dis = 1 – sim) ; sa valeur est comprise entre 100% (différence totale) et 0 (similitude totale).

Attention : dans certains cas, un indice de dissimilarité peut varier de 0 à +\(\infty\)**. Il n’existe alors pas d’indice de similarité complémentaire. C’est le cas précisément de la distance euclidienne que nous avons exploré jusqu’ici.

Tous les indices de similarité / dissimilarité peuvent servir à construire des matrices de distances.

6.2.2.1 Indice de Bray-Curtis

L’indice de dissimilarité de Bray-Curtis, aussi appelé coefficient de Czecanowski est calculé comme suit :

\[\mathrm{D_{Bray-Curtis}}_{j,k}=\frac{\sum_{i=1}^{n}\left|y_{ij}-y_{ik}\right|}{\sum_{i=1}^{n}(y_{ij}+y_{ik})}\]

Dans SciViews R nous utiliserons dissimilarity(DF, method = "bray"). cet indice s’utilise pour mesurer, entre autres, la similitude entre échantillons sur base du dénombrement d’espèces. Si le nombre d’espèces est très variable (espèces dominantes versus espèces rares), nous devons transformer les données pour éviter de donner trop de poids aux espèces les plus abondantes (par exemple à l’aide de transformations \(log(x+1)\), double racine carrée …).

Une caractéristique essentielle de cet indice (contrairement à la distance euclidienne) est que toute double absence n’est pas prise en compte dans le calcul. C’est souvent pertinent dans le cadre de son utilisation comme le dénombrement d’espèces. En effet, quelle information utile retire-t-on de doubles zéros dans un tableau répertoriant la faune belge pour le crocodile du Nil et le tigre de Sibérie par exemple ? Aucune ! Ils sont tous deux systématiquement absents des dénombrements, mais cette double absence n’apporte aucune information utile pour caractériser la faune belge par ailleurs.

L’indice de similarité de Bray-Curtis (sim) est complémentaire à l’indice de dissimilarité correspondant (dis tel que calculé ci-dessus) :

\[sim = 1 – dis\]

6.2.2.2 Indice de Canberra

L’indice de dissimilarité de Canberra est apparenté à l’indice de Bray-Curtis mais il pondère les individus en fonction du nombre d’occurrences afin de donner le même poids à chacun. Il se calcule comme suit :

\[\mathrm{D_{Canberra}}_{j,k}=\frac{1}{nz}\sum_{i'=1}^{nz}\frac{\left|y_{i'j}-y_{i'k}\right|}{\left|y_{i'j}\right|+\left|y_{i'k}\right|}\]

\(nz\) est le nombre de valeurs non nulles simultanément dans le tableau de départ. Tous les cas contribuent ici de manière égale. C’est un point positif, mais il faut faire attention à ce que cet indice a souvent tendance à donner, au contraire, trop d’importance aux dénombrements très rares observés une seule fois ou un petit nombre de fois ! Dans SciViews R, nous utiliserons dissimilarity(DF, method = "canberra").

Toute double absence n’est pas prise en compte ici également. Seuls les indices ne dépendant pas des doubles zéros sont utilisables pour des dénombrements d’espèces ou des présences-absences. Ainsi, pour ce type de données, notre choix se portera sur :

  • Bray-Curtis si l’on souhaite que le résultat soit dominé par les espèces les plus abondantes.

  • Canberra si notre souhait est de donner la même importance à toutes les espèces, mais avec un risque d’importance exagérée des espèces rares par rapport au nombre relatif d’observations pour l’ensemble du jeu de données.

  • Bray-Curtis sur données transformées (\(log(x+1)\) ou double racine carrée) pour un compromis entre les deux avec prise en compte de toutes les espèces, mais dépondération partielle des espèces les plus abondantes. C’est souvent un bon compromis.

Attention : Si les volumes échantillonnés entre stations ne sont pas comparables, il faut standardiser (moyenne nulle et écart type un) les données selon les échantillons avant de faire les calculs de distances. L’argument scale = TRUE pourra être ajouté à l’appel de dissimilarity() pour ce faire.

De même que pour Bray-Curtis, l’indice de similarité sim se calcule à partir de l’indice de dissimilarité dis tel que ci-dessus comme \(sim = 1 - dis\).

6.2.2.3 Distance Euclidienne

Nous savons déjà que c’est la distance géométrique entre les points dans un espace à n dimensions :

\[\mathrm{D_{Euclidean}}_{j,k}=\sqrt{\sum_{i=1}^{n}(y_{ij}-y_{ik})^2}\]

Dans SciViews R, cette distance peut être calculée avec dissimilarity(DF, method = "euclidean"). Attention : en anglais, c’est euclidean avec un “e”, pas euclidian ! Cet indice de dissimilarité est utile pour des mesures quantitatives, pour des données environnementales, etc. Il faut que les mesures soient toutes effectuées dans les mêmes unités. Si ce n’est pas le cas, nous devons les standardiser avant le calcul, avec scale = TRUE. Il n’existe pas d’indice de similarité complémentaire.

6.2.2.4 Distance de Manhattan

La distance de Manhattan, encore appelée “city-block distance” en anglais, est un indice de dissimilarité qui, contrairement à la distance euclidienne ne mesure pas la distance géométrique entre les points en ligne droite, mais via un trajet composé de segments de droites parallèles aux axes. C’est comme si la distance euclidenne reliait les points à vol d’oiseau, alors qu’avec la distance de Manhattan, nous devions contourner les blocs de maisons du quartier pour aller d’un point A à un point B (d’où le nom de cette métrique). Elle se calcule comme suit :

\[\mathrm{D_{Manhattan}}_{j,k}=\sum_{i=1}^{n}|y_{ij}-y_{ik}|\]

Dans SciViews R, nous utiliserons dissimilarity(DF, method = "manhattan"). Ici aussi, seul l’indice de dissimilarité est défini. L’indice de similarité complémentaire n’existe pas, car la valeur de l’indice de dissimlarité n’est pas bornée à droite et peut varier de zéro (dissimilarité nulle, cela signifie alors que les deux individus sont identiques) à l’infini pour une différence maximale.

6.2.3 Utilisation des indices

  • Les distances euclidienne ou de Manhattan sont à préférer pour les mesures environnementales ou de manière générale pour les variables quantitatives continues.

  • Les distances de Bray-Curtis ou Canberra sont meilleures pour les dénombrements d’espèces (ou à chaque fois que les doubles zéros ne portent aucune information utile). Il s’agit souvent de variables quantitatives discrètes prenant des valeurs nulles ou positives.

À vous de jouer !
h5p

6.2.4 Propriétés des indices

Les indices varient en 0 et 1 (0 et 100%), mais les distances sont utilisées aussi comme indices de dissimilarité et peuvent varier entre 0 et \(+\infty\).

Un indice est dit métrique si il a les propriétés suivantes :

  • Minimum 0 : \(I_{j, k} = 0\) si \(j = k\)

  • Positif : \(I_{j, k}>0\) si \(j \neq k\)

  • Symétrique : \(I_{j, k}=I_{k, j}\)

  • Inégalité triangulaire : \(I_{j, k} + I_{k, l} >= I_{j, l}\)

La dernière propriété d’inégalité triangulaire est la plus difficile à obtenir, et n’est pas toujours nécessaire. Nous pouvons montrer que certains indices qui ne respectent pas cette dernière propriété sont pourtant utiles dans le contexte. Nous dirons alors d’un indice que c’est une semi-métrique s’il répond à toutes les conditions sauf la quatrième. Enfin, un indice est dit non métrique dans tous les autres cas. Le tableau suivant reprend les métriques que nous avons vues jusqu’ici, et rajoute d’autres candidats potentiels (la distance Chi carré, l’indice de corrélation ou de variance/covariance) en indiquant leur type :

Distance Type
Bray-Curtis semi-métrique
Canberra métrique
Euclidienne métrique
Manhattan métrique
Chi carré métrique
(corrélation) (non métrique)
(variance/covariance) (non métrique)
Pour en savoir plus
  • Vous pouvez aussi transposer le tableau pour calculer la distance entre ses colonnes en utilisant l’argument transpose = TRUE dans dissimilarity(). Par exemple, dans le cas d’un tableau “espèces - station” (dénombrement d’espèces en différentes stations), nous pouvons comparer les stations du point de vue de la composition en espèces, mais nous pouvons aussi comparer les espèces du point de vue de leur répartition entre les stations. Pour passer d’un calcul à l’autre, nous transposerons donc le tableau (les colonnes deviennent les lignes et inversement).

  • Pour bien comprendre la logique derrière ces indices, il est utile de comprendre leurs équations. Si elles sont pour vous trop abstraites, une façon efficace de comprendre consiste à faire le calcul à la main. Par exemple dans le cas de l’indice de Canberra, la notion de nombre de données non nulles \(nz\) n’est pas évident. Effectuons un calcul à la main détaillé sur le tableau fictif suivant concernant trois espèces A, B, et C dénombré en trois stations sta1, sta2 et sta3 dans le tableau nommé ex1 :

A B C
sta1 4 0 2
sta2 3 0 10
sta3 1 8 0

Pour rappel, la dissimilarité de Canberra se calcule comme suit :

\[\mathrm{D_{Canberra}}_{j,k}=\frac{1}{nz}\sum_{i'=1}^{nz}\frac{\left|y_{i'j}-y_{i'k}\right|}{y_{i'j}+y_{i'k}}\]

où :

  • nz est le nombre d’observations non nulles simultanément dans les deux vecteurs comparés (les doubles zéros ne sont pas pris en compte)
  • i’ est l’itérateur sur toutes les valeurs non doubles zéros

Voici le détail du calcul (notez bien comment le double zéro pour l’espèce B entre les stations sta1 et sta2 est pris en compte dans le calcul) :

sta1_2 <- (1/2) * ((abs(4 - 3)) / (4 + 3) + (abs(2 - 10)) / (2 + 10))
round(sta1_2, 3)
# [1] 0.405
sta1_3 <- (1/3) * (abs(4 - 1) / (4 + 1) + abs(0 - 8) / (0 + 8) +
  abs(2 - 0) / (2 + 0))
round(sta1_3, 3)
# [1] 0.867
sta2_3 <- (1/3) * (abs(3 - 1) / (3 + 1) + abs(0 - 8) / (0 + 8) +
  abs(10 - 0) / (10 + 0))
round(sta2_3, 3)
# [1] 0.833

La matrice finale est la suivante :

sta1 sta2
sta2 0.405
sta3 0.867 0.833

Vérifions en laissant R faire le calcul :

(dissimilarity(ex1, method = "canberra"))
# Dissimilarity matrix with metric: canberra
# # A data.frame: 3 x 3
#   labels sta1  sta2 
#   <chr>  <dst> <dst>
# 1 sta1              
# 2 sta2   0.405      
# 3 sta3   0.867 0.833

  1. Le nombre de paires uniques et distinctes (donc dans les cas où jk) possibles parmi n items est \(n(n-1)/2\). Donc pour nos 1262 individus, nous avons 795.691 paires à calculer.↩︎