Metrics Spectral

spectral

FUNCTION	DESCRIPTION
spectral_gap	Compute spectral gap between leading eigenvalues.
spectral_gap_preservation	Compare spectral gaps between source and target eigenvalue spectra.
eigenvalue_reconstruction_error	Compute reconstruction error between eigenvalues and target eigenvalues.
eigenvector_alignment	Measure alignment between dominant eigenspaces from two methods.

Functions

spectral_gap(eigenvalues: torch.Tensor, n_states: int | None = None, relative: bool = False, eps: float = torch.finfo(torch.float32).eps) -> torch.Tensor

Compute spectral gap between leading eigenvalues.

The spectral gap is the difference between the n_states-th eigenvalue and the (n_states + 1)-th eigenvalue. A larger spectral gap indicates better separation between dominant and subdominant eigenspaces, which is desirable for identifying metastable states or dominant structures.

PARAMETER	DESCRIPTION
eigenvales	Eigenvalues in descending order, shape (n_eigenvalues, ). TYPE: Tensor
n_states	Index for gap computation. If None, finds maximum gap. TYPE: int | None, optional, by default None DEFAULT: None
relative	If True, return relative gap (gap / eigval[n_states - 1]). TYPE: bool, optional, by default False DEFAULT: False
eps	Small value to avoid division by zero. TYPE: float, optional, by default torch.finfo(torch.float32).eps DEFAULT: eps

RETURNS	DESCRIPTION
Tensor	Spectral gap value.

Examples:

>>> eigval = torch.tensor([1.0, 0.9, 0.7, 0.1, 0.05])
>>> spectral_gap(eigval, n_states=3)
0.6
>>> spectral_gap(eigval, n_states=3, relative=True)
0.857...

Find maximum spectral gap:

>>> spectral_gap(eigval, n_states=None)
0.6

Source code in spectre/metrics/spectral.py

def spectral_gap(
    eigenvalues: torch.Tensor,
    n_states: int | None = None,
    relative: bool = False,
    eps: float = torch.finfo(torch.float32).eps,
) -> torch.Tensor:
    """
    Compute spectral gap between leading eigenvalues.

    The spectral gap is the difference between the `n_states`-th eigenvalue
    and the (`n_states` + 1)-th eigenvalue. A larger spectral gap indicates
    better separation between dominant and subdominant eigenspaces, which is
    desirable for identifying metastable states or dominant structures.

    Parameters
    ----------
    eigenvales : torch.Tensor
        Eigenvalues in descending order, shape (n_eigenvalues, ).

    n_states : int | None, optional, by default None
        Index for gap computation. If None, finds maximum gap.

    relative : bool, optional, by default False
        If True, return relative gap (gap / eigval[n_states - 1]).

    eps : float, optional, by default torch.finfo(torch.float32).eps
        Small value to avoid division by zero.

    Returns
    -------
    torch.Tensor
        Spectral gap value.

    Examples
    --------
    >>> eigval = torch.tensor([1.0, 0.9, 0.7, 0.1, 0.05])
    >>> spectral_gap(eigval, n_states=3)
    0.6
    >>> spectral_gap(eigval, n_states=3, relative=True)
    0.857...

    Find maximum spectral gap:

    >>> spectral_gap(eigval, n_states=None)
    0.6
    """
    check_1d(eigenvalues)

    if not isinstance(eps, float):
        raise TypeError(f"Expected `eps` to be a float, got {type(eps).__name__}.")
        check_in_interval(eps, "[0, inf)")

    if eigenvalues.numel() < 2:
        raise ValueError(f"Expected at least 2 eigenvalues, got {eigenvalues.numel()}.")

    if n_states is None:
        # Find index with maximum gap
        gaps = eigenvalues[:-1] - eigenvalues[1:]
        n_states = torch.argmax(gaps).item() + 1
        gap = gaps[n_states - 1]
    else:
        if not isinstance(n_states, int):
            raise TypeError(
                f"Expected `n_states` to be an integer, got {type(n_states).__name__}."
            )
        check_in_interval(n_states, f"[1, {len(eigenvalues)})")

        gap = eigenvalues[n_states - 1] - eigenvalues[n_states]

    if relative:
        if eigenvalues[n_states - 1].abs() < eps:
            raise ValueError(
                "Cannot compute relative gap: leading eigenvalue is near zero."
            )
        gap = gap / eigenvalues[n_states - 1]

    return gap

spectral_gap_preservation(eigenvalues: torch.Tensor, eigenvalues_target: torch.Tensor, n_states: int | None = None, reduce: str = 'relative_error', eps: float = torch.finfo(torch.float32).eps) -> torch.Tensor

Compare spectral gaps between source and target eigenvalue spectra.

Evaluates how well target preserves the spectral gap structure of source.

PARAMETER	DESCRIPTION
eigenvalues	Source eigenvalues, shape (n_eigen, ). TYPE: Tensor
eigenvalues_target	Target eigenvalues, shape (n_eigen, ). TYPE: Tensor
n_states	Index for gap computation. If None, uses maximum gap from source. TYPE: int | None, optionalby default None DEFAULT: None
reduce	Comparison reduce: "relative_error", "absolute_error", "correlation". TYPE: str, optional, by default "relative_error" DEFAULT: 'relative_error'
eps	Small value to avoid division by zero. TYPE: float, optional, by default torch.finfo(torch.float32).eps DEFAULT: eps

RETURNS	DESCRIPTION
Tensor	Preservation reduce value.

Examples:

>>> eigenvalues = torch.tensor([1.0, 0.9, 0.7, 0.1])
>>> eigenvalues_target = torch.tensor([0.98, 0.88, 0.68, 0.12])
>>> spectral_gap_preservation(eigenvalues, eigenvalues_target, n_states=3)
0.033...

Using correlation reduce:

>>> spectral_gap_preservation(
...     eigenvalues,
...     eigenvalues_target,
...     n_states=3,
...     reduce="correlation",
... )
0.999...

Source code in spectre/metrics/spectral.py

def spectral_gap_preservation(
    eigenvalues: torch.Tensor,
    eigenvalues_target: torch.Tensor,
    n_states: int | None = None,
    reduce: str = "relative_error",
    eps: float = torch.finfo(torch.float32).eps,
) -> torch.Tensor:
    """
    Compare spectral gaps between source and target eigenvalue spectra.

    Evaluates how well target preserves the spectral gap structure of source.

    Parameters
    ----------
    eigenvalues : torch.Tensor
        Source eigenvalues, shape (n_eigen, ).

    eigenvalues_target : torch.Tensor
        Target eigenvalues, shape (n_eigen, ).

    n_states : int | None, optionalby default None
        Index for gap computation. If None, uses maximum gap from source.

    reduce : str, optional, by default "relative_error"
        Comparison reduce: "relative_error", "absolute_error", "correlation".

    eps : float, optional, by default torch.finfo(torch.float32).eps
        Small value to avoid division by zero.

    Returns
    -------
    torch.Tensor
        Preservation reduce value.

    Examples
    --------
    >>> eigenvalues = torch.tensor([1.0, 0.9, 0.7, 0.1])
    >>> eigenvalues_target = torch.tensor([0.98, 0.88, 0.68, 0.12])
    >>> spectral_gap_preservation(eigenvalues, eigenvalues_target, n_states=3)
    0.033...

    Using correlation reduce:

    >>> spectral_gap_preservation(
    ...     eigenvalues,
    ...     eigenvalues_target,
    ...     n_states=3,
    ...     reduce="correlation",
    ... )
    0.999...
    """
    check_1d(eigenvalues)
    check_1d(eigenvalues_target)
    check_same_len(eigenvalues, eigenvalues_target)
    check_same_device(eigenvalues, eigenvalues_target)

    gap = spectral_gap(eigenvalues, n_states=n_states, relative=False)
    gap_target = spectral_gap(eigenvalues_target, n_states=n_states, relative=False)

    if reduce == "relative_error":
        if abs(gap) < eps:
            raise ValueError("Cannot compute relative error: source gap is near zero.")
        return abs(gap_target - gap) / abs(gap)

    elif reduce == "absolute_error":
        return abs(gap_target - gap)

    elif reduce == "correlation":
        # Compute correlation of all spectral gaps
        gaps_source = eigenvalues[:-1] - eigenvalues[1:]
        gaps_target = eigenvalues_target[:-1] - eigenvalues_target[1:]
        return torch.corrcoef(torch.stack([gaps_source, gaps_target]))[0, 1]

    else:
        raise ValueError(
            f"Unknown reduce '{reduce}'. Options: 'relative_error', "
            f"'absolute_error', 'correlation'."
        )

eigenvalue_reconstruction_error(eigenvalues: torch.Tensor, eigenvalues_target: torch.Tensor, reduce: str = 'mse', weights: torch.Tensor | None = None, eps: float = torch.finfo(torch.float32).eps) -> torch.Tensor

Compute reconstruction error between eigenvalues and target eigenvalues.

PARAMETER	DESCRIPTION
eigenvalues	Predicted eigenvalues, shape (n_eigen, ). TYPE: Tensor
eigenvalues_target	True/reference eigenvalues, shape (n_eigen, ). TYPE: Tensor
reduce	Reduction operation: "mse", "mae", "relative", "correlation". TYPE: str, optional, by default "mse" DEFAULT: 'mse'
weights	Optional weights for eigenvalues (higher weight for leading eigenvalues), shape (n_eigen, ). TYPE: torch.Tensor | None, optional, by default None DEFAULT: None
eps	Small value to avoid division by zero. TYPE: float, optional, by default torch.finfo(torch.float32).eps DEFAULT: eps

RETURNS	DESCRIPTION
Tensor	Reconstruction error value.

Examples:

>>> eig = torch.tensor([0.95, 0.85, 0.65, 0.15])
>>> eig_target = torch.tensor([1.0, 0.9, 0.7, 0.1])
>>> eigenvalue_reconstruction_error(eig, eig_target, reduce="mse")
0.00625

With weights emphasizing leading eigenvalues:

>>> weights = torch.tensor([4.0, 2.0, 1.0, 0.5])
>>> eigenvalue_reconstruction_error(eig, eig_target, weights=weights)
0.00833...

Source code in spectre/metrics/spectral.py

def eigenvalue_reconstruction_error(
    eigenvalues: torch.Tensor,
    eigenvalues_target: torch.Tensor,
    reduce: str = "mse",
    weights: torch.Tensor | None = None,
    eps: float = torch.finfo(torch.float32).eps,
) -> torch.Tensor:
    """
    Compute reconstruction error between eigenvalues and target eigenvalues.

    Parameters
    ----------
    eigenvalues : torch.Tensor
        Predicted eigenvalues, shape (n_eigen, ).

    eigenvalues_target : torch.Tensor
        True/reference eigenvalues, shape (n_eigen, ).

    reduce : str, optional, by default "mse"
        Reduction operation: "mse", "mae", "relative", "correlation".

    weights : torch.Tensor | None, optional, by default None
        Optional weights for eigenvalues (higher weight for leading eigenvalues),
        shape (n_eigen, ).

    eps : float, optional, by default torch.finfo(torch.float32).eps
        Small value to avoid division by zero.

    Returns
    -------
    torch.Tensor
        Reconstruction error value.

    Examples
    --------
    >>> eig = torch.tensor([0.95, 0.85, 0.65, 0.15])
    >>> eig_target = torch.tensor([1.0, 0.9, 0.7, 0.1])
    >>> eigenvalue_reconstruction_error(eig, eig_target, reduce="mse")
    0.00625

    With weights emphasizing leading eigenvalues:

    >>> weights = torch.tensor([4.0, 2.0, 1.0, 0.5])
    >>> eigenvalue_reconstruction_error(eig, eig_target, weights=weights)
    0.00833...
    """
    check_1d(eigenvalues)
    check_1d(eigenvalues_target)
    check_same_len(eigenvalues, eigenvalues_target)
    check_same_device(eigenvalues, eigenvalues_target)

    if weights is not None:
        check_1d(weights)
        check_same_len(weights, eigenvalues)
        check_sample_weights(weights)
        weights = weights / weights.sum()

    diff = eigenvalues - eigenvalues_target

    if reduce == "mse":
        error = diff**2
        if weights is not None:
            error = error * weights
        return error.mean()

    elif reduce == "mae":
        error = diff.abs()
        if weights is not None:
            error = error * weights
        return error.mean()

    elif reduce == "relative":
        # Relative error: mean(|target - true| / |true|)
        denom = eigenvalues_target.abs().clamp(min=eps)
        error = diff.abs() / denom
        if weights is not None:
            error = error * weights
        return error.mean()

    elif reduce == "correlation":
        # Pearson correlation coefficient
        return torch.corrcoef(torch.stack([eigenvalues, eigenvalues_target]))[0, 1]

    else:
        raise ValueError(
            f"Unknown reduce operation '{reduce}'. Options: 'mse', 'mae', "
            f"'relative', 'correlation'."
        )

eigenvector_alignment(eigenvectors: torch.Tensor, eigenvectors_target: torch.Tensor, n_eigen: int | None = None, reduce: str = 'subspace_overlap') -> torch.Tensor

Measure alignment between dominant eigenspaces from two methods.

Computes how well the leading eigenvectors from a target method align with those from a source method.

PARAMETER	DESCRIPTION
eigenvectors	Source eigenvectors, shape (n_samples, n_eigenvectors). TYPE: Tensor
eigenvectors_target	Target eigenvectors, shape (n_samples, n_eigenvectors). TYPE: Tensor
n_eigen	Number of leading components to compare. If None, uses all components. TYPE: int | None, optional, by default None DEFAULT: None
reduce	Alignment reduce: "subspace_overlap", "grassmann", "canonical_angles". TYPE: str, optional, by default "subspace_overlap" DEFAULT: 'subspace_overlap'

RETURNS	DESCRIPTION
Tensor	Alignment reduce value. Higher values indicate better alignment.

Examples:

>>> eigenvectors = torch.randn(100, 5)
>>> eigenvectors_target = eigenvectors + 0.1 * torch.randn(100, 5)
>>> eigenvector_alignment(eigenvectors, eigenvectors_target, n_eigen=3)
0.95...

Using Grassmann distance:

>>> eigenvector_alignment(
...     eigenvectors,
...     eigenvectors_target,
...     n_eigen=3,
...     reduce="grassmann",
... )
0.12...

Source code in spectre/metrics/spectral.py

def eigenvector_alignment(
    eigenvectors: torch.Tensor,
    eigenvectors_target: torch.Tensor,
    n_eigen: int | None = None,
    reduce: str = "subspace_overlap",
) -> torch.Tensor:
    """
    Measure alignment between dominant eigenspaces from two methods.

    Computes how well the leading eigenvectors from a target method align
    with those from a source method.

    Parameters
    ----------
    eigenvectors : torch.Tensor
        Source eigenvectors, shape (n_samples, n_eigenvectors).

    eigenvectors_target : torch.Tensor
        Target eigenvectors, shape (n_samples, n_eigenvectors).

    n_eigen : int | None, optional, by default None
        Number of leading components to compare. If None, uses all components.

    reduce : str, optional, by default "subspace_overlap"
        Alignment reduce: "subspace_overlap", "grassmann", "canonical_angles".

    Returns
    -------
    torch.Tensor
        Alignment reduce value. Higher values indicate better alignment.

    Examples
    --------
    >>> eigenvectors = torch.randn(100, 5)
    >>> eigenvectors_target = eigenvectors + 0.1 * torch.randn(100, 5)
    >>> eigenvector_alignment(eigenvectors, eigenvectors_target, n_eigen=3)
    0.95...

    Using Grassmann distance:

    >>> eigenvector_alignment(
    ...     eigenvectors,
    ...     eigenvectors_target,
    ...     n_eigen=3,
    ...     reduce="grassmann",
    ... )
    0.12...
    """
    check_2d(eigenvectors)
    check_2d(eigenvectors_target)
    check_same_len(eigenvectors, eigenvectors_target)
    check_same_device(eigenvectors, eigenvectors_target)

    n_samples, n_eigenvectors = eigenvectors.shape

    if n_eigen is None:
        n_eigen = n_eigenvectors
    else:
        if not isinstance(n_eigen, int):
            raise TypeError(
                f"Expected `n_eigen` to be an integer, got {type(n_eigen).__name__}."
            )
        check_in_interval(n_eigen, f"[1, {n_eigenvectors}]")

    # Select leading components
    U_source = eigenvectors[:, :n_eigen]
    U_target = eigenvectors_target[:, :n_eigen]

    # Orthonormalize the subspaces using QR decomposition
    U_source, _ = torch.linalg.qr(U_source)
    U_target, _ = torch.linalg.qr(U_target)

    # Compute U^T V (cross-correlation matrix)
    M = U_source.T @ U_target  # Shape: (n_eigen, n_eigen)

    if reduce == "subspace_overlap":
        # Frobenius norm of M (higher = better alignment)
        # For orthonormal bases: Perfect alignment: ||M||_F = sqrt(n_eigen)
        overlap = torch.norm(M, p="fro") / (n_eigen**0.5)
        return overlap

    elif reduce == "grassmann":
        # Grassmann distance based on canonical angles
        # Compute singular values of M
        try:
            s = torch.linalg.svdvals(M)
            # Clamp to [0, 1] to avoid numerical issues with arccos
            s = torch.clamp(s, 0.0, 1.0)
            # Grassmann distance: sqrt(sum(arccos(s_i)^2))
            distance = torch.sqrt(torch.sum(torch.acos(s) ** 2))
            return distance
        except RuntimeError:
            raise ValueError("Failed to compute SVD for Grassmann distance.")

    elif reduce == "canonical_angles":
        # Return average canonical angle (in radians)
        try:
            s = torch.linalg.svdvals(M)
            # Clamp to [0, 1] to avoid numerical issues with arccos
            s = torch.clamp(s, 0.0, 1.0)
            return torch.acos(s).mean()
        except RuntimeError:
            raise ValueError("Failed to compute SVD for canonical angles.")

    else:
        raise ValueError(
            f"Unknown reduce '{reduce}'. Options: 'subspace_overlap', "
            f"'grassmann', 'canonical_angles'."
        )