The optimal transport (OT) problem has gained significant traction in modern machine learning for its ability to: (1) provide versatile metrics, such as Wasserstein distances and their variants, and (2) determine optimal couplings between probability measures. To reduce the computational complexity of OT solvers, methods like entropic regularization and sliced optimal transport have been proposed. The sliced OT framework improves efficiency by comparing one-dimensional projections (slices) of high-dimensional distributions. However, despite their computational efficiency, sliced-Wasserstein approaches lack a transportation plan between the input measures, limiting their use in scenarios requiring explicit coupling. In this paper, we address two key questions: Can a transportation plan be constructed between two probability measures using the sliced transport framework? If so, can this plan be used to define a metric between the measures? We propose a ‘lifting’ operation to extend one-dimensional optimal transport plans back to the original space of the measures. By computing the expectation of these lifted plans, we derive a new transportation plan, termed expected sliced transport (EST) plans. We further prove that using the EST plan to weight the sum of the individual Euclidean costs $\|x - y\|^p$ for moving from $x$ to $y$ results in a valid metric between the input discrete probability measures. Finally, we demonstrate the connection between our approach and the recently proposed min-SWGG, along with illustrative numerical examples that support our theoretical findings.

This paper introduces a clustering algorithm, called principal coordinate clustering. It takes in a similarity matrix SWof a data matrix W and computes the singular value decomposition of SWto determine the principal coordinates to convert the clustering problem to a simpler domain. It is a relative of spectral clustering, however, principal coordinate clustering is easier to interpret, and gives a clear understanding of why it performs well. In a fashion, this gives intuition behind why spectral clustering works from a more simple, linear algebra perspective, beyond the typical explanations via graph cuts, or other techniques. Moreover, it was demonstrated through experimentation on real and synthetic data that the proposed method performs equally well on average as spectral clustering, and that the method has the ability to scale quite easily to truly large data.

This paper provides a comprehensive analysis of skeleton decomposition used for segmentation of data W = [w1···WN] ⊂ ℝddrawn from a union U = ∪i=1MSiof linearly independent subspaces {Si}i=1Mof dimensions of {di}i=1M. Our previous work developed a generalized theoretical framework for computing similarity matrices by matrix factorization. Skeleton decomposition is a special case of this general theory. First, a square sub-matrix A ϵ ℝr×rof W with the same rank r as W is found. Then, the corresponding row restriction R of W is constructed. This leads to P = A-1ℝ and corresponding similarity matrix SW= (pTp)dmax, where dmaxis the maximum subspace dimension. Since most of the data matrices are low-rank in many subspace segmentation problems, this is computationally efficient compared to the other constructions of similarity matrices. It is also shown (with some limitations) that center-of-mass based sorting of data columns in SWcan be used to quickly assess clustering performance while algorithm development in both noisy or noise-free cases.

The motion segmentation problem is a special case of the general subspace segmentation problem that clusters data drawn from an unknown union of subspaces. This paper provides a nonlinear model for general subspace segmentation problem and presents an algorithm to compute the optimal solution for noiseless data. We also provide a combined algorithm that addresses issues with noise to some extent. Furthermore, a devised algorithm that specifically targets motion segmentation has been developed and applied to the Hopkins 155 Dataset. It generates the best segmentation rate to the date.

This letter presents a clustering algorithm for high dimensional data that comes from a union of lower dimensional subspaces of equal and known dimensions. The algorithm estimates a local subspace for each data point, and computes the distances between the local subspaces and the points to convert the problem to a one-dimensional data clustering problem. The algorithm is reliable in the presence of noise, and applied to the Hopkins 155 Dataset, it generates the best results to date for motion segmentation. The two motion, three motion, and overall segmentation rates for the video sequences are 99.43%, 98.69%, and 99.24%, respectively.

We consider the problem of registering (aligning) two images to sub-pixel accuracy by optimization of objective functions constructed from the images' intensity values. We show that some widely used interpolation methods can introduce multiple local optima in the energy of the interpolated image which, if not counter-balanced by other terms, can cause local optima in registration objective functions including the sum of squared differences, cross correlation, and mutual information. We discuss different solutions to address the problem based on high degree B-spline interpolation, low pass filtering the images, and stochastic integration. Numerical examples using synthetic and real signals and images are shown.

We look at the problem of extracting geometric functions such as spatial transformations and curves delineating edges to subpixel accuracy from noisy, sampled data. We analyze the stochastic properties of continuous B-spline interpolation to show that, in general noisy circumstances, sub-pixel accuracy is not obtainable when using low degree B-splines. Results using magnetic resonance image (MRI) data are shown

Nonrigid registration of medical images is important for a number of applications such as the creation of population averages, atlas-based segmentation, or geometric correction of functional magnetic resonance imaging (fMRI) images to name a few. In recent years, a number of methods have been proposed to solve this problem, one class of which involves maximizing a mutual information (MI)-based objective function over a regular grid of splines. This approach has produced good results but its computational complexity is proportional to the compliance of the transformation required to register the smallest structures in the image. Here, we propose a method that permits the spatial adaptation of the transformation's compliance. This spatial adaptation allows us to reduce the number of degrees of freedom in the overall transformation, thus speeding up the process and improving its convergence properties. To develop this method, we introduce several novelties: 1) we rely on radially symmetric basis functions rather than B-splines traditionally used to model the deformation field; 2) we propose a metric to identify regions that are poorly registered and over which the transformation needs to be improved; 3) we partition the global registration problem into several smaller ones; and 4) we introduce a new constraint scheme that allows us to produce transformations that are topologically correct. We compare the approach we propose to more traditional ones and show that our new algorithm compares favorably to those in current use.

I. Introduction Wavelets are the result of collective efforts that recognized common threads between ideas and concepts that had been independently developed and investigated by distinct research communities. They provide a unifying framework for decomposing images, volumes, and time-series data into their elementary constituents across scale. Although a relatively recent construct, wavelets have become a tool of choice for engineers, physicists, and mathematicians, leading to efficient solutions in time and space frequency analysis problems, as well as a multitude of other applications. One of the consequences is that wavelet methods of analysis and representation are presently having a significant impact on the science of medical imaging and the diagnosis of disease and screening protocols. Because of a powerful underlying mathematical theory, they offer exciting opportunities for the design of new multiresolution image processing algorithms, and novel acquisition methods such as wavelet-encoded magnetic resonance imaging (MRI). This special issue of the IEEE Transactions on Medical Imaging focuses on these recent developments and highlights progress that has been accomplished in the areas related to medical imaging. II. Sizing the Wave: Some Facts and Figures Rarely has a mathematical concept generated so much response and enthusiasm within and between the engineering and mathematical research communities at large. To give a rough idea of the phenomenon, we provide a brief chronology. While wavelets have been traced all the way back to Alfred Haar in 1910 [1], for many, the starting point of their modern history coincides with two publications in the late 1980s by S. Mallat [2] and I. Daubechies [3]. These groundbreaking papers established a solid mathematical footing which would both shape and define the field. In a nutshell, S. Mallat identified the important concept of multiresolution analysis which is the corner stone of modern wavelet theory, while I. Daubechies constructed the first orthogonal wavelet bases that were compactly supported. These two contributions count among the most cited papers in the scientific literature (over 1500 SCI citations each). From that point on, the number of contributions relating to wavelet-applications and theory has increased steadily on the order of 9000 journal papers published to date. This trend is likely to continue as suggested by the strong response to the call for papers for this special issue (over 30 submissions). Wavelets have become so popular that distinct communities continue to have conferences and scientific journals entirely devoted to them. There are already historical anecdotes and folklore associated with them; an entertaining account of which can be found in the book of B. Burke Hubbard [4]. Readers who want to dive deeper into the subject have the daunting task of choosing among over 200 books written on wavelets. Our only advice in this regard is: in case of doubt, stick with the classics. Given the size of the phenomenon, it is no surprise that wavelets have had an impact on a number of disciplines, medical imaging being no exception. A first record of activity in this particular area is the workshop on wavelets in medicine and biology that took place at the annual IEEE-EMBS meeting, Baltimore, MD, 1992. The first journal paper describing a wavelet application in medical imaging—noise reduction in MRI by soft-thresholding in the wavelet domain—also appeared in 1992 [5]. Note that this work, which is often overlooked, provided the earliest description of a wavelet denoising method that has become extremely popular through the impulsion of Donoho et al. A large palette of wavelet applications in medical imaging is provided in [6]. Two complementary review articles are also available; the first gives a complete account of the activity taking place from the beginning to 1996 [7], while the second covers the more recent papers until 2000 [8]. So far, the primary applications of wavelets in medical imaging have been the following: Compression of medical images. CT reconstruction; local tomography. Wavelet denoising (MRI, ultrasound). Wavelet-based feature extraction; texture and statistical descriptors Medical image enhancement (e.g. fluoroscopy, mammography). Analysis of functional images of the brain [positron emission tomography (PET), functional MRI (fMRI)]. Wavelet-encoded MRI: Most of these topics are still active areas of research, as illustrated by the papers that are published in this special issue. III. Scanning Through the Issue Functional imaging is an area were methods of wavelet processing hold great promise. This particular line of research was initiated by U. Ruttimann, a creative researcher and good friend, who sadly passed away shortly before the publication of his paper in this very journal [9]. Another first rate statistician who was also active in this area at an early stage is J. Raz. By a sad coincidence, he also suffered a sudden death about a week before he was to present his latest results on wavelet analysis of fMRI [10]. Despite the tragic loss of these two pioneers, research in this area is alive and well as exemplified by the first three papers of this issue. Turkheimer et al. [11] consider the problem of the analysis of dynamic PET data; in particular, they advocate the use of a linear (James-Stein) wavelet estimator as an alternative to the more classical wavelet shrinkage or threshold detectors. Hossein-Zadeh et al. propose a wavelet-technique for the detection of activation in fMRI data [12]. Their contribution is twofold: first, the use of a redundant wavelet transform for better translation invariance, and second, a nonparametric detection method based on a randomization procedure. F. Meyer also considers fMRI time series but applies wavelets differently, within the context of a generalized linear model, to detrend the data [13]; that is, to get rid of signal drifts and disturbances that are not related to the stimulus. Two other areas where wavelets have achieved great success is signal denoising (typically, by simple thresholding in the wavelet domain) and tomographic reconstruction, mainly, because the Radon operator is well localized in a wavelet basis. Pizurica et al. [14] propose a novel wavelet denoising method that uses a local statistical model for improved signal estimation and noise suppression. Willett and Nowak [15] introduce a new multiscale image model, using piecewise planar basis functions, and apply their method to the reconstruction of photon-limited data (with Poisson noise). They develop penalized maximum likelihood methods for image denoising, deconvolution, and tomographic reconstruction. Kalifa et al. present a direct method for the efficient reconstruction of PET and single photom emission computed tomography data [16]. Their approach includes a nonlinear noise reduction step that is implemented by thresholding in the transformed domain; the key here is to select a transform (wavelet packet) that is optimized for the problem and data at hand (sparse representation of the signal and near diagonalization of the Radon operator). Bonnet et al. [17] also develop a direct approach for the reconstruction of cone-beam data which is known to be challenging. In essence, their approach is a wavelet adaptation of the Feldkamp algorithm. The special issue also features two contributions relating to ultrasound imaging. Michailovich and Adam consider the problem of the estimation of the spectrum of a ultrasound pulse [18]. Specifically, they develop a modified (outlier-resistant) wavelet estimator that they apply to the log-spectrum of the radio-freqency sequence. Lee et al. present a pattern recognition system that uses an M-band wavelet filterbank to extract fractal and texture features from ultrasonic images of the liver [19]. They report promising classification results, differientiating normal liver, cirrhosis, and hepatoma using a hierarchical classifier. The recent development of commercial digital mammography imaging systems not only provides a significant improvement in image quality for traditional screening, but translates into a wealth of information for the analysis and detection of mammographic features by computer. The papers by Lemaur et al. [20]. and Heinlein et al. [21] focus on the goals of early detection and visual enhancement of microcalcifications, respectively. In the former, the regularity of a wavelet basis is used to identify microcalcification in clusters. The identification of microcalcifaction in clusters as opposed to individual occurrences is of clinical significance as clusters may suggest the likelihood of malignancy. In the later paper, a discretization of the continuous wavelet transform is developed which allows a filterbank to be adapted for the enhancement of mammographic features. This implementation allows for the reconstruction of modified wavelet coefficients at arbitrary scales and orientations without the introduction of artifacts or loss of completeness. The integration of such an interactive enhancement tool into digital mammographic screening systems will be of great importance as the wealth of dynamic range (contrast) provided by digital detectors become generally available to radiologist through the introduction of lower cost softcopy display systems. The paper of Davatzikos et al. [22] offers another illustration of the versatility of wavelets. Their proposal is to represent the contours of a shape in a wavelet bases and to use this sparse representation to derive active shape models. Their results are promising and significant in terms of providing an automated solution to problems in volume quantification. The amount of data generated by modern imaging devices is often very large, and ever increasing. Thus, an important problem is to find efficient ways of compressing and encoding this information

We prove that the exact reconstruction of a function s from its samples s(x/sub i/) on any "sufficiently dense" sampling set {x/sub i/}/sub i/spl isin/I//spl sub/R/sup n/, where I is a countable indexing set, can be obtained for a large class of spline-like spaces that belong to LP(R/sup n/). Moreover, the reconstruction can be implemented using fast algorithms. Since, a special case is the space of bandlimited functions, our result generalizes the classical Shannon-Whittacker (1949) sampling theorem on regular sampling and the Paley-Wiener (1934) theorem on nonuniform sampling.

We introduce a new method for initializing the multi-wavelet decomposition algorithm. The approach assumes that the input signal is contained within some well-defined subspace of L/sub 2/ (e.g. space of bandlimited functions). The initialization algorithm is the orthogonal projection of the input signal into the space defined by the multi-scaling function. Unlike an interpolation approach, the projection method will always have a solution. We provide examples and implementation details.

Describes an image processing system for the automatic assessment of ejection fraction (EF) from noisy SPECT left ventricular myocardial perfusion images. A segmentation scheme consisting of an elliptical coordinate transformation, matched filtering and dynamic contour tracking detects the endo- and epicardial boundaries. Computation of EF is performed based on the epi- rather than the endocardial boundary, which proves to be more robust for images with low signal-to-noise ratios. The computation incorporates anatomical constraints of constant myocardial mass and smooth cardiac variation. The algorithm was tested on different image modalities and shows good linear agreement with EFs obtained from conventional approaches based on planar gated blood pool imaging (PET: y=8.7+1.07x, r=0.84, Technetium-99m MIBI SPECT: y=5.0+0.80x, r=0.90, Thallium-201 SPECT: y=23.5+0.82x, r=0.77, RMS error=10.0). The proposed scheme may be an alternative for EF computation without the need for additional image acquisitions.

We present an overview of the various uses of the wavelet transform (WT) in medicine and biology. We start by describing the wavelet properties that are the most important for biomedical applications. In particular we provide an interpretation of the the continuous wavelet transform (CWT) as a prewhitening multiscale matched filter. We also briefly indicate the analogy between the WT and some of the the biological processing that occurs in the early components of the auditory and visual system. We then review the uses of the WT for the analysis of 1-D physiological signals obtained by phonocardiography, electrocardiography (ECG), mid electroencephalography (EEG), including evoked response potentials. Next, we provide a survey of wavelet developments in medical imaging. These include biomedical image processing algorithms (e.g., noise reduction, image enhancement, and detection of microcalcifications in mammograms), image reconstruction and acquisition schemes (tomography, and magnetic resonance imaging (MRI)), and multiresolution methods for the registration and statistical analysis of functional images of the brain (positron emission tomography (PET) and functional MRI (fMRI)). In each case, we provide the reader with same general background information and a brief explanation of how the methods work.

We present examples of a new type of wavelet basis functions that are orthogonal across shifts but not across scales. The analysis functions are piecewise linear while the synthesis functions are polynomial splines of degree n (odd). The approximation power of these representations is essentially as good as that of the corresponding Battle-Lemarie orthogonal wavelet transform, with the difference that the present wavelet synthesis filters have a much faster decay. This last property, together with the fact that these transformations are almost orthogonal, may be useful for image coding applications.

The purpose of this paper is to derive optimal spline algorithms for the enlargement or reduction of digital images by arbitrary (noninteger) scaling factors. In our formulation, the original and rescaled signals are each represented by an interpolating polynomial spline of degree n with step size one and Delta, respectively. The change of scale is achieved by determining the spline with step size Delta that provides the closest approximation of the original signal in the L(2)-norm. We show that this approximation can be computed in three steps: (i) a digital prefilter that provides the B-spline coefficients of the input signal, (ii) a resampling using an expansion formula with a modified sampling kernel that depends explicitly on Delta, and (iii) a digital postfilter that maps the result back into the signal domain. We provide explicit formulas for n=0, 1, and 3 and propose solutions for the efficient implementation of these algorithms. We consider image processing examples and show that the present method compares favorably with standard interpolation techniques. Finally, we discuss some properties of this approach and its connection with the classical technique of bandlimiting a signal, which provides the asymptotic limit of our algorithm as the order of the spline tends to infinity.

The authors are concerned with the derivation of general methods for the L/sub 2/ approximation of signals by polynomial splines. The main result is that the expansion coefficients of the approximation are obtained by linear filtering and sampling. The authors apply those results to construct a L/sub 2/ polynomial spline pyramid that is a parametric multiresolution representation of a signal. This hierarchical data structure is generated by repeated application of a REDUCE function (prefilter and down-sampler). A complementary EXPAND function (up-sampler and post-filter) allows a finer resolution mapping of any coarser level of the pyramid. Four equivalent representations of this pyramid are considered, and the corresponding REDUCE and EXPAND filters are determined explicitly for polynomial splines of any order n (odd). Some image processing examples are presented. It is demonstrated that the performance of the Laplacian pyramid can be improved significantly by using a modified EXPAND function associated with the dual representation of a cubic spline pyramid.>

The authors describe a novel digital filtering algorithms for the processing and representation of signals using polynomial splines. The classical polynomial spline interpolation problem is considered. It is found that it can be solved efficiently by recursive digital filtering. This result also yields a simple procedure for signal differentiation. Filters that efficiently solve the problem of smoothing spline approximations are derived. This technique is a regularized version of spline interpolation and is therefore less sensitive to noise. It is applied to the design of a robust edge detection algorithm with an adjustable scale parameter. A filtering/sampling algorithm for least squares spline approximation is described. This data reduction technique is applied to the generation of a cubic spline image pyramid that is found to compare favorably with the Gauss/Laplace pyramid.>

The least-squares polynomial spline approximation of a signal g(t) in L/sub 2/(R) is obtained by projecting g(t) on S/sup n/(R) (the space of polynomial splines of order n). It is shown that this process can be linked to the classical problem of cardinal spline interpolation by first convolving g(t) with a B-spline of order n. More specifically, the coefficients of the B-spline interpolation of order 2n+1 of the sampled filtered sequence are identical to the coefficients of the least-squares approximation of g(t) of order n. It is shown that this approximation can be obtained from a succession of three basic operations: prefiltering, sampling, and postfiltering, which confirms the parallel with the classical sampling/reconstruction procedure for bandlimited signals. The frequency responses of these filters are determined for three equivalent spline representations using alternative sets of shift-invariant basis functions of S/sup n/(R): the standard expansion in terms of B-spline coefficients, a representation in terms of sampled signal values, and a representation using orthogonal basis functions.>

A family of nonorthogonal polynomial spline wavelet transforms is considered. These transforms are fully reversible and can be implemented efficiently. The corresponding wavelet functions have a compact support. It is proven that these B-spline wavelets converge to Gabor functions (modulated Gaussian) pointwise and in all L/sub p/-norms with 1>

Least squares approximation problems that are regularized with specified highpass stabilizing kernels are discussed. For each problem, there is a family of discrete regularization filters (R-filters) which allow an efficient determination of the solutions. These operators are stable symmetric lowpass filters with an adjustable scale factor. Two decomposition theorems for the z-transform of such systems are presented. One facilitates the determination of their impulse response, while the other allows an efficient implementation through successive causal and anticausal recursive filtering. A case of special interest is the design of R-filters for the first- and second-order difference operators. These results are extended for two-dimensional signals and, for illustration purposes, are applied to the problem of edge detection. This leads to a very efficient implementation (8 multiplies+10 adds per pixel) of the optimal Canny edge detector based on the use of a separable second-order R-filter.>

Efficient algorithms for the continuous representation of a discrete signal in terms of B-splines (direct B-spline transform) and for interpolative signal reconstruction (indirect B-spline transform) with an expansion factor m are described. Expressions for the z-transforms of the sampled B-spline functions are determined and a convolution property of these kernels is established. It is shown that both the direct and indirect spline transforms involve linear operators that are space invariant and are implemented efficiently by linear filtering. Fast computational algorithms based on the recursive implementations of these filters are proposed. A B-spline interpolator can also be characterized in terms of its transfer function and its global impulse response (cardinal spline of order n). The case of the cubic spline is treated in greater detail. The present approach is compared with previous methods that are reexamined from a critical point of view. It is concluded that B-spline interpolation correctly applied does not result in a loss of image resolution and that this type of interpolation can be performed in a very efficient manner.>