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

Probleemstelling

Het artikel behandelt het inverse probleem van het reconstrueren van een potentieel coëfficiënt $q$ in een elliptische partiële differentiaalvergelijking (PDV) op basis van gedistribueerde, ruisgevoelige observaties van de oplossing $u$. De onderliggende vergelijking is:
$$\begin{cases} -\Delta u + qu = f & \text{in } \Omega, \\ u = 0 & \text{op } \partial\Omega, \end{cases}$$
waarbij $\Omega \subset \mathbb{R}^s$ ($s=1,2$) een begrensd domein is, $f$ een bekende bron is, en het potentieel $q$ behoort tot een toegestane verzameling van begrensde functies. De observaties worden gegeven door $u^\eta(x) = u(x) + \eta(x)$, waarbij $\eta$ de meetruis representeert.

De kernuitdaging ligt in de aard van de onbekende functie $q$. In veel praktische toepassingen (bijv. beeldvorming, elastografie) kan $q$ scherpe discontinuïteiten of niet-gladde kenmerken vertonen. Conventionele Bayesiaanse benaderingen die vertrouwen op Gaussische priors falen vaak in deze scenario's omdat ze een excessieve gladheid afdwingen, wat leidt tot vervaagde reconstructies die abrupte sprongen of stapdiscontinuïteiten niet kunnen vastleggen. Hoewel Total Variation (TV) priors dit aanpakken door stuksgewijs constante oplossingen te bevorderen, hebben zij beperkingen bij het definiëren van distributies op algemene topologische ruimtes en het vastleggen van complexe structurele informatie.

Methodologie

De auteurs stellen een nieuwe Persistent Homology (PH)-gebaseerde Bayesiaanse prior voor binnen een hybride raamwerk. De methodologie bestaat uit de volgende componenten:

1. Bayesiaans Raamwerk:
   Het probleem wordt geformuleerd met behulp van de stelling van Bayes. De posterieure distributie $\mu_{post}$ wordt afgeleid van een likelihood-functie (gebaseerd op Gaussische ruis) en een prior-maat $\mu_{pr}$. Om te garanderen dat de prior goed gedefinieerd is in oneindig-dimensionale ruimtes, construeren de auteurs de prior met betrekking tot een Gaussische referentiemaat $\mu_0$ (meestal gedefinieerd op $L^2(\Omega)$). De posterieure dichtheid wordt uitgedrukt via de Radon-Nikodym-afgeleide:
   $$\frac{d\mu_{post}}{d\mu_0}(q) \propto \exp(-\Phi(q; u^\eta) - R(q))$$
   waarbij $\Phi$ de data-misfit is en $R(q)$ een regularisatieterm.

2. Constructie van Persistent Homology (PH):

   • Topologische Kenmerken: De auteurs maken gebruik van persistente homologie om de topologische kenmerken van kandidaatfuncties te kwantificeren. Door het domein te discretiseren in een simpliciaal complex en sublevel-sets te definiëren, volgen zij de "geboorte" en "dood" van homologieklassen (verbonden componenten in 1D, lussen in 2D) naarmate de drempelwaarde van de functiewaarde toeneemt.
   • Persistentieparen: Deze kenmerken worden gecodeerd als persistentieparen $(\sigma, \tau)$, waarbij $\sigma$ het geboorte-simplex is en $\tau$ het dood-simplex. De persistentie wordt gedefinieerd als $|q(\tau) - q(\sigma)|$.
   • Persistentieafstand: Een metriek $\|q\|_{per}$ wordt gedefinieerd als de som van de persistenties van alle paren. In 1D is aangetoond dat dit nauw verwant is aan de Total Variation (TV), specifiek $\|q\|_{per} + \text{globale range} = TV(q)$.

3. De PH-Gaussische Prior:
   De auteurs introduceren een regularisatieterm $R_{PH}(q)$ gebaseerd op de persistentieafstand. De prior-dichtheid wordt gedefinieerd als:
   $$\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)$$

   • Wegingsstrategie: De term $\alpha_j(q)$ is een wegeingsfunctie die afhankelijk is van de persistentie-magnitude. Kleine persistenties (waarschijnlijk ruis) worden zwaar bestraft (grote $\alpha_j$), terwijl grote persistenties (significante structurele kenmerken) minder worden bestraft (kleine $\alpha_j$).
   • Generalisatie: Indien $\alpha_j \equiv 1$, reduceert de prior tot de standaard TV-prior. Zo wordt de PH-prior gepresenteerd als een generalisatie van TV die topologische weging incorporeert.
   • 2D Uitbreiding: Voor 2D vierkante domeinen wordt de prior geconstrueerd door de 1D PH-prior toe te passen op elke rij en kolom van de gediscretiseerde functiematrix en de resultaten op te tellen.

4. Sampling Algoritme:
   Om te samplen uit de posterieure verdeling gebruiken de auteurs een Delayed Acceptance Metropolis-Hastings (MH) algoritme. De potentiële functional $V(q)$ wordt gedecomposeerd in twee delen: een computationeel dure data-misfit term en een goedkopere regularisatieterm. Het algoritme gebruikt een tweestaps-acceptatiemechanisme:

   • Fase 1: Een pre-screening stap met behulp van de regularisatieterm $R(q)$ om voorstellen te filteren.
   • Fase 2: Een correctiestap met de volledige posterieure potentiaal indien het voorstel de eerste fase passeert.
     Deze aanpak verbetert de computationele efficiëntie, met name wanneer het evalueren van de forward operator kostbaar is.

Belangrijkste Bijdragen

• Nieuwe Prior Formulering: Het artikel introduceert een hybride prior die de probabilistische raamwerken van Gaussische maten combineert met topologische restricties afgeleid van persistente homologie.
• Topologische Regularisatie: In tegen tegenover standaard TV-priors, codeert de PH-prior expliciet de topologische structuur van de onbekende functie, waarbij onderscheid wordt gemaakt tussen ruis (kortlevende topologische kenmerken) en significante kenmerken (langlevende kenmerken) via een gewogen persistentieafstand.
• Generalisatie van TV: De auteurs demonstreren dat de PH-prior een generalisatie is van de TV-prior, die in staat is om op algemenere topologische ruimtes te opereren en rijkere structurele informatie te bewaren.
• Efficiënte Sampling: De ontwikkeling van een delayed acceptance algoritme, afgestemd op de decompositie van de PH-Gaussische potentiaal, maakt efficiënte sampling mogelijk in hoog-dimensionale inverse problemen.

Numerieke Resultaten

De auteurs valideren hun methode via numerieke experimenten in zowel 1D als 2D settings, waarbij de PH-Gaussische prior vergelijkt met:

1. Standaard Gaussische priors (met diverse covariantie-operatoren).
2. TV-Gaussische priors.

Bevindingen:

• Niet-gladde Functies: In gevallen die te maken hebben met stuksgewijs constante functies of functies met scherpe sprongen (Voorbeelden 3, 4, 5, 6, 7), falen standaard Gaussische priors om discontinuïteiten te vangen, wat resulteert in oscillerende of overmatig gladde resultaten. Zowel TV-Gaussische als PH-Gaussische priors presteren significant beter dan de standaard Gaussische prior.
• Oscillerende Functies: In Voorbeeld 1 (een Weierstrass-type functie) vertoonde de PH-Gaussische prior een superieure nauwkeurigheid vergeleken met de TV-Gaussische prior, met zichtbaar lagere foutsporen. De PH-prior onderdrukte effectief hoogfrequente ruis terwijl de complexe oscillerende structuur behouden bleef.
• 2D Reconstructie: In 2D voorbeelden met circulaire discontinuïteiten leverde de PH-Gaussische prior reconstructies die vergelijkbaar met of licht beter waren dan de TV-Gaussische prior, met name bij het afhandelen van complexe topologische structuren (Voorbeeld 7).
• Computationele Kosten: Het delayed acceptance algoritme bleek iets sneller dan standaard pCN-schema's. Hoewel de PH-prior zelf computationeel duurder is om te evalueren dan de TV-prior, rechtvaardigt de verbeterde nauwkeurigheid deze kosten in specifieke niet-gladde scenario's.

Betekenis en Claims

Het artikel claimt dat de voorgestelde PH-gebaseerde prior een bescheiden maar consistente verbetering biedt ten opzichte van de klassieke TV-prior, met name voor functies met complexe topologische kenmerken of hoogfrequente oscillaties.

De auteurs benadrukken dat het primaire voordeel van hun aanpak de toepasbaarheid is: door de prior te definiëren met betrekking tot een Gaussische referentiemaat en enkel te vereisen dat onbekende functies tot een geschikte topologische ruimte behoren, vermijdt de methode de strikte regulariteitsrestricties van klassieke benaderingen. Deze flexibiliteit maakt de PH-prior tot een robuust instrument voor inverse problemen waarbij de doelfunctie niet-glad is en complexe structurele informatie bezit die standaard gladheids-priors niet kunnen vangen. Het werk positioneert de PH-prior niet als een vervanging voor alle bestaande methoden, maar als een gespecialiseerde, topologisch geïnformeerde uitbreiding van de TV-prior voor uitdagende inverse potentieelproblemen.