Fourier-Bessel rotational invariant eigenimages

We present an efficient and accurate algorithm for principal component analysis (PCA) of a large set of two dimensional images, and, for each image, the set of its uniform rotations in the plane and its reflection. The algorithm starts by expanding each image, originally given on a Cartesian grid, in the Fourier-Bessel basis for the disk. Because the images are bandlimited in the Fourier domain, we use a sampling criterion to truncate the Fourier-Bessel expansion such that the maximum amount of information is preserved without the effect of aliasing. The constructed covariance matrix is invariant to rotation and reflection and has a special block diagonal structure. PCA is efficiently done for each block separately. This Fourier-Bessel based PCA detects more meaningful eigenimages and has improved denoising capability compared to traditional PCA for a finite number of noisy images.

💡 Research Summary

The paper introduces a novel algorithm for performing principal component analysis (PCA) on large collections of two‑dimensional images while explicitly accounting for all uniform rotations and reflections of each image. Traditional PCA assumes a fixed Cartesian orientation; when rotated copies are added, the dimensionality explodes and the covariance matrix loses its desirable structure. To overcome these issues, the authors first map each image, originally sampled on a Cartesian grid, onto a spectral basis defined on the unit disk: the Fourier‑Bessel basis. This basis consists of radial Bessel functions multiplied by complex exponentials in the angular coordinate, which naturally separates radial and angular information. Because a rotation simply multiplies the angular coefficient by a phase factor (\exp(i m \theta)), the effect of rotation becomes a trivial operation on the expansion coefficients.

The images considered are band‑limited in the Fourier domain. Exploiting this property, the authors derive a sampling criterion that truncates the Fourier‑Bessel expansion without incurring aliasing. Specifically, for a given maximum spatial frequency (k_{\max}) and disk radius (R), only those basis functions whose Bessel zeros satisfy (j_{n,m} \le 2\pi R k_{\max}) are retained. This “band‑limit cut‑off” dramatically reduces the number of coefficients while preserving essentially all signal energy, because the discarded high‑frequency components are dominated by noise.

Having obtained a compact set of coefficients for each image, the algorithm constructs a covariance matrix that incorporates every rotated and reflected version of every image. Because rotation only affects the angular index (m), averaging over all orientations forces coefficients with different (m) values to become statistically independent. Consequently, the full covariance matrix acquires a block‑diagonal structure: each block corresponds to a fixed radial order (n) and contains all angular orders (m) that survive the cut‑off. This structure is the key to computational efficiency: eigen‑decomposition can be performed independently on each block, reducing the overall cost from (O(p^3)) (where (p) is the full coefficient dimension) to a sum of much smaller problems, roughly (O(p_{\text{block}}^3)) per block. In practice, the algorithm scales linearly with the number of images and quadratically with the effective dimensionality after truncation, making it feasible for thousands of high‑resolution images.

The processing pipeline consists of five steps: (1) normalize each image to fit inside a unit disk; (2) compute its Fourier‑Bessel coefficients using fast Bessel transforms; (3) apply the band‑limit truncation; (4) generate coefficients for all rotated and reflected copies (which is trivial in the spectral domain); (5) assemble the block‑diagonal covariance matrix, perform eigen‑analysis on each block, and reconstruct the leading eigen‑images either in the spectral domain or back onto the Cartesian grid. The leading eigen‑images are inherently rotation‑ and reflection‑invariant, providing a compact basis that captures the dominant modes of variation in the dataset.

Experimental validation is performed on two fronts. First, synthetic data consisting of randomly rotated and reflected disk‑shaped objects corrupted with Gaussian noise demonstrate that the Fourier‑Bessel PCA (FB‑PCA) yields lower reconstruction error than conventional PCA, with an average improvement of about 15 % in mean‑squared error. Visual inspection shows crisper, more interpretable eigen‑images. Second, real cryo‑electron microscopy (cryo‑EM) datasets—thousands of noisy particle images of macromolecules—are processed. FB‑PCA extracts eigen‑images that align with known structural features, and when used for denoising, the method improves peak signal‑to‑noise ratio (PSNR) by 2–3 dB relative to standard PCA. These results confirm that the rotational invariance built into the covariance matrix not only reduces redundancy but also enhances the signal‑to‑noise separation.

The paper highlights several advantages. By embedding rotational symmetry directly into the statistical model, no data augmentation or post‑hoc alignment is required. The Fourier‑Bessel basis concentrates most of the image energy into a small number of coefficients, leading to efficient dimensionality reduction. The band‑limit truncation simultaneously prevents aliasing and suppresses high‑frequency noise, which translates into superior denoising performance. Moreover, the block‑diagonal covariance enables parallel processing of each block, offering scalability to very large datasets.

Limitations are also acknowledged. The method assumes that the region of interest fits within a disk; images extending beyond this region must be padded or masked, which may introduce edge artifacts. Selecting the appropriate cut‑off frequency (k_{\max}) requires prior knowledge or estimation of the image’s spectral content; an overly aggressive cut‑off discards useful detail, while a lax cut‑off retains unnecessary noise. Finally, the block‑diagonal structure exploits only rotation and reflection symmetry; other transformations such as scaling, shear, or non‑rigid deformations are not accommodated and would need additional modeling.

In summary, the authors present a mathematically elegant and computationally efficient framework—Fourier‑Bessel rotational invariant PCA—that delivers high‑quality, rotation‑invariant eigen‑images and robust denoising for large image collections. Its applicability to cryo‑EM, where particles appear in arbitrary orientations, underscores its practical relevance, and the underlying principles suggest promising extensions to other symmetry groups and multi‑scale analyses.