Levenberg-Marquardt

Nelle sezioni precedenti, gli algoritmi di risoluzione di sistemi non lineari sono stati distinti tra metodi di discesa del gradiente e metodi di Gauss-Newton. Per una trattazione più approfondita si rimanda a (MBT04).

Nei metodi Gauss-Newton, quando $\mathbf{J}^{\top}\mathbf{J}$ è definita positiva, il metodo fornisce sempre una direzione di discesa del costo. Tuttavia, quando $\mathbf{J}^{\top}\mathbf{J}$ diventa singolare, il metodo può diventare numericamente instabile.

La tecnica proposta da Levenberg-Marquardt cerca di combinare i punti di forza di Gauss-Newton e della discesa del gradiente, traendone vantaggio da entrambi.

L'algoritmo di Levenberg-Marquardt (LM) è una tecnica iterativa ormai considerata standard per la risoluzione di problemi non lineari multivariabili. Una descrizione dettagliata dell'algoritmo è disponibile in (Lou05).

LM può essere visto come composto da una fase iniziale di discesa del gradiente, più lenta ma stabile, seguita da un risolutore di tipo Gauss-Newton, più veloce ma meno robusto.

L'algoritmo risolve una versione modificata dell'equazione (3.49), caso particolare dell'equazione (3.35), nota come augmented normal equations:

$$\mathbf{N} \boldsymbol\delta_\beta = -\mathbf{J}_{r}^{\top} \mathbf{r}$$ (3.53)

dove $\mathbf{N} = \mathbf{H}_{S} + \mu \mathbf{I}$ e $\mu > 0$ è un fattore di attenuazione (damping factor). Quando $\mu$ è elevato, $\mathbf{N}$ è quasi diagonale e l'algoritmo si comporta come una discesa del gradiente. Quando $\mu$ è vicino a zero, LM approssima il metodo di Newton.

A differenza della ricerca lungo una linea (line search), LM implementa il concetto di trust region, adattando dinamicamente la regione entro cui si assume valida la linearizzazione del modello.

Come nel metodo di Gauss-Newton, LM sfrutta l'approssimazione dell'Hessiana:

$$\mathbf{H}_S(\boldsymbol\beta) \approx \mathbf{J}_{r}^{\top}\mathbf{J}_{r}$$ (3.54)

valida quando la loss function è quadratica.

La scelta iniziale e l'aggiornamento del parametro $\mu$ tra le iterazioni è lasciata al risolutore, e diverse strategie sono proposte in letteratura.

Una delle implementazioni più diffuse (Nie99) propone di inizializzare $\mu$ come:

$$\mu_0 = \tau \max \trace \mathbf{H}$$ (3.55)

dove $\tau$ è scelto liberamente dall'utente in base alla fiducia nella stima iniziale di $\boldsymbol\beta$.

L'aggiornamento di $\mu$ è guidato dal gain ratio $\rho$:

$$\rho = \frac{ S(\boldsymbol\beta) - S(\boldsymbol\beta + \bo... ...\mu \boldsymbol\delta_\beta + \mathbf{J}^{\top} \mathbf{r} ) }$$ (3.56)

Un valore elevato di $\rho$ indica che la linearizzazione del modello è efficace, e si può ridurre $\mu$. Viceversa, se $\rho$ è basso o negativo, $\mu$ va aumentato per avvicinarsi a un comportamento da discesa del gradiente. Quando $\rho \approx 1$, si ha una buona corrispondenza tra modello e dati.

L'aggiornamento di $\mu$ può essere gestito secondo la seguente regola: