hofa.rgz

This module contains the Numpy implementation for performing higher-order Fourier analysis using the spectral approach developed by Candela, González-Sánchez, and Szegedy.

Conventions

The following conventions are used for functions and data structures in this file:

• Any finite abelian group \(Z\) is isomorphic to \(Z/n_1Z \times Z/n_2Z \times ... \times Z/n_kZ\). Therefore, a function \(f: Z \to \mathbb{C}\) is represented as a numpy.ndarray tensor with shape (n_1, n_2, ..., n_k).

• A \(Z\)-matrix, where \(Z = Z/n_1Z \times ... \times Z/n_kZ\), is represented as a numpy.ndarray tensor with shape (n_1, ..., n_k, n_1, ..., n_k).

Classes

LayerRegularizer	Abstract base class for an abstract cumulative distribution function.
StandardLayerRegularizer	An implementation of the LayerRegularizer class with common functionality.
RegularizationResult	A class containing the output of the regularition algorithm

Functions

regularize(→ RegularizationResult)	Perform the regularization of a function using certain LayerRegularizer s.
regularize_after_setup(f, per_layer_regularizers)	Internal process that implements the main steps of the regularization.
regularize_and_decompose(f, K, layer_regularizer)	Compute the eigenvalues and eigenvectors of a regularized Z-matrix.
weighed_projection(f, eigvals, eigvects, layer_regularizer)	Projects a function to certain eigenspaces with weights depending on eigenvalues.
default_reg_list(k[, rng])	Returns a list of pre-defined regularizers.
move_diag_to_rows(M)	Creates a new Z-matrix where each row is a Z-diagonal of M.
move_rows_to_diag(M)	Reconstructs the Z-matrix from its row representation.
apply_on_diag(M, K)	Apply a regularizer to all Z-diagonals of a Z-matrix.
t_prod_itself(g)	Compute the outer product of a Z-function with itself conjugate.

Module Contents

class hofa.rgz.LayerRegularizer

Bases: abc.ABC

Abstract base class for an abstract cumulative distribution function.

This class will contain all necessary information to project a function onto the eigenspace of a certain matrix where the projection to the different eigenspaces is weighted according to an LayerRegularizer.alpha() method.

The user must instantiate this class (or use the class StandardLayerRegularizer provided by this module).

abstractmethod setup(f: numpy.ndarray, depth: int, layer: int)

Method to be called once at the beginning of the regularization process.

As input, it takes the function f to regularize, the total order depth of regularization, (1 for linear/classical Fourier analysis, 2 for quadratic, etc.), and the position layer of this regularizer in the algorithm. That is, layer will be an integer between 0 and depth-1, representing which part of the regularization process this instance of LayerRegularizer will perform. For example, 0 means Fourier regularization (acting on the squared norms of the Fourier coefficients of some multiplicative derivative), 1 means quadratic, etc.

Parameters:

• f (np.ndarray) – A anrray containing the function f to regularize.

• depth (int) – An integer representing the total order of regularization. (1 for Fourier regularization, 2 for quadratic Fourier, etc.)

• layer (int) – An integer specifying the position of this regularizer in the process. It will be an integer between 0 and depth-1. For example, 0 means Fourier regularization, 1 means quadratic regularization, etc.

Returns:

None

Raises:

NotImplementedError – If the method is not overridden in a subclass.

abstractmethod update(f: numpy.ndarray, height: int)

Update the object with the provided parameters at each step of the regularization.

This is an abstract method that should be implemented by subclasses to update the state of the LayerRegularizer object at each step of the regularization process. Specifically, for an instance of LayerRegularizer where the setup() method was called with depth = k and layer = l, the update() method will be called k-l times, with height between k and l+1. At each call, the implementation of the update() method can modify the internal state of the LayerRegularizer class using the provided function f.

For example, when doing cubic regularization (k=3), the setup() method will be called with layer=1 and the original function to regularize. In this case, the update() method will be called |Z| + 1 times, once with height = 3 (with f being the original function to regularize \(f_0\)), and |Z| times with height=2 (with f being \(\Delta_t f_0\) for \(t \in Z\)). The purpose of calling the method with these different heights is to keep track of all multiplicative derivatives made during the regularization process.

Parameters:

• f (np.ndarray) – A array representing a multiplicative derivative \(\Delta_{t_1}\cdots\Delta_{t_{k-h}} f_0\), where \(f_0\) is the original function to regularize. At each step, f changes according to the value of height and the regularization depth.

• height (int) – An integer representing the number of multiplicative derivatives taken.

Returns:

None

Raises:

NotImplementedError – If the method is not overridden in a subclass.

abstractmethod alpha(eigs: numpy.ndarray)

This is the regularizer that will be used to weight the different eigenspaces.

The implementation of any class can use the information about the specific function to regularize, see StandardLayerRegularizer

Parameters:

eigs (np.ndarray) – An array of eigenvalues to obtain the weights of the different eigenspaces.

Returns:

An array of the same shape as eigs with the weights of the different eigenspaces.

Return type:

np.ndarray

Raises:

NotImplementedError – If the method is not overridden in a subclass.

abstractmethod eigh(M: numpy.ndarray)

A method to obtain a number of eigenvalues and eigenvectors.

As the LayerRegularizer.alpha() function typically cancels the contribution of eigenspaces with small eigenvalues, there is usually no need to compute the full eigendecomposition of the matrix M. Thus, instead of using a typical method such as numpy.linalg.eigh, which would compute the full eigendecomposition, the user may define a custom method which is more efficient.

Parameters:

M (np.ndarray) – A self-adjoint matrix to perform regularization on.

Returns:

A tuple containing the eigenvalues and eigenvectors of the matrix M. The format of this tuple follows the numpy standard of ordering the eigenvalues increasingly and returning the corresponding i-th eigenvector indexed by [...,i].

Return type:

tuple of (np.ndarray, np.ndarray)

Raises:

NotImplementedError – If the method is not overridden in a subclass.

__deepcopy__(memo)

General-purpose deepcopy implementation for ABC.

Handles both __dict__-based and __slots__-based attributes, and preserves cycles and shared references.

Parameters:

memo (dict) – Dictionary of already-copied objects to handle recursion.

Returns:

A deep copy of this object.

Return type:

LayerRegularizer

copy()

Create a deep copy of this object.

This method returns a new instance of the concrete subclass, with all attributes deep-copied. Cycles and shared references are handled correctly.

Returns:

A new independent deep copy of this object.

Return type:

Base

Note:

Relies on the object’s __deepcopy__ implementation. Concrete subclasses with special internal state (e.g., RNGs) will automatically have their custom deepcopy logic executed.

class hofa.rgz.StandardLayerRegularizer(alpha_method: Callable = hofa.cdf.relu, param: float | int = 1.2, mode: str = 'dynamic-original', lin_alg_method: str = 'sparse', num_eigen: str | int = 'dynamic', rng: int | numpy.random.Generator | None = None)

Bases: LayerRegularizer

An implementation of the LayerRegularizer class with common functionality.

This class contains the implementation of a LayerRegularizer with the functionality described in the paper.

Parameters:

• alpha_method (Callable) –

  A function that defines the alpha calculation method. By default, it uses hofa.cdf.relu() as the method. Here we can use any user-defined method that has a signature of the form function( f : numpy.ndarray , coefficient : float | int) -> numpy.ndarray. Examples of admisible functions are:

  • hofa.cdf.relu()

  • hofa.cdf.cutoff()

  • hofa.cdf.avg_cutoff()

  • hofa.cdf.top_eig()

  • hofa.cdf.u2_dual()

• param (float | int, optional) – A floating-point or integer parameter used for various purposes depending on the mode below. Default value is 1.2.

• mode (str, optional) –

  A string representing the mode of operation. This determines how the alpha_method will interpret parameters and perform regularization. For the rest of this explanation \(k\) will denote the order or regularization (1 classical Fourier, 2 quadratic, etc.), \(f_0\) will be the original function to regularize, \(f_1\) will denote a multiplicative derivative, \(f_2\) two multiplicative derivatives, etc. Note that and an StandardLayerRegularizer object of height \(0 \le j < k\) will perform regularization of height \(j\), which means that it will act on a function \(f_{k-j-1}\), i.e., on the \(k-j-1\) multiplicative derivatives of the original function. In order to choose an appropriate parameter, the StandardLayerRegularizer object will keep track of all the typical deviations of the multiplicative derivatives made to compute \(f_{k-j-1}\). Those typical deviations will be denoted \(\sigma_0\) (for \(f_0\)), \(\sigma_1\) (for \(f_1\)) etc. Possible values:

  • 'dynamic-original': In this mode, we use the variace of the original function but modified using the fact that we are using it on a multiplicative derivative of ome order. In this mode, the parameter passed to the alpha_method will be param multiplied by \(\sigma_0^{2^{k-j-1}}\sqrt{|Z|/\log|Z|}\).

  • 'dynamic-strict': In this mode, we use the variace of the function we are regularizing. In this mode, the parameter passed to the alpha_method will be: param multiplied by \(\sigma_{k-j-1}\sqrt{|Z|/\log|Z|}\).

  • 'literal': Uses the parameter param directly.

  Default is ‘dynamic-original’.

• lin_alg_method (str, optional) –

  The method for computing the eigendecomposition of a matrix. Possible values:

  • 'sparse': Uses scipy.sparse.linalg.eigsh to compute eigenvalues.

  • 'full': Uses numpy.linalg.eigh to compute eigenvalues.

  Default is 'sparse'.

• num_eigen (str | int, optional) –

  The number of eigenvalues and eigenvectors to compute when using 'sparse' method in lin_alg_method. Possible values:

  • 'dynamic': Computed dynamically according to a heuristic.

  • an int: Specifies the number of eigenvalues to compute.

  Default is 'dynamic'.

• rng (int | np.random.Generator | None) – Either an integer seed, a np.random.Generator, or None. It is used for deterministic behaviour of scipy.sparse.linalg.eigsh

alpha_method

param = 1.2

mode = 'dynamic-original'

lin_alg_method = 'sparse'

num_eigen = 'dynamic'

total_depth = 0

layer = 0

array_sigma_sq

group = 'no group'

group_size = 0

_initialized = False

setup(f: numpy.ndarray, depth: int, layer: int)

Setup method of all StandardLayerRegularizer in a list to perform regularization

This method will records some important information that this object will need in order to perform the regularization. Namely, it records the size of the group, the total order of regularization, the position of this regularizer in the process, and initialize an array to keep track of the variances of the different multiplicative derivatives involved in the process (to be used for computing the parameter for the StandardLayerRegularizer.alpha()).

Parameters:

• f (np.ndarray) – An array containing the original function f to regularize.

• depth (int) – An integer representing the total order of regularization, i.e., 1 for Fourier regularization, 2 for quadratic Fourier, etc.

• layer (int) – An integer specifying the position of this regularizer in the process. It will be an integer between 0 and depth-1, with 0 meaning that it will do the Fourier regularization part, 1 the quadratic, etc.

Returns:

None

_require_initialized()

Ensures that the StandardLayerRegularizer has been properly initialized before use.

This method checks whether the internal state of the instance of StandardLayerRegularizer has been initialized via the StandardLayerRegularizer.setup() method. It is intended to be called at the beginning of any method that requires access to the runtime state, such as StandardLayerRegularizer.update() or other internal operations.

Raises:

RuntimeError – If the LayerRegularizer has not been initialized with StandardLayerRegularizer.setup().

This method does not modify any internal state. Its sole purpose is to enforce the correct usage pattern: construction → setup → algorithm execution.

update(f: numpy.ndarray, height: int)

Update the object with the provided parameters at each step of the regularization.

This method records the variance of the different multiplicative derivatives to be used later in StandardLayerRegularizer.alpha(). More precisely, depth=k and layer=l be the parameters that were used to call the StandardLayerRegularizer.setup() function for this object. Let height=h and \(Z\) be the group where f is defined. Recall that this object has to perform regularization of functions of the form \(\Delta_{t_1}\cdots\Delta_{t_{k-l-1}} f_0\). Then this method will keep track of the variances of \(f_0\), \(\Delta_{t_{k-l-1}} f_0\), …, \(\Delta_{t_1}\cdots\Delta_{t_{k-l-1}} f_0\) so that this information is used when StandardLayerRegularizer.alpha() is called.

Parameters:

• f (np.ndarray) – A multiplicative derivative of the form \(\Delta_{t_1}\cdots\Delta_{t_{k-h}} f_0\).

• height (int) – An integer between l+1 and k representing the number of multiplicative derivatives taken.

Returns:

None

Return type:

None

alpha(eigs: numpy.ndarray)

This method applies the alpha_method with a parameter computed according to mode and param as explained in StandardLayerRegularizer.

Parameters:

eigs (np.ndarray) – An array of eigenvalues to obtain the weights of the different eigenspaces.

Returns:

An array of the same shape as eigs with the weights of the different eigenspaces.

Return type:

np.ndarray

eigh(M: numpy.ndarray)

A method to obtain a number of eigenvalues and eigenvectors.

Depending on lin_alg_method given in the constructor, this method will use either numpy.linalg.eigh or scipy.sparse.linalg.eigsh to compute (part of) the eigendecomposition of a matrix M.

Parameters:

M (np.ndarray) – A self-adjoint matrix to perform regularization.

Returns:

A tuple containing the eigenvalues and eigenvectors.

Return type:

tuple(np.ndarray, np.ndarray)

Note:

The returned tuple of eigenvalues and eigenvectors follow the same format as numpy.linalg.eigh.

__deepcopy__(memo)

Deepcopy implementation that creates a new RNG for the copy.

Parameters:

memo (dict) – Dictionary of already-copied objects to handle recursion.

Returns:

A new copy of this object with an independent RNG.

Return type:

StandardLayerRegularizer

class hofa.rgz.RegularizationResult

Bases: NamedTuple

A class containing the output of the regularition algorithm

This class contains the 3 outputs of the regularization algorithm. That is, the regularization of the function, the eigenvalues of the matrix \(\mathcal{K}(f \otimes \overline{f})\) in increasing order, and their corresponding eigenvectors.

Parameters:

• regularization (np.ndarray) – The result of the regularization.

• eigenvalues – A 1-dimensional array of shape (N,) with the eigenvalues in increasing order.

• eigenvectors (np.ndarray) – The corresponding eigenvectors. It has shape (f.shape, N) where N is the number of eigenvalues returned.

regularization: numpy.ndarray

eigenvalues: numpy.ndarray

eigenvectors: numpy.ndarray

hofa.rgz.regularize(f: numpy.ndarray, per_layer_regularizers: List[LayerRegularizer] | int = 2, rng: int | numpy.random.Generator | None = None) → RegularizationResult

Perform the regularization of a function using certain LayerRegularizer s.

The behaviour of this function depends on the number and type of per_layer_regularizers. If per_layer_regularizers is an integer, then it is used a default list of LayerRegularizer s given by default_reg_list() is used. Otherwise it is used the list provided by the user. It returns the result of the regularization algorithm and returns a RegularizationResult object containing the regularization of f, the eigenvalues of the matrix \(\mathcal{K}(f \otimes \overline{f})\) and its corresponding eigenvectors. The eigenvalues are a N dimensional numpy.andarray where N is the number of eigenvalues returned by the RegularizationResult.eigh() method of per_layer_regularizers[-1] and are returned in increasing order. The eigenvectors are returned as a numpy.andarray of shape (f.shape,N). The eigenvalue in position i corresponds to the eigenvector indexed at [...,i]. There is some special behaviour in the following cases.

• If len(per_layer_regularizers)=0,

  this function returns as a regularization an numpy.andarray array with shape equal to that of f, the average of f as the sole eigenvalue and the constant 1 numpy.andarray of shape f.shape as the sole eigenvector.

• If len(per_layer_regularizers)=1,

  this function returns the squares of the Fourier coefficients along with the full Fourier basis of the group \(Z\) as the eigenvalues (ordered in increasing order) and the full set of Fourier characters of \(Z\) as the eigenvectors.

Parameters:

• f (numpy.ndarray) – An array representing the function to regularize.

• per_layer_regularizers (list of LayerRegularizer objects or an integer.) – A list of LayerRegularizer objects that will be set up using the array f or an integer to use a default list of such length. The LayerRegularizer object at index i will perform regularization of order i+1. That is, per_layer_regularizers[0] will act on multiplicative derivatives of the form \(\Delta_{t_1}\cdots\Delta_{t_{k-1}} f\) where len(per_layer_regularizers)=k. Similarly, per_layer_regularizers[1] will act on multiplicative derivatives of the form \(\Delta_{t_1}\cdots\Delta_{t_{k-2}} f\), etc.

• rng (int | np.random.Generator | None) – To have predictable behaviour in case the user uses an integer in the list of regularizers.

Returns:

The regularization of f according to per_layer_regularizers and the eigenvalues and eigenvectors corresponding to the last step in the process.

Return type:

RegularizationResult

hofa.rgz.regularize_after_setup(f: numpy.ndarray, per_layer_regularizers: List[LayerRegularizer])

Internal process that implements the main steps of the regularization.

This function, which requires being called with LayerRegularizer s already set up, performs the main recursive step in the regularization algorithm.

Parameters:

• f (np.ndarray) – An array representing the function to regularize.

• per_layer_regularizers (List[LayerRegularizer]) – A list of LayerRegularizer objects.

Returns:

The result of the regularization of f, the eigenvalues, and the eigenvectors of the last regularization step.

Return type:

tuple(np.ndarray,np.ndarray,np.ndarray)

hofa.rgz.regularize_and_decompose(f: numpy.ndarray, K: Callable, layer_regularizer: LayerRegularizer)

Compute the eigenvalues and eigenvectors of a regularized Z-matrix.

This function takes a Z-function f (represented as a numpy.ndarray), computes the Z-matrix \(M = f \otimes f^*\), applies a regularizer K to the Z-diagonals of M to form the Z-matrix \(\mathcal{K}(M)\), and then returns the eigenvalues in increasing order and eigenvectors of K(M) ordered according to its corresponding eigenvalue.

Parameters:

• f (np.ndarray) – A Z-function, represented as a numpy.ndarray.

• K (Callable) – A regularizer, i.e., an operator that acts on Z-functions. The operator K must be “invariant,” meaning \(K(f(\cdot+t)) = K(f)(\cdot+t)\) and \(K(\overline{f}) = \overline{K(f)}\) If these conditions are not met, the behavior of the function is unpredictable.

• layer_regularizer (LayerRegularizer) – The regularizer that contains the method for computing the eigendecomposition of the matrix once the diagonals are regularized.

Returns:

• A numpy.ndarray of shape (total_dim,) , where total_dim is the

  product of the elements in f.shape, i.e., the size of the group Z containing the eigenvalues of the Z-matrix \(\mathcal{K}(M)\), appropriately normalized and ordered increasingly.

• A numpy.ndarray of shape (f.shape, total_dim), where total_dim

  is as before and eigenvectors[...,i] is the eigenvector corresponding to eigenvalue[i].

Return type:

tuple (np.ndarray, np.ndarray)

Note:

The operator K must satisfy the invariance conditions for predictable behavior.

hofa.rgz.weighed_projection(f: numpy.ndarray, eigvals: numpy.ndarray, eigvects: numpy.ndarray, layer_regularizer: LayerRegularizer)

Projects a function to certain eigenspaces with weights depending on eigenvalues.

This function first computes a weight for each eigenvalue, represented by eigvals. Then, it projects f to all eigenspaces. Such projections are summed in a weighted manner depending on the weights computer by the layer_regularizer variable according to its alpha method.

Parameters:

• f (numpy.ndarray) – A \(Z\)-function, represented as a numpy.ndarray. The function \(f: Z \to \mathbb{C}\) is assumed to be provided in this form.

• eigvals (np.ndarray) – A dimensional array of eigenvalues. It must have shape (N,).

• eigvects (np.ndarray) – An array of orthonormal eigenvectors. It must have shape (f.shape, N) where N is the dimension of eigvals.

• layer_regularizer (LayerRegularizer) – A regularizer to be used to weight the corresponding eigenspaces.

Returns:

A numpy.ndarray of shape equal to f.shape with the weighted projection of f to the eigenspaces weighted according to layer_regularizer.

Return type:

np.ndarray

Note:

This function is meant to be used in conjunction with eig(f,K). Undefined behaviour otherwise.

hofa.rgz.default_reg_list(k: int, rng: int | numpy.random.Generator | None = None)

Returns a list of pre-defined regularizers.

This function gives back a list with pre-defined regularizers that have been tested in practice. Namely, it returns a list of StandardLayerRegularizer where all but the two last ones are initialized with num_eigen = 6, the previous to the last one using the empty constructor, and the last one uses alpha_method = hofa.cdf.cutoff.

Parameters:

• k (int) – The length of the list. If smaller or equal to 0 is the same as passing 0.

• rng (int | np.random.Generator | None) – A random number generator for reproducibilty purposes. It is used to fix the seed of the StandardLayerRegularizer

Returns:

A list of default LayerRegularizers of length k

Return type:

List[LayerRegularizer]

hofa.rgz.move_diag_to_rows(M: numpy.ndarray)

Creates a new Z-matrix where each row is a Z-diagonal of M.

Inverse of move_rows_to_diag(). This function takes a Z-matrix M of shape (n_1, ..., n_k, n_1, ..., n_k), and returns a new matrix whose rows are the Z-diagonals of M. The output matrix has the same shape as M. The returned matrix does not share memory with the input M.

Parameters:

M (numpy.ndarray) – A Z-matrix of shape (n_1, ..., n_k, n_1, ..., n_k), where \(n_i \ge 1\). Undefined behavior occurs if this condition is not met.

Returns:

A Z-matrix of the same shape as M, where the elements of M[i, i+z] are transferred to the output position [z, i].

Return type:

numpy.ndarray

Note:

This function reverses the operation performed by move_rows_to_diag().

hofa.rgz.move_rows_to_diag(M: numpy.ndarray)

Reconstructs the Z-matrix from its row representation.

Inverse of move_diag_to_rows(). This function takes a Z-matrix M of shape (n_1, ..., n_k, n_1, ..., n_k), and returns a new matrix where the rows of M become the Z-diagonals of the returned matrix. The output matrix has the same shape as M, and the returned values do not share memory with the input matrix.

Parameters:

M (numpy.ndarray) – A Z-matrix of shape (n_1, ..., n_k, n_1, ..., n_k), where \(n_i \ge 1\). Undefined behavior occurs if this condition is not met.

Returns:

A Z-matrix of the same shape as M, such that for each pair [i, z], output[i, i+z] equals M[z, i].

Return type:

numpy.ndarray

Note:

This function reverses the operation performed by move_diag_to_rows().

hofa.rgz.apply_on_diag(M: numpy.ndarray, K: Callable)

Apply a regularizer to all Z-diagonals of a Z-matrix.

This function applies the operator K to each Z-diagonal of the input Z-matrix M. The operation is performed on the diagonals indexed by [i, i+z], where the function K is applied to the corresponding values of M[i, i+z].

Parameters:

• M (numpy.ndarray) – A Z-matrix of shape (n_1, ..., n_k, n_1, ..., n_k), where \(n_i \ge 1\). Undefined behavior occurs if this condition is not met.

• K (Callable) – A regularizer operator that acts on Z-functions. The operator K must be “invariant”.

Returns:

A Z-matrix of the same shape as M, where for each pair [i, z], output[i, i+z] equals K(M[i, i+z]).

Return type:

np.ndarray

Note:

The operator K must satisfy the invariance conditions for predictable behavior.

hofa.rgz.t_prod_itself(g: numpy.ndarray)

Compute the outer product of a Z-function with itself conjugate.

This function takes a Z-function g and returns the Z-matrix formed by computing the outer product of g with its conjugate transpose, denoted as \(g\otimes g^*\).

Parameters:

g (numpy.ndarray) – A Z-function. The shape of g must be compatible for the outer product operation.

Returns:

A Z-matrix of shape (g.shape, g.shape), representing the outer product of g and its conjugate transpose.

Return type:

numpy.ndarray

Note:

The outer product is performed element-wise, and the resulting matrix has a shape of (g.shape, g.shape).