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

Enunciado do Problema

O artigo aborda o problema inverso de reconstruir um coeficiente de potencial $q$ em uma EDP elíptica a partir de observações distribuídas e ruidosas da solução $u$. A equação governante é:
$$\begin{cases} -\Delta u + qu = f & \text{em } \Omega, \\ u = 0 & \text{em } \partial\Omega, \end{cases}$$
onde $\Omega \subset \mathbb{R}^s$ ($s=1,2$) é um domínio limitado, $f$ é uma fonte conhecida e o potencial $q$ pertence a um conjunto admissível de funções limitadas. As observações são dadas por $u^\eta(x) = u(x) + \eta(x)$, onde $\eta$ representa o ruído de medição.

O desafio central reside na natureza da função desconhecida $q$. Em muitas aplicações práticas (ex: imagem, elastografia), $q$ pode exibir descontinuidades agudas ou características não suaves. Abordagens Bayesianas convencionais que dependem de priors Gaussianos frequentemente falham nesses cenários porque impõem uma suavidade excessiva, levando a reconstruções borradas que não conseguem capturar saltos abruptos ou descontinuidades de degrau. Embora os priors de Variação Total (TV) abordem isso promovendo soluções constantes por partes, eles possuem limitações na definição de distribuições em espaços topológicos gerais e na captura de informações estruturais complexas.

Metodologia

Os autores propõem um novo prior Bayesiano baseado em Homologia Persistente (PH) dentro de um framework híbrido. A metodologia consiste nos seguintes componentes:

1. Framework Bayesiano:
   O problema é formulado utilizando o teorema de Bayes. A distribuição posterior $\mu_{post}$ é derivada de uma função de verossimilhança (baseada em ruído Gaussiano) e de uma medida prior $\mu_{pr}$. Para garantir que o prior seja bem definido em espaços de dimensão infinita, os autores constroem o prior em relação a uma medida de referência Gaussiana $\mu_0$ (tipicamente definida em $L^2(\Omega)$). A densidade posterior é expressa via a derivada de Radon-Nikodym:
   $$\frac{d\mu_{post}}{d\mu_0}(q) \propto \exp(-\Phi(q; u^\eta) - R(q))$$
   onde $\Phi$ é o desajuste de dados (data misfit) e $R(q)$ é um termo de regularização.

2. Construção da Homologia Persistente (PH):

   • Características Topológicas: Os autores utilizam a homologia persistente para quantificar as características topológicas das funções candidatas. Ao discretizar o domínio em um complexo simplicial e definir conjuntos de subnível, eles rastreiam o "nascimento" e a "morte" de classes de homologia (componentes conexos em 1D, loops em 2D) conforme o limiar do valor da função aumenta.
   • Pares de Persistência: Essas características são codificadas como pares de persistência $(\sigma, \tau)$, onde $\sigma$ é o simplex de nascimento e $\tau$ é o simplex de morte. A persistência é definida como $|q(\tau) - q(\sigma)|$.
   • Distância de Persistência: Uma métrica $\|q\|_{per}$ é definida como a soma das persistências de todos os pares. Em 1D, demonstra-se que esta está intimamente relacionada à Variação Total (TV), especificamente $\|q\|_{per} + \text{amplitude global} = TV(q)$.

3. O Prior Gaussiano-PH:
   Os autores introduzem um termo de regularização $R_{PH}(q)$ baseado na distância de persistência. A densidade do prior é definida como:
   $$\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)$$

   • Estratégia de Ponderação: O termo $\alpha_j(q)$ é uma função de peso dependente da magnitude da persistência. Pequenas persistências (prováveis ruídos) são fortemente penalizadas (grande $\alpha_j$), enquanto grandes persistências (características estruturais significativas) são penalizadas menos (pequeno $\alpha_j$).
   • Generalização: Se $\alpha_j \equiv 1$, o prior reduz-se ao prior TV padrão. Assim, o prior PH é apresentado como uma generalização de TV que incorpora ponderação topológica.
   • Extensão para 2D: Para domínios quadrados 2D, o prior é construído aplicando o prior PH de 1D a cada linha e coluna da matriz da função discretizada e somando os resultados.

4. Algoritmo de Amostragem:
   Para amostrar da posterior, os autores empregam um algoritmo de Metropolis-Hastings (MH) de Aceitação Retardada (Delayed Acceptance). O potencial funcional $V(q)$ é decomposto em duas partes: um termo de desajuste de dados computacionalmente caro e um termo de regularização mais barato. O algoritmo utiliza um mecanismo de aceitação de dois estágios:

   • Estágio 1: Um passo de pré-triagem usando o termo de regularização $R(q)$ para filtrar propostas.
   • Estágio 2: Um passo de correção usando o potencial completo da posterior se a proposta passar no primeiro estágio.
     Esta abordagem melhora a eficiência computacional, particularmente quando a avaliação do operador direto é custosa.

Principais Contribuições

• Nova Formulação de Prior: O artigo introduz um prior híbrido que combina o framework probabilístico de medidas Gaussianas com restrições topológicas derivadas da homologia persistente.
• Regularização Topológica: Diferente dos priors de TV padrão, o prior PH codifica explicitamente a estrutura topológica da função desconhecida, distinguindo entre ruído (características topológicas de vida curta) e características significativas (características de vida longa) através de uma distância de persistência ponderada.
• Generalização de TV: Os autores demonstram que o prior PH é uma generalização do prior TV, capaz de operar em espaços topológicos mais gerais e preservar informações estruturais mais ricas.
• Amostragem Eficiente: O desenvolvimento de um algoritmo de aceitação retardada adaptado para a decomposição do potencial PH-Gaussiano permite a amostragem eficiente em problemas inversos de alta dimensão.

Resultados Numéricos

Os autores validam o método através de experimentos numéricos em configurações 1D e 2D, comparando o prior PH-Gaussiano contra:

1. Priors Gaussianos padrão (com vários operadores de covariância).
2. Priors TV-Gaussianos.

Descobertas:

• Funções Não Suaves: Em casos envolvendo funções constantes por partes ou funções com saltos agudos (Exemplos 3, 4, 5, 6, 7), os priors Gaussianos padrão falham em capturar descontinuidades, produzindo resultados oscilatórios ou excessivamente suaves. Os priors TV-Gaussiano e PH-Gaussiano superam significativamente o prior Gaussiano padrão.
• Funções Oscilatórias: No Exemplo 1 (uma função do tipo Weierstrass), o prior PH-Gaussiano demonstrou precisão superior em comparação ao TV-Gaussiano, com traços de erro visivelmente menores. O prior PH efetivamente suprimiu o ruído de alta frequência enquanto preservou a complexa estrutura oscilatória.
• Reconstrução 2D: Em exemplos 2D envolvendo descontinuidades circulares, o prior PH-Gaussiano forneceu reconstruções comparáveis ou ligeiramente melhores que o TV-Gaussiano, particularmente no tratamento de estruturas topológicas complexas (Exemplo 7).
• Custo Computacional: O algoritmo de aceitação retardada mostrou-se ligeiramente mais rápido que os esquemas pCN padrão. Embora o prior PH seja computacionalmente mais caro de avaliar do que o prior TV, a melhoria na precisão justifica o custo em cenários não suaves específicos.

Significância e Alegações

O artigo afirma que o prior proposto baseado em PH oferece uma melhoria modesta, porém consistente, sobre o clássico prior TV, particularmente para funções com características topológicas complexas ou oscilações de alta frequência.

Os autores enfatizam que a principal vantagem de sua abordagem é a sua aplicabilidade: ao definir o prior em relação a uma medida de referência Gaussiana e exigir apenas que as funções desconhecidas pertençam a um espaço topológico adequado, o método evita as restrições rigorosas de regularidade das abordagens clássicas. Esta flexibilidade torna o prior PH uma ferramenta robusta para problemas inversos onde a função alvo é não suave e possui informações estruturais intrincadas que os priors de suavidade padrão não conseguem capturar. O trabalho posiciona o prior PH não como um substituto para todos os métodos existentes, mas como uma extensão especializada e topologicamente informada do prior TV para problemas de potencial inverso desafiadores.