Regularization Methods

Regularization Methods#

Overview#

Pyralysis supports a variety of regularization (penalty) functions to encode prior knowledge or promote desirable properties in reconstructed images. Regularizers can be combined with data fidelity terms in the ObjectiveFunction for flexible, composite objectives.

Note

All regularizers support masks, multi-channel images, and have efficient gradient and proximal operator implementations for use in modern optimization algorithms.

Available Regularizers#

L1 Norm Regularizer#

Promotes sparsity in the image (many pixels exactly zero).

Mathematical Formulation:

\[\ell_1(I) = \sum_i \|I_i\|\]

Typical Use Case: Sparse solutions, e.g., point sources, compressed sensing.

Implementation:

from pyralysis.optimization.terms.regularizers import L1Norm

l1_term = L1Norm(penalization_factor=0.01)

# Apply to specific channel in multi-channel image
l1_channel = L1Norm(parameter_index=0, penalization_factor=0.01)

Entropy Regularizer#

Favors similarity to a prior image (maximum entropy method).

Mathematical Formulation:

\[S = \sum_i I_i \log(I_i / G)\]

Typical Use Case: Biases solution toward a reference/prior, useful in MEM.

Implementation:

from pyralysis.optimization.terms.regularizers import Entropy

entropy_term = Entropy(
    prior_image=prior_model,
    penalization_factor=0.005
)

Tikhonov Regularizer#

Penalizes large values (L2 norm squared).

Mathematical Formulation:

\[f(I) = \frac{1}{2}\|I\|_2^2\]

Typical Use Case: Smooths the solution, discourages large pixel values.

Implementation:

from pyralysis.optimization.terms.regularizers import Tikhonov

tikhonov_term = Tikhonov(penalization_factor=0.001)

Total Flux Regularizer#

Controls the total sum of pixel values (flux).

Mathematical Formulation:

\[f(I) = \sum_k I_k\]

Typical Use Case: Enforces or penalizes total flux, e.g., for normalization.

Implementation:

from pyralysis.optimization.terms.regularizers import TotalFlux

flux_term = TotalFlux(penalization_factor=0.1)

TSV (Total Squared Variation) Regularizer#

Promotes smoothness while preserving edges (quadratic TV).

Mathematical Formulation:

\[TSV(I) = \sum_i [ (\partial_x I_i)^2 + (\partial_y I_i)^2 ]\]

Typical Use Case: Edge-preserving smoothing, denoising, super-resolution.

Implementation:

from pyralysis.optimization.terms.regularizers import TSV

tsv_term = TSV(penalization_factor=0.01)

L2 Norm Regularizer#

Penalizes the L2 norm (Euclidean norm) of the image.

Mathematical Formulation:

\[\|I\|_2 = \sqrt{\sum_i I_i^2}\]

Typical Use Case: General smoothness, different from Tikhonov (which uses L2 squared).

Implementation:

from pyralysis.optimization.terms.regularizers import L2Norm

l2_term = L2Norm(penalization_factor=0.001)

Note: This is different from Tikhonov regularization, which uses the L2 squared norm \(\frac{1}{2}\|I\|_2^2\). The L2 norm is non-differentiable at zero, while Tikhonov is differentiable everywhere.

Huber Loss Regularizer#

Robust loss function that combines L1 and L2 behavior, providing a smooth transition.

Mathematical Formulation:

\[\begin{split}L_\delta(x) = \begin{cases} \frac{1}{2}x^2 & \text{if } |x| \leq \delta \\ \delta(|x| - \frac{\delta}{2}) & \text{if } |x| > \delta \end{cases}\end{split}\]

Typical Use Case: Robust optimization in the presence of outliers, handling calibration errors, RFI mitigation in radio interferometry.

Implementation:

from pyralysis.optimization.terms.regularizers import Huber

huber_term = Huber(delta=1.0, penalization_factor=0.01)

Parameters:

  • delta: Threshold parameter determining L1/L2 transition (default: 1.0)

  • Smaller values: More L1-like (robust to outliers)

  • Larger values: More L2-like (smooth)

Isotropic Total Variation (TV) Regularizer#

True Total Variation using Euclidean norm of gradients (rotationally invariant).

Mathematical Formulation:

\[TV_{iso}(I) = \sum_{i} \sqrt{ \left(\frac{\partial I}{\partial x}\right)_i^2 + \left(\frac{\partial I}{\partial y}\right)_i^2 }\]

Typical Use Case: Edge-preserving smoothing, image denoising, superior to TSV for preserving sharp features.

Implementation:

from pyralysis.optimization.terms.regularizers import TVIsotropic

tv_iso_term = TVIsotropic(penalization_factor=0.01, epsilon=1e-8)

Parameters:

  • epsilon: Regularization parameter to avoid division by zero in gradient (default: 1e-8)

Note: This is the true non-smooth TV, different from TSV (smooth approximation). Requires specialized optimization methods (FISTA, SDMM) due to non-differentiability.

Anisotropic Total Variation (TV) Regularizer#

Total Variation using L1 norm of gradients (separable, faster than isotropic).

Mathematical Formulation:

\[TV_{aniso}(I) = \sum_{i} \left| \frac{\partial I}{\partial x} \right|_i + \left| \frac{\partial I}{\partial y} \right|_i\]

Typical Use Case: Edge-preserving smoothing, faster than isotropic TV, good for images with axis-aligned features.

Implementation:

from pyralysis.optimization.terms.regularizers import TVAnisotropic

tv_aniso_term = TVAnisotropic(penalization_factor=0.01)

Note: Anisotropic TV is separable and typically faster to compute than isotropic TV, but is not rotationally invariant. Both TV variants are non-smooth and require proximal optimization methods (FISTA, SDMM).

Regularizer Compatibility with Optimizers#

Gradient-Based Optimizers (CG, L-BFGS): These optimizers work with smooth (differentiable) regularizers:

  • Entropy, Tikhonov, L2Norm, Huber, TotalFlux, TSV

For gradient-based optimizers and line search, see Optimizers and Line Searchers (overview: Optimization in Pyralysis).

Proximal Optimizers (FISTA, SDMM): These optimizers work with regularizers that have proximal operators (both smooth and non-smooth):

  • All regularizers implement proximal operators and can be used with FISTA/SDMM

  • Non-differentiable regularizers (L1Norm, TVIsotropic, TVAnisotropic, L2Norm) require proximal methods

  • Smooth regularizers can also be used with proximal methods for flexibility

For details on proximal optimization, see Compressed Sensing & Proximal Methods.

Key Features#

  • Multi-channel Support: All regularizers support multi-channel images with optional channel selection via parameter_index

  • Masking: Full support for boolean masks to restrict operations to specific regions

  • Efficient Gradients: Optimized gradient computations using Dask arrays for scalability

  • Proximal Operators: All regularizers implement proximal operators for FISTA/SDMM compatibility

  • Flexible Penalization: Configurable penalization factors for each regularizer

  • Optimized Operations: All mask-based operations use efficient arithmetic instead of conditional branching

TSV Regularizer Details#

The TSV (Total Squared Variation) regularizer is a smooth approximation of Total Variation that promotes edge-preserving smoothness. It computes finite differences along spatial dimensions and uses FFT-based proximal operators for efficient optimization.

Advantages over Tikhonov: - Preserves sharp features while smoothing noise - Better edge preservation - More suitable for images with complex structures

Mathematical Details: TSV computes finite differences along spatial dimensions and uses FFT-based proximal operators:

\[\nabla TSV(I) = D^T D\, I\]

where \(D\) denotes the finite difference operator along the spatial dimensions.

Using Regularizers in Objective Functions#

Combine multiple regularizers with data fidelity terms:

from pyralysis.optimization import ObjectiveFunction
from pyralysis.optimization.terms import Chi2, L1Norm, TSV, Entropy

# Create objective function with multiple terms
terms = [
    Chi2(model_visibility=mv, normalize=True),           # Data fidelity
    L1Norm(penalization_factor=0.01),                    # Sparsity
    TSV(penalization_factor=0.005),                      # Edge-preserving smoothness
    Entropy(prior_image=prior_model, penalization_factor=0.001)  # Prior bias
]

obj_func = ObjectiveFunction(terms=terms, parameter=image)

Advanced Usage with Masks#

Apply regularization only to specific regions:

from pyralysis.reconstruction.mask import Mask

# Create mask for specific region
mask = Mask(data=region_of_interest)

# Apply regularization with mask
l1_term = L1Norm(penalization_factor=0.01)
l1_term.function(mask=mask)  # Only considers masked pixels

Multi-channel Support#

Handle spectral and polarization data efficiently:

# Apply to specific channel in multi-channel image
l1_channel = L1Norm(parameter_index=0, penalization_factor=0.01)  # Channel 0 only

# Apply to all channels
l1_all = L1Norm(parameter_index=-1, penalization_factor=0.01)     # All channels

Performance Optimization#

  • Gradient Persistence: Enable persist_gradient=True in ObjectiveFunction for caching gradients

  • Proximal Persistence: Enable persist_proximal=True in regularizer terms for caching proximal operator results (useful for FISTA/SDMM)

  • Chunking: Configure Dask chunking for optimal memory usage

  • Masking: Use masks to limit computation to relevant regions

  • Channel Selection: Use parameter_index to avoid unnecessary computations

Proximal Operator Persistence:

For proximal optimization methods (FISTA, SDMM), enable persistence to maintain consistent iteration times:

from pyralysis.optimization.terms.regularizers import L1Norm, Tikhonov

# Enable proximal persistence
l1_term = L1Norm(penalization_factor=0.01, persist_proximal=True)
tikh_term = Tikhonov(penalization_factor=0.001, persist_proximal=True)

# Proximal results are cached in memory, similar to persist_gradient

Best Practices#

  1. Start Simple: Begin with basic regularizers (L1, Tikhonov) before adding complexity

  2. Balance Terms: Ensure data fidelity and regularization terms are properly balanced

  3. Use Masks: Apply regularization only to regions where it’s meaningful

  4. Monitor Convergence: Watch for convergence issues that may indicate regularization conflicts

  5. Validate Results: Always validate that regularization improves image quality

Example: Progressive Regularization

# Start with basic regularization
basic_terms = [Chi2(mv), L1Norm(0.01)]

# Add complexity gradually
advanced_terms = [Chi2(mv), L1Norm(0.01), TSV(0.005)]

# Full regularization suite
full_terms = [Chi2(mv), L1Norm(0.01), TSV(0.005), Entropy(prior, 0.001)]

See Also#

Optimizers | Compressed Sensing & Proximal Methods | Objective Function Framework