Search

Optimization over the set of matrices $X$ that satisfy $X^\top B X = I_p$, referred to as the generalized Stiefel manifold, appears in many applications involving sampled covariance matrices such as the canonical correlation analysis (CCA), independent component analysis (ICA), and the generalized eigenvalue problem (GEVP). Solving these problems is typically done by iterative methods that require a fully formed $B$. We propose a cheap stochastic iterative method that solves the optimization problem while having access only to a random estimates of $B$. Our method does not enforce the constraint in every iteration; instead, it produces iterations that converge to critical points on the generalized Stiefel manifold defined in expectation. The method has lower per-iteration cost, requires only matrix multiplications, and has the same convergence rates as its Riemannian optimization counterparts that require the full matrix $B$. Experiments demonstrate its effectiveness in various machine learning applications involving generalized orthogonality constraints, including CCA, ICA, and the GEVP., Comment: This v2 is the camera-ready version for ICML 2024

We propose to use low-rank matrix approximation using the component-wise L1-norm for direct imaging of exoplanets. Exoplanet detection by direct imaging is a challenging task for three main reasons: (1) the host star is several orders of magnitude brighter than exoplanets, (2) the angular distance between exoplanets and star is usually very small, and (3) the images are affected by the noises called speckles that are very similar to the exoplanet signal both in shape and intensity. We first empirically examine the statistical noise assumptions of the L1 and L2 models, and then we evaluate the performance of the proposed L1 low-rank approximation (L1-LRA) algorithm based on visual comparisons and receiver operating characteristic (ROC) curves. We compare the results of the L1-LRA with the widely used truncated singular value decomposition (SVD) based on the L2 norm in two different annuli, one close to the star and one far away., Comment: 13 pages, 4 figures, BNAIC/BeNeLearn 2023

Orthogonality constraints naturally appear in many machine learning problems, from Principal Components Analysis to robust neural network training. They are usually solved using Riemannian optimization algorithms, which minimize the objective function while enforcing the constraint. However, enforcing the orthogonality constraint can be the most time-consuming operation in such algorithms. Recently, Ablin & Peyr\'e (2022) proposed the Landing algorithm, a method with cheap iterations that does not enforce the orthogonality constraint but is attracted towards the manifold in a smooth manner. In this article, we provide new practical and theoretical developments for the landing algorithm. First, the method is extended to the Stiefel manifold, the set of rectangular orthogonal matrices. We also consider stochastic and variance reduction algorithms when the cost function is an average of many functions. We demonstrate that all these methods have the same rate of convergence as their Riemannian counterparts that exactly enforce the constraint. Finally, our experiments demonstrate the promise of our approach to an array of machine-learning problems that involve orthogonality constraints.

We consider the problem of minimizing a differentiable function with locally Lipschitz continuous gradient on a stratified set and present a first-order algorithm designed to find a stationary point of that problem. Our assumptions on the stratified set are satisfied notably by the determinantal variety (i.e., matrices of bounded rank), its intersection with the cone of positive-semidefinite matrices, and the set of nonnegative sparse vectors. The iteration map of the proposed algorithm applies a step of projected-projected gradient descent with backtracking line search, as proposed by Schneider and Uschmajew (2015), to its input but also to a projection of the input onto each of the lower strata to which it is considered close, and outputs a point among those thereby produced that maximally reduces the cost function. Under our assumptions on the stratified set, we prove that this algorithm produces a sequence whose accumulation points are stationary, and therefore does not follow the so-called apocalypses described by Levin, Kileel, and Boumal (2022). We illustrate the apocalypse-free property of our method through a numerical experiment on the determinantal variety.

Low-rank matrix completion is the task of recovering unknown entries of a matrix by assuming that the true matrix admits a good low-rank approximation. Sometimes additional information about the variables is known, and incorporating this information into a matrix completion model can lead to a better completion quality. We consider the situation where information between the column/row entities of the matrix is available as a weighted graph. In this framework, we address the problem of completing missing entries in air temperature data recorded by weather stations. We construct test sets by holding back data at locations that mimic real-life gaps in weather data. On such test sets, we show that adequate spatial and temporal graphs can significantly improve the accuracy of the completion obtained by graph-regularized low-rank matrix completion methods., Comment: Accepted by EUSIPCO2023

We consider the problem of minimizing a differentiable function with locally Lipschitz continuous gradient over the real determinantal variety, and present a first-order algorithm designed to find stationary points of that problem. This algorithm applies steps of a retraction-free descent method proposed by Schneider and Uschmajew (2015), while taking the numerical rank into account to attempt rank reductions. We prove that this algorithm produces a sequence of iterates the accumulation points of which are stationary, and therefore does not follow the so-called apocalypses described by Levin, Kileel, and Boumal (2022). Moreover, the rank reduction mechanism of this algorithm requires at most one rank reduction attempt per iteration, in contrast with the one of the $\mathrm{P}^2\mathrm{GDR}$ algorithm introduced by Olikier, Gallivan, and Absil (2022) which can require a number of rank reduction attempts equal to the rank of the iterate in the worst-case scenario., Comment: arXiv admin note: text overlap with arXiv:2201.03962

We compare two first-order low-rank optimization algorithms, namely $\text{P}^2\text{GD}$ (Schneider and Uschmajew, 2015), which has been proven to be apocalypse-prone (Levin et al., 2021), and its apocalypse-free version $\text{P}^2\text{GDR}$ obtained by equipping $\text{P}^2\text{GD}$ with a suitable rank reduction mechanism (Olikier et al., 2022). Here an apocalypse refers to the situation where the stationarity measure goes to zero along a convergent sequence whereas it is nonzero at the limit. The comparison is conducted on two simple examples of apocalypses, the original one (Levin et al., 2021) and a new one. We also present a potential side effect of the rank reduction mechanism of $\text{P}^2\text{GDR}$ and discuss the choice of the rank reduction parameter.

Tangent and normal cones play an important role in constrained optimization to describe admissible search directions and, in particular, to formulate optimality conditions. They notably appear in various recent algorithms for both smooth and nonsmooth low-rank optimization where the feasible set is the set $\mathbb{R}_{\leq r}^{m \times n}$ of all $m \times n$ real matrices of rank at most $r$. In this paper, motivated by the convergence analysis of such algorithms, we study, by computing inner and outer limits, the continuity of the correspondence that maps each $X \in \mathbb{R}_{\leq r}^{m \times n}$ to the tangent cone to $\mathbb{R}_{\leq r}^{m \times n}$ at $X$. We also deduce results about the continuity of the corresponding normal cone correspondence. Finally, we show that our results include as a particular case the $a$-regularity of the Whitney stratification of $\mathbb{R}_{\leq r}^{m \times n}$ following from the fact that this set is a real algebraic variety, called the real determinantal variety.

We consider the problem of minimizing a differentiable function with locally Lipschitz continuous gradient on the real determinantal variety, and present a first-order algorithm designed to find stationary points of that problem. This algorithm applies steps of steepest descent with backtracking line search on the variety, as proposed by Schneider and Uschmajew (2015), but by taking the numerical rank into account to perform suitable rank reductions. We prove that this algorithm produces sequences of iterates the accumulation points of which are stationary, and therefore does not follow the so-called apocalypses described by Levin, Kileel, and Boumal (2021).

Direct imaging of exoplanets is a challenging task due to the small angular distance and high contrast relative to their host star, and the presence of quasi-static noise. We propose a new statistical method for direct imaging of exoplanets based on a likelihood ratio detection map, which assumes that the noise after the background subtraction step obeys a Laplacian distribution. We compare the method with two detection approaches based on signal-to-noise ratio (SNR) map after performing the background subtraction by the widely used Annular Principal Component Analysis (AnnPCA). The experimental results on the Beta Pictoris data set show the method outperforms SNR maps in terms of achieving the highest true positive rate (TPR) at zero false positive rate (FPR).

We propose a differential geometric approach for building families of low-rank covariance matrices, via interpolation on low-rank matrix manifolds. In contrast with standard parametric covariance classes, these families offer significant flexibility for problem-specific tailoring via the choice of “anchor” matrices for interpolation, for instance over a grid of relevant conditions describing the underlying stochastic process. The interpolation is computationally tractable in high dimensions, as it only involves manipulations of low-rank matrix factors. We also consider the problem of covariance identification, i.e., selecting the most representative member of the covariance family given a data set. In this setting, standard procedures such as maximum likelihood estimation are nontrivial because the covariance family is rank-deficient; we resolve this issue by casting the identification problem as distance minimization. We demonstrate the utility of these differential geometric families for interpolation and identification in a practical application: wind field covariance approximation for unmanned aerial vehicle navigation.

We propose a differential geometric approach for building families of low-rank covariance matrices, via interpolation on low-rank matrix manifolds. In contrast with standard parametric covariance classes, these families offer significant flexibility for problem-specific tailoring via the choice of “anchor” matrices for interpolation, for instance over a grid of relevant conditions describing the underlying stochastic process. The interpolation is computationally tractable in high dimensions, as it only involves manipulations of low-rank matrix factors. We also consider the problem of covariance identification, i.e., selecting the most representative member of the covariance family given a data set. In this setting, standard procedures such as maximum likelihood estimation are nontrivial because the covariance family is rank-deficient; we resolve this issue by casting the identification problem as distance minimization. We demonstrate the utility of these differential geometric families for interpolation and identification in a practical application: wind field covariance approximation for unmanned aerial vehicle navigation.

We propose the SH model, a simplified version of the well-known SIR compartmental model of infectious diseases. With optimized parameters and initial conditions, this time-invariant two-parameter two-dimensional model is able to fit COVID-19 hospitalization data over several months with high accuracy (e.g., the root relative squared error is below 10% for Belgium over the period from 2020-03-15 to 2020-07-15). Moreover, we observed that, when the model is trained on a suitable three-week period around the first hospitalization peak for Belgium, it forecasts the subsequent two months with mean absolute percentage error (MAPE) under 4%. We repeated the experiment for each French department and found 14 of them where the MAPE was below 20%. However, when the model is trained in the increase phase, it is less successful at forecasting the subsequent evolution., Comment: Paper home page: https://sites.uclouvain.be/absil/2020.05

This paper concerns the minimax center of a collection of linear subspaces. When the subspaces are $k$-dimensional subspaces of $\mathbb{R}^n$, this can be cast as finding the center of a minimum enclosing ball on a Grassmann manifold, Gr$(k,n)$. For subspaces of different dimension, the setting becomes a disjoint union of Grassmannians rather than a single manifold, and the problem is no longer well-defined. However, natural geometric maps exist between these manifolds with a well-defined notion of distance for the images of the subspaces under the mappings. Solving the initial problem in this context leads to a candidate minimax center on each of the constituent manifolds, but does not inherently provide intuition about which candidate is the best representation of the data. Additionally, the solutions of different rank are generally not nested so a deflationary approach will not suffice, and the problem must be solved independently on each manifold. We propose and solve an optimization problem parametrized by the rank of the minimax center. The solution is computed using a subgradient algorithm on the dual. By scaling the objective and penalizing the information lost by the rank-$k$ minimax center, we jointly recover an optimal dimension, $k^*$, and a central subspace, $U^* \in$ Gr$(k^*,n)$ at the center of the minimum enclosing ball, that best represents the data., Comment: 26 pages

The Grassmann manifold of linear subspaces is important for the mathematical modelling of a multitude of applications, ranging from problems in machine learning, computer vision and image processing to low-rank matrix optimization problems, dynamic low-rank decompositions and model reduction. With this mostly expository work, we aim to provide a collection of the essential facts and formulae on the geometry of the Grassmann manifold in a fashion that is fit for tackling the aforementioned problems with matrix-based algorithms. Moreover, we expose the Grassmann geometry both from the approach of representing subspaces with orthogonal projectors and when viewed as a quotient space of the orthogonal group, where subspaces are identified as equivalence classes of (orthogonal) bases. This bridges the associated research tracks and allows for an easy transition between these two approaches. Original contributions include a modified algorithm for computing the Riemannian logarithm map on the Grassmannian that is advantageous numerically but also allows for a more elementary, yet more complete description of the cut locus and the conjugate points. We also derive a formula for parallel transport along geodesics in the orthogonal projector perspective, formulae for the derivative of the exponential map, as well as a formula for Jacobi fields vanishing at one point., Comment: 38 pages, 4 figures

We consider a Canonical Polyadic (CP) decomposition approach to low-rank tensor completion (LRTC) by incorporating external pairwise similarity relations through graph Laplacian regularization on the CP factor matrices. The usage of graph regularization entails benefits in the learning accuracy of LRTC, but at the same time, induces coupling graph Laplacian terms that hinder the optimization of the tensor completion model. In order to solve graph-regularized LRTC, we propose efficient alternating minimization algorithms by leveraging the block structure of the underlying CP decomposition-based model. For the subproblems of alternating minimization, a linear conjugate gradient subroutine is specifically adapted to graph-regularized LRTC. Alternatively, we circumvent the complicating coupling effects of graph Laplacian terms by using an alternating directions method of multipliers. Based on the Kurdyka-{\L}ojasiewicz property, we show that the sequence generated by the proposed algorithms globally converges to a critical point of the objective function. Moreover, the complexity and convergence rate are also derived. In addition, numerical experiments including synthetic data and real data show that the graph regularized tensor completion model has improved recovery results compared to those without graph regularization, and that the proposed algorithms achieve gains in time efficiency over existing algorithms.

Invariance transformations of polyadic decompositions of matrix multiplication tensors define an equivalence relation on the set of such decompositions. In this paper, we present an algorithm to efficiently decide whether two polyadic decompositions of a given matrix multiplication tensor are equivalent. With this algorithm, we analyze the equivalence classes of decompositions of several matrix multiplication tensors. This analysis is relevant for the study of fast matrix multiplication as it relates to the question of how many essentially different fast matrix multiplication algorithms there exist. This question has been first studied by de~Groote, who showed that for the multiplication of $2\times2$ matrices with $7$ active multiplications, all algorithms are essentially equivalent to Strassen's algorithm. In contrast, the results of our analysis show that for the multiplication of larger matrices, (e.g., $2\times3$ by $3\times2$ or $3\times3$ by $3\times3$ matrices), two decompositions are very likely to be essentially different. We further provide a necessary criterion for a polyadic decomposition to be equivalent to a polyadic decomposition with integer entries. Decompositions with specific integer entries, e.g., powers of two, provide fast matrix multiplication algorithms with better efficiency and stability properties. This condition can be tested algorithmically and we present the conclusions obtained for the decompositions of small/medium matrix multiplication tensors.

We present a method to compute a fitting curve B to a set of data points d0,...,dm lying on a manifold M. That curve is obtained by blending together Euclidean B\'ezier curves obtained on different tangent spaces. The method guarantees several properties among which B is C1 and is the natural cubic smoothing spline when M is the Euclidean space. We show examples on the sphere S2 as a proof of concept., Comment: in Proceedings of iTWIST'18, Paper-ID: 18, Marseille, France, November, 21-23, 2018

We demonstrate 100-Gbps silicon-contacted germanium waveguide p-i-n photodetectors integrated on imec's silicon photonics platform. The performance of 14 and 20 μm long devices is compared. The responsivity of the devices is 0.74 and 0.92 A/W at 1550 nm, respectively.

We present the Vortex Image Processing (VIP) library, a python package dedicated to astronomical high-contrast imaging. Our package relies on the extensive python stack of scientific libraries and aims to provide a flexible framework for high-contrast data and image processing. In this paper, we describe the capabilities of VIP related to processing image sequences acquired using the angular differential imaging (ADI) observing technique. VIP implements functionalities for building high-contrast data processing pipelines, encompass- ing pre- and post-processing algorithms, potential sources position and flux estimation, and sensitivity curves generation. Among the reference point-spread function subtraction techniques for ADI post-processing, VIP includes several flavors of principal component analysis (PCA) based algorithms, such as annular PCA and incremental PCA algorithm capable of processing big datacubes (of several gigabytes) on a computer with limited memory. Also, we present a novel ADI algorithm based on non-negative matrix factorization (NMF), which comes from the same family of low-rank matrix approximations as PCA and provides fairly similar results. We showcase the ADI capabilities of the VIP library using a deep sequence on HR8799 taken with the LBTI/LMIRCam and its recently commissioned L-band vortex coronagraph. Using VIP we investigated the presence of additional companions around HR8799 and did not find any significant additional point source beyond the four known planets. VIP is available at http://github.com/vortex-exoplanet/VIP and is accompanied with Jupyter notebook tutorials illustrating the main functionalities of the library.

We present the Vortex Image Processing (VIP) library, a python package dedicated to astronomical high-contrast imaging. Our package relies on the extensive python stack of scientific libraries and aims to provide a flexible framework for high-contrast data and image processing. In this paper, we describe the capabilities of VIP related to processing image sequences acquired using the angular differential imaging (ADI) observing technique. VIP implements functionalities for building high-contrast data processing pipelines, encompass- ing pre- and post-processing algorithms, potential sources position and flux estimation, and sensitivity curves generation. Among the reference point-spread function subtraction techniques for ADI post-processing, VIP includes several flavors of principal component analysis (PCA) based algorithms, such as annular PCA and incremental PCA algorithm capable of processing big datacubes (of several gigabytes) on a computer with limited memory. Also, we present a novel ADI algorithm based on non-negative matrix factorization (NMF), which comes from the same family of low-rank matrix approximations as PCA and provides fairly similar results. We showcase the ADI capabilities of the VIP library using a deep sequence on HR8799 taken with the LBTI/LMIRCam and its recently commissioned L-band vortex coronagraph. Using VIP we investigated the presence of additional companions around HR8799 and did not find any significant additional point source beyond the four known planets. VIP is available at http://github.com/vortex-exoplanet/VIP and is accompanied with Jupyter notebook tutorials illustrating the main functionalities of the library.

This paper investigates the influence of TSV noise coupling on nearby devices based on an extended 3D TSV circuit model. This model not only takes into account the complex RF field distributions in bulk Si, but also incorporates the anomalous TSV capacitance inversion behavior, which has been found to occur due to the presence of fixed charges in the backside passivation layer after wafer thinning. The extended 3D TSV circuit model is validated by the excellent agreement between the simulation results and experimental data. It demonstrates that the inversion behavior of the TSV capacitance increases the noise coupling to adjacent devices mainly in the low frequency range. Furthermore, we show that noise mitigation techniques can be easily implemented in this 3D circuit model to predict the extent of noise coupling alleviation.

Data processing constitutes a critical component of high-contrast exoplanet imaging. Its role is almost as important as the choice of a coronagraph or a wavefront control system, and it is intertwined with the chosen observing strategy. Among the data processing techniques for angular differential imaging (ADI), the most recent is the family of principal component analysis (PCA) based algorithms. PCA serves, in this case, as a subspace projection technique for constructing a reference point spread function (PSF) that can be subtracted from the science data for boosting the detectability of potential companions present in the data. Unfortunately, when building this reference PSF from the science data itself, PCA comes with certain limitations such as the sensitivity of the lower dimensional orthogonal subspace to non-Gaussian noise. Inspired by recent advances in machine learning algorithms such as robust PCA, we aim to propose a localized subspace projection technique that surpasses current PCA-based post-processing algorithms in terms of the detectability of companions at near real-time speed, a quality that will be useful for future direct imaging surveys. We used randomized low-rank approximation methods recently proposed in the machine learning literature, coupled with entry-wise thresholding to decompose an ADI image sequence locally into low-rank, sparse, and Gaussian noise components (LLSG). This local three-term decomposition separates the starlight and the associated speckle noise from the planetary signal, which mostly remains in the sparse term. We tested the performance of our new algorithm on a long ADI sequence obtained on beta Pictoris with VLT/NACO. Compared to a standard PCA approach, LLSG decomposition reaches a higher signal-to-noise ratio and has an overall better performance in the receiver operating characteristic space. (abridged).

Data processing constitutes a critical component of high-contrast exoplanet imaging. Its role is almost as important as the choice of a coronagraph or a wavefront control system, and it is intertwined with the chosen observing strategy. Among the data processing techniques for angular differential imaging (ADI), the most recent is the family of principal component analysis (PCA) based algorithms. PCA serves, in this case, as a subspace projection technique for constructing a reference point spread function (PSF) that can be subtracted from the science data for boosting the detectability of potential companions present in the data. Unfortunately, when building this reference PSF from the science data itself, PCA comes with certain limitations such as the sensitivity of the lower dimensional orthogonal subspace to non-Gaussian noise. Inspired by recent advances in machine learning algorithms such as robust PCA, we aim to propose a localized subspace projection technique that surpasses current PCA-based post-processing algorithms in terms of the detectability of companions at near real-time speed, a quality that will be useful for future direct imaging surveys. We used randomized low-rank approximation methods recently proposed in the machine learning literature, coupled with entry-wise thresholding to decompose an ADI image sequence locally into low-rank, sparse, and Gaussian noise components (LLSG). This local three-term decomposition separates the starlight and the associated speckle noise from the planetary signal, which mostly remains in the sparse term. We tested the performance of our new algorithm on a long ADI sequence obtained on beta Pictoris with VLT/NACO. Compared to a standard PCA approach, LLSG decomposition reaches a higher signal-to-noise ratio and has an overall better performance in the receiver operating characteristic space. (abridged).

The third edition of the "international - Traveling Workshop on Interactions between Sparse models and Technology" (iTWIST) took place in Aalborg, the 4th largest city in Denmark situated beautifully in the northern part of the country, from the 24th to 26th of August 2016. The workshop venue was at the Aalborg University campus. One implicit objective of this biennial workshop is to foster collaboration between international scientific teams by disseminating ideas through both specific oral/poster presentations and free discussions. For this third edition, iTWIST'16 gathered about 50 international participants and features 8 invited talks, 12 oral presentations, and 12 posters on the following themes, all related to the theory, application and generalization of the "sparsity paradigm": Sparsity-driven data sensing and processing (e.g., optics, computer vision, genomics, biomedical, digital communication, channel estimation, astronomy); Application of sparse models in non-convex/non-linear inverse problems (e.g., phase retrieval, blind deconvolution, self calibration); Approximate probabilistic inference for sparse problems; Sparse machine learning and inference; "Blind" inverse problems and dictionary learning; Optimization for sparse modelling; Information theory, geometry and randomness; Sparsity? What's next? (Discrete-valued signals; Union of low-dimensional spaces, Cosparsity, mixed/group norm, model-based, low-complexity models, ...); Matrix/manifold sensing/processing (graph, low-rank approximation, ...); Complexity/accuracy tradeoffs in numerical methods/optimization; Electronic/optical compressive sensors (hardware)., Comment: 69 pages, 22 extended abstracts, iTWIST'16 website: http://www.itwist16.es.aau.dk

Many problems in the sciences and engineering can be rephrased as optimization problems on matrix search spaces endowed with a so-called manifold structure. This book shows how to exploit the special structure of such problems to develop efficient numerical algorithms. It places careful emphasis on both the numerical formulation of the algorithm and its differential geometric abstraction--illustrating how good algorithms draw equally from the insights of differential geometry, optimization, and numerical analysis. Two more theoretical chapters provide readers with the background in differential geometry necessary to algorithmic development. In the other chapters, several well-known optimization methods such as steepest descent and conjugate gradients are generalized to abstract manifolds. The book provides a generic development of each of these methods, building upon the material of the geometric chapters. It then guides readers through the calculations that turn these geometrically formulated methods into concrete numerical algorithms. The state-of-the-art algorithms given as examples are competitive with the best existing algorithms for a selection of eigenspace problems in numerical linear algebra. Optimization Algorithms on Matrix Manifolds offers techniques with broad applications in linear algebra, signal processing, data mining, computer vision, and statistical analysis. It can serve as a graduate-level textbook and will be of interest to applied mathematicians, engineers, and computer scientists.

Many problems in the sciences and engineering can be rephrased as optimization problems on matrix search spaces endowed with a so-called manifold structure. This book shows how to exploit the special structure of such problems to develop efficient numerical algorithms. It places careful emphasis on both the numerical formulation of the algorithm and its differential geometric abstraction--illustrating how good algorithms draw equally from the insights of differential geometry, optimization, and numerical analysis. Two more theoretical chapters provide readers with the background in differential geometry necessary to algorithmic development. In the other chapters, several well-known optimization methods such as steepest descent and conjugate gradients are generalized to abstract manifolds. The book provides a generic development of each of these methods, building upon the material of the geometric chapters. It then guides readers through the calculations that turn these geometrically formulated methods into concrete numerical algorithms. The state-of-the-art algorithms given as examples are competitive with the best existing algorithms for a selection of eigenspace problems in numerical linear algebra. Optimization Algorithms on Matrix Manifolds offers techniques with broad applications in linear algebra, signal processing, data mining, computer vision, and statistical analysis. It can serve as a graduate-level textbook and will be of interest to applied mathematicians, engineers, and computer scientists.

Many problems in the sciences and engineering can be rephrased as optimization problems on matrix search spaces endowed with a so-called manifold structure. This book shows how to exploit the special structure of such problems to develop efficient numerical algorithms. It places careful emphasis on both the numerical formulation of the algorithm and its differential geometric abstraction--illustrating how good algorithms draw equally from the insights of differential geometry, optimization, and numerical analysis. Two more theoretical chapters provide readers with the background in differential geometry necessary to algorithmic development. In the other chapters, several well-known optimization methods such as steepest descent and conjugate gradients are generalized to abstract manifolds. The book provides a generic development of each of these methods, building upon the material of the geometric chapters. It then guides readers through the calculations that turn these geometrically formulated methods into concrete numerical algorithms. The state-of-the-art algorithms given as examples are competitive with the best existing algorithms for a selection of eigenspace problems in numerical linear algebra. Optimization Algorithms on Matrix Manifolds offers techniques with broad applications in linear algebra, signal processing, data mining, computer vision, and statistical analysis. It can serve as a graduate-level textbook and will be of interest to applied mathematicians, engineers, and computer scientists.

The implicit objective of the biennial "international - Traveling Workshop on Interactions between Sparse models and Technology" (iTWIST) is to foster collaboration between international scientific teams by disseminating ideas through both specific oral/poster presentations and free discussions. For its second edition, the iTWIST workshop took place in the medieval and picturesque town of Namur in Belgium, from Wednesday August 27th till Friday August 29th, 2014. The workshop was conveniently located in "The Arsenal" building within walking distance of both hotels and town center. iTWIST'14 has gathered about 70 international participants and has featured 9 invited talks, 10 oral presentations, and 14 posters on the following themes, all related to the theory, application and generalization of the "sparsity paradigm": Sparsity-driven data sensing and processing; Union of low dimensional subspaces; Beyond linear and convex inverse problem; Matrix/manifold/graph sensing/processing; Blind inverse problems and dictionary learning; Sparsity and computational neuroscience; Information theory, geometry and randomness; Complexity/accuracy tradeoffs in numerical methods; Sparsity? What's next?; Sparse machine learning and inference., Comment: 69 pages, 24 extended abstracts, iTWIST'14 website: http://sites.google.com/site/itwist14

Let f be a real-valued function on a Riemannian submanifold of a Euclidean space, and let f̄ be a local extension of f. We show that the Riemannian Hessian of f can be conveniently obtained from the Euclidean gradient and Hessian of f̄ by means of two m

Optimization on manifolds is a rapidly developing branch of nonlinear optimization. Its focus is on problems where the smooth geometry of the search space can be leveraged to design efficient numerical algorithms. In particular, optimization on manifolds is well-suited to deal with rank and orthogonality constraints. Such structured constraints appear pervasively in machine learning applications, including low-rank matrix completion, sensor network localization, camera network registration, independent component analysis, metric learning, dimensionality reduction and so on. The Manopt toolbox, available at www.manopt.org, is a user-friendly, documented piece of software dedicated to simplify experimenting with state of the art Riemannian optimization algorithms. We aim particularly at reaching practitioners outside our field.

Communities play a crucial role to describe and analyse modern networks. However, the size of those networks has grown tremendously with the increase of computational power and data storage. While various methods have been developed to extract community structures, their computational cost or the difficulty to parallelize existing algorithms make partitioning real networks into communities a challenging problem. In this paper, we propose to alter an efficient algorithm, the Louvain method, such that communities are defined as the connected components of a tree-like assignment graph. Within this framework, we precisely describe the different steps of our algorithm and demonstrate its highly parallelizable nature. We then show that despite its simplicity, our algorithm has a partitioning quality similar to the original method on benchmark graphs and even outperforms other algorithms. We also show that, even on a single processor, our method is much faster and allows the analysis of very large networks.

We study the problem of separating audio sources from a single linear mixture. The goal is to find a decomposition of the single channel spectrogram into a sum of individual contributions associated to a certain number of sources. In this paper, we consider an informed source separation problem in which the input spectrogram is partly annotated. We propose a convex formulation that relies on a nuclear norm penalty to induce low rank for the contributions. We show experimentally that solving this model with a simple subgradient method outperforms a previously introduced nonnegative matrix factorization (NMF) technique, both in terms of source separation quality and computation time., Comment: Submitted to ESANN 2013 conference

Synchronization of rotations is the problem of estimating a set of rotations R_i in SO(n), i = 1, ..., N, based on noisy measurements of relative rotations R_i R_j^T. This fundamental problem has found many recent applications, most importantly in structural biology. We provide a framework to study synchronization as estimation on Riemannian manifolds for arbitrary n under a large family of noise models. The noise models we address encompass zero-mean isotropic noise, and we develop tools for Gaussian-like as well as heavy-tail types of noise in particular. As a main contribution, we derive the Cram\'er-Rao bounds of synchronization, that is, lower-bounds on the variance of unbiased estimators. We find that these bounds are structured by the pseudoinverse of the measurement graph Laplacian, where edge weights are proportional to measurement quality. We leverage this to provide interpretation in terms of random walks and visualization tools for these bounds in both the anchored and anchor-free scenarios. Similar bounds previously established were limited to rotations in the plane and Gaussian-like noise.

Approximate matrix factorization techniques with both nonnegativity and orthogonality constraints, referred to as orthogonal nonnegative matrix factorization (ONMF), have been recently introduced and shown to work remarkably well for clustering tasks such as document classification. In this paper, we introduce two new methods to solve ONMF. First, we show athematical equivalence between ONMF and a weighted variant of spherical k-means, from which we derive our first method, a simple EM-like algorithm. This also allows us to determine when ONMF should be preferred to k-means and spherical k-means. Our second method is based on an augmented Lagrangian approach. Standard ONMF algorithms typically enforce nonnegativity for their iterates while trying to achieve orthogonality at the limit (e.g., using a proper penalization term or a suitably chosen search direction). Our method works the opposite way: orthogonality is strictly imposed at each step while nonnegativity is asymptotically obtained, using a quadratic penalty. Finally, we show that the two proposed approaches compare favorably with standard ONMF algorithms on synthetic, text and image data sets., Comment: 17 pages, 8 figures. New numerical experiments (document and synthetic data sets)

It is well-known that the H2-norm and H∞-norm of a transfer function can differ arbitrarily since both norms reflect fundamentally different properties. However, if the pole structure of the transfer function is known it is possible to bound the H∞-no

In this paper we derive a canonical lower bound for the autocovariance function of any unbiased transfer-function estimator. As a generalization of the Cramér-Rao bound, the Cramér-Rao kernel that we define can be derived without parametrizing the model

Various optimal control strategies exist in the literature. Prominent approaches are Robust Control and Linear Quadratic Regulators, the first one being related to the H∞ norm of a system, the second one to the H2 norm. In 1994, F. De Bruyne et al [1] s

Many problems in the sciences and engineering can be rephrased as optimization problems on matrix search spaces endowed with a so-called manifold structure. This book shows how to exploit the special structure of such problems to develop efficient numeric

The classes of continuous-time flows on Rn×p that induce the same flow on the set of p- dimensional subspaces of Rn×p are described. The power flow is briefly reviewed in this framework, and a subspace generalization of the Rayleigh quotient flow [Linea

Consider a space M endowed with two or more Lie group actions. Under a certain condition on the orbits of the Lie group actions, we show how to construct a flow on M that projects to prescribed flows on the orbit spaces of the group actions. Hence, in order to design a flow that converges to the intersection of given orbits, it suffices to design flows on the various orbit spaces that display convergence to the desired orbits, and then to lift these flows to M using the proposed procedure. We illustrate the technique by creating a flow for principal component analysis. The flow projects to a flow on the Grassmann manifold that achieves principal subspace analysis and to a flow on the "shape" manifold that converges to the set of orthonormal matrices.

In the early eighties Lojasiewicz [in Seminari di Geometric, 1982-1983, Università di Bologna, Istituto di Geometria, Dipartimento di Matematica, 1984, pp. 115-117] proved that a bounded solution of a gradient flow for an analytic cost function converges

We propose a Newton-like iteration that evolves on the set of fixed dimensional subspaces of ℝn and converges locally cubically to the invariant subspaces of a symmetric matrix. This iteration is compared in terms of numerical cost and global behavior w

We give simple formulas for the canonical metric, gradient, Lie derivative, Riemannian connection, parallel translation, geodesics and distance on the Grassmann manifold of p-planes in ℝn. In these formulas, p-planes are represented as the column space

A continuous-time differential equation inspired by the Rayleigh quotient iteration for a symmetric matrix A is studied. The set of all continuous solutions, termed the Rayleigh quotient flow, is shown to be a time-scaled version of the Newton flow for Rayleigh quotient cost functional. The scaling factor ensures that the rate of variation of the Rayleigh quotient is constant and positive along solutions. This interpretation leads to a precise phase portrait for the Rayleigh quotient flow. It is shown that a complete (non-degenerate) solution of the Rayleigh quotient flow visits each of the eigenvectors of A in ascending order.

We study the global behaviour of a Newton algorithm on the Grassmann manifold for invariant subspace computation. It is shown that the basins of attraction of the invariant subspaces may collapse in case of small eigenvalue gaps. A Levenberg-Marquardt-like modification of the algorithm with low numerical cost is proposed. A simple strategy for choosing the parameter is shown to dramatically enlarge the basins of attraction of the invariant subspaces while preserving the fast local convergence.

The classical Rayleigh quotient iteration (RQI) allows one to compute a one-dimensional invariant subspace of a symmetric matrix A. Here we propose a generalization of the RQI which computes a p-dimensional invariant subspace of A. Cubic convergence is preserved and the cost per iteration is low compared to other methods proposed in the literature.