A persistent-homology-based Bayesian prior for potential coefficient reconstruction in an elliptic PDE

Énoncé du problème

L'article traite du problème inverse de la reconstruction d'un coefficient de potentiel $q$ dans une équation aux dérivées partielles (EDP) elliptique à partir d'observations distribuées et bruitées de la solution $u$. L'équation direct est :
$$\begin{cases} -\Delta u + qu = f & \text{dans } \Omega, \\ u = 0 & \text{sur } \partial\Omega, \end{cases}$$
où $\Omega \subset \mathbb{R}^s$ ($s=1,2$) est un domaine borné, $f$ est une source connue, et le potentiel $q$ appartient à un ensemble admissible de fonctions bornées. Les observations sont données par $u^\eta(x) = u(x) + \eta(x)$, où $\eta$ représente le bruit de mesure.

Le défi central réside dans la nature de la fonction inconnue $q$. Dans de nombreuses applications pratiques (par exemple, l'imagerie, l'élastographie), $q$ peut présenter des discontinuités nettes ou des caractéristiques non lisses. Les approches bayésiennes conventionnelles reposant sur des priors gausiens échouent souvent dans ces scénarios car elles imposent une lissité excessive, conduisant à des reconstructions floues qui ne peuvent capturer les sauts brusques ou les discontinuités de type échelon. Bien que les priors de Variation Totale (TV) traitent cela en favorisant des solutions par morceaux constantes, ils présentent des limites pour définir des distributions sur des espaces topologiques généraux et pour capturer des informations structurelles complexes.

Méthodologie

Les auteurs proposent un nouveau prior bayésien basé sur l'homologie persistante (PH) au sein d'un cadre hybride. La méthodologie se compose des éléments suivants :

1. Cadre Bayésien :
   Le problème est formulé à l'aide du théorème de Bayes. La distribution a posteriori $\mu_{post}$ est dérivée d'une fonction de vraisemblance (basée sur un bruit gaussien) et d'une mesure a priori $\mu_{pr}$. Pour garantir que le prior soit bien défini dans des espaces de dimension infinie, les auteurs construisent le prior par rapport à une mesure de référence gaussienne $\mu_0$ (typiquement définie sur $L^2(\Omega)$). La densité a posteriori est exprimée via la dérivée de Radon-Nikodym :
   $$\frac{d\mu_{post}}{d\mu_0}(q) \propto \exp(-\Phi(q; u^\eta) - R(q))$$
   où $\Phi$ est l'écart aux données (misfit) et $R(q)$ est un terme de régularisation.

2. Construction de l'Homologie Persistante (PH) :

   • Caractéristiques Topologiques : Les auteurs utilisent l'homologie persistante pour quantifier les caractéristiques topologiques des fonctions candidates. En discrétisant le domaine en un complexe simplicial et en définissant des ensembles de sous-niveaux, ils suivent la « naissance » et la « mort » des classes d'homologie (composantes connexes en 1D, boucles en 2D) à mesure que le seuil de la valeur de la fonction augmente.
   • Paires de Persistance : Ces caractéristiques sont encodées sous forme de paires de persistance $(\sigma, \tau)$, où $\sigma$ est le simplexe de naissance et $\tau$ est le simplexe de mort. La persistance est définie comme $|q(\tau) - q(\sigma)|$.
   • Distance de Persistance : Une métrique $\|q\|_{per}$ est définie comme la somme des persistances de toutes les paires. En 1D, il est démontré que cela est étroitement lié à la Variation Totale (TV), spécifiquement $\|q\|_{per} + \text{étendue globale} = TV(q)$.

3. Le Prior PH-Gaussien :
   Les auteurs introduisent un terme de régularisation $R_{PH}(q)$ basé sur la distance de persistance. La densité du prior est définie par :
   $$\frac{d\mu_{pr}}{d\mu_0}(q) \propto \exp\left(-\lambda \sum_{(x_j, \tilde{x}_j) \in P(q)} \alpha_j(q) |q(x_j) - q(\tilde{x}_j)|\right)$$

   • Stratégie de Pondération : Le terme $\alpha_j(q)$ est une fonction de poids dépendant de la magnitude de la persistance. Les petites persistances (probablement du bruit) sont lourdement pénalisées (grand $\alpha_j$), tandis que les grandes persistances (caractéristiques structurelles significatives) sont moins pénalisées (petit $\alpha_j$).
   • Généralisation : Si $\alpha_j \equiv 1$, le prior se réduit au prior TV standard. Ainsi, le prior PH est présenté comme une généralisation de la TV qui incorpore une pondération topologique.
   • Extension 2D : Pour les domaines carrés en 2D, le prior est construit en appliquant le prior PH 1D à chaque ligne et colonne de la matrice de la fonction discrétisée, puis en sommant les résultats.

4. Algorithme d'Échantillonnage :
   Pour échantillonner de l'a posteriori, les auteurs emploient un algorithme de Metropolis-Hastings (MH) à acceptation différée. Le potentiel de la fonction $V(q)$ est décomposé en deux parties : un terme d'écart aux données coûteux en calcul et un terme de régularisation moins coûteux. L'algorithme utilise un mécanisme d'acceptation à deux étapes :

   • Étape 1 : Une étape de pré-filtrage utilisant le terme de régularisation $R(q)$ pour filtrer les propositions.
   • Étape 2 : Une étape de correction utilisant le potentiel complet de l'a posteriori si la proposition passe la première étape.
     Cette approche améliore l'efficacité computationnelle, particulièrement lorsque l'évaluation de l'opérateur direct est coûteuse.

Contributions Clés

• Nouvelle Formulation de Prior : L'article introduit un prior hybride qui combine le cadre probabiliste des mesures gaussiennes avec des contraintes topologiques dérivées de l'homologie persistante.
• Régularisation Topologique : Contrairement aux priors TV standards, le prior PH encode explicitement la structure topologique de la fonction inconnue, distinguant le bruit (caractéristiques topologiques de courte durée) des caractéristiques significatives (caractéristiques de longue durée) grâce à une distance de persistance pondérée.
• Généralisation de la TV : Les auteurs démontrent que le prior PH est une généralisation du prior TV, capable d'opérer sur des espaces topologiques plus généraux et de préserver des informations structurelles plus riches.
• Échantillonnage Efficace : Le développement d'un algorithme d'acceptation différée adapté à la décomposition du potentiel PH-gaussien permet un échantillonnage efficace dans des problèmes inverses de haute dimension.

Résultats Numériques

Les auteurs valident la méthode via des expériences numériques dans des contextes 1D et 2D, comparant le prior PH-gaussien aux suivants :

1. Priors gausiens standards (avec divers opérateurs de covariance).
2. Priors TV-gausiens.

Résultats :

• Fonctions Non-Lisses : Dans les cas impliquant des fonctions par morceaux constantes ou des fonctions avec des sauts brusques (Exemples 3, 4, 5, 6, 7), les priors gausiens standards échouent à capturer les discontinuités, produisant des résultats oscillatoires ou excessivement lisses. Les priors TV-gausiens et PH-gausiens surpassent significativement le prior gaussien standard.
• Fonctions Oscillatoires : Dans l'Exemple 1 (une fonction de type Weierstrass), le prior PH-gaussien a démontré une précision supérieure par rapport au prior TV-gaussien, avec des traces d'erreur visiblement plus faibles. Le prior PH supprime efficacement le bruit à haute fréquence tout en préservant la structure oscillatoire complexe.
• Reconstruction 2D : Dans les exemples 2D impliquant des discontinuités circulaires, le prior PH-gaussien a fourni des reconstructions comparables ou légèrement meilleures que le prior TV-gaussien, particulièrement pour la gestion des structures topologiques complexes (Exemple 7).
• Coût Computationnel : L'algorithme d'acceptation différée s'est avéré légèrement plus rapide que les schémas pCN standards. Bien que le calcul du prior PH soit plus coûteux que celui du prior TV, l'amélioration de la précision justifie ce coût dans des scénarios non lisses spécifiques.

Signification et Revendications

L'article affirme que le prior proposé basé sur la PH offre une amélioration modeste mais cohérente par rapport au la classe de la TV classique, particulièrement pour les fonctions possédant des caractéristiques topologiques complexes ou des oscillations à haute fréquence.

Les auteurs soulignent que le principal avantage de leur approche est son applicabilité : en définissant le prior par rapport à une mesure de référence gaussienne et en exigeant seulement que les fonctions inconnues appartiennent à un espace topologique approprié, la méthode évite les contraintes de régularité strictes des approches classiques. Cette flexibilité fait du prior PH un outil robuste pour les problèmes inverses où la fonction cible est non lisse et possède des informations structurelles complexes que les priors de lissité standards ne peuvent capturer. Ce travail positionne le prior PH non pas comme un remplacement de toutes les méthodes existantes, mais comme une extension spécialisée et topologiquement informée du prior TV pour les problèmes complexes de reconstruction de potentiel.