Homologie persistante

SMACUD1L - Master IMD


Calcul (homologie standard)

  • Complexe simplicial \( K \) ≃ liste de tuples
  • liste de tuples ≃ matrices de bord,
    \( \partial_{q+1}(\langle v_1, \cdots, v_q \rangle) = \sum_{i=1}^q (-1)^{i+1} \cdot \langle v_1, \cdots, \hat{v}_i, \cdots, v_q \rangle \)
  • \( H_q(K) = \ker(\partial_q) / \text{im}(\partial_{q+1}) = (\mathbb{Z}/2 \mathbb{Z})^{\beta_q} \)
  • \( \beta_q = \dim(\ker(\partial_q)) - \dim(\text{im}(\partial_{q+1})) \)
  • \( \dim(\ker(\partial_q)) \), \( \dim(\text{im}(\partial_{q+1})) \) ? ⟶ Diagonalisation
Exemple 8

\( \beta_0 = 4-3, \beta_1 = 2-1, \beta_2 = 0-0 \)

Pour calculer les groupes d'homologie avec coefficients dans une anneau ( \( \mathbb{Z} \) ), il faut calculer la forme normale de Smith des matrices de bord.

\( H_q = \mathbb{Z}^{4-3} \oplus \mathbb{Z}/2\mathbb{Z} \oplus \mathbb{Z}/6\mathbb{Z} \)

Calcul (homologie persistante)

Algorithme incrémental pour l'homologie

Algorithme pour calculer les nombres de Betti d'un complexe simplicial Delfinado + Edelsbrunner, 1992

Il faut trier les simplexes de façon à avoir une filtration :

  • \(K = \{\sigma^1, \ldots, \sigma^n\}\)
  • Pour tout j, \(K^j = \{\sigma^1, \ldots, \sigma^j\}\) est un complexe simplicial
β0, β1, β2 ← 0
for j = 1...n
    q ← dim(σj)
    if σj appartient à un cycle de Kj
        βq ← βq + 1
    else
        βq-1 ← βq-1 - 1
return (β0, β1, β2)

Idée de la preuve :

  • Ajouter un q-simplexe incrémente \(\beta_q\) ou décrémente \(\beta_{q-1}\) (Théorème d'Euler-Poincaré)
  • Si \(\sigma^j\) appartient à un cycle de \(K^j\) alors il y a un nouveau cycle par rapport à \(K^{j-1}\), et ce n'est pas un bord.
Exercice 7

Calculer les nombres de Betti

14

Deux types de simplexes :

  1. Positif : il appartient à un cycle dans \(K^j\)
    ⇒ « il ajoute un trou »
  2. Négatif : il n'appartient pas à un cycle dans \(K^j\)
    ⇒ son bord est un cycle non trivial dans \(K^{j-1}\) (exercice)
    ⇒ « il supprime un trou »

Algorithme

L'algorithme incrémental peut être adapté pour calculer les paires de persistance. Il nous faudra :

  1. Bon choix des bases d'homologie
  2. Coupler les simplexes négatifs (un trou disparaît) avec des simplexes positifs (un trou apparaît)

Pour simplifier, on considère des filtrations avec un seul simplexe entrant à la fois.

Quand on arrive à un simplexe positif \( \sigma^i \) :

  • nouveau cycle : \( g^i = \sigma^i + \sum_{l \in L} \sigma^l \)
  • notons que pour tout \( l \in L, l < i \)
  • on peut se restreindre à des \(\sigma^l\) négatifs (exercice)
  • base de \( H(K^i) \) = base de \( H(K^{i-1}) \cup \{ g^i \}\)*
  • Correspondance élément d'une base ↔ simplexe positif ↔ \( i \in [1, n] \)

Quand on arrive à un simplexe négatif \( \sigma^j \) :

  • \( \partial(\sigma^j) \) cycle non trivial dans \( K^{j-1} \)
  • \( \partial(\sigma^j) = \sum_{i \in I} g^i \) dans \( H(K^{j-1}) \)
  • \( \sum_{i \in I} g^i \sim 0 \) dans \( H(K^{j}) \) ⇒ ils ne sont pas linéairement indépendants
  • base de \( H(K^j) \) = base de \( H(K^{j-1}) \setminus \{ g^k \}, k \in I \) (lequel ?)
  • Correspondance simplexe négatif ↔ ensemble de simplexes positifs précédents ↔ \( I \subset [ 1, j-1 ] \)

Pseudocode :

L0, L1, L2 ← ∅ // listes des paires
B0, B1, B2 ← ∅ // bases d'homologie
for j = 1...n
    q ← dim(σj)
    if σj est positif
        Bq ← Bq ∪ {j}
    else
        I ← ∂(σj)
        i ← max(I)  // c'est lui, le plus jeune !
        Bq-1 ← Bq-1 - {i}
        Lq-1 ← Lq-1 ∪ {(i,j)}
foreach σi non couplé
    q ← dim(σi)
    Lq ← Lq ∪ {(i,n+1)}
return (L0, L1, L2)

\( D_q = \{ (a_i, a_j) \vert (i,j) \in L_q \} \)

Exercice 8

Calculer les paires de persistance \( (i, j) \)

1

k-Triangle Lemma

Soit \( Q^{i,j} = ]-\infty, i] \times ]j, \infty[ \), alors
\( \beta_q^{i,j} = \#(L_q \cap Q^{i,j}) \) (exercice)

C'est la clé de la preuve de l'algorithme.

Preuve de l'algorithme

On fixe \( q \) et \( i \), et on le démontre par induction sur \( k = j-i \geq 0 \).

Cas de base : Si \( i = j \) alors \( \beta_q^{i,i} = \beta_q(K^i) \).

D'autre part, \( \#(L_q \cap Q^{i,i}) \) =

  • = nombre de paires \( (a,b) \in L_q, a \leq i < b\)
  • = nombre \( q \)-simplexes positifs moins nombre \( (q+1) \)-simplexes négatifs à l'étape i
  • = \( \beta_q(K^i) \) (comme l'algorithme incrémental)

Induction : Vrai pour \( \beta_q^{i,j} \). Deux possibilités pour \( \beta_q^{i,j+1} \) :
\( \beta_q^{i,j+1} = \beta_q^{i,j} \) ou \( \beta_q^{i,j+1} = \beta_q^{i,j} -1 \)

Si \( \sigma^{j+1} \) est un \(q\)-simplexe positif :

  • \( \beta_q^{i,j+1} = \beta_q^{i,j} \)
  • il y une paire \( (j+1, b) \in L_q, b > j+1 \)
  • \( (j+1, b) \notin Q^{i,j} \Rightarrow \#(L_q \cap Q^{i,j+1}) = \#(L_q \cap Q^{i,j}) \)
  • \( \beta_q^{i,j+1} = \#(L_q \cap Q^{i,j+1}) \)

Si \( \sigma^{j+1} \) est un \( (q+1) \)-simplexe négatif :

  • \( \partial(\sigma^{j+1}) = \sum_{i \in I} g^i \)
  • \( \beta_q^{i,j+1} = \beta_q^{i,j} -1 \Leftrightarrow \max(I) \leq i \) (exercice)
  • \( \max(I) \leq i \Leftrightarrow \#(L_q \cap Q^{i,j+1}) = \#(L_q \cap Q^{i,j}) - 1\)
  • \( \beta_q^{i,j+1} = \#(L_q \cap Q^{i,j+1}) \) □

Ce qui reste :

  • Démontrer que \( D_q = \{ (a_i, a_j) \vert (i,j) \in L_q \} \).
  • Filtration générale, avec plusieurs simplexes entrant en même temps.
  • L'algorithme explicite pour savoir si un simplexe est positif, négatif, et son bord dans la base de \( H_q^{j} \)

L'algorithme explicite est similaire à celui du calcul de l'homologie avec la Forme Normale de Smith et sa complexité en temps est \( \mathcal{O}(n^3) \)

Morale
  • Algorithme aussi basé sur la Forme Normale de Smith.
  • Quand plusieurs trous fusionnent, ils tuent le plus jeune.

Applications

Quelques applications :

  1. Traitement d'images : reconstruction de contours
  2. Robotique : comment attraper un objet
  3. Visualisation scientifique : distinguer les parties importantes
  4. Analyse de formes : mesurer les trous

1. Traitement d'images

Comment reconnaître les contours ?

  • Entrée : ensemble fini de points \( P \subset \mathbb{R}^2 \)
    Il devrait ressembler à un échantillon sur un graphe planaire
  • Sortie : un graphe planaire \( G \)

Si l'ensemble de points représente assez bien la forme, on devrait la retrouver.

Algorithme :

  1. Construis une filtration avec des boules croissantes centrées sur les points (α-complexe)
  2. Calcule le diagramme de persistance de dimension 1 (composantes connexes ⇒ très rapide)
  3. Ignore les paires les plus proches de la diagonale à partir d'un seuil (widest gap) ⟵ clé de l'algorithme
  4. Le graphe \( G \) est formée par les cycles d'homologie des paires qui restent

Alpha complexe : filtration sur un ensemble de points \( P \subset R^2 \)

  • Complexe : triangulation de Delaunay d'un ensemble de points
  • Fonction de filtration \( f \) :
    • sur les 2-simplexes : rayon du cercle circonscrit du triangle
    • sur les 1-simplexes : rayon du cercle circonscrit de l'arête (s'il ne contient pas d'autre point), minimum des valeurs des cofaces sinon
    • sur les 0-simplexes : 0
  • Propriété : \( f^{-1}(-\infty, a] \cong \) union des boules de rayon \( a \) centrées sur les points de \( P \)
Exemple 9
Exercice 9

Quelle est la sortie de l'algorithme pour ces données ?

2. Robotique

Question : comment attraper un objet à l'aide d'une main robotique et de la gravité ?

Observation : il y a plusieurs façons de l'attraper ; certaines sont meilleurs qu'autres

Filtration : espace de translations invalides + un rideau qui descend

\( (x,y) \in D_1 \) ⇒ il y a une façon d'attraper l'objet et la quantité d’énergie nécessaire pour s'échapper est \( y - x \)

En pratique :

  • Translations et rotations ⇒ espace de configurations sans collision \( F \subset \mathbb{R}^2 \times \mathbb{S}^1 \)
  • Complexe : échantillon sur les configurations invalides et triangulation de Delaunay
  • Filtration définie sur les arêtes : \( -\infty \) si longueur sous un seuil, \( \min(-y) \) sinon. On regarde \( D_2 \)
Exemple 10

3. Visualisation scientifique

Problème : visualiser la formation de doigts dans un écoulement d'eau salée.

Données :

  • Pour chaque pas de temps \( k \), un ensemble fini \( P_k \subset \mathbb{R}^4 \)
  • Un élément \( (x, y, z, c) \) correspond à un point avec coordonnées spatiales \( (x, y, z) \) avec concentration de sel \( c \)

Étapes :

  1. Interpoler la concentration sur l'espace
  2. Séparer l'eau salée de l'eau non-salée (seuillage)
  3. Détecter les différents doigts (segmentation)
  4. Relier les différents doigts entre les pas de temps

Détails de l'étape 3 (segmentation) :

  • Fonction de filtration : distance géodésique (négative) vers la base
  • Seuillage du diagramme de persistance de dimension 0 proportionnel à l'image de la fonction (10%)
  • Chaque sommet est associé avec le maximum le plus proche
Exemple 11

4. Analyse de formes

Question : comment distinguer les trous d'un objet ?

Érosions + dilatations = filtration

  • \( D_q \cap Q^{0,0} \) = \( \beta_q \) paires \( (- a_i, b_i) \) avec \( a_i, b_i \geq 0 \)
  • Une paire \( (a_i, b_i) \) pour chaque trou
  • \( a_i \) : de combien il faut éroder l'objet pour casser ce trou
    = épaisseur
  • \( b_i \) : de combien il faut dilater l'objet pour remplir ce trou
    = largeur

De plus, le couplage de l'algorithme de l'homologie persistante permet de définir les deux centres de chaque trou.

Stabilité : la distance entre les paires de deux objets \( X \) et \( Y \) est bornée par \[ d_H(X, Y) + d_H(X^c, Y^c) \]

En français : si on perturbe la frontière de l'objet, l'épaisseur et la largeur de chaque trou reste à peu près la même.

Morale
  • La stabilité des diagrammes permet de dégager le bruit
  • En même temps, on peut distinguer les features les plus importantes (potentiellement)
  • On utilise souvent l'homologie persistante en dimension 0 (calcul beaucoup plus rapide)

Conclusion

  • Homologie persistante = extension de l'homologie à une suite de complexes
  • Stabilité ⇒ correction, classement, comparaison
  • Complexité worst-case cubique, mais presque linéaire en pratique
  • Très à la mode car elle permet de faire de petites erreurs

Quelques ressources :

  • Applied Algebraic Topology Network [YouTube]
  • Cours de Topological data analysis avec exercices [cours]
  • Geometric and Topological Inference [livre]
  • Computational topology: an introduction [livre]
  • Topological Data Analysis for Scientific Visualization [livre]