Publications

We propose new variants of the sketch-and-project method for solving large scale ridge regression problems. First, we propose a new momentum alternative and provide a theorem showing it can speed up the convergence of sketch-and-project, through a fast sublinear convergence rate. We carefully delimit under what settings this new sublinear rate is faster than the previously known linear rate of convergence of sketch-and-project without momentum. Second, we consider combining the sketch-and-project method with new modern sketching methods such as Count sketch, SubCount sketch (a new method we propose), and subsampled Hadamard transforms. We show experimentally that when combined with the sketch-and-project method, the (Sub)Count sketch is very effective on sparse data and the standard Subsample sketch is effective on dense data. Indeed, we show that these sketching methods, combined with our new momentum scheme, result in methods that are competitive even when compared to the conjugate gradient method on real large scale data. On the contrary, we show the subsampled Hadamard transform does not perform well in this setting, despite the use of fast Hadamard transforms, and nor do recently proposed acceleration schemes work well in practice. To support all of our experimental findings, and invite the community to validate and extend our results, with this paper we are also releasing an open source software package: RidgeSketch. We designed this object-oriented package in Python for testing sketch-and-project methods and benchmarking ridge regression solvers. RidgeSketch is highly modular, and new sketching methods can easily be added as subclasses. We provide code snippets of our package in the appendix.

Recovering the phase of the Fourier transform is a ubiquitous problem in imaging applications from astronomy to nanoscale X-ray diffraction imaging. Despite the efforts of a multitude of scientists, from astronomers to mathematicians, there is, as yet, no satisfactory theoretical or algorithmic solution to this class of problems. Written for mathematicians, physicists and engineers working in image analysis and reconstruction, this book introduces a conceptual, geometric framework for the analysis of these problems, leading to a deeper understanding of the essential, algorithmically independent, difficulty of their solutions. Using this framework, the book studies standard algorithms and a range of theoretical issues in phase retrieval and provides several new algorithms and approaches to this problem with the potential to improve the reconstructed images. The book is lavishly illustrated with the results of numerous numerical experiments that motivate the theoretical development and place it in the context of practical applications.

An efficient, reliable, and interpretable global solution method, the Deep learning-based algorithm for Heterogeneous Agent Models (DeepHAM), is proposed for solving high dimensional heterogeneous agent models with aggregate shocks. The state distribution is approximately represented by a set of optimal generalized moments. Deep neural networks are used to approximate the value and policy functions, and the objective is optimized over directly simulated paths. In addition to being an accurate global solver, this method has three additional features. First, it is computationally efficient in solving complex heterogeneous agent models, and it does not suffer from the curse of dimensionality. Second, it provides a general and interpretable representation of the distribution over individual states, which is crucial in addressing the classical question of whether and how heterogeneity matters in macroeconomics. Third, it solves the constrained efficiency problem as easily as it solves the competitive equilibrium, which opens up new possibilities for studying optimal monetary and fiscal policies in heterogeneous agent models with aggregate shocks.

In this paper, we consider the reconstruction of the shape and the impedance function of an obstacle from measurements of the scattered field at a collection of receivers outside the object. The data is assumed to be generated by plane waves impinging on the unknown obstacle from multiple directions and at multiple frequencies. This inverse problem can be reformulated as an optimization problem: that of finding band-limited shape and impedance functions which minimize the L2 distance between the computed value of the scattered field at the receivers and the given measurement data. The optimization problem is highly non-linear, non-convex, and ill-posed. Moreover, the objective function is computationally expensive to evaluate (since a large number of Helmholtz boundary value problems need to be solved at every iteration in the optimization loop). The recursive linearization approach (RLA) proposed by Chen has been successful in addressing these issues in the context of recovering the sound speed of an inhomogeneous object or the shape of a sound-soft obstacle. We present an extension of the RLA for the recovery of both the shape and impedance functions of the object. The RLA is, in essence, a continuation method in frequency where a sequence of single frequency inverse problems is solved. At each higher frequency, one attempts to recover incrementally higher resolution features using a step assumed to be small enough to ensure that the initial guess obtained at the preceding frequency lies in the basin of attraction for Newton’s method at the new frequency. We demonstrate the effectiveness of this approach with several numerical examples. Surprisingly, we find that one can recover the shape with high accuracy even when the measurements are generated by sound-hard or sound-soft objects, eliminating the need to know the precise boundary conditions appropriate for modeling the object under consideration. While the method is effective in obtaining high-quality reconstructions for many complicated geometries and impedance functions, a number of interesting open questions remain regarding the convergence behavior of the approach. We present numerical experiments that suggest underlying mechanisms of success and failure, pointing out areas where improvements could help lead to robust and automatic tools for the solution of inverse obstacle scattering problems.

In recent years, tremendous progress has been made on numerical algorithms for solving partial differential equations (PDEs) in a very high dimension, using ideas from either nonlinear (multilevel) Monte Carlo or deep learning. They are potentially free of the curse of dimensionality for many different applications and have been proven to be so in the case of some nonlinear Monte Carlo methods for nonlinear parabolic PDEs. In this paper, we review these numerical and theoretical advances. In addition to algorithms based on stochastic reformulations of the original problem, such as the multilevel Picard iteration and the deep backward stochastic differential equations method, we also discuss algorithms based on the more traditional Ritz, Galerkin, and least square formulations. We hope to demonstrate to the reader that studying PDEs as well as control and variational problems in very high dimensions might very well be among the most promising new directions in mathematics and scientific computing in the near future.

A major problem in motor control is understanding how the brain plans and executes proper movements in the face of delayed and noisy stimuli. A prominent framework for addressing such control problems is Optimal Feedback Control (OFC). OFC generates control actions that optimize behaviorally relevant criteria by integrating noisy sensory stimuli and the predictions of an internal model using the Kalman filter or its extensions. However, a satisfactory neural model of Kalman filtering and control is lacking because existing proposals have the following limitations: not considering the delay of sensory feedback, training in alternating phases, and requiring knowledge of the noise covariance matrices, as well as that of systems dynamics. Moreover, the majority of these studies considered Kalman filtering in isolation, and not jointly with control. To address these shortcomings, we introduce a novel online algorithm which combines adaptive Kalman filtering with a model free control approach (i.e., policy gradient algorithm). We implement this algorithm in a biologically plausible neural network with local synaptic plasticity rules. This network performs system identification and Kalman filtering, without the need for multiple phases with distinct update rules or the knowledge of the noise covariances. It can perform state estimation with delayed sensory feedback, with the help of an internal model. It learns the control policy without requiring any knowledge of the dynamics, thus avoiding the need for weight transport. In this way, our implementation of OFC solves the credit assignment problem needed to produce the appropriate sensory-motor control in the presence of stimulus delay.

Gauge symmetries play a key role in physics appearing in areas such as quantum field theories of the fundamental particles and emergent degrees of freedom in quantum materials. Motivated by the desire to efficiently simulate many-body quantum systems with exact local gauge invariance, gauge equivariant neural-network quantum states are introduced, which exactly satisfy the local Hilbert space constraints necessary for the description of quantum lattice gauge theory with Zd gauge group on different geometries. Focusing on the special case of Z2 gauge group on a periodically identified square lattice, the equivariant architecture is analytically shown to contain the loop-gas solution as a special case. Gauge equivariant neural-network quantum states are used in combination with variational quantum Monte Carlo to obtain compact descriptions of the ground state wavefunction for the Z2 theory away from the exactly solvable limit, and to demonstrate the confining/deconfining phase transition of the Wilson loop order parameter.

The Laplace-Beltrami problem on closed surfaces embedded in three dimensions arises in many areas of physics, including molecular dynamics (surface diffusion), electromagnetics (harmonic vector fields), and fluid dynamics (vesicle deformation). Using classical potential theory,the Laplace-Beltrami operator can be pre-/post-conditioned with integral operators whose kernel is translation invariant, resulting in well-conditioned Fredholm integral equations of the second-kind. These equations have the standard Laplace kernel from potential theory, and therefore the equations can be solved rapidly and accurately using a combination of fast multipole methods (FMMs) and high-order quadrature corrections. In this work we detail such a scheme, presenting two alternative integral formulations of the Laplace-Beltrami problem, each of whose solution can be obtained via FMM acceleration. We then present several applications of the solvers, focusing on the computation of what are known as harmonic vector fields, relevant for many applications in electromagnetics. A battery of numerical results are presented for each application, detailing the performance of the solver in various geometries.

The brain effortlessly solves blind source separation (BSS) problems, but the algorithm it uses remains elusive. In signal processing, linear BSS problems are often solved by Independent Component Analysis (ICA). To serve as a model of a biological circuit, the ICA neural network (NN) must satisfy at least the following requirements: 1. The algorithm must operate in the online setting where data samples are streamed one at a time, and the NN computes the sources on the fly without storing any significant fraction of the data in memory. 2. The synaptic weight update is local, i.e., it depends only on the biophysical variables present in the vicinity of a synapse. Here, we propose a novel objective function for ICA from which we derive a biologically plausible NN, including both the neural architecture and the synaptic learning rules. Interestingly, our algorithm relies on modulating synaptic plasticity by the total activity of the output neurons. In the brain, this could be accomplished by neuromodulators, extracellular calcium, local field potential, or nitric oxide.

Recently introduced shadow tomography protocols use “classical shadows” of quantum states to predict many target functions of an unknown quantum state. Unlike full quantum state tomography, shadow tomography does not insist on accurate recovery of the density matrix for high rank mixed states. Yet, such a protocol makes multiple accurate predictions with high confidence, based on a moderate number of quantum measurements. One particular influential algorithm, proposed by Huang et al. [Huang, Kueng, and Preskill, Nat. Phys. 16, 1050 (2020)], requires additional circuits for performing certain random unitary transformations. In this paper, we avoid these transformations but employ an arbitrary informationally complete positive operator-valued measure and show that such a procedure can compute k-bit correlation functions for quantum states reliably. We also show that, for this application, we do not need the median of means procedure of Huang et al. Finally, we discuss the contrast between the computation of correlation functions and fidelity of reconstruction of low rank density matrices.