# pylops.utils.estimators.trace_hutchpp#

Trace of linear operator using the Hutch++ method.

Returns an estimate of the trace of a linear operator using the Hutch++ method [1].

Parameters
Oppylops.LinearOperator

Linear operator to compute trace on.

nevalint, optional

Maximum number of matrix-vector products compute. Defaults to 10% of shape[1].

samplerstr, optional

Sample sketching matrices from the following distributions:

• “gaussian”: Mean zero, unit variance Gaussian.

• “rayleigh”: Sample from mean zero, unit variance Gaussian and normalize the columns.

backendstr, optional

Backend used to densify matrix (numpy or cupy). Note that this must be consistent with how the operator has been created.

Returns
tracefloat

Operator trace.

Raises
ValueError

If neval is smaller than 3.

NotImplementedError

If the sampler is not one of the available samplers.

Notes

This function follows Algorithm 1 of [1]. Let $$m$$ = shape[1] and $$k$$ = neval.

1. Sample sketching matrices $$\mathbf{S} \in \mathbb{R}^{m \times \lfloor k/3\rfloor}$$, and $$\mathbf{G} \in \mathbb{R}^{m \times \lfloor k/3\rfloor}$$, from sub-Gaussian distributions.

2. Compute reduced QR decomposition of $$\mathbf{Op}\,\mathbf{S}$$, retaining only $$\mathbf{Q}$$.

3. Return $$\operatorname{tr}(\mathbf{Q}^T\,\mathbf{Op}\,\mathbf{Q}) + \frac{1}{\lfloor k/3\rfloor}\operatorname{tr}\left(\mathbf{G}^T(\mathbf{I} - \mathbf{Q}\mathbf{Q}^T)\,\mathbf{Op}\,(\mathbf{I} - \mathbf{Q}\mathbf{Q}^T)\mathbf{G}\right)$$

Use the Rademacher sampler unless you know what you are doing.

1(1,2)

Meyer, R. A., Musco, C., Musco, C., & Woodruff, D. P. (2021). Hutch++: Optimal Stochastic Trace Estimation. In Symposium on Simplicity in Algorithms (SOSA) (pp. 142–155). Philadelphia, PA: Society for Industrial and Applied Mathematics. link