GKB Methods for X-Ray Computed Tomography with an Unmatched Back Projector

Sintesi Tecnica: Metodi GKB per la Tomografia Computerizzata a Raggi X con un Back Projector Non Corrispondente

Enunciato del Problema
Nella tomografia computerizzata (CT) su larga scala, la ricostruzione di un'immagine $x$ da dati misurati $b$ è formulata come un problema inverso lineare $Ax \approx b$. Nelle implementazioni pratiche, il proiettore diretto ($A$) e il back proiettore ($B$) sono spesso generati utilizzando schemi di discretizzazione differenti basati su approssimazioni fisiche. Di conseguenza, $B$ non è equivalente alla trasposta di $A$ ($B \neq A^\top$), creando una "coppia di proiettori non corrispondenti" (unmatched projector pair).

Questa discrepanza porta a equazioni normali inconsistenti. I solutori iterativi standard come il Gradiente Coniugato (CG) e il Metodo dei Residui Minimi Generalizzato (GMRES), quando applicati alla formulazione dei minimi quadrati di questi sistemi non corrispondenti, spesso non riescono a produrre ricostruzioni significative. Sebbene siano state proposte varianti precondizionate di GMRES (AB-GMRES e BA-GMRES) per gestire questi sistemi, esse soffrono di "semiconvergenza". In questo fenomeno, la soluzione converge inizialmente verso l'immagine reale, ma successivamente diverge con l'aumentare delle iterazioni, amplificando il rumore a causa dell'inclusione di componenti associate a piccoli valori singolari. Inoltre, i sistemi proiettati in GMRES coinvolgono matrici di Hessenberg superiori, che possono essere computazionalmente costose da decomporre man mano che la dimensione dello spazio di Krylov cresce.

Metodologia
Il documento propone due nuovi algoritmi, AB-GKB e BA-GKB, che adattano il metodo di bidiagonalizzazione di Golub–Kahan (GKB) per gestire coppie di proiettori non corrispondenti. A differenza di GMRES, che proietta il problema su un sottospazio che genera una matrice di Hessenberg superiore, GKB proietta il problema su un sottospazio a dimensione inferiore che genera una matrice bidiagonale inferiore sparsa ($C_k$).

I metodi proposti funzionano come segue:

• AB-GKB: Applica il processo GKB al sistema $ABx = Bb$ (trattando $B$ come un precondizionatore a destra). La soluzione viene recuperata tramite $x = BQ_k u_k$, dove $Q_k$ è una base ortogonale e $u_k$ è la soluzione del sistema bidiagonale proiettato.
• BA-GKB: Applica il processo GKB al sistema $BAx = Bb$ (trattando $B$ come un precondizionatore a sinistra). La soluzione viene recuperata tramite $x = Q_k u_k$.

Entrambi gli algoritmi evitano la formazione esplicita degli operatori composti $AB$o$BA$, affidandosi invece a successivi prodotti matrice-vettore con $A$, $B$ e le loro trasposte. I problemi proiettati di piccola scala vengono risolti utilizzando la Decomposizione ai Valori Singolari (SVD) della matrice bidiagonale sparsa $C_k$. Gli autori investigano due criteri di arresto per mitigare la semiconvergenza: il Principio della Discrepanza (DP) e lo Stagnazione della Norma del Residuo (RNS).

Contributi Chiave

1. Proposta Algoritmica: L'introduzione di AB-GKB e BA-GKB come forme precondizionate di GKB specificamente progettate per problemi di CT con back proiettori non corrispondenti.
2. Analisi della Complessità Teorica: Un confronto dettagliato delle operazioni in virgola mobile (FLOPS) tra i proposti metodi GKB e gli esistenti metodi AB/BA-GMRES. L'analisi rivela che, sebbene i metodi GMRES siano computazionalmente più veloci per iterazione (richiedendo meno operazioni), i metodi GKB generano una matrice proiettata più sparsa ($C_k$) rispetto alla matrice di Hessenberg $H_k$ in GMRES. La SVD della matrice bidiagonale $C_k$ è significativamente più economica ($O(k^3)$) rispetto a quella della matrice di Hessenberg $H_k$ ($O(k^4)$).
3. Validazione Numerica: Test estesi su tre esempi di fantocci CT (128x128, 180x180 e 256x256) con diversi livelli di rumore e diverse dimensioni del sistema (sottodeterminati e sovradeterminati).
4. Valutazione dei Criteri di Arresto: Una valutazione comparativa tra DP e RNS, dimostrando che DP termina costantemente le iterazioni vicino al minimo errore di ricostruzione relativo (RRE).

Risultati
Gli esperimenti numerici attraverso vari scenari di CT hanno prodotto le seguenti osservazioni:

• Comportamento di Semiconvergenza: I metodi GKB proposti hanno dimostrato una maggiore resilienza contro la semiconvergenza rispetto ad AB-GMRES e BA-GMRES. Mentre i metodi GMRES spesso divergevano rapidamente dopo poche iterazioni (ad esempio, 7–10 iterazioni), i metodi GKB si sono stabilizzati per una durata maggiore (ad esempio, 40–100 iterazioni) prima che l'amplificazione del rumore diventasse dominante.
• Qualità della Ricostruzione: In termini di Errore di Ricostruzione Relativo (RRE), i metodi GKB hanno generalmente ottenuto tassi di errore inferiori rispetto ai metodi GMRES, in particolare nelle iterazioni avanzate, a condizione che l'iterazione venisse interrotta correttamente.
• Costo Computazionale: Coerentemente con l'analisi teorica, i metodi GMRES sono stati più veloci in termini di tempo di esecuzione totale. Per un numero fisso di iterazioni, GMRES richiedeva meno FLOPS. Nello specifico, AB-GMRES richiedeva circa i 3/5 del costo di AB-GKB, e BA-GMRES richiedeva circa la metà del costo di BA-GKB.
• Criteri di Arresto: Il Principio della Discrepanza (DP) si è rivelato superiore a RNS. DP ha interrotto costantemente le iterazioni vicino all'RRE ottimale per tutti i metodi, mentre RNS spesso terminava troppo presto o troppo tardi a seconda dell'algoritmo, portando a una qualità di ricostruzione subottimale.
• Dimensionalità del Sistema: L'efficienza relativa delle formulazioni AB rispetto a BA dipendeva dal fatto che il sistema fosse sovradeterminato ($m > n$) o sottodeterminato ($m < n$). Nei casi sovradeterminati, i metodi basati su BA erano generalmente più veloci; nei casi sottodeterminati, i metodi basati su AB erano più costosi.

Significatività e Rivendicazioni
Il documento sostiene che, sebbene gli algoritmi AB-GKB e BA-GKB comportino un costo computazionale per iterazione più elevato rispetto agli approcci basati su GMRES, essi offrono un vantaggio distinto in termini di stabilità e qualità della soluzione per problemi di proiettore non corrispondente. Sfruttando la scarsità della matrice bidiagonale generata da GKB, i metodi forniscono un'alternativa robusta che ritarda l'insorgere della semiconvergenza, permettendo ricostruzioni più stabili senza sacrificare il potenziale di accuratezza dei metodi di sottospazio di Krylov.

Gli autori sottolineano che l'efficacia di questi metodi dipende fortemente dalla scelta di un criterio di arresto; il Principio della Discrepanza è identificato come il meccanismo più affidabile per terminare le iterazioni prima che la corruzione da rumore degradi la soluzione. Il lavoro conclude che questi algoritmi sono efficaci per fornire soluzioni stabili in problemi inversi di CT su larga scala dove proiettori corrispondenti non sono disponibili, sebbene il lavoro futuro possa concentrarsi sull'integrazione di tecniche di regolarizzazione avanzate per migliorare ulteriormente l'efficienza.