Low-Rank Hessian Approximation

hessian

Wrapper for low-rank Hessian approximation in Hippylib.

This module provides a wrapper for the corresponding functionality in Hippylib. Configuration is done via dataclasses, all input and output vectors are numpy arrays.

Internally, the module uses the double-pass algorithm in Hippylib, to compute a randomized truncated SVD of the Hessian operator of the provided Hippylib inference model.

Classes:

Name	Description
LowRankHessianSettings	Configuration for the low-rank computation

Functions:

Name	Description
compute_low_rank_hessian	Compute the low-rank Hessian approximation of the inference model posterior at a given evaluation point.

spin.hippylib.hessian.LowRankHessianSettings dataclass

Configuration for computation of the low-rank Hessian approximation.

Attributes:

Name	Type	Description
inference_model	hl.Model	Hippylib inference model to use for computations.
num_eigenvalues	int	Number of eigenvalues to compute.
num_oversampling	int	Number of oversampling eigenvalues (for numerical stability).
gauss_newton_approximation	bool	Whether to use the Gauss-Newton approximation.

spin.hippylib.hessian.compute_low_rank_hessian

compute_low_rank_hessian(settings: LowRankHessianSettings, evaluation_point: Iterable[npt.NDArray[np.floating], npt.NDArray[np.floating], npt.NDArray[np.floating]]) -> tuple[npt.NDArray[np.floating], Iterable[npt.NDArray[np.floating]]]

Compute a low-rank Hessian approximation of the inference model posterior.

The Hessian is computed about the given evaluation point, based on randomized truncated SVD. More precisely, we consider the case where the Hessian is the second parametric derivative of a negative posterior functional, defined through a Gaussian prior and noise model,

\[ \pi_{post}(\mathbf{m}) \propto \exp\Big( -\frac{1}{2}||\mathbf{F}(m)-\mathbf{d}_{obs}||_{R_{\text{noise}}}^2 -\frac{1}{2}||\mathbf{m}-\mathbf{m}_{pr}||_{R_{\text{prior}}}^2\Big). \]

Consequently, the Hessian can be divided into a misfit and a prior term,

\[ \mathbf{H}(\mathbf{m}) = \mathbf{H}_{\text{misfit}}(\mathbf{m}) + \mathbf{R}_{\text{prior}} \in \mathbb{R}^{N\times N}. \]

To obtain an approximation of \(\mathbf{H}\) at the point \(\mathbf{m}\), we solve a generalized eigenvalue problem for the first \(r\) eigenvalues and eigenvectors of the misfit Hessian,

\[ \mathbf{H}_{\text{misfit}}\mathbf{v}_i = \lambda_i\mathbf{R}_{\text{prior}}\mathbf{v}_i, \quad \lambda_1 \geq \ldots \geq \lambda_r. \]

The low-rank approximation of the Hessian can then be constructed as

\[ \begin{gather*} \mathbf{H} = \mathbf{R}_{\text{prior}} + \mathbf{R}_{\text{prior}}\mathbf{V}_r\mathbf{D}_r\mathbf{V}_r^T\mathbf{R}_{\text{prior}} + \mathcal{O}\big(\sum_{i=r+1}^N\lambda_i\big), \\ \mathbf{V}_r = [\mathbf{v}_1,\ldots,\mathbf{v}_r], \\ \mathbf{D}_r = \text{diag}(\lambda_1,\ldots,\lambda_r). \end{gather*} \]

The inverse can be appoximated with the Sherman-Morrison-Woodbury formula. Apparently, the approximation is good if the trailing eigenvalues \(\lambda_i, i=r+1,\ldots,N\) are small, corresponding to compactness of the Hessian operator.

Parameters:

Name	Type	Description	Default
settings	LowRankHessianSettings	Configuration for the low-rank computation	required
evaluation_point	tuple[npt.NDArray, npt.NDArray, npt.NDArray]	The evaluation point for the Hessian.	required

Returns:

Type	Description
tuple[npt.NDArray[np.floating], Iterable[npt.NDArray[np.floating]]]	tuple[npt.NDArray[np.floating], Iterable[npt.NDArray[np.floating]]]: Arrays for the computed eigenvalues and eigenvectors