Alexander V. Mamonov

Alexander V. Mamonov, Yen-Hsi Richard Tsai

We consider a problem of identification of point sources in time dependent advection-diffusion systems with a non-linear reaction term. The linear counterpart of the problem in question can be reduced to solving a system of non-linear algebraic equations via the use of adjoint equations. We extend this approach by constructing an algorithm that solves the problem iteratively to account for the non-linearity of the reaction term. We study the question of improving the quality of source identification by adding more measurements adaptively using the solution obtained previously with a smaller number of measurements.

Liliana Borcea, Vladimir Druskin, Alexander V. Mamonov et al.

The motivation of this work is an inverse problem for the acoustic wave equation, where an array of sensors probes an unknown medium with pulses and measures the scattered waves. The goal of the inversion is to determine from these measurements the structure of the scattering medium, modeled by a spatially varying acoustic impedance function. Many inversion algorithms assume that the mapping from the unknown impedance to the scattered waves is approximately linear. The linearization, known as the Born approximation, is not accurate in strongly scattering media, where the waves undergo multiple reflections before they reach the sensors in the array. Thus, the reconstructions of the impedance have numerous artifacts. The main result of the paper is a novel, linear-algebraic algorithm that uses a reduced order model (ROM) to map the data to those corresponding to the single scattering (Born) model. The ROM construction is based only on the measurements at the sensors in the array. The ROM is a proxy for the wave propagator operator, that propagates the wave in the unknown medium over the duration of the time sampling interval. The output of the algorithm can be input into any off-the-shelf inversion software that incorporates state of the art linear inversion algorithms to reconstruct the unknown acoustic impedance.

Liliana Borcea, Vladimir Druskin, Alexander V. Mamonov et al.

We propose a novel numerical inversion algorithm for the coefficients of parabolic partial differential equations, based on model reduction. The study is motivated by the application of controlled source electromagnetic exploration, where the unknown is the subsurface electrical resistivity and the data are time resolved surface measurements of the magnetic field. The algorithm presented in this paper considers inversion in one and two dimensions. The reduced model is obtained with rational interpolation in the frequency (Laplace) domain and a rational Krylov subspace projection method. It amounts to a nonlinear mapping from the function space of the unknown resistivity to the small dimensional space of the parameters of the reduced model. We use this mapping as a nonlinear preconditioner for the Gauss-Newton iterative solution of the inverse problem. The advantage of the inversion algorithm is twofold. First, the nonlinear preconditioner resolves most of the nonlinearity of the problem. Thus the iterations are less likely to get stuck in local minima and the convergence is fast. Second, the inversion is computationally efficient because it avoids repeated accurate simulations of the time-domain response. We study the stability of the inversion algorithm for various rational Krylov subspaces, and assess its performance with numerical experiments.

Liliana Borcea, Vladimir Druskin, Alexander V. Mamonov et al.

We introduce a novel algorithm for nonlinear processing of data gathered by an active array of sensors which probes a medium with pulses and measures the resulting waves. The algorithm is motivated by the application of array imaging. We describe it for a generic hyperbolic system that applies to acoustic, electromagnetic or elastic waves in a scattering medium modeled by an unknown coefficient called the reflectivity. The goal of imaging is to invert the nonlinear mapping from the reflectivity to the array data. Many existing imaging methodologies ignore the nonlinearity i.e., operate under the assumption that the Born (single scattering) approximation is accurate. This leads to image artifacts when multiple scattering is significant. Our algorithm seeks to transform the array data to those corresponding to the Born approximation, so it can be used as a pre-processing step for any linear inversion method. The nonlinear data transformation algorithm is based on a reduced order model defined by a proxy wave propagator operator that has four important properties. First, it is data driven, meaning that it is constructed from the data alone, with no knowledge of the medium. Second, it can be factorized in two operators that have an approximately affine dependence on the unknown reflectivity. This allows the computation of the Fréchet derivative of the reflectivity to the data mapping which gives the Born approximation. Third, the algorithm involves regularization which balances numerical stability and data fitting with accuracy of the order of the standard deviation of additive data noise. Fourth, the algebraic nature of the algorithm makes it applicable to scalar (acoustic) and vectorial (elastic, electromagnetic) wave data without any specific modifications.

Vladimir Druskin, Alexander V. Mamonov, Mikhail Zaslavsky

We developed a novel reduced-order multi-scale method for solving large time-domain wavefield simulation problems. Our algorithm consists of two main stages. During the first "off-line" stage the fine-grid operator (of the graph Laplacian type} is partitioned on coarse cells (subdomains). Then projection-type multi-scale reduced order models (ROMs) are computed for the coarse cell operators. The off-line stage is embarrassingly parallel as ROM computations for the subdomains are independent of each other. It also does not depend on the number of simulated sources (inputs) and it is performed just once before the entire time-domain simulation. At the second "on-line" stage the time-domain simulation is performed within the obtained multi-scale ROM framework. The crucial feature of our formulation is the representation of the ROMs in terms of matrix Stieltjes continued fractions (S-fractions). The layered structure of the S-fraction introduces several hidden layers in the ROM representation, that results in the block-tridiagonal dynamic system within each coarse cell. This allows us to sparsify the obtained multi-scale subdomain operator ROMs and to reduce the communications between the adjacent subdomains which is highly beneficial for a parallel implementation of the on-line stage. Our approach suits perfectly the high performance computing architectures, however in this paper we present rather promising numerical results for a serial computing implementation only. These results include 3D acoustic and multi-phase anisotropic elastic problems.

Liliana Borcea, Fernando Guevara Vasquez, Alexander V. Mamonov

We propose a discrete approach for solving an inverse problem for Schrödinger's equation in two dimensions, where the unknown potential is to be determined from boundary measurements of the Dirichlet to Neumann map. For absorptive potentials, and in the continuum, it is known that by using the Liouville identity we obtain an inverse conductivity problem. Its discrete analogue is to find a resistor network that matches the measurements, and is well understood. Here we show how to use a discrete Liouville identity to transform its solution to that of Schrödinger's problem. The discrete Schrödinger potential given by the discrete Liouville identity can be used to reconstruct the potential in the continuum in two ways. First, we can obtain a direct but coarse reconstruction by interpreting the values of the discrete Schrödinger potential as averages of the continuum Schrödinger potential on a special sensitivity grid. Second, the discrete Schrödinger potential may be used to reformulate the conventional nonlinear output least squares optimization formulation of the inverse Schrödinger problem. Instead of minimizing the boundary measurement misfit, we minimize the misfit between the discrete Schrödinger potentials. This results in a better behaved optimization problem that converges in a single Gauss-Newton iteration, and gives good quality reconstructions of the potential, as illustrated by the numerical results.

Liliana Borcea, Wei Li, Alexander V. Mamonov et al.

We consider the problem of optical imaging of small nonlinear scatterers in random media. We propose an extension of coherent interferometric imaging (CINT) that applies to scatterers that emit second-harmonic light. We compare this method to a nonlinear version of migration imaging and find that the images obtained by CINT are more robust to statistical fluctuations. This finding is supported by a resolution analysis that is carried out in the setting of geometrical optics in random media. It is also consistent with numerical simulations for which the assumptions of the geometrical optics model do not hold.

Malihehsadat Chavooshi, Alexander V. Mamonov

We propose a novel approach to unsupervised learning by constructing a non-linear embedding of the data into a low-dimensional space followed by any conventional clustering algorithm. The embedding promotes clusterability of the data and is comprised of two mappings: the encoder of an autoencoder neural network and the output of UMAP algorithm. The autoencoder is trained with a composite loss function that incorporates both a conventional data reconstruction as a regularization component and a clustering-promoting component built using the spectral graph theory. The two embeddings and the subsequent clustering are integrated into a three-stage unsupervised learning framework, referred to as Autoencoded UMAP-Enhanced Clustering (AUEC). When applied to MNIST data, AUEC significantly outperforms the state-of-the-art techniques in terms of clustering accuracy.

Vladimir Druskin, Alexander V. Mamonov, Mikhail Zaslavsky

Graph-Laplacians and their spectral embeddings play an important role in multiple areas of machine learning. This paper is focused on graph-Laplacian dimension reduction for the spectral clustering of data as a primary application. Spectral embedding provides a low-dimensional parametrization of the data manifold which makes the subsequent task (e.g., clustering) much easier. However, despite reducing the dimensionality of data, the overall computational cost may still be prohibitive for large data sets due to two factors. First, computing the partial eigendecomposition of the graph-Laplacian typically requires a large Krylov subspace. Second, after the spectral embedding is complete, one still has to operate with the same number of data points. For example, clustering of the embedded data is typically performed with various relaxations of k-means which computational cost scales poorly with respect to the size of data set. In this work, we switch the focus from the entire data set to a subset of graph vertices (target subset). We develop two novel algorithms for such low-dimensional representation of the original graph that preserves important global distances between the nodes of the target subset. In particular, it allows to ensure that target subset clustering is consistent with the spectral clustering of the full data set if one would perform such. That is achieved by a properly parametrized reduced-order model (ROM) of the graph-Laplacian that approximates accurately the diffusion transfer function of the original graph for inputs and outputs restricted to the target subset. Working with a small target subset reduces greatly the required dimension of Krylov subspace and allows to exploit the conventional algorithms (like approximations of k-means) in the regimes when they are most robust and efficient.

Vladimir Druskin, Alexander V. Mamonov, Mikhail Zaslavsky

We introduce a novel nonlinear imaging method for the acoustic wave equation based on data-driven model order reduction. The objective is to image the discontinuities of the acoustic velocity, a coefficient of the scalar wave equation from the discretely sampled time domain data measured at an array of transducers that can act as both sources and receivers. We treat the wave equation along with transducer functionals as a dynamical system. A reduced order model (ROM) for the propagator of such system can be computed so that it interpolates exactly the measured time domain data. The resulting ROM is an orthogonal projection of the propagator on the subspace of the snapshots of solutions of the acoustic wave equation. While the wavefield snapshots are unknown, the projection ROM can be computed entirely from the measured data, thus we refer to such ROM as data-driven. The image is obtained by backprojecting the ROM. Since the basis functions for the projection subspace are not known, we replace them with the ones computed for a known smooth kinematic velocity model. A crucial step of ROM construction is an implicit orthogonalization of solution snapshots. It is a nonlinear procedure that differentiates our approach from the conventional linear imaging methods (Kirchhoff migration and reverse time migration - RTM). It resolves all dynamical behavior captured by the data, so the error from the imperfect knowledge of the velocity model is purely kinematic. This allows for almost complete removal of multiple reflection artifacts, while simultaneously improving the resolution in the range direction compared to conventional RTM.

Alexander V. Mamonov, Vladimir Druskin, Mikhail Zaslavsky

We introduce a novel nonlinear seismic imaging method based on model order reduction. The reduced order model (ROM) is an orthogonal projection of the wave equation propagator operator on the subspace of the snapshots of the solutions of the wave equation. It can be computed entirely from the knowledge of the measured time domain seismic data. The image is a backprojection of the ROM using the subspace basis for the known smooth kinematic velocity model. The implicit orthogonalization of solution snapshots is a nonlinear procedure that differentiates our approach from the conventional linear methods (Kirchhoff, RTM). It allows for the removal of multiple reflection artifacts. It also enables us to estimate the magnitude of the reflectors similarly to the true amplitude migration algorithms.

Alexander V. Mamonov, Isabel N. Figueiredo, Pedro N. Figueiredo et al.

Colorectal polyps are important precursors to colon cancer, a major health problem. Colon capsule endoscopy (CCE) is a safe and minimally invasive examination procedure, in which the images of the intestine are obtained via digital cameras on board of a small capsule ingested by a patient. The video sequence is then analyzed for the presence of polyps. We propose an algorithm that relieves the labor of a human operator analyzing the frames in the video sequence. The algorithm acts as a binary classifier, which labels the frame as either containing polyps or not, based on the geometrical analysis and the texture content of the frame. The geometrical analysis is based on a segmentation of an image with the help of a mid-pass filter. The features extracted by the segmentation procedure are classified according to an assumption that the polyps are characterized as protrusions that are mostly round in shape. Thus, we use a best fit ball radius as a decision parameter of a binary classifier. We present a statistical study of the performance of our approach on a data set containing over 18,900 frames from the endoscopic video sequences of five adult patients. The algorithm demonstrates a solid performance, achieving 47% sensitivity per frame and over 81% sensitivity per polyp at a specificity level of 90%. On average, with a video sequence length of 3747 frames, only 367 false positive frames need to be inspected by a human operator.