Stochastic Zeroth-Order Method for Computing Generalized Rayleigh Quotients

Technische Samenvatting: Stochastische Zeroth-Order Methode voor het Berekenen van Gegeneraliseerde Rayleigh-Quotiënten

Probleemstelling
Het artikel behandelt de maximalisatie van de gegeneraliseerde Rayleigh-quotiënt, gedefinieerd als:
$$R(A, B) = \max_{v \in \mathbb{R}^d \setminus \{0\}} \frac{\langle v, Av \rangle}{\langle v, Bv \rangle}$$
waarbij $A \in \mathbb{R}^{d \times d}$ en $B \in \text{Sym}^d_{\succ 0}$ (symmetrisch positief definiet). Dit probleem is fundamenteel in de numerieke lineaire algebra en komt voor bij spectrale equivalentie, gegeneraliseerde enkelvoudige waardeontbindingen en stabiliteitsanalyse van partiële differentiaalvergelijkingen (via de numerieke abscissa).

Conventionele methoden voor het oplossen van dit probleem vertrouwen doorgaans op matrix-adjunct producten ($A^T$) en de expliciete berekening van $B^{-1}$. De auteurs identificeren twee primaire beperkingen van deze benaderingen:

1. Computationele Kosten en Stabiliteit: Het berekenen van $B^{-1}$ of het uitvoeren van matrix-vector producten met deze matrix is kostbaar en gevoelig voor numerieke onnauwkeurigheden.
2. Adjoint Mismatch: In toepassingen zoals computertomografie wordt de transponster $A^T$ vaak vervangen door een benaderde, computationeel hanteerbare operator. Deze "adjoint mismatch" leidt tot significante reconstructiefouten, waardoor methoden die exacte adjoints vereisen ongeschikt zijn.

Methodologie
De auteurs stellen een stochastische zeroth-order Riemannische algoritme voor dat direct opereert op de variëteit $S^{d-1}_B = \{v \in \mathbb{R}^d : \langle v, Bv \rangle = 1\}$. De methode vermijdt het volledig berekenen van $A^T$ of $B^{-1}$, en vertrouwt uitsluitend op matrix-vector producten met $A$ en $B$.

De kern van het algoritme (Algoritme 1 voor een enkele steekproef, Algoritme 2 voor meerdere steekproeven) verloopt als volgt:

1. Initialisatie: Start met een willekeurige vector $v_0$ op de variëteit.
2. Sampling: Trek in elke iteratie $k$ een willekeurige richting $x_k$ uit de raakruimte $T_{v_k} S^{d-1}_B$. Dit wordt bereikt door een standaard Gaussische vector te projecteren op de raakruimte en te normaliseren.
3. Optimale Stapgrootte: In plaats van een vaste of backtracking stapgrootte te gebruiken, lost het algoritme het eendimensionale maximalisatieprobleem langs de gesamplede richting analytisch op:
   $$\tau_k = \arg\max_{\tau \in \mathbb{R}} r(A, B, v_k + \tau x_k)$$
   Stelling 3.1 biedt een gesloten vorm oplossing voor $\tau_k$ die enkel binnenproducten van $v_k$ en $x_k$ met $A$ en $B$ bevat.
4. Update: De iteratie wordt bijgewerkt via de retractie-afbeelding: $v_{k+1} = R_{v_k}(\tau_k x_k)$.

De multi-sample variant (Algoritme 2) middelt meerdere gesamplede richtingen om de variantie van de gradiëntschatting te verminderen voordat de stapgrootte wordt berekend.

Theoretische Bijdragen
Het artikel stelt strikte convergentiegaranties vast onder de voorwaarde dat de dimensie van de eigensubruimte die overeenkomt met de maximale eigenwaarde, $\dim(G_{\max})$, kleiner is dan $d-1$.

• Niet-terminatie: De auteurs bewijzen dat als $\dim(G_{\max}) < d-1$, het algoritme met waarschijnlijkheid één niet voortijdig stopt (d.w.z. de gradiënt-proxy $b_k$ verdwijnt niet).
• Convergentie naar Globale Maximalisatoren: De sequentie van objectieve waarden $r(A, B, v_k)$ convergeert bijna zeker naar het globale maximum $R(A, B)$. Bovendien convergeert de sequentie van iteraties naar de verzameling van globale maximalisatoren $G_{\max}$.
• Convergentiesnelheden:
  • De norm van de Riemannische gradiënt verdwijnt met een sublineaire snelheid van $O(1/n)$ met waarschijnlijkheid één.
  • Specifiek, $\min_{k=0,\dots,n} \| \text{grad } f(v_k) \|^2 \leq \frac{C}{n+1}$, waarbij $C$ afhankelijk is van de conditiegetal van $B$ en de norm van $A$.
• Zeroth-Order Interpretatie: De methode wordt geherinterpreteerd als een stochastische zeroth-order Riemannische gradiëntascens. De auteurs tonen aan dat de update-richting $b_k x_k$ een onbevooroordeelde schatter is van de Riemannische gradiënt (geschaald met $1/(d-1)$), waarmee de voorgestelde methode wordt gekoppeld aan bestaande stochastische optimalisatiekaders op variëteiten.

Numerieke Resultaten
Uitgebreide numerieke experimenten valideren de theoretische bevindingen en vergelijken de voorgestelde methode met de meest geavanceerde alternatieven:

• Vergelijking met ZO-RGA: De voorgestelde methode presteert aanzienlijk beter dan de stochastische zeroth-order Riemannische gradiëntascens (ZO-RGA) uit [32], met name in hoge dimensies ($d=100, 500$). De analytische selectie van de stapgrootte $\tau_k$ blijkt superieur aan constante stapgroottes of Armijo backtracking gebruikt in ZO-RGA.
• Vergelijking met Gen-Oja: In een Karhunen-Loève probleem convergeert de voorgestelde methode sneller en levert het nauwkeurigere eigenfunctie-approximaties op dan zowel deterministische als stochastische Gen-Oja methoden.
• Robuustheid: Het algoritme vertoont robuustheid tegen slecht geconditioneerde matrices $B$, hoewel de convergentietijd toeneemt met het conditiegetal $\kappa(B)$.
• Steekproefgrootte ($m$): Het verhogen van het aantal steekproeven $m$ (van 1 naar 100) verbetert de convergentiesnelheid aanzienlijk en vermindert de variantie, met minimale extra runtime kosten dankzij vectorisatie.
• Rayleigh-Ritz Subprobleem: De auteurs vergelijken hun methode met een Rayleigh-Ritz benadering op de subruimte gespannen door de huidige iteratie en $m$ steekproeven. Hoewel de Rayleigh-Ritz methode sneller kan zijn voor specifieke combinaties van grote $m$ en $d$, biedt de voorgestelde methode een betere afweging tussen computationele complexiteit ($O(dm)$ versus $O(d^2m)$) en convergentiesnelheid.

Betekenis en Claims
Het artikel beweert een "eenvoudig en efficiënt" algoritme te introduceren dat de maximalisatie van de gegeneraliseerde Rayleigh-quotiënt oplost zonder dat daarvoor adjoint- of inversieberekeningen nodig zijn. De primaire betekenis ligt in:

1. Adjoint-Vrije Werking: Het biedt een bewezen convergerende oplossing voor problemen waar $A^T$ niet beschikbaar of mismatch is, een veelvoorkomend scenario in beeldvorming en inverse problemen.
2. Theoretische Rigor: Het biedt de eerste convergentieanalyse voor een stochastische zeroth-order Riemannische methode specifiek toegepast op de gegeneraliseerde Rayleigh-quotiënt, waarbij bijna zekere convergentie naar globale maximalisatoren en sublineaire snelheden worden vastgesteld.
3. Prestaties: Empirische resultaten tonen aan dat de methode bestaande stochastische en deterministische baselines overtreft, met name in hoog-dimensionale settings en bij het werken met slecht geconditioneerde matrices.

Het werk breidt een eerdere preprint [41] uit door het kader te generaliseren naar Riemannische meetkunde en gegeneraliseerde Rayleigh-quotiënten, inclusief aanvullende verduidelijkingen en substantiële verbeteringen in de convergentieanalyse.