A key property of neural networks is their capacity of adapting to data during training. Yet, our current mathematical understanding of feature learning and its relationship to generalization remain limited. In this work, we provide a random matrix analysis of how fully-connected two-layer neural networks adapt to the target function after a single, but aggressive, gradient descent step. We rigorously establish the equivalence between the updated features and an isotropic spiked random feature model, in the limit of large batch size. For the latter model, we derive a deterministic equivalent description of the feature empirical covariance matrix in terms of certain low-dimensional operators. This allows us to sharply characterize the impact of training in the asymptotic feature spectrum, and in particular, provides a theoretical grounding for how the tails of the feature spectrum modify with training. The deterministic equivalent further yields the exact asymptotic generalization error, shedding light on the mechanisms behind its improvement in the presence of feature learning. Our result goes beyond standard random matrix ensembles, and therefore we believe it is of independent technical interest. Different from previous work, our result holds in the challenging maximal learning rate regime, is fully rigorous and allows for finitely supported second layer initialization, which turns out to be crucial for studying the functional expressivity of the learned features. This provides a sharp description of the impact of feature learning in the generalization of two-layer neural networks, beyond the random features and lazy training regimes.

In this manuscript, we analyze a solvable model of flow or diffusion-based generative model. We consider the problem of learning a model parametrized by a two-layer auto-encoder, trained with online stochastic gradient descent, on a high-dimensional target density with an underlying low-dimensional manifold structure. We derive a tight asymptotic characterization of low-dimensional projections of the distribution of samples generated by the learned model, ascertaining in particular its dependence on the number of training samples. Building on this analysis, we discuss how mode collapse can arise, and lead to model collapse when the generative model is re-trained on generated synthetic data.

In this manuscript, we investigate the problem of how two-layer neural networks learn features from data, and improve over the kernel regime, after being trained with a single gradient descent step. Leveraging the insight from (Ba et al., 2022), we model the trained network by a spiked Random Features (sRF) model. Further building on recent progress on Gaussian universality (Dandi et al., 2023), we provide an exact asymptotic description of the generalization error of the sRF in the high-dimensional limit where the number of samples, the width, and the input dimension grow at a proportional rate. The resulting characterization for sRFs also captures closely the learning curves of the original network model. This enables us to understand how adapting to the data is crucial for the network to efficiently learn non-linear functions in the direction of the gradient - where at initialization it can only express linear functions in this regime.

Motivated by the recent application of approximate message passing (AMP) to the analysis of convex optimizations in multi-class classifications [Loureiro, et. al., 2021], we present a convergence analysis of AMP dynamics with non-separable multivariate nonlinearities. As an application, we present a complete (and independent) analysis of the motivated convex optimization problem.

It has been observed that the performances of many high-dimensional estimation problems are universal with respect to underlying sensing (or design) matrices. Specifically, matrices with markedly different constructions seem to achieve identical performance if they share the same spectral distribution and have “generic” singular vectors. We prove this universality phenomenon for the case of convex regularized least squares (RLS) estimators under a linear regression model with additive Gaussian noise. Our main contributions are two-fold: (1) We introduce a notion of universality classes for sensing matrices, defined through a set of deterministic conditions that fix the spectrum of the sensing matrix and precisely capture the notion of generic singular vectors; (2) We show that for all sensing matrices that lie in the same universality class, the dynamics of the proximal gradient descent algorithm for solving the regression problem, as well as the performance of RLS estimators themselves (under additional strong convexity conditions) are asymptotically identical. In addition to including i.i.d. Gaussian and rotational invariant matrices as special cases, our universality class also contains highly structured, strongly dependent, and even (nearly) deterministic matrices. Examples of the latter include randomly signed versions of incoherent tight frames and randomly subsampled Hadamard transforms. As a consequence of this universality principle, the asymptotic performance of regularized linear regression on many structured matrices constructed with limited randomness can be characterized by using the rotationally invariant ensemble as an equivalent yet mathematically more tractable surrogate.

We prove a universality theorem for learning with random features. Our result shows that, in terms of training and generalization errors, a random feature model with a nonlinear activation function is asymptotically equivalent to a surrogate linear Gaussian model with a matching covariance matrix. This settles a so-called Gaussian equivalence conjecture based on which several recent papers develop their results. Our method for proving the universality theorem builds on the classical Lindeberg approach. Major ingredients of the proof include a leave-one-out analysis for the optimization problem associated with the training process and a central limit theorem, obtained via Stein’s method, for weakly correlated random variables.

In sparse linear regression, the SLOPE estimator generalizes LASSO by penalizing different coordinates of the estimate according to their magnitudes. In this paper, we present a precise performance characterization of SLOPE under the i.i.d. Gaussian design. We focus on the asymptotic regime where the number of unknown parameters grows in proportion to the number of observations. Our asymptotic characterization enables us to derive the fundamental limits of SLOPE in both estimation and variable selection settings. We also provide a computationally feasible way to optimally design the regularizing sequences to reach the aforementioned fundamental limits. In both settings, we show that the optimal design problem can be formulated as certain infinite-dimensional convex optimization problems that have efficient and accurate finite-dimensional approximations. Numerical simulations verify all our asymptotic predictions. They demonstrate the superiority of our optimal regularizing sequences over other designs used in the existing literature.

Randomly perturbing networks during the training process is a commonly used approach to improving generalization performance. In this paper, we present a theoretical study of one particular way of random perturbation, which corresponds to injecting artificial noise to the training data. We provide a precise asymptotic characterization of the training and generalization errors of such randomly perturbed learning problems on a random feature model. Our analysis shows that Gaussian noise injection in the training process is equivalent to introducing a weighted ridge regularization, when the number of noise injections tends to infinity. The explicit form of the regularization is also given. Numerical results corroborate our asymptotic predictions, showing that they are accurate even in moderate problem dimensions. Our theoretical predictions are based on a new correlated Gaussian equivalence conjecture that generalizes recent results in the study of random feature models.

This paper proposes a new algorithm, named Householder Dice (HD), for simulating dynamics on dense random matrix ensembles with rotational invariance. Examples include the Gaussian ensemble, the Haar-distributed random orthogonal ensemble, and their complex-valued counterparts. A “direct” approach to the simulation, where one first generates a dense$n \times n$matrix from the ensemble, requires at least$\mathcal {O}(n^{2})$resource in space and time. The HD algorithm overcomes this$\mathcal {O}(n^{2})$bottleneck by using the principle of deferred decisions: rather than fixing the entire random matrix in advance, it lets the randomness unfold with the dynamics. At the heart of this matrix-free algorithm is an adaptive and recursive construction of (random) Householder reflectors. These orthogonal transformations exploit the group symmetry of the matrix ensembles, while simultaneously maintaining the statistical correlations induced by the dynamics. The memory and computation costs of the HD algorithm are$\mathcal {O}(nT)$and$\mathcal {O}(nT^{2})$, respectively, with$T$being the number of iterations. When$T \ll n$, which is nearly always the case in practice, the new algorithm leads to significant reductions in runtime and memory footprint. Numerical results demonstrate the promise of the HD algorithm as a new computational tool in the study of high-dimensional random systems.

We consider a high-dimensional mixture of two Gaussians in the noisy regime where even an oracle knowing the centers of the clusters misclassifies a small but finite fraction of the points. We provide a rigorous analysis of the generalization error of regularized convex classifiers, including ridge, hinge and logistic regression, in the high-dimensional limit where the number $n$ of samples and their dimension $d$ go to infinity while their ratio is fixed to $\alpha=n/d$. We discuss surprising effects of the regularization that in some cases allows to reach the Bayes-optimal performances. We also illustrate the interpolation peak at low regularization, and analyze the role of the respective sizes of the two clusters.

We consider a commonly studied supervised classification of a synthetic dataset whose labels are generated by feeding a one-layer non-linear neural network with random iid inputs. We study the generalization performances of standard classifiers in the high-dimensional regime where $\alpha=\frac{n}{d}$ is kept finite in the limit of a high dimension $d$ and number of samples $n$. Our contribution is three-fold: First, we prove a formula for the generalization error achieved by $\ell_2$ regularized classifiers that minimize a convex loss. This formula was first obtained by the heuristic replica method of statistical physics. Secondly, focussing on commonly used loss functions and optimizing the $\ell_2$ regularization strength, we observe that while ridge regression performance is poor, logistic and hinge regression are surprisingly able to approach the Bayes-optimal generalization error extremely closely. As $\alpha \to \infty$ they lead to Bayes-optimal rates, a fact that does not follow from predictions of margin-based generalization error bounds. Third, we design an optimal loss and regularizer that provably leads to Bayes-optimal generalization error.

In Generalized Linear Estimation (GLE) problems, we seek to estimate a signal that is observed through a linear transform followed by a component-wise, possibly nonlinear and noisy, channel. In the Bayesian optimal setting, Generalized Approximate Message Passing (GAMP) is known to achieve optimal performance for GLE. However, its performance can significantly deteriorate whenever there is a mismatch between the assumed and the true generative model, a situation frequently encountered in practice. In this paper, we propose a new algorithm, named Generalized Approximate Survey Propagation (GASP), for solving GLE in the presence of prior or model misspecifications. As a prototypical example, we consider the phase retrieval problem, where we show that GASP outperforms the corresponding GAMP, reducing the reconstruction threshold and, for certain choices of its parameters, approaching Bayesian optimal performance. Furthermore, we present a set of state evolution equations that can precisely characterize the performance of GASP in the high-dimensional limit.

In sparse linear regression, the SLOPE estimator generalizes LASSO by assigning magnitude-dependent regularizations to different coordinates of the estimate. In this paper, we present an asymptotically exact characterization of the performance of SLOPE in the high-dimensional regime where the number of unknown parameters grows in proportion to the number of observations. Our asymptotic characterization enables us to derive optimal regularization sequences to either minimize the MSE or to maximize the power in variable selection under any given level of Type-I error. In both cases, we show that the optimal design can be recast as certain infinite-dimensional convex optimization problems, which have efficient and accurate finite-dimensional approximations. Numerical simulations verify our asymptotic predictions. They also demonstrate the superiority of our optimal design over LASSO and a regularization sequence previously proposed in the literature.

We present a theoretical analysis of the training process for a single-layer GAN fed by high-dimensional input data. The training dynamics of the proposed model at both microscopic and macroscopic scales can be exactly analyzed in the high-dimensional limit. In particular, we prove that the macroscopic quantities measuring the quality of the training process converge to a deterministic process characterized by an ordinary differential equation (ODE), whereas the microscopic states containing all the detailed weights remain stochastic, whose dynamics can be described by a stochastic differential equation (SDE). This analysis provides a new perspective different from recent analyses in the limit of small learning rate, where the microscopic state is always considered deterministic, and the contribution of noise is ignored. From our analysis, we show that the level of the background noise is essential to the convergence of the training process: setting the noise level too strong leads to failure of feature recovery, whereas setting the noise too weak causes oscillation. Although this work focuses on a simple copy model of GAN, we believe the analysis methods and insights developed here would prove useful in the theoretical understanding of other variants of GANs with more advanced training algorithms.

For many modern applications in science and engineering, data are collected in a streaming fashion carrying time-varying information, and practitioners need to process them with a limited amount of memory and computational resources in a timely manner for decision making. This often is coupled with the missing data problem, such that only a small fraction of data attributes are observed. These complications impose significant, and unconventional, constraints on the problem of streaming Principal Component Analysis (PCA) and subspace tracking, which is an essential building block for many inference tasks in signal processing and machine learning. This survey article reviews a variety of classical and recent algorithms for solving this problem with low computational and memory complexities, particularly those applicable in the big data regime with missing data. We illustrate that streaming PCA and subspace tracking algorithms can be understood through algebraic and geometric perspectives, and they need to be adjusted carefully to handle missing data. Both asymptotic and non-asymptotic convergence guarantees are reviewed. Finally, we benchmark the performance of several competitive algorithms in the presence of missing data for both well-conditioned and ill-conditioned systems.

Finding the sparse representation of a signal in an overcomplete dictionary has attracted a lot of attention over the past years. This paper studies ProSparse, a new polynomial complexity algorithm that solves the sparse representation problem when the underlying dictionary is the union of a Vandermonde matrix and a banded matrix. Unlike our previous work, which establishes deterministic (worst-case) sparsity bounds for ProSparse to succeed, this paper presents a probabilistic average-case analysis of the algorithm. Based on a generating-function approach, closed-form expressions for the exact success probabilities of ProSparse are given. The success probabilities are also analyzed in the high-dimensional regime. This asymptotic analysis characterizes a sharp phase transition phenomenon regarding the performance of the algorithm.

Solving a large-scale regularized linear inverse problem using multiple processors is important in various real-world applications due to the limitations of individual processors and constraints on data sharing policies. This paper focuses on the setting where the matrix is partitioned column-wise. We extend the algorithmic framework and the theoretical analysis of approximate message passing (AMP), an iterative algorithm for solving linear inverse problems, whose asymptotic dynamics are characterized by state evolution (SE). In particular, we show that column-wise multiprocessor AMP (C-MP-AMP) obeys an SE under the same assumptions when the SE for AMP holds. The SE results imply that (i) the SE of C-MP-AMP converges to a state that is no worse than that of AMP and (ii) the asymptotic dynamics of C-MP-AMP and AMP can be identical. Moreover, for a setting that is not covered by SE, numerical results show that damping can improve the convergence performance of C-MP-AMP.

A novel data-driven reconstruction algorithm for quantum image sensors is proposed. Binary observations are efficiently decoded by modeling the reconstruction structure as a two-layer neural network, where optimal coefficients are obtained via error backpropagation. Such a model encapsulates the structure of state-of-the-art algorithms, yet it presents a considerably faster alternative which adapts to input examples without a priori statistical information. Simulations on natural and synthetic datasets show accurate reconstructions with structural similarities consistent with the state of the art, while requiring approximately 5 times less computational cost.

We study a simple spectral method that serves as a key ingredient in a growing line of work using efficient iterative algorithms for estimating signals in nonconvex settings. Unlike previous work, which focuses on the phase retrieval setting and provides only bounds on the performance, we consider arbitrary generalized linear sensing models and provide an exact characterization of the performance of the spectral method in the high-dimensional regime. Our analysis reveals a phase transition phenomenon that depends on the sampling ratio. When the ratio is below a critical threshold, the estimates given by the spectral method are no better than random guesses drawn uniformly from the hypersphere; above the threshold, however, the estimates become increasingly aligned with the underlying signal. Worked examples and numerical simulations are provided to illustrate and verify the analytical predictions.

We analyze the dynamics of an online algorithm for independent component analysis in the high-dimensional scaling limit. As the ambient dimension tends to infinity, and with proper time scaling, we show that the time-varying joint empirical measure of the target feature vector and the estimates provided by the algorithm will converge weakly to a deterministic measured-valued process that can be characterized as the unique solution of a nonlinear PDE. Numerical solutions of this PDE, which involves two spatial variables and one time variable, can be efficiently obtained. These solutions provide detailed information about the performance of the ICA algorithm, as many practical performance metrics are functionals of the joint empirical measures. Numerical simulations show that our asymptotic analysis is accurate even for moderate dimensions. In addition to providing a tool for understanding the performance of the algorithm, our PDE analysis also provides useful insight. In particular, in the high-dimensional limit, the original coupled dynamics associated with the algorithm will be asymptotically “decoupled”, with each coordinate independently solving a 1-D effective minimization problem via stochastic gradient descent. Exploiting this insight to design new algorithms for achieving optimal trade-offs between computational and statistical efficiency may prove an interesting line of future research.

Many patch-based image denoising algorithms can be formulated as applying a smoothing filter to the noisy image. Expressed as matrices, the smoothing filters must be row normalized, so that each row sums to unity. Surprisingly, if we apply a column normalization before the row normalization, the performance of the smoothing filter can often be significantly improved. Prior works showed that such performance gain is related to the Sinkhorn-Knopp balancing algorithm, an iterative procedure that symmetrizes a row-stochastic matrix to a doubly stochastic matrix. However, a complete understanding of the performance gain phenomenon is still lacking. In this paper, we study the performance gain phenomenon from a statistical learning perspective. We show that Sinkhorn-Knopp is equivalent to an expectation-maximization (EM) algorithm of learning a Gaussian mixture model of the image patches. By establishing the correspondence between the steps of Sinkhorn-Knopp and the EM algorithm, we provide a geometrical interpretation of the symmetrization process. This observation allows us to develop a new denoising algorithm called Gaussian mixture model symmetric smoothing filter (GSF). GSF is an extension of the Sinkhorn-Knopp and is a generalization of the original smoothing filters. Despite its simple formulation, GSF outperforms many existing smoothing filters and has a similar performance compared with several state-of-the-art denoising algorithms.

We present a novel algorithm - ProSparse Denoise - that can solve the sparsity recovery problem in the presence of noise when the dictionary is the union of Fourier and identity matrices. The algorithm is based on a proper use of Cadzow routine and Prony's method and exploits the duality of Fourier and identity matrices. The algorithm has low complexity compared to state of the art algorithms for sparse recovery since it relies on the Fast Fourier Transform (FFT) algorithm. We provide conditions on the noise that guarantees the correct recovery of the sparsity pattern. Our approach outperforms state of the art algorithms such as Basis Pursuit De-noise and Subspace Pursuit when the dictionary is the union of Fourier and identity matrices.

We study the dynamics of an online algorithm for learning a sparse leading eigenvector from samples generated from a spiked covariance model. This algorithm combines the classical Oja's method for online PCA with an element-wise nonlinearity at each iteration to promote sparsity. In the high-dimensional limit, the joint empirical measure of the underlying sparse eigenvector and its estimate provided by the algorithm is shown to converge weakly to a deterministic, measure-valued process. This scaling limit is characterized as the unique solution of a nonlinear PDE, and it provides exact information regarding the asymptotic performance of the algorithm. For example, performance metrics such as the cosine similarity and the misclassification rate in sparse support recovery can be obtained by examining the limiting dynamics. A steady-state analysis of the nonlinear PDE also reveals an interesting phase transition phenomenon. Although our analysis is asymptotic in nature, numerical simulations show that the theoretical predictions are accurate for moderate signal dimensions.

Estimating low-rank positive-semidefinite (PSD) matrices from symmetric rank-one measurements is of great importance in many applications, such as high-dimensional data processing, quantum state tomography, and phase retrieval. When the rank is known a priori, this problem can be regarded as solving a system of quadratic equations of a low-dimensional subspace. The authors develop a fast iterative algorithm based on an adaptation of the Kaczmarz method, which is traditionally used for solving overdetermined linear systems. In particular, the authors characterize the dynamics of the algorithm when the measurement vectors are composed of standard Gaussian entries in the online setting. Numerical simulations demonstrate the compelling performance of the proposed algorithm.

Standard convolution as a model of radiometric degradation is in majority of cases inaccurate as the blur varies in space and we are thus required to work with a computationally demanding space-variant model. Space-variant degradation can be approximately decomposed to a set of standard convolutions. We explain in detail the properties of the space-variant degradation operator and show two possible decomposition models and two approximation approaches. Our target application is space-variant image deconvolution, on which we illustrate theoretical differences between these models. We propose a computationally efficient restoration algorithm that belongs to a category of alternating direction methods of multipliers, which consists of four update steps with closed-form solutions. Depending on the used decomposition, two variations of the algorithm exist with distinct properties. We test the effectiveness of the decomposition models under different levels of approximation on synthetic and real examples, and conclude the letter by drawing several practical observations.

We propose a sampling scheme that can perfectly reconstruct a collection of spikes on the sphere from samples of their low-pass filtered observations. The proposed algorithm can reconstruct K spikes from (K + √K)2spatial samples, thus improving over previously known FRI sampling schemes on the sphere by a factor of up to four. Further, we show how multiple sound source localization (SSL) by a spherical microphone array can be transformed into a spherical FRI sampling problem. We certify the effectiveness of the proposed algorithm by using it to solve the SSL problem.

The problem of finding the sparse representation of a signal has attracted a lot of attention over the past years. In particular, uniqueness conditions and reconstruction algorithms have been established by relaxing a non-convex optimisation problem. The finite rate of innovation (FRI) theory is an alternative approach that solves the sparsity problem using algebraic methods based around Prony's algorithm. Recent extensions to this framework have shown that it is possible to recover sparse representations beyond the uniqueness limits, that is, finding all the possible sparse representations that fit the observation for the case of signals which are sparse in the union of Fourier and canonical bases. In this paper, we show the application of such methods to the case of the union of DCT and Haar basis. We present an extension that takes advantage of the even symmetry of the cosine functions to build an algorithm that can operate over the observed vector and in a dual domain. We also analyse the case of the union of frames. Simulation results confirm the validity of this new approach and show that it outperforms state of the art algorithms in a number scenarios.

We study a class of smoothing filters for image denoising. Expressed as matrices, these smoothing filters must be row normalized so that each row sums to unity. Surprisingly, if one applies a column normalization to the matrix before the row normalization, the denoising quality can often be significantly improved. This column-row normalization corresponds to one iteration of a symmetrization process called the Sinkhorn-Knopp balancing algorithm. However, a complete understanding of the performance gain phenomenon is lacking. In this paper, we analyze the performance gain from a Gaussian mixture model (GMM) perspective. We show that the symmetrization is equivalent to an expectation-maximization (EM) algorithm for learning the GMM. Moreover, we make modifications to the symmetrization procedure and present a new denoising algorithm. Experimental results show that the new algorithm achieves comparable denoising results to some state-of-the-art methods.

Traditional Finite Rate of Innovation (FRI) theory has considered the problem of sampling continuous-time signals. This framework can be naturally extended to the case where the input is a discrete-time signal. Here we present a novel approach which uses both the traditional FRI sampling scheme, based on the annihilating filter method, and the fact that in this new setup the null space of the problem to be solved is finite dimensional. In the noiseless scenario, we show that this new approach is able to perfectly recover the original signal at the critical sampling rate. We also present simulation results in the noisy scenario where this new approach improves performances in terms of the mean squared error (MSE) of the reconstructed signal when compared to the canonical FRI algorithms and compressed sensing (CS).

We propose a randomized version of the nonlocal means (NLM) algorithm for large-scale image filtering. The new algorithm, called Monte Carlo nonlocal means (MCNLM), speeds up the classical NLM by computing a small subset of image patch distances, which are randomly selected according to a designed sampling pattern. We make two contributions. First, we analyze the performance of the MCNLM algorithm and show that, for large images or large external image databases, the random outcomes of MCNLM are tightly concentrated around the deterministic full NLM result. In particular, our error probability bounds show that, at any given sampling ratio, the probability for MCNLM to have a large deviation from the original NLM solution decays exponentially as the size of the image or database grows. Second, we derive explicit formulas for optimal sampling patterns that minimize the error probability bound by exploiting partial knowledge of the pairwise similarity weights. Numerical experiments show that MCNLM is competitive with other state-of-the-art fast NLM algorithms for single-image denoising. When applied to denoising images using an external database containing ten billion patches, MCNLM returns a randomized solution that is within 0.2 dB of the full NLM solution while reducing the runtime by three orders of magnitude.

We consider the classical problem of finding the sparse representation of a signal in a pair of bases. When both bases are orthogonal, it is known that the sparse representation is unique when the sparsity K of the signal satisfies K <; 1/μ(D), where μ(D) is the mutual coherence of the dictionary. Furthermore, the sparse representation can be obtained in polynomial time by basis pursuit (BP), when K <; 0.91/μ(D). Therefore, there is a gap between the unicity condition and the one required to use the polynomial-complexity BP formulation. For the case of general dictionaries, it is also well known that finding the sparse representation under the only constraint of unicity is NP-hard. In this paper, we introduce, for the case of Fourier and canonical bases, a polynomial complexity algorithm that finds all the possible K-sparse representations of a signal under the weaker condition that K <; √2/μ(D). Consequently, when K <; 1/μ(D), the proposed algorithm solves the unique sparse representation problem for this structured dictionary in polynomial time. We further show that the same method can be extended to many other pairs of bases, one of which must have local atoms. Examples include the union of Fourier and local Fourier bases, the union of discrete cosine transform and canonical bases, and the union of random Gaussian and canonical bases.

We consider the problem of distinguishing between two hypotheses: that a sequence of signals on a large graph consists entirely of noise, or that it contains a realization of a random walk buried in the noise. The problem of computing the error exponent of the optimal detector is simple to formulate, but exhibits deep connections to problems known to be difficult, such as computing Lyapunov exponents of products of random matrices and the free entropy density of statistical mechanical systems. We describe these connections, and define an algorithm that efficiently computes the error exponent of the Neyman-Pearson detector. We also derive a closed-form formula, derived from a statistical mechanics-based approximation, for the error exponent on an arbitrary graph of large size. The derivation of this formula is not entirely rigorous, but it closely matches the empirical results in all our experiments. This formula explains a phase transition phenomenon in the error exponent: below a threshold SNR, the error exponent is nearly constant and near zero, indicating poor performance; above the threshold, there is rapid improvement in performance as the SNR increases. The location of the phase transition depends on the entropy rate of the random walk.

Non-local means (NLM) is a popular denoising scheme. Conceptually simple, the algorithm is computationally intensive for large images. We propose to speed up NLM by using random sampling. Our algorithm picks, uniformly at random, a small number of columns of the weight matrix, and uses these “representatives” to compute an approximate result. It also incorporates an extra column-normalization of the sampled columns, a form of symmetrization that often boosts the denoising performance on real images. Using statistical large deviation theory, we analyze the proposed algorithm and provide guarantees on its performance. We show that the probability of having a large approximation error decays exponentially as the image size increases. Thus, for large images, the random estimates generated by the algorithm are tightly concentrated around their limit values, even if the sampling ratio is small. Numerical results confirm our theoretical analysis: the proposed algorithm reduces the run time of NLM, and thanks to the symmetrization step, actually provides some improvement in peak signal-to-noise ratios.

Sensors of most digital cameras are made of silicon that is inherently sensitive to both the visible and near-infrared parts of the electromagnetic spectrum. In this paper, we address the problem of color and NIR joint acquisition. We propose a framework for the joint acquisition that uses only a single silicon sensor and a slightly modified version of the Bayer color-filter array that is already mounted in most color cameras. Implementing such a design for an RGB and NIR joint acquisition system requires minor changes to the hardware of commercial color cameras. One of the important differences between this design and the conventional color camera is the post-processing applied to the captured values to reconstruct full resolution images. By using a CFA similar to Bayer, the sensor records a mixture of NIR and one color channel in each pixel. In this case, separating NIR and color channels in different pixels is equivalent to solving an under-determined system of linear equations. To solve this problem, we propose a novel algorithm that uses the tools developed in the field of compressive sensing. Our method results in high-quality RGB and NIR images (the average PSNR of more than 30 dB for the reconstructed images) and shows a promising path towards RGB and NIR cameras.

The spectral theory of graphs provides a bridge between classical signal processing and the nascent field of graph signal processing. In this paper, a spectral graph analogy to Heisenberg's celebrated uncertainty principle is developed. Just as the classical result provides a tradeoff between signal localization in time and frequency, this result provides a fundamental tradeoff between a signal's localization on a graph and in its spectral domain. Using the eigenvectors of the graph Laplacian as a surrogate Fourier basis, quantitative definitions of graph and spectral “spreads” are given, and a complete characterization of the feasibility region of these two quantities is developed. In particular, the lower boundary of the region, referred to as the uncertainty curve, is shown to be achieved by eigenvectors associated with the smallest eigenvalues of an affine family of matrices. The convexity of the uncertainty curve allows it to be found to within ε by a fast approximation algorithm requiringO(ε-1/2) typically sparse eigenvalue evaluations. Closed-form expressions for the uncertainty curves for some special classes of graphs are derived, and an accurate analytical approximation for the expected uncertainty curve of Erd-s-Rényi random graphs is developed. These theoretical results are validated by numerical experiments, which also reveal an intriguing connection between diffusion processes on graphs and the uncertainty bounds.

The classical uncertainty principle provides a fundamental tradeoff in the localization of a signal in the time and frequency domains. In this paper we describe a similar tradeoff for signals defined on graphs. We describe the notions of “spread” in the graph and spectral domains, using the eigenvectors of the graph Laplacian as a surrogate Fourier basis. We then describe how to find signals that, among all signals with the same spectral spread, have the smallest graph spread about a given vertex. For every possible spectral spread, the desired signal is the solution to an eigenvalue problem. Since localization in graph and spectral domains is a desirable property of the elements of wavelet frames on graphs, we compare the performance of some existing wavelet transforms to the obtained bound.

We present a blind estimation algorithm for multi-input and multi-output (MIMO) systems with sparse common support. Key to the proposed algorithm is a matrix generalization of the classical annihilating filter technique, which allows us to estimate the nonlinear parameters of the channels through an efficient and noniterative procedure. An attractive property of the proposed algorithm is that it only needs the sensor measurements at a narrow frequency band. By exploiting this feature, we can derive efficient sub-Nyquist sampling schemes which significantly reduce the number of samples that need to be retained at each sensor. Numerical simulations verify the accuracy of the proposed estimation algorithm and its robustness in the presence of noise.

We present a novel regularization framework for demosaicking by viewing image as smooth signal on a weighted graph. The restoration problem is formed as a minimization of variation of the signal on graph. As an initial step, we build a weight matrix which measures the similarity between every pair of pixels, from an estimate of the full color image. After that, a two-stage optimization is carried out: first, we assume that graph Laplacian is signal dependent and solve a non-quadratic problem by gradient descent; then, we pose a variational problem on graph with a static Laplacian, under the constraint of consistency with the available samples in each color component. Performance evaluation shows that our approach can improve the previous demosaicking methods both quantitively and visually, by alleviating the artifical effect. Moreover, the mapping from image to signal on graph provides a general method for image processing.

We study a new image sensor that is reminiscent of a traditional photographic film. Each pixel in the sensor has a binary response, giving only a 1-bit quantized measurement of the local light intensity. To analyze its performance, we formulate the oversampled binary sensing scheme as a parameter estimation problem based on quantized Poisson statistics. We show that, with a single-photon quantization threshold and large oversampling factors, the Cramér-Rao lower bound (CRLB) of the estimation variance approaches that of an ideal unquantized sensor, i.e., as if there were no quantization in the sensor measurements. Furthermore, the CRLB is shown to be asymptotically achievable by the maximum-likelihood estimator (MLE). By showing that the log-likelihood function of our problem is concave, we guarantee the global optimality of iterative algorithms in finding the MLE. Numerical results on both synthetic data and images taken by a prototype sensor verify our theoretical analysis and demonstrate the effectiveness of our image reconstruction algorithm. They also suggest the potential application of the oversampled binary sensing scheme in high dynamic range photography.

We consider the problem of estimating room geometry from the acoustic room impulse response (RIR). Existing approaches addressing this problem exploit the knowledge of multiple RIRs. In contrast, we are interested in reconstructing the room geometry from a single RIR — a 1-D function of time. We discuss the uniqueness of the mapping between the geometry of a planar polygonal room and a single RIR. In addition to this theoretical analysis, we also propose an algorithm that performs the “blindfolded” room estimation. Furthermore, the derived results are used to construct an algorithm for localization in a known room using only a single RIR. Verification of the theoretical developments with numerical simulations is given before concluding the paper.

Spectral factorization is a classical tool in signal processing and communications. It also plays a critical role in X-ray crystallography, in the context of phase retrieval. In this work, we study the problem of sparse spectral factorization, aiming to recover a one-dimensional sparse signal from its autocorrelation. We present a sufficient condition for the recovery to be unique, and propose an iterative algorithm that can obtain the original signal (up to a sign change, time-shift and time-reversal). Numerical simulations verify the effectiveness of the proposed algorithm.

We study the spatiotemporal sampling of a diffusion field generated by K point sources, aiming to fully reconstruct the unknown initial field distribution from the sample measurements. The sampling operator in our problem can be described by a matrix derived from the diffusion model. We analyze the important properties of the sampling matrices, leading to precise bounds on the spatial and temporal sampling densities under which perfect field reconstruction is feasible. Moreover, our analysis indicates that it is possible to compensate linearly for insufficient spatial sampling densities by oversampling in time. Numerical simulations on initial field reconstruction under different spatiotemporal sampling densities confirm our theoretical results.

Joint processing of visible (RGB) and near-infrared (NIR) images has recently found some appealing applications, which make joint capturing a pair of visible and NIR images an important problem. In this paper, we propose a new method to design color filter arrays (CFA) and demosaicing matrices for acquiring NIR and visible images using a single sensor. The proposed method modifies the optimum CFA algorithm proposed in [1] by taking advantage of the NIR/visible correlation in the design process. Simulation results show that by applying the proposed method, the quality of demo-saiced NIR and visible images is increased by about 1 dB in peak signal-to-noise ratio over the results of the optimum CFA algorithm. It is also shown that better visual quality can be obtained by using the proposed algorithm.

We propose a novel algorithm for sparse system identification in the frequency domain. Key to our result is the observation that the Fourier transform of the sparse impulse response is a simple sum of complex exponentials, whose parameters can be efficiently determined from only a narrow frequency band. From this perspective, we present a sub-Nyquist sampling scheme, and show that the original continuous-time system can be learned by considering an equivalent low-rate discrete system. The impulse response of that discrete system can then be adaptively obtained by a novel frequency-domain LMS filter, which exploits the parametric structure of the model. Numerical experiments confirm the effectiveness of the proposed scheme for sparse system identification tasks.

Color image demosaicking is a key process in the digital imaging pipeline. In this paper, we study a well-known and influential demosaicking algorithm based upon alternating projections (AP), proposed by Gunturk, Altunbasak and Mersereau in 2002. Since its publication, the AP algorithm has been widely cited and compared against in a series of more recent papers in the demosaicking literature. Despite good performances, a limitation of the AP algorithm is its high computational complexity. We provide three main contributions in this paper. First, we present a rigorous analysis of the convergence property of the AP demosaicking algorithm, showing that it is a contraction mapping, with a unique fixed point. Second, we show that this fixed point is in fact the solution to a constrained quadratic minimization problem, thus, establishing the optimality of the AP algorithm. Finally, using the tool of polyphase representation, we show how to obtain the results of the AP algorithm in a single step, implemented as linear filtering in the polyphase domain. Replacing the original iterative procedure by the proposed one-step solution leads to substantial computational savings, by about an order of magnitude in our experiments.

We study the spatial-temporal sampling of a linear diffusion field, and show that it is possible to compensate for insufficient spatial sampling densities by oversampling in time. Our work is motivated by the following issue often encountered in sensor network sampling, namely increasing the temporal sampling density is often easier and less expensive than increasing the spatial sampling density of the network. For the case of sampling a diffusion field, we show that, to achieve trade-off between spatial and temporal sampling, the spatial arrangement of the sensors must satisfy certain conditions. We provide in this paper the precise relationships between the achievable reduction of spatial sampling density, the required temporal oversampling rate, the spatial arrangement of the sensors, and the bound for the condition numbers of the resulting sampling and reconstruction procedures.

We consider the task of recovering correlated vectors at a central decoder based on fixed linear measurements obtained by distributed sensors. A general formulation of the problem is proposed, under both a universal and an almost sure reconstruction requirement. We then study a specific correlation model which involves a filter that is sparse in the time domain. While this sparsity assumption does not allow reducing the description cost in the universal case, we show that large gains can be achieved in the almost sure scenario by means of a novel distributed scheme based on annihilating filters. The robustness of the proposed method is also investigated.

Digital camera sensors are inherently sensitive to the near-infrared (NIR) part of the light spectrum. In this paper, we propose a general design for color filter arrays that allow the joint capture of visible/NIR images using a single sensor. We pose the CFA design as a novel spatial domain optimization problem, and provide an efficient iterative procedure that finds (locally) optimal solutions. Numerical experiments confirm the effectiveness of the proposed CFA design, which can simultaneously capture high quality visible and NIR image pairs.

SensorScope is a collaborative project between network, signal processing, and environmental researchers that aims at providing a cheap and out-of-the-box environmental monitoring system based on a wireless sensor network. It has been successfully used in a number of deployments to gather hundreds of megabytes of environmental data. With data gathering techniques well mastered, the efficient processing of the huge amounts of the acquired information to allow for useful exploitation has become an increasingly important issue. In this paper, we present a number of challenging and relevant signal processing tasks that arise from the SensorScope project. We believe the resolution of these problems will benefit from a better understanding of the underlying physical processes. We show an example to demonstrate how physical correlations between different sensing modalities can help reduce the sampling rate.

The search for alias-free sampling lattices for a given frequency support, in particular those lattices achieving minimum sampling densities, is a fundamental issue in various applications of signal and image processing. In this paper, we propose an efficient computational procedure to find all alias-free integral sampling lattices for a given frequency support with minimum sampling density. Central to this algorithm is a novel condition linking the alias-free sampling with the Fourier transform of the indicator function defined on the frequency support. We study the computation of these Fourier transforms based on the divergence theorem, and propose a simple closed-form formula for a fairly general class of support regions consisting of arbitrary N-dimensional polytopes, with polygons in 2D and polyhedra in 3D as special cases. The proposed algorithm can be useful in a variety of applications involving the design of efficient acquisition schemes for multidimensional bandlimited signals.

In 1992, Bamberger and Smith proposed the directional filter bank (DFB) for an efficient directional decomposition of 2-D signals. Due to the nonseparable nature of the system, extending the DFB to higher dimensions while still retaining its attractive features is a challenging and previously unsolved problem. We propose a new family of filter banks, named NDFB, that can achieve the directional decomposition of arbitrary N-dimensional (N > or =2) signals with a simple and efficient tree-structured construction. In 3-D, the ideal passbands of the proposed NDFB are rectangular-based pyramids radiating out from the origin at different orientations and tiling the entire frequency space. The proposed NDFB achieves perfect reconstruction via an iterated filter bank with a redundancy factor of N in N-D. The angular resolution of the proposed NDFB can be iteratively refined by invoking more levels of decomposition through a simple expansion rule. By combining the NDFB with a new multiscale pyramid, we propose the surfacelet transform, which can be used to efficiently capture and represent surface-like singularities in multidimensional data.

The contourlet transform was proposed as a directional multiresolution image representation that can efficiently capture and represent singularities along smooth object boundaries in natural images. Its efficient filter bank construction as well as low redundancy make it an attractive computational framework for various image processing applications. However, a major drawback of the original contourlet construction is that its basis images are not localized in the frequency domain. In this paper, we analyze the cause of this problem, and propose a new contourlet construction as a solution. Instead of using the Laplacian pyramid, we employ a new multiscale decomposition defined in the frequency domain. The resulting basis images are sharply localized in the frequency domain and exhibit smoothness along their main ridges in the spatial domain. Numerical experiments on image denoising show that the proposed new contourlet transform can significantly outperform the original transform both in terms of PSNR (by several dB 's) and in visual quality, while with similar computational complexity.

Directional information is an important and unique feature of multidimensional signals. As a result of a separable extension from 1D bases, the multidimensional wavelet transform has very limited directionality. Furthermore, different directions are mixed in certain wavelet subbands. In this paper, we propose a new transform that fixes this frequency mixing problem by using a simple "add-on" to the wavelet transform. In the 2D case, it provides one lowpass subband and six directional highpass subbands at each scale. Just like the wavelet transform, the proposed transform is nonredundant, and can be easily extended to higher dimensions. Though nonseparable in essence, the proposed transform has an efficient implementation based on 1D operations only.

Many signals of interest can be characterized by a finite number of parameters per unit of time. Instead of spanning a single linear space, these signals often lie on a union of spaces. Under this setting, traditional sampling schemes are either inapplicable or very inefficient. We present a framework for sampling these signals based on an injective projection operator, which "flattens" the signals down to a common low dimensional representation space while still preserving all the information. Standard sampling procedures can then be applied on that space. We show the necessary and sufficient conditions for such operators to exist and provide the minimum sampling rate for the representation space, which indicates the efficiency of this framework. These results provide a new perspective on the sampling of signals with finite rate of innovation and can serve as a guideline for designing new algorithms for a class of problems in signal processing and communications.