# Linear Operators¶

We provide a collection of linear operators with efficient JAX based implementations that are relevant in standard signal/image processing problems. We also provide a bunch of utilities to combine and convert linear operators.

This module is inspired by pylops although the implementation approach is different.

A linear operator $$T : X \to Y$$ connects a model space $$X$$ to a data space $$Y$$.

A linear operator satisfies following laws:

(1)$T (x + y) = T (x) + T (y)$
(2)$T (\alpha x) = \alpha T(x)$

Thus, for a general linear combination:

(3)$T (\alpha x + \beta y) = \alpha T (x) + \beta T (y)$

We are concerned with linear operators $$T : \mathbb{F}^n \to \mathbb{F}^m$$ where $$\mathbb{F}$$ is either the field of real numbers or complex numbers. $$X = \mathbb{F}^n$$ is the model space and $$Y = \mathbb{F}^m$$ is the data space.

Such a linear operator can be represented by a two dimensional matrix $$A$$.

The forward operation is given by:

(4)$y = A x.$

The corresponding adjoint operation is given by:

(5)$\hat{x} = A^H y$

We represent a linear operator by a pair of functions times and trans. The times function implements the forward operation while the trans function implements the adjoint operation.

An inverse problem consists of computing $$x$$ given $$y$$ and $$A$$.

1D, 2D, ND operators

• A 1D operator takes a 1D array as input and returns a 1D array as output. E.g. identity, pad, etc.

• A 2D operator takes a 2D array as input and returns a 2D array as output. E.g. conv2D, dwt2D, etc.

• An ND operator takes an ND array as input and returns an ND array as output.

The vectors may be stored using multi-dimensional arrays in memory. E.g., images are usually stored in 2-3 dimensional arrays. The operators themselves may work directly on multi-dimensions arrays. E.g. a 2D convolution operator can be applied directly to an image to result in another image.

In other words, the vectors from the model space as well as data space may be stored in memory using 1D,2D,…,ND array. They should still be treated as vectors for the purposes of this module.

This is a departure from pylops convention where input to a linear operator must be flattened into a 1D array and output needs to be reshaped again. In this library, the input and output to a 2D linear operator would be a 2D array.

axis parameter in a 1D linear operator

• A 1D linear operator may get an ND array as input.

• In this case, the axis parameter to the operator specifies the axis along which the linear operator is to be applied.

• The input ND array will be broken into slices of 1D arrays along the axis and the linear operator will be applied separately to each slice.

• Then the slices will be combined to generate the output ND array.

• E.g. if the input is a matrix then:

• axis=0 means apply the linear operator over each column (along axis=0)

• axis=1 means apply the linear operator over each row (along axis=1)

This is based on the convention followed by numpy.apply_along_axis.

## Data types¶

 Operator Represents a finite linear operator $$T : A -> B$$ where $$A$$ and $$B$$ are finite vector spaces.

## Basic operators¶

 identity(in_dim[, out_dim, axis]) Returns an identity linear operator from model space to data space matrix(A[, axis]) Converts a matrix into a linear operator Converts a real matrix into a linear operator sparse_real_matrix(A[, axis]) Converts a sparse real matrix into a linear operator scalar_mult(alpha, n) Returns a linear operator T such that $$T v = \alpha v$$ diagonal(d[, axis]) Returns a linear operator which mimics multiplication by a diagonal matrix zero(in_dim[, out_dim, axis]) Returns a linear operator which maps everything to 0 vector in data space Returns an operator which flips the order of entries in input upside down Returns an operator which computes the sum of a vector dot(v[, adjoint, axis]) Returns a linear operator T such that $$T x = \langle v , x \rangle = v^H x$$ pad_zeros(n, before, after) Adds zeros before and after a vector. An operator which constructs a symmetric vector by pre-pending the input in reversed order restriction(n, I[, axis]) An operator which computes y = x[I] over an index set I reshape(in_shape, out_shape[, order]) Returns a linear operator which reshapes vectors from model space to data space arr2vec(shape) Returns a linear operator which reshapes arrays to vectors heaviside(n[, axis, normalized]) Returns a linear operator implements the Heaviside step function inv_heaviside(n[, axis, normalized]) Returns a linear operator that computes the inverse of Heaviside/cumsum on input

## Operator calculus¶

It is possible to combine one or more linear operators to create new linear operators. The functions in this section provide different ways to combine linear operators.

 Returns the negative of a linear operator $$T = -A$$ scale(A, alpha) Returns the linear operator $$T = \alpha A$$ for the operator $$A$$ add(A, B) Returns the sum of two linear operators $$T = A + B$$ subtract(A, B) Returns a linear operator $$T = A - B$$ compose(A, B[, ignore_compatibility, shape]) Returns the composite linear operator $$T = AB$$ such that $$T(x)= A(B(x))$$ Returns the transpose of a given operator $$T = A^T$$ Returns the adjoint of a given operator $$T = A^H$$ hcat(A, B) Returns the linear operator $$T = [A \, B]$$ power(A, p) Returns the linear operator $$T = A^p$$ block_diag(operators[, axis]) Returns a block diagonal operator from 2 or more operators windowed_op(m, T[, overlap, axis]) A wrapper to convert an operator into an overcomplete windowed operator

## Signal processing operators¶

 running_average(n, length[, axis]) Computes a running average of entries in x fir_filter(n, h[, axis]) Implements an FIR filter defined by coeffs convolve(n, h[, offset, axis]) Implements a convolution operator with the filter h convolve2D(shape, h[, offset, axes]) Performs 2 dimensional convolution on the input array convolveND(shape, h[, offset, axes]) Performs N dimensional convolution on input array

## Orthonormal transforms and bases¶

 fft(m[, n, axis, mode]) 1D FFT operator dwt(n[, wavelet, level, axis, basis]) Returns a 1D Discrete Wavelet Transform operator dwt2D(shape[, wavelet, level, axes, basis]) Returns a 2D Discrete Wavelet Transform operator Returns an operator which represents the DFT orthonormal basis Returns an operator which represents the DCT-II orthonormal basis Returns an operator which represents the Walsh Hadamard Transform Basis

## Unions of bases¶

 Returns an operator for a two-ortho basis dictionary consisting of Dirac basis and Fourier basis

## Random compressive sensing operators¶

 gaussian_dict(key, m[, n, normalize_atoms, axis]) An operator which represents a Gaussian sensix matrix (with normalized columns) rademacher_dict(key, m[, n, …]) An operator which represents a Rademacher sensing matrix random_onb_dict(key, n[, axis]) An operator representing a random orthonormal basis random_orthonormal_rows_dict(key, m[, n, axis]) An operator whose rows are orthonormal (sampled from a random orthonormal basis) sparse_binary_dict(key, m[, n, d, …]) An operator which represents a sparse binary sensing matrix

## Operators for special matrices¶

 circulant(n, c[, axis]) Circulant matrix operator

## Derivatives (finite differences)¶

 first_derivative(n[, dx, kind]) Computes the first derivative second_derivative(n[, dx]) tv(n[, kind, axis]) Returns a total variation linear operator for 1D signals tv2D(shape[, kind]) Returns a total variation linear operator for 2D images

## Convenience operators¶

These operators are technically not linear on $$\mathbb{F}^n \to \mathbb{F}^m$$

 Returns the real parts of a vector of complex numbers

## Operator parts¶

 column(T, i) Returns the i-th column of the operator T columns(T, indices) Returns the i-th column of the operator T

## Properties of a linear operator¶

These are still experimental and not efficient.

 normest(operator[, max_iters, error_tolerance]) Estimates the norm of a linear operator by power method normest_jit(operator[, max_iters, …]) Estimates the norm of a linear operator by power method Computes the upper frame bound for a linear operator using full SVD

## Utilities¶

 jit(operator) Returns the same linear operator with compiled times and trans functions Converts a linear operator to a matrix Converts the adjoint of a linear operator to a matrix Converts a linear operator to a matrix in complex numbers