Next Article in Journal
Periodic DFTB for Supported Clusters: Implementation and Application on Benzene Dimers Deposited on Graphene
Next Article in Special Issue
Inverse Modeling of Hydrologic Parameters in CLM4 via Generalized Polynomial Chaos in the Bayesian Framework
Previous Article in Journal
Clustering Analysis for the Pareto Optimal Front in Multi-Objective Optimization
Previous Article in Special Issue
Bayesian Instability of Optical Imaging: Ill Conditioning of Inverse Linear and Nonlinear Radiative Transfer Equation in the Fluid Regime
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

DIAS: A Data-Informed Active Subspace Regularization Framework for Inverse Problems

1
Department of Aerospace Engineering and Engineering Mechanics, UT Austin, Austin, TX 78712, USA
2
The Oden Institute of Computational Engineering and Sciences, UT Austin, Austin, TX 78712, USA
3
Department of Aerospace Engineering and Engineering Mechanics, The Oden Institute for Computational Engineering and Sciences, UT Austin, Austin, TX 78712, USA
*
Author to whom correspondence should be addressed.
Submission received: 25 January 2022 / Revised: 28 February 2022 / Accepted: 8 March 2022 / Published: 11 March 2022
(This article belongs to the Special Issue Inverse Problems with Partial Data)

Abstract

:
This paper presents a regularization framework that aims to improve the fidelity of Tikhonov inverse solutions. At the heart of the framework is the data-informed regularization idea that only data-uninformed parameters need to be regularized, while the data-informed parameters, on which data and forward model are integrated, should remain untouched. We propose to employ the active subspace method to determine the data-informativeness of a parameter. The resulting framework is thus called a data-informed (DI) active subspace (DIAS) regularization. Four proposed DIAS variants are rigorously analyzed, shown to be robust with the regularization parameter and capable of avoiding polluting solution features informed by the data. They are thus well suited for problems with small or reasonably small noise corruptions in the data. Furthermore, the DIAS approaches can effectively reuse any Tikhonov regularization codes/libraries. Though they are readily applicable for nonlinear inverse problems, we focus on linear problems in this paper in order to gain insights into the framework. Various numerical results for linear inverse problems are presented to verify theoretical findings and to demonstrate advantages of the DIAS framework over the Tikhonov, truncated SVD, and the TSVD-based DI approaches.

1. Introduction with Related Work and Novelties

Over the past few decades, the development of new regularization techniques has played an important role in addressing the ill-posedness of inverse problems. A few example applications include image reconstruction in X-ray tomography [1,2] and various partial differential equation constrained inversion problems, for instance, inverse scattering and parameter identification [3,4,5,6,7]. Popular techniques, such as Tikhonov regularization and truncated singular value decomposition (truncated SVD), are ubiquitous in practical inverse problems [8,9,10]. However, one particular challenge of Tikhonov-based regularization strategies is the need to specify a suitable choice of the regularization parameter. Suboptimal choices can lead to excessively smooth or unstable reconstructions. While methods such as Morozov’s discrepancy principle, the L-curve criterion, and cross validation are popular for choosing a satisfactory regularization parameter, these are often computationally expensive and may require the computation of multiple inverse solutions. Even when the “right” regularization parameter is chosen, the smoothing effect leading to smeared reconstructions is unavoidable. On the other hand, spectral-decomposition-based methods, such as truncated SVD, aim to avoid regularizing data-informative features while applying infinite regularization on the rest. It is, however, not trivial to determine how many dominant modes should be retained. Furthermore, these SVD-typed approaches are typically suitable only for linear inverse problems.
Other regularization techniques have been proposed to combat the smoothing effect of Tikhonov regularization. Total variation (TV) regularization is designed with an anisotropic diffusion mechanism to ensure discontinuities and sharp interfaces in the inverse solution to be maximally preserved [9,10,11]. One problem with TV regularization is that, due to the non-differentiability of the TV functional, it could produce a staircasing effect [12]. To overcome this issue, smooth approximations and sophisticated optimization methods have been developed [9,13]. One of the reasons why Tikhonov and TV regularizations are popular is their convexity. This is a particularly appealing feature for inversion techniques using optimization. Though less popular, non-convex regularization strategies [14,15,16] are viable alternatives. However, inverse solutions with non-convex regularizations also require advanced optimization methods, such as alternating direction method of multipliers (ADMM) [17,18,19] or iteratively reweighted least squares (IRLS) [20].
In our previous work [21], inspired by the truncated SVD method, we put forward a new regularization approach, called data-informed (DI) regularization. DI was designed to overcome the ill-posed nature of inverse problems by placing regularization only where it is necessary in order to preserve the fidelity of the reconstructions. This is accomplished by a two-step process: first, find data-informed directions in the solution space, and then apply Tikhonov regularization only in the data-uninformed directions. We theoretically showed that DI is a valid regularization strategy. Numerically, we demonstrated that the reconstruction accuracy is robust for a wide range of regularization parameter values and DI outperforms both Tikhonov and truncated SVD (TSVD) for various inverse and imaging problems. Since DI, as well as TSVD, exploits the SVD of the forward operator, extension is necessary for nonlinear inverse problems. One straightforward extension for the Newton-based optimizer is to apply the DI approach at each Newton iteration by linearizing the forward operator around the current parameter estimate. However, this approach can require significant additional computation, especially in high dimensions.
Meanwhile, one recent tool for studying sensitivity is the active subspace (AS) method introduced in [22]. The active subspace is designed to study the sensitivity of a given nonlinear function with respect to the input parameters. The key idea is to identify a set of the directions in the input space that, on average, contribute most to the variation of the function. The beauty of the approach is that it requires the computation of the active subspace only once for any nonlinear input–output map. An application where AS is particularly useful is dimension reduction [23,24,25,26,27]. Due to the computational cost involved in dealing with high-dimensional parameter spaces in traditional Markov chain Monte Carlo simulations (MCMC), the active subspace method was used to accelerate MCMC [28]. The work in [29] combined the active subspace method and proper orthogonal decomposition with interpolation to obtain a better reconstruction of the modal coefficients when a small number of snapshots are used. More recently, the active subspace method was adopted to develop neural network based surrogate models [30].
In this work, we equip the DI idea with the active space approach and we call this combination the data-informed active subspace (DIAS) regularization strategy. The DIAS method retains all of the DI advantages, including being robust with respect to the regularization parameter, and avoids polluting solution features informed by the data while estimating the data-informed subspace better than the original DI approach. Unlike DI or Tikhonov, DIAS takes into account the uncertain nature of the inverse problem via the active subspace using the prior distribution of the unknown parameter. More importantly, it is applicable not only for linear, but also seamlessly for nonlinear inverse problems. Indeed, the active subspace, thanks to its global nature, is computed once at the beginning regardless of the linear or nonlinear nature of the inverse problem under consideration. Furthermore, the DIAS approach can effectively reuse any Tikhonov regularization codes/libraries. In fact, DIAS regularization can be considered a special Tikhonov regularization once the DIAS regularization operator is constructed. Table 1 summarizes the main advantageous features of DIAS approaches over the classical Tikhonov regularization method and the DI method [21].
The paper is structured as follows. Section 2 briefly reviews the standard (uncentered) active subspace method and discusses its centered variation. The relationship between uncentered and centered active subspaces is then rigorously investigated. Each of these active subspace approaches can be incorporated into the DI framework with either the original data misfit or an approximate one. For linear inverse problems, Section 3 proposes two DIAS regularization methods with approximate data misfit: one using centered active subspaces (cDIAS-A) and the other using uncentered active subspaces (DIAS-A). An important result is that the truncated SVD approach is a special case of cDIAS-A. Similarly, two DIAS approaches with the original (full) data misfit—one using centered active subspaces (cDIAS-F) and the other using uncentered active subspaces (DIAS-F)—are presented for linear inverse problems in Section 4. It is in this section that the practical aspect of the four DIAS variants is discussed. The full data misfit variants, DIAS-F and cDIAS-F, are non-intrusive, while the approximate misfit counterparts are intrusive. We also show that the DI approach in [21] is a special case of cDIAS-F. Various numerical results for linear inverse problems are presented in Section 5 to expose the pros and cons of each of the DIAS variants and to compare them with Tikhonov regularization. Section 6 concludes the paper with future work.

2. Centered versus Uncentered Active Subspace Methods for Linear Inverse Problems

In this section, we first recall some results on the active subspaces (AS) method that are useful for our developments. A detailed derivation of the active subspace can be consulted in [28]. The main issue that we investigate is the difference between the uncentered AS and the centered AS. To begin, let f ( x ) be a differentiable function from R p to R and its gradient vector be denoted as f ( x ) R p . The key object in the AS method is the uncentered covariance matrix C R p × p defined by
C = f ( x ) f ( x ) T ρ ( x ) d x = W D W T ,
where ρ ( x ) is a probability density in R p (we assume it to be a centered Gaussian N 0 , Γ in this paper, where Γ is the covariance matrix), and W , D are matrices containing eigenpairs of C . The eigenvalue and eigenvector matrices can be partitioned as
D = D 1 D 2 W = W 1 W 2 ,
where D 1 is the diagonal matrix containing the r largest eigenvalues of C , and W 1 the corresponding eigenvectors. The active subspace —the subspace on which f ( x ) is most sensitive on average—is defined to be the span of the columns of W 1 . This is in turn determined by the r largest eigenvalues. Likewise, the inactive subspace—the subspace on which f ( x ) is almost invariant on average—is defined by the span of the columns of W 2 . It is thus sensible to “ignore” the inactive subspace without compromising much of the accuracy in computing f x (see [28] for a rigorous proof). If r p , significant computational gain can be achieved by focusing on only the active variables. One way to eliminate inactive variables is to integrate them out, as we now describe. Any x R p can be written as
x = W 1 W 1 T x + W 2 W 2 T x = W 1 y + W 2 z ,
where y = W 1 T x are called the active variables, and z = W 2 T x the inactive variables. The density function ρ ( x ) can be therefore considered as the joint density between the active and inactive variables:
ρ ( x ) = ρ W 1 y + W 2 z = ρ y , z .
In the AS approach, one typically approximates f ( x ) using the active variables by integrating out the inactive ones
f ( x ) g ( y ) = f W 1 y + W 2 z ρ z | y d z ,
where ρ ( z | y ) is the conditional density of z given y . Note that if ρ x = N 0 , Γ , which is the case in this paper, ρ z | y is trivial (see Section 3.1). The integral evaluation is less straightforward, in fact, computationally prohibited, for nonlinear f x , and it is typically approximated using Monte Carlo sampling.
From a statistical perspective, the uncentered covariance matrix of the form (1) is not common, though it has been investigated. Comparisons between centered and uncentered covariance matrices have been carried out for general covariance matrices in the context of principal component analysis [31,32,33,34]. We now address the difference between uncentered and centered covariance matrices in the context of active subspaces. As we will show, more can be said in this context by exploiting the structure of the inverse problems. To that end, let us introduce the centered version of the AS covariance matrix
C ˜ = f ( x ) f ¯ ( x ) f ( x ) f ¯ ( x ) T ρ ( x ) d x ,
where the gradient mean is given by
f ¯ ( x ) = f ( x ) ρ ( x ) d x .
For notational convenience, we will refer to the active subspace based on (1) as the uncentered active subspace, and the one based on (3) as the centered active subspace.
To understand the advantages/disadvantages of uncentered AS and centered AS, we restrict ourselves in the linear inversion setting. Consider the following additive noise observational model
d = A x + e ,
with A R n × p , e N 0 , Λ , and the data d R n . The inverse problem is to determine x R p given d . Posing this problem as a least squares problem in a weighted Euclidean norm, we minimize the data misfit:
min x 1 2 A x d Λ 1 2 .
To overcome an ill-conditioning issue due to the ill-posed nature of the inverse problem, we can employ the classical Tikhonov regularization approach:
min x 1 2 A x d Λ 1 2 + 1 2 x Γ 1 2 ,
where Γ 1 is a given weight matrix. In the context of Bayesian inverse problems, Γ is typically chosen as the covariance of the prior distribution of x . Following [22], we determine the active subspaces based on the data misfit, i.e.,
f ( x ) : = 1 2 A x d Λ 1 2 ,
whose gradient is
f ( x ) = A T Λ 1 A x d .
Since we choose ρ x = N 0 , Γ , it is easy to see that the uncentered covariance matrix C in (1) can be written as [22]
C = f ( x ) f ( x ) T ρ ( x ) d x = A T Λ 1 A Γ A T + d d T Λ 1 A .
Similarly, the mean of the gradient with respect to ρ x is
f ¯ ( x ) = A T Λ 1 A x d ρ ( x ) d x = A T Λ 1 d ,
and thus the centered covariance matrix C ˜ in (3) becomes
C ˜ = f ( x ) f ¯ ( x ) f ( x ) f ¯ ( x ) T ρ ( x ) d x = C f ¯ ( x ) f ¯ ( x ) T = A T Λ 1 A Γ A T Λ 1 A .
In order to gain insights into the difference between centered and uncentered active subspaces, for the rest of Section 2, we assume Γ = I . Let the full SVD of the noise covariance whitened forward operator be Λ 1 2 A = U Σ V T . Since Λ 1 2 d resides in the column space of U , there exists β R n such that
Λ 1 2 d = U β .
Lemma 1.
Let D : = diag λ 1 , , λ p and D ˜ : = diag λ ˜ 1 , , λ ˜ p be eigenvalue matrices of C and C ˜ , respectively, such that λ i λ i + 1 and λ ˜ i λ ˜ i + 1 .
  • For k min n , p λ ˜ k = σ k 4 , and λ ˜ k = 0 for n < k p .
  • Define γ = Σ T β . For 2 k p , λ ˜ k λ k λ ˜ k 1 , and λ ˜ 1 λ 1 λ ˜ 1 + γ T γ .
Proof. 
The first assertion is obvious as
C ˜ = V Σ T Σ Σ T Σ V T
using the SVD of Λ 1 2 A . For the second assertion, we have
C = A T Λ 1 A A T + d d T Λ 1 A = V Σ T U T U Σ Σ T U T + U β β T U T U Σ V T = V D ˜ + γ γ T Θ V T = V Θ V T .
Since C is a similarity transformation of Θ , we seek the relationship between the spectra of Θ and C ˜ . Now, since Θ is a rank-one perturbation of the diagonal matrix D ˜ , using a standard interlacing eigenvalue perturbation result (see, for example, ([35], Section 5) and [36]) concludes the proof. □
Lemma 1 shows that the eigenvalues of C are not smaller than those of C ˜ , but this does not reveal insights on how d could make a difference going from C ˜ to C . To wit, let us consider a special case of d = β σ k u k , where u k is the k-th column of U and β is some number, then it is straightforward to see that
C = V diag λ ˜ 1 , , λ ˜ i , , λ ˜ k 1 + β 2 , , λ ˜ p V T .
Now if β 2 is sufficiently large such that λ ˜ k 1 + β 2 > λ ˜ i where i r < k , then λ i = λ ˜ k 1 + β 2 . As a direct consequence, while v k is not part of the centered active subspace, it is for the uncentered one. In other words, it is important to see that the uncentered AS approach takes both the forward operator A and the data into account when constructing the active subspace, while the centered AS approach is purely determined by the spectrum of the forward operator. When k r , both approaches are similar. However, the uncentered AS is expected to outperform the centered counterpart when k > r since eigendirections, classified as inactive in the centered approach, are in fact active when taking the data into account; we verify this fact later in Section 5. The proof of Lemma 1 also implies that, due to the symmetry of Θ , the eigenvectors of C are not only a reordering, but also a rotation of the eigenvectors of C ˜ in general since γ γ T is not necessarily diagonal.
Remark 1.
Note that we can alternatively first perform the whitening to transform the inverse problem to the standard setting with Λ = I and Γ = I , and then compute the active subspaces. This simplifies the exposition. However, as shown in Appendix A, the active subspaces for the whitened problem change, and the corresponding DIAS solutions are less accurate.

3. Data-Informed Active Subspace Regularization for Linear Inverse Problems: Approximate Data Misfit

In our previous work [21], we proposed a data-informed (DI) regularization approach for linear inverse problems. The workhorse behind the DI approach is the generalized singular value decomposition taking advantage of the linear nature of the forward operator. In Section 2, we compare and contrast the centered and uncentered active subspace approaches. In this section, inspired by the DI idea, we construct a data-informed active subspace regularization approach for both centered and uncentered cases. We take advantage of the linear structure to approximate the data misfit. We consider the full data misfit in Section 4.

3.1. DIAS with Uncentered Active Subspace (DIAS-A)

The goal of this section is to explicitly derive the DIAS regularization for linear inverse problems using the uncentered active subspace method. Suppose that the active subspace W 1 and the inactive subspace W 2 are identified. The joint prior probability density ρ x can be written as
ρ ( x ) = ρ ( y , z ) = c x exp 1 2 x Γ 1 2 = c x exp 1 2 y T z T W 1 T Γ 1 W 1 W 1 T Γ 1 W 2 W 2 T Γ 1 W 1 W 2 T Γ 1 W 2 y z ,
where c x is a normalization constant. It follows that the conditional probability ρ ( z | y ) is
ρ ( z | y ) = c z , y exp 1 2 z z ¯ W 2 T Γ 1 W 2 2 ,
where z ¯ = W 2 T Γ 1 W 2 1 W 2 T Γ 1 W 1 y , and c z , y is a normalization constant. From (2), the AS approximation g y of the data misfit term f ( x ) in (5) is now simplified using (6)
g ( y ) = 1 2 A W 1 y + W 2 z d Λ 1 2 c z , y exp 1 2 z z ¯ W 2 T Γ 1 W 2 2 d z , = 1 2 A I W 2 W 2 T Γ 1 W 2 1 W 2 T Γ 1 W 1 y d Λ 1 2 + c 2
where
c 2 = c z , y 2 A W 2 z z ¯ Λ 1 2 exp 1 2 z z ¯ W 2 T Γ 1 W 2 2 d z ,
which is further simplified as
c 2 = c z , y 2 Trace A W 2 z z ¯ Λ 1 2 exp 1 2 z z ¯ W 2 T Γ 1 W 2 2 d z = c z , y 2 Trace A W 2 z z ¯ Λ 1 2 exp 1 2 z z ¯ W 2 T Γ 1 W 2 2 d z = c z , y 2 Trace z z ¯ T W 2 T A T Λ 1 A W 2 z z ¯ exp 1 2 z z ¯ W 2 T Γ 1 W 2 2 d z = c z , y 2 Trace W 2 T A T Λ 1 A W 2 z z ¯ z z ¯ T exp 1 2 z z ¯ W 2 T Γ 1 W 2 2 d z = c z , y 2 Trace W 2 T A T Λ 1 A W 2 z z ¯ z z ¯ T exp 1 2 z z ¯ W 2 T Γ 1 W 2 d z = c ˜ z , y 2 Trace W 2 T A T Λ 1 A W 2 W 2 T Γ 1 W 2 1 ,
where we have used the commutativity of matrix trace and integral, and the invariance of matrix trace under cyclic permutations. Here, c ˜ z , y is the product of c z , y and the constant from the integration.
Since c 2 is a constant, the inverse problem (4) can be equivalently written as
min x 1 2 M x d Λ 1 2 + 1 2 x Γ 1 2 ,
where, for brevity, we have defined
M : = A I W 2 W 2 T Γ 1 W 2 1 W 2 T Γ 1 W 1 W 1 T .
It should be pointed out that we have used the relationship y = W 1 T x to write the inverse problem in terms of the original parameter x .
Proposition 1.
The term M can be simplified further as
M = A I W 2 W 2 T Γ 1 W 2 1 W 2 T Γ 1 W 1 W 1 T = A W 1 W 1 T .
Proof. 
Let the singular value decomposition of Γ 1 2 be Γ 1 2 = R K R T . We have
Γ 1 2 W 1 = R K R T W I 0 , Γ 1 2 W 2 = R K R T W 0 I .
We can rewrite M as
M = A I W 2 W 2 T Γ 1 W 2 1 W 2 T Γ 1 W 1 W 1 T = A W 1 W 1 T A Γ 1 2 Γ 1 2 W 2 W 2 T Γ 1 2 Γ 1 2 W 2 1 W 2 T Γ 1 2 Γ 1 2 W 1 W 1 T Γ 1 2 T Γ 1 2 .
Since W 2 is a full column rank matrix, so is Γ 1 2 W 2 . Thus
Γ 1 2 W 2 T Γ 1 2 W 2 1 Γ 1 2 W 2 T = W 2 T Γ 1 2 = 0 I W T R K 1 R T .
The second equality is simply the pseudo-inverse using the singular value decomposition. Substituting (8) into the T in the Equation (9) yields
T = R K R T W 0 I 0 I W T R K 1 R T R K R T W I 0 I 0 W T R K R T = 0 .
It is well known that the Tikhonov approach (7) applies regularization on all parameters and this is the reason why a Tikhonov solution is typically smeared out. As in the DI approach [21], we argue that only data-uninformed parameters should be regularized while the data-informed ones should remain untouched to preserve the fidelity of inverse solutions.
In this paper, we treat the active parameters as the data-informed parameters and the inactive parameters as the data-uninformed parameters. Recall that ρ x = N 0 , Γ , thus z = W 2 T x N 0 , W 2 T Γ W 2 . In particular, we replace the Tikhonov regularization in the inverse problem (7) with the DIAS regularization, in which we regularize only the inactive modes,
min x 1 2 M x d Λ 1 2 + 1 2 W 2 T x W 2 T Γ W 2 1 2 ,
which yields the following inverse solution:
x DIAS - A = M T Λ 1 M + W 2 W 2 T Γ W 2 1 W 2 T 1 M T Λ 1 d .

3.2. Dias with Centered Active Subspace (cDIAS-A)

Let the eigenvalue decomposition of C ˜ be
C ˜ = V E V T ,
and similar to the uncentered AS, we partition the eigenvalue and eigenvector matrices as
E = E 1 E 2 , V = V 1 V 2
in which E 1 is the diagonal matrix containing the r largest eigenvalues and V 1 consists of the corresponding eigenvectors. For notational convenience, we still denote the active and inactive parameters as y and z , respectively. Whether they correspond to the uncentered or the centered AS should be clear from the context.
It is easy to see that any parameter x R p can be uniquely expressed as
x = f ¯ ( x ) + V 1 y + V 2 z = V 1 y y 0 + V 2 z z 0 ,
where y 0 : = V 1 T f ¯ and z 0 : = V 2 T f ¯ . Note that the relationship between the active and original parameters is given by y = V 1 T x + y 0 , and similarly z = V 2 T x + z 0 for the relationship between the inactive and original parameters.
Now following the same procedure as in Section 3.1, we can obtain the DIAS inverse solution with centered AS:
x cDIAS - A = G T Λ 1 G + V 2 V 2 T Γ V 2 1 V 2 T 1 G T Λ 1 d ,
where
G : = A I V 2 V 2 T Γ 1 V 2 1 V 2 T Γ 1 V 1 V 1 T = A V 1 V 1 T .
The second equality can be achieved in the same manner as Proposition 1.
Proposition 2.
If Γ = β 2 I , then the cDIAS-A solution x cDIAS - A reduces to the truncated SVD solution.
Proof. 
Given Γ = β 2 I , then V in (12) is the right eigenvectors of Λ 1 2 A . In this case, the solution (13) reads
x cDIAS - A = V 1 V 1 T A T Λ 1 A V 1 V 1 T + 1 β 2 I V 1 V 1 T V 1 V 1 T A T Λ 1 d = V diag 1 σ i 2 i = 1 r 0 0 0 + 0 0 0 diag β 2 r + 1 p V T V 1 V 1 T A T Λ 1 d = V diag 1 σ i 2 i = 1 r 0 0 diag β 2 r + 1 p diag σ i i = 1 r 0 0 0 U T Λ 1 2 d = V diag 1 σ i i = 1 r 0 0 0 U T Λ 1 2 d = i = 1 r u i T Λ 1 2 d σ i v i ,
where r is the dimension of the active subspace. This is exactly the solution of the truncated SVD method of (10). □

4. Data-Informed Active Subspace Regularization for Linear Inverse Problems: Non-Approximate Data Misfit (DIAS-F)

Following the standard practice in the AS method, we intrusively replaced the original data misfit with an approximate one in Section 3. Though this approach is clean and significantly reduces the number of optimization variables to manage for linear inverse problems, its intrusive nature requires a complete overhaul of the original inverse codes. This section avoids this issue by considering a DIAS formulation without changing the data misfit at the expense of working with the original high-dimensional parameter x . As this approach does not approximate the data misfit, it could provide more robust solutions as demonstrated numerically in Section 5.
The uncentered and centered DIAS formulations in this case simply read
min x 1 2 A x d Λ 1 2 + 1 2 W 2 T x W 2 T Γ W 2 1 2 , and min x 1 2 A x d Λ 1 2 + 1 2 V 2 T x V 2 T Γ V 2 1 2 ,
respectively. Thus, the inverse solution using the uncentered AS, analogous to (11), is
x DIAS - F = A T Λ 1 A + W 2 W 2 T Γ W 2 1 W 2 T 1 A T Λ 1 d ,
and the centered AS solution, similar to (13), reads
x cDIAS - F = A T Λ 1 A + V 2 V 2 T Γ V 2 1 V 2 T 1 A T Λ 1 d .
Proposition 3.
If Γ = β 2 I , then the cDIAS-F solution (15) reduces to the DI solution proposed in [21].
Proof. 
Let A ¯ = Λ 1 2 A and d ¯ = Λ 1 2 d , then the cDIAS-F solution (15) becomes
x cDIAS - F = A ¯ T A ¯ + 1 β 2 I V 1 V 1 T 1 A ¯ T d ¯ .
which is exactly the DI solution in [21]. □
Theorem 1.
Let λ i be the i-th eigenvalue of the active subspace matrix C or C ˜ , and
r ε : = max i : 1 i p and λ i ε , ε > 0 .
Set L T L = W 2 W 2 T Γ W 2 1 W 2 T and L T L = V 2 V 2 T Γ V 2 1 V 2 T for the uncentered active subspace and centered active subspace, respectively. Suppose that the nullspace of A is trivial, i.e., N A = 0 . Consider the inverse problem
min x 1 2 Λ 1 2 B x Λ 1 2 d 2 2 + 1 2 x 2 2 ,
using the DIAS regularization approaches with rank- r ε active subspace, where B = M , or G , or A , which correspond to the DIAS-A, cDIAS-A, or (c)DIAS-F approaches, respectively. Define the “reconstruction operator” [7,37], approximating the map from observations d to parameter x , as
R ε : = H r ε 1 B T Λ 1 2 ,
where
H r ε = B T Λ 1 B + L T L .
The following hold:
(i) 
The inverse problem with rank- r ε DIAS approaches, i.e., the optimization problem (16), is well-posed in the Hadamard sense.
(ii) 
DIAS techniques are regularization strategies [7,37] in the following sense:
lim ε 0 R ε Λ 1 2 A x = x .
(iii) 
The rank- r ε DIAS technique is an admissible regularization method.
Proof. 
The fact that H r ε is invertible for all approaches is clear.
(i)
It is sufficient to show that the solution is continuous with respect to the observational data. Indeed, we have
x DIAS 2 H r ε 1 2 A 2 Λ 1 2 d 2 1 ϕ A 2 Λ 1 2 d 2
where ϕ is the smallest singular value of H r ε . Note that, owing to the invertibility of H r ε , ϕ is bounded away from 0 and approaches the smallest singular value of A T Λ 1 A , which is larger than 0 due to the triviality of the nullspace of A .
(ii)
The matrix inversion is a continuous mapping (which is the composition of continuous maps, i.e., the matrix determinant and the matrix adjoint), thus
lim ε 0 R ε Λ 1 2 A x = lim ε 0 H r ε 1 A T Λ 1 A x = lim ε 0 H r ε 1 A T Λ 1 A x = A T Λ 1 A 1 A T Λ 1 Λ 1 2 A x = x ,
where we have used the following facts
  • The rank of active subspace r ε p as ε 0 , and thus L T L = 0 .
  • As r ε p , we have B = A .
(iii)
It is sufficient to show that
sup y R ε Λ 1 2 y x 2 : Λ 1 2 A x y 2 ε 0 as ε 0 ,
for any x . We have
R ε Λ 1 2 y x 2 R ε Λ 1 2 A x x 2 + R ε Λ 1 2 A x y 2 H r ε 1 B T Λ 1 A I 2 x 2 + R ε 2 ε H r ε 1 2 B T Λ 1 A B T Λ 1 B 2 + L T L 2 x 2 + ε ϕ B 2 Λ 1 2 2 1 ϕ σ ( ε ) + γ ( ε ) x 2 + ε ϕ B 2 Λ 1 2 2 ,
where σ ( ε ) is the maximum singular value of B T Λ 1 A B T Λ 1 B , γ ( ε ) is the maximum eigenvalue of L T L . We note that, as ε 0 , σ ( ε ) 0 and γ ( ε ) 0 .

5. Numerical Results

We now test the proposed DIAS approaches against the Tikhonov method on a variety of popular linear inverse problems. In particular, we consider one-dimensional (1D) deconvolution, various benchmark problems from [38], and X-ray tomography. In all linear inverse problems, we assume Λ = δ 2 I and Γ = β 2 I . Under this assumption, we have α = β 2 δ 2 . Furthermore, by Propositions 2 and 3, the cDIAS-A and cDIAS-F methods become the truncated SVD (TSVD) and the DI approaches, respectively. Thus, for clarity we use TSVD and DI in the places of cDIAS-A and cDIAS-F for all examples. Additionally, by a % additive white noise we mean additive Gaussian noise with the standard deviation of a % of the maximum value of synthetic data. In the 1D deconvolution problem, we numerically investigate the difference in the inverse solutions, using the uncentered AS approaches and the centered AS approaches. In particular, we highlight the reordering and rotation effect of the uncentered AS induced by the data d (see Section 2 for a theoretical discussion).

5.1. 1D Deconvolution Problem

We consider the one-dimensional deconvolution problem (see, for example, [9,39])
d ( s j ) = 0 1 a ( s j , t ) x ( t ) d t + e ( s j ) , j = 1 , , n ,
where a ( s , t ) = 1 2 π μ 2 exp 1 2 μ 2 t s 2 is a Gaussian convolution kernel with μ = 0.05 , s j = j n ( 0 j < n ), n = 1000 , and e is 5% additive white noise. The exact solution, which we aim to reconstruct from the noisy data d , is given as
x ( t ) = 10 t 0.5 exp 50 t 0.5 2 0.8 + 1.6 t , t 0 , 1 .
Figure 1a shows the projected components of the true solution x t in the uncentered active eigenspace W 1 and the centered active eigenspace V 1 . In the uncentered eigenspace, the true solution lies almost entirely in the first eigenmode, while it predominantly lies in the second and sixth modes of the centered eigenspace. The relative error between projection x r with r dimensional AS and the true solution x t is shown in Figure 1b. It can be seen that the uncentered AS provides more accurate projection, even with one active direction ( r = 1 ). The result also shows that the centered AS needs at least 10 active directions to start being comparable to the uncentered counterpart in terms of accuracy. These numerical observations verify the reordering and rotating results of Lemma 1.
The uncentered AS eigenvector reordering and rotating effects induced by the data d are presented in Table 2, where we compute the cosine of the angle between the centered and uncentered active modes. As can be seen, w 1 is slightly rotated from v 2 . That is, compared to the centered AS method, the uncentered AS, under the influence of the data, reorders the eigenvectors of the forward operator so that v 2 is in the first position and slightly rotates it to align better with the directions in which the data are informative. The most significant shortcoming of the DI method and others using the basis V is that they misclassify v 1 as the most data-informed direction, while for the data under consideration, v 2 is much more informative. Indeed, the relative error in Figure 1b shows that the true parameter is almost orthogonal to v 1 .
Figure 2 presents solutions for the 1D deconvolution problem with various combinations of active subspace dimension r and regularization parameter α . For the one-dimensional active subspace r = 1 , the TSVD method is not able to reconstruct the solution since the first mode of V 1 contributes very little to the true solution. On the contrary, DIAS-A yields a reasonable solution, as its first mode is the most data-informative one. For large dimensional active subspace, the DIAS-A and TSVD solutions are almost the same. These observations are consistent with the discussion following Lemma 1. Recall that the uncentered eigenvectors are a reordering and rotation of the centered ones. By considering a sufficiently large number of modes to be data informed, the subspaces spanned by V 1 and W 1 become more similar. This can also be clearly seen in Table 2 along the diagonal. Notice that v i w i except for the first two modes. This is because the first two eigenvectors are swapped and slightly rotated. Furthermore, inherited from the DI approach [21], DIAS-A, and TSVD (in fact cDIAS-A) solutions are almost invariant with respect to the regularization parameter. This is not surprising since the DIAS approach, by design, regularizes only the data-uninformed modes. Its solution should remain unaffected if sufficient data-informed modes are taken into account.
Figure 2 also shows that the DIAS-F, DI, and Tikhonov solutions are indistinguishable for small regularization parameter regardless of whether r = 1 or r = 10 . This is due to two reasons: (1) the DIAS approaches and Tikhonov regularize inactive modes in the same way, and (2) Tikhonov regularization has little effect on the active modes when the regularization parameter is small, thus having little impact on the inverse solution. The situation changes for larger regularization parameters, especially for r = 10 . Both DIAS-F and DI, by leaving the active parameters untouched, are insensitive to the regularization while Tikhonov oversmooths the solution by heavily regularizing the active parameters.
Another important observation that Figure 2 conveys is the difference among DIAS approaches with centered and uncentered active subspaces. For r = 1 , the uncentered approaches DIAS-A and DIAS-F outperform the centered counterparts TSVD and DI, especially for large regularization parameters. In other words, uncentered approaches are more robust to the regularization parameter. The reason is that uncentered methods do not penalize the data-informed direction w 1 while the centered ones—without taking data into account—regularize v 2 , the most data-informative direction in the basis V . As discussed in Section 2, for sufficiently large active subspace dimension r = 10 , all the DIAS solutions are similar since all methods end up spanning the same subspace. However, at the optimal regularization parameter α = 1 , determined by the L-curve method [38], the DIAS-F solution is visibly closer to the ground truth than TSVD and Tikhonov.

5.2. Benchmark Problems

In this section, we apply the DIAS regularization approach with centered AS and uncentered AS methods to six benchmark problems from regularization tools [38]. We briefly describe the Shaw benchmark [40] here, and we refer the readers to [38] for the descriptions of other benchmark problems. It is in fact a deconvolution problem in which the kernel is given as
a ( s , t ) = cos s + cos t sin u u 2 ,
where u = π sin ( s ) + sin ( t ) , and the true solution is given by
x ( t ) = a 1 exp c 1 t t 1 2 + a 2 exp c 2 t t 2 2 ,
where a 1 = 2 , c 1 = 6 , t 1 = 0.8 , a 2 = 1 , c 2 = 2 , t 2 = 0.5 . For all benchmark problems, the domain 0 , 1 is divided into 1000 subdomains. These problems provide good examples to study the robustness and accuracy properties of the DIAS approaches. Observational data for each problem are corrupted with 1% additive white noise. Figure 3 measures the relative error between the ground truth x t and its projection x r on both centered and uncentered active subspaces with various values of the active subspace dimension r. For Shaw, heat, gravity, and Phillips problems, projecting the exact solution on the uncentered AS results in a lower error than on the centered AS. For the Deriv2 and Foxgood problems, the results are identical for both centered and uncentered AS. The reason is that the data for these two problems do not provide new information and thus the active subspaces are entirely determined by the forward operator A .
Figure 4, Figure 5, Figure 6, Figure 7, Figure 8 and Figure 9 show the following: (1) top row—relative errors in the inverse solutions for different regularization methods for two values of regularization parameter ( α o p t i m a l and 500 × α o p t i m a l ); (2) middle row—inverse solutions with Tikhonov regularization and DIAS regularizations with one dimensional active subspace ( r = 1 ) for two values of regularization parameter; and (3) bottom row—inverse solutions with Tikhonov regularization and DIAS regularizations with optimal active subspace dimension for two values of the regularization parameter. Here, the optimal regularization parameters are chosen based on the L-curve method with Tikhonov regularization [38], and the optimal AS dimension is found experimentally for each method. It turns out that the optimal AS dimension is the same for all methods.
Around the optimal AS dimension (top rows of all figures), regardless of what the regularization parameter is, all DIAS regularizations have similar accuracy and they are more accurate than Tikhonov, as expected. As can be seen in the middle rows, when the active subspace dimension is from r = 1 to r = 6 , less than the optimal dimension r optimal = 7 for the Shaw problem, the full data misfit methods outperform the approximated misfit counterparts. We provide the reason only for the Shaw problem, as it is the same for others. When taking r = 1 , the approximate misfit approaches completely remove six other active modes, which are misclassified as inactive, in addition to removing the truly inactive modes. The inverse solutions thus lack the important contribution from these modes, leading to inaccurate reconstructions. Even when the active subspace is chosen to be too small, the full misfit methods only regularize the misclassified modes rather than eliminating them from the solution entirely.
The results also show that the DIAS regularizations with full data misfit are at least as good as Tikhonov regularization with α o p t i m u m and are more accurate with 500 × α o p t i m u m . Again, the reason is that for a reasonably small regularization parameter, the smoothing effect from Tikhonov regularization does not change the solution significantly. On the bottom rows where the AS dimensions are optimal, DIAS solutions are similar and outperform the Tikhonov counterparts for both values of regularization parameters.
From middle to bottom rows, the DIAS solutions with approximate data misfit change from worse to comparable to the DIAS solutions with original data misfit, and thus from worse to more accurate than Tikhonov solutions. This is not surprising: Equation (2) shows that as the active subspace dimension increases, the error due to the data misfit approximation decreases. On the other hand, for the same reasons as in the deconvolution problem, the uncentered AS methods (DIAS-A and DIAS-F) are more accurate than the corresponding centered AS counterparts (TSVD and DI) for all problems with any active subspace dimension and with any regularization parameter.
Note that for these benchmarks, we take 500 × α optimal as a large regularization parameter case to show that DIAS regularization is robust with respect to the choice of regularization parameter while Tikhonov is not. The DIAS solutions are much more accurate than the Tikhonov solutions, as the latter overregularizes all modes in this case. Additionally, we observe that the full data misfit approaches become similar to the approximate misfit ones as the regularization parameter increases. To see why, we recall from Equation (7) that the only difference between approximate and full misfit approaches is the removal of the inactive variables in the former. When the regularization parameter approaches infinity, it annihilates the contribution of the inactive subspace in the inverse solution. In fact, they behave like TSVD in the limit and we know that TSVD is equivalent to applying infinite regularization on the truncated modes [21]. Thus both approaches would yield identical solutions in the limit. For the case we consider here, the regularization is sufficiently large for us to already see this asymptotic behavior.

5.3. X-ray Tomography

X-ray tomographic reconstruction is another well-known linear inverse problem. A more detailed description of this problem can be found in [9]. Synthetic data are generated with 1%, 3% and 5% white noise added to be realistic. Tikhonov solutions using the optimal regularization parameter α = 1.3 × 10 5 , obtained by the L-curve method [9], are compared with DIAS approaches.
Figure 10 depicts the eigenvalues of C and C ˜ . While the first eigenvalue of the uncentered active subspace matrix is significantly larger than the first eigenvalue of the centered active subspace matrix, this only hints that w 1 is more important than any of the vectors in V . Figure 11 shows that there is indeed a striking difference between the first eigenvector of the two active subspaces. Visual inspection of the eigenvectors w 1 and v 1 makes it obvious that w 1 contributes significantly to the solution, while the contribution of v 1 to the solution is much less pronounced.
The Tikhonov solutions with different regularization parameters and the original image are shown in the Figure 12 for the case of 1 % noise. Clearly, underregularization results in noisy solutions and overregularization yields overly smoothed and blurry solutions.
Table 3 and Figure 13 show the relative error and the reconstructed images for the underregularization case. The approximated misfit methods (DIAS-A and TSVD) remove the inactive variables so they are not prone to noise pollution amplified by inverting the small singular values corresponding to these data-uninformed modes. This has a regularization effect known as regularization by truncation. Since the full misfit approaches and Tikhonov capture all modes, their solutions are much more vulnerable to noise in the underregularization regime. For optimal regularization case in Table 4 and Figure 14, the full misfit methods DIAS-F and DI are more accurate than their approximate misfit counterparts for small active subspace dimensions, while all methods give similar results when the active subspace dimension is sufficiently large. For the overregularization case in Table 5 and Figure 15, the approximate and full misfit approaches provide comparable solutions, as they both behave like the TSVD method. Moreover, DIAS regularization solutions are robust to the regularization parameter as opposed to the Tikhonov solution. Another observation is that DIAS methods perform poorly when the active subspace dimension is taken to be too large ( r = 12,000). This is not surprising since as r p , all regularization is removed and the problem becomes ill-posed again.
To better understand the robustness of the DIAS approaches to various noise levels, we perform a sensitivity analysis, increasing the noise levels to 3% and 5%. One might expect that the DIAS method would be especially sensitive to noise since it is a data-driven approach. The following discussion shows, however, that the DIAS approaches maintain their advantages, even in the presence of significant noise. The results with higher noise levels are presented in Table 6, Table 7 and Table 8 for 3%, and Table 9, Table 10 and Table 11 for 5%. It can be observed that approximately the first 500 modes of the active subspace are not perturbed by noise (1%, 3% and 5%). For example, at r = 500 , the relative error of DIAS-A and TSVD are approximately 49% and 56%, respectively, for all noise levels. Meanwhile, when r = 5000 , the corresponding figures for both approaches in the cases of 1%, 3% and 5% are (27.39%, 29.10%), (33.95%, 35.10%) and (43.49 % and 44.39%), respectively. This is consistent with the well-known fact that higher modes of the active subspace (corresponding to smaller eigenvalues) contain more noise. In other words, as the noise level increases, the 5000-th mode moves to the inactive side of the space because it becomes dominated by noise rather than useful information.

6. Conclusions

We have presented four different data-informed regularization approaches for linear inverse problems using active subspaces (DIAS): two with approximate data misfit, and two with full data misfit. We rigorously showed the connection between the centered and uncentered active subspaces and the consequences on the performance of the corresponding DIAS approaches. For linear inverse problems, we showed that the TSVD and the DI regularization methods are members of DIAS regularization. Regularizing only the inactive directions is fundamental to the success of the DIAS approaches. All four DIAS regularization methods are robust to a wide range of regularization parameter values and outperform Tikhonov regularization for many choices of regularization parameter. The uncentered DIAS approaches are more robust and more accurate than their centered counterparts (the TSVD and DI approaches). Among DIAS regularizations methods, DIAS-F (uncentered active subspace with the original data misfit) has the best compromise: it is a data-informed non-intrusive approach.
By data-informed, we mean that the method balances the information encoded in the forward operator and information gained from the particular realization of the data. By non-intrusive, we mean that the method provides the ability to reuse existing inverse codes with minor modification only on the regularization. Various numerical results have demonstrated that the DIAS-F approach is the most effective method presented. In particular, excellent results can be obtained from DIAS-F with only a one-dimensional data informed subspace.
For problems with significant noise in the data, DIAS regularization methods could result in noisy reconstructions unless a small number of active directions are taken, since the less data-informed directions, reflecting the noise in the data, still need some regularization to smooth out the noise. Ongoing work is to equip the DIAS regularization approach with a mechanism to combat high noise scenarios while ensuring the fidelity of inverse solutions. Another appealing feature of the DIAS framework is that it is applicable not only for linear, but also seamlessly for nonlinear inverse problems. Detailed numerical investigations and rigorous analysis of the DIAS approach for nonlinear inverse problems are part of ongoing work.

Author Contributions

Conceptualization, H.N., J.W. and T.B.-T.; Methodology, H.N., J.W. and T.B.-T.; Project administration, H.N., J.W. and T.B.-T.; Supervision, T.B.-T.; Validation, H.N., J.W. and T.B.-T.; Writing—original draft, H.N., J.W. and T.B.-T.; Writing—review and editing, H.N., J.W. and T.B.-T. All authors have read and agreed to the published version of the manuscript.

Funding

This work was partially funded by the National Science Foundation awards NSF-2108320, NSF-1808576 and NSF-CAREER-1845799; by the Department of Energy award DE-SC0018147 and DE-26083989; and by 2020 ConTex award; and by 2021 UT-Portugal CoLab award.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Data are contained within the article.

Acknowledgments

The authors would like to thank the anonymous reviewers for their critical and thoughtful comments that helped improve the paper substantially.

Conflicts of Interest

The authors declare no conflict of interest.

Appendix A. Dias Regularization with Whitened Active Subspaces

One can alternatively whiten the inverse problem to reduce to the standard form. Let us define the following whitening transformation: x = Γ 1 2 x ˜ , x ˜ N 0 , I , A ˜ = Λ 1 2 A Γ 1 2 , and d ˜ = Λ 1 2 d .
The whitened inverse problem reads
min x ˜ 1 2 A ˜ x ˜ d ˜ 2 + 1 2 x ˜ 2 ,
The corresponding misfit function and its gradient are given as
f ( x ˜ ) = 1 2 A ˜ x ˜ d ˜ 2 , x ˜ f ( x ˜ ) = A ˜ T A ˜ x ˜ d ˜ .
By simple algebraic manipulations, we obtain
C x ˜ = x ˜ f ( x ˜ ) x ˜ f ( x ˜ ) T ρ ( x ˜ ) d x ˜ = A ˜ T A ˜ A ˜ T + d ˜ d ˜ T A ˜ = Γ 1 2 C Γ 1 2 , C ˜ x ˜ = Γ 1 2 C ˜ Γ 1 2 .
We only consider the DIAS regularization with full misfit (similar derivation can be done for approximate misfit). The uncentered and centered DIAS solutions in this case are
x ˜ DIAS - F = A ˜ T A ˜ + α I W ˜ 1 W ˜ 1 T 1 A ˜ T d ˜ , x ˜ cDIAS - F = A ˜ T A ˜ + α I V ˜ 1 V ˜ 1 T 1 A ˜ T d ˜ .
The 1D deconvolution DIAS inverse solutions with non-identity regularization, inspired by the relax boundary prior covariance in [41], show that the whitened active subspace is less accurate (Figure A1).
Figure A1. The relative error of original active subspace and whitened active subspace in the overregularization case.
Figure A1. The relative error of original active subspace and whitened active subspace in the overregularization case.
Computation 10 00038 g0a1

References

  1. Natterer, F. The Mathematics of Computerized Tomography; SIAM: Philadelphia, PA, USA, 2001. [Google Scholar]
  2. Natterer, F.; Wübbeling, F. Mathematical Methods in Image Reconstruction; SIAM: Philadelphia, PA, USA, 2001. [Google Scholar]
  3. Kravaris, C.; Seinfeld, J.H. Identification of parameters in distributed parameter systems by regularization. SIAM J. Control. Optim. 1985, 23, 217–241. [Google Scholar] [CrossRef] [Green Version]
  4. Banks, H.T.; Kunisch, K. Estimation Techniques for Distributed Parameter Systems; Springer: Berlin/Heidelberg, Germany, 2012. [Google Scholar]
  5. Bui-Thanh, T.; Ghattas, O.; Martin, J.; Stadler, G. A computational framework for infinite-dimensional Bayesian inverse problems Part I: The linearized case, with application to global seismic inversion. SIAM J. Sci. Comput. 2013, 35, A2494–A2523. [Google Scholar] [CrossRef] [Green Version]
  6. Bui-Thanh, T.; Burstedde, C.; Ghattas, O.; Martin, J.; Stadler, G.; Wilcox, L.C. Extreme-scale UQ for Bayesian inverse problems governed by PDEs. In Proceedings of the SC’12: International Conference on High Performance Computing, Networking, Storage and Analysis, Salt Lake City, UT, USA, 10–16 November 2012; pp. 1–11. [Google Scholar]
  7. Colton, D.; Kress, R. Inverse Acoustic and Electromagnetic Scattering, 2nd ed.; Applied Mathematical Sciences; Springer: Berlin/Heidelberg, Germany; New York, NY, USA; Tokyo, Japan, 1998; Volume 93. [Google Scholar]
  8. Hansen, P.C. Discrete Inverse Problems: Insight and Algorithms; SIAM: Philadelphia, PA, USA, 2010. [Google Scholar]
  9. Mueller, J.L.; Siltanen, S. Linear and Nonlinear Inverse Problems with Practical Applications; Society for Industrial and Applied Mathematics: Philadelphia, PA, USA, 2012. [Google Scholar]
  10. Rudin, L.; Osher, S.; Fatemi, E. Nonlinear total variation based noise removal algorithms. Physica D 1992, 60, 259–268. [Google Scholar] [CrossRef]
  11. Beck, A.; Teboulle, M. Fast Gradient-Based Algorithms for Constrained Total Variation Image Denoising and Deblurring Problems. IEEE Trans. Image Process. 2009, 18, 2419–2434. [Google Scholar] [CrossRef] [Green Version]
  12. Nikolova, M. Weakly Constrained Minimization: Application to the Estimation of Images and Signals Involving Constant Regions. J. Math. Imaging Vis. 2004, 21, 155–175. [Google Scholar] [CrossRef]
  13. Goldstein, T.; Osher, S. The slit Bregman method for L1-regularized problems. SIAM J. Imaging Sci. 2009, 2, 323–343. [Google Scholar] [CrossRef]
  14. Ramirez-Giraldo, J.; Trzasko, J.; Leng, S.; Yu, L.; Manduca, A.; McCollough, C.H. Nonconvex prior image constrained compressed sensing (NCPICCS): Theory and simulations on perfusion CT. Med. Phys. 2011, 38, 2157–2167. [Google Scholar] [CrossRef] [Green Version]
  15. Babacan, S.D.; Mancera, L.; Molina, R.; Katsaggelos, A.K. Non-convex priors in Bayesian compressed sensing. In Proceedings of the 2009 17th European Signal Processing Conference, Glasgow, UK, 24–28 August 2009; pp. 110–114. [Google Scholar]
  16. Nikolova, M. Analysis of the recovery of edges in images and signals by minimizing nonconvex regularized least-squares. Multiscale Model. Simul. 2005, 4, 960–991. [Google Scholar] [CrossRef]
  17. Chartrand, R.; Wohlberg, B. A Nonconvex ADMM Algorithm for Group Sparsity with Sparse Groups. In Proceedings of the 2013 IEEE International Conference on Acoustics, Speech and Signal Processing, Vancouver, BC, Canada, 26–31 May 2013; pp. 6009–6013. [Google Scholar]
  18. Boley, D. Local linear convergence of the alternating direction method of multipliers on quadratic or linear programs. SIAM J. Optim. 2013, 23, 2183–2207. [Google Scholar] [CrossRef] [Green Version]
  19. Boyd, S.; Parikh, N.; Chu, E.; Peleato, B.; Eckstein, J. Distributed Optimization and Statistical Learning via the Alternating Direction Method of Multipliers. Found. Trends Mach. Learn. 2010, 3, 1–122. [Google Scholar] [CrossRef]
  20. Chartrand, R.; Yin, W. Iteratively reweighted algorithms for compressive sensing. In Proceedings of the 2008 IEEE International Conference on Acoustics, Speech and Signal Processing, Las Vegas, NV, USA, 31 March–4 April 2008; pp. 3869–3872. [Google Scholar]
  21. Wittmer, J.; Bui-Thanh, T. Data-Informed Regularization for Inverse and Imaging Problems. In Handbook of Mathematical Models and Algorithms in Computer Vision and Imaging; Springer International Publishing: Cham, Switzerland, 2021. [Google Scholar]
  22. Constantine, P.G.; Dow, E.; Wang, Q. Active subspace methods in theory and practice: Applications to kriging surfaces. SIAM J. Sci. Comput. 2014, 36, A1500–A1524. [Google Scholar] [CrossRef]
  23. Constantine, P.G.; Emory, M.; Larsson, J.; Iaccarino, G. Exploiting active subspaces to quantify uncertainty in the numerical simulation of the HyShot II scramjet. J. Comput. Phys. 2015, 302, 1–20. [Google Scholar] [CrossRef] [Green Version]
  24. Diaz, P.; Constantine, P.; Kalmbach, K.; Jones, E.; Pankavich, S. A modified SEIR model for the spread of Ebola in Western Africa and metrics for resource allocation. Appl. Math. Comput. 2018, 324, 141–155. [Google Scholar] [CrossRef] [Green Version]
  25. Constantine, P.G.; Zaharatos, B.; Campanelli, M. Discovering an active subspace in a single-diode solar cell model. Stat. Anal. Data Mining Asa Data Sci. J. 2015, 8, 264–273. [Google Scholar] [CrossRef] [Green Version]
  26. Cui, C.; Zhang, K.; Daulbaev, T.; Gusak, J.; Oseledets, I.; Zhang, Z. Active subspace of neural networks: Structural analysis and universal attacks. SIAM J. Math. Data Sci. 2020, 2, 1096–1122. [Google Scholar] [CrossRef]
  27. Lam, R.R.; Zahm, O.; Marzouk, Y.M.; Willcox, K.E. Multifidelity dimension reduction via active subspaces. SIAM J. Sci. Comput. 2020, 42, A929–A956. [Google Scholar] [CrossRef] [Green Version]
  28. Constantine, P.G.; Kent, C.; Bui-Thanh, T. Accelerating Markov chain Monte Carlo with active subspaces. SIAM J. Sci. Comput. 2016, 38, A2779–A2805. [Google Scholar] [CrossRef] [Green Version]
  29. Demo, N.; Tezzele, M.; Rozza, G. A non-intrusive approach for the reconstruction of POD modal coefficients through active subspaces. Comptes Rendus Mécanique 2019, 347, 873–881. [Google Scholar] [CrossRef]
  30. O’Leary-Roseberry, T.; Villa, U.; Chen, P.; Ghattas, O. Derivative-informed projected neural networks for high-dimensional parametric maps governed by PDEs. Comput. Methods Appl. Mech. Eng. 2022, 388, 114199. [Google Scholar] [CrossRef]
  31. Cadima, J.; Jolliffe, I.T. On Relationships between Uncentred and Column-Centred Principal Component Analysis. Pak. J. Stat. 2009, 25, 473–503. [Google Scholar]
  32. Honeine, P. An eigenanalysis of data centering in machine learning. arXiv 2014, arXiv:1407.2904. [Google Scholar]
  33. Jolliffe, I.T.; Cadima, J. Principal component analysis: A review and recent developments. Philos. Trans. R. Soc. Math. Phys. Eng. Sci. 2016, 374, 20150202. [Google Scholar] [CrossRef] [PubMed]
  34. Alexandris, N.; Gupta, S.; Koutsias, N. Remote sensing of burned areas via PCA, Part 1; centering, scaling and EVD vs. SVD. Open Geospat. Data Softw. Stand. 2017, 2, 17. [Google Scholar] [CrossRef] [Green Version]
  35. Golub, G.H. Some Modified Matrix Eigenvalue Problems. SIAM Rev. 1973, 15, 318–334. [Google Scholar] [CrossRef] [Green Version]
  36. Wilkinson, J. The Algebraic Eigenvalue Problem; Claredon Press: Oxford, UK, 1965. [Google Scholar]
  37. Kirsch, A. An Introduction to the Mathematical Theory of Inverse Problems, 2nd ed.; Applied Mathematical Sciences; Springer: New York, NY, USA, 2011; Volume 120. [Google Scholar]
  38. Hansen, P.C. Regularization tools version 4.0 for Matlab 7.3. Numer. Algorithms 2007, 46, 189–194. [Google Scholar] [CrossRef]
  39. Calvetti, D.; Somersalo, E. An Introduction to Bayesian Scientific Computing: Ten Lectures on Subjective Computing; Springer: Berlin/Heidelberg, Germany, 2007; Volume 2. [Google Scholar]
  40. Shaw, C., Jr. Improvement of the resolution of an instrument by numerical solution of an integral equation. J. Math. Anal. Appl. 1972, 37, 83–112. [Google Scholar] [CrossRef] [Green Version]
  41. Calvetti, D. Preconditioned iterative methods for linear discrete ill-posed problems from a Bayesian inversion perspective. J. Comp. Appl. Math. 2007, 2, 378–395. [Google Scholar] [CrossRef] [Green Version]
Figure 1. 1D deconvolution problem. (a): uncentered AS subspace versus centered AS subspace in representing the true solution. (b): relative error of DIAS-A and TSVD solutions as a function of active subspace dimension r.
Figure 1. 1D deconvolution problem. (a): uncentered AS subspace versus centered AS subspace in representing the true solution. (b): relative error of DIAS-A and TSVD solutions as a function of active subspace dimension r.
Computation 10 00038 g001
Figure 2. Solutions for 1D deconvolution problem with different active subspace dimensions r = 1 , 10 and different regularization parameter values α = 10 4 , 1 , 10 4 .
Figure 2. Solutions for 1D deconvolution problem with different active subspace dimensions r = 1 , 10 and different regularization parameter values α = 10 4 , 1 , 10 4 .
Computation 10 00038 g002
Figure 3. The comparison between centered AS and uncentered AS in recovering the true solution x t . x r is the projection of x t in the active subspaces.
Figure 3. The comparison between centered AS and uncentered AS in recovering the true solution x t . x r is the projection of x t in the active subspaces.
Computation 10 00038 g003
Figure 4. Shaw problem, α optimal = 1.35 × 10 2 . Comparison of different DIAS regularizations and Tikhonov regularization. Top row: relative errors in the inverse solutions for different regularization methods for two values of regularization parameter. Middle row: inverse solutions with DIAS regularizations with one-dimensional active subspace and Tikhonov regularization for two values of regularization parameter. Bottom row: inverse solutions with DIAS regularizations with the optimal dimension active subspace ( r optimal = 7 ) and Tikhonov regularization for two values of regularization parameter.
Figure 4. Shaw problem, α optimal = 1.35 × 10 2 . Comparison of different DIAS regularizations and Tikhonov regularization. Top row: relative errors in the inverse solutions for different regularization methods for two values of regularization parameter. Middle row: inverse solutions with DIAS regularizations with one-dimensional active subspace and Tikhonov regularization for two values of regularization parameter. Bottom row: inverse solutions with DIAS regularizations with the optimal dimension active subspace ( r optimal = 7 ) and Tikhonov regularization for two values of regularization parameter.
Computation 10 00038 g004
Figure 5. Gravity, α optimal = 0.138 . Comparison of different DIAS regularizations and Tikhonov regularization. Top row: relative errors in the inverse solutions for different regularization methods for two values of regularization parameter. Middle row: inverse solutions with DIAS regularizations with one-dimensional active subspace and Tikhonov regularization for two values of regularization parameter. Bottom row: inverse solutions with DIAS regularizations with the optimal dimension active subspace ( r optimal = 8 ) and Tikhonov regularization for two values of regularization parameter.
Figure 5. Gravity, α optimal = 0.138 . Comparison of different DIAS regularizations and Tikhonov regularization. Top row: relative errors in the inverse solutions for different regularization methods for two values of regularization parameter. Middle row: inverse solutions with DIAS regularizations with one-dimensional active subspace and Tikhonov regularization for two values of regularization parameter. Bottom row: inverse solutions with DIAS regularizations with the optimal dimension active subspace ( r optimal = 8 ) and Tikhonov regularization for two values of regularization parameter.
Computation 10 00038 g005
Figure 6. Deriv2, α optimal = 6.6 × 10 4 . Comparison of different DIAS regularizations and Tikhonov regularization. Top row: relative errors in the inverse solutions for different regularization methods for two values of regularization parameter. Middle row: inverse solutions with DIAS regularizations with one-dimensional active subspace and Tikhonov regularization for two values of regularization parameter. Bottom row: inverse solutions with DIAS regularizations with the optimal dimension active subspace ( r optimal = 15 ) and Tikhonov regularization for two values of regularization parameter.
Figure 6. Deriv2, α optimal = 6.6 × 10 4 . Comparison of different DIAS regularizations and Tikhonov regularization. Top row: relative errors in the inverse solutions for different regularization methods for two values of regularization parameter. Middle row: inverse solutions with DIAS regularizations with one-dimensional active subspace and Tikhonov regularization for two values of regularization parameter. Bottom row: inverse solutions with DIAS regularizations with the optimal dimension active subspace ( r optimal = 15 ) and Tikhonov regularization for two values of regularization parameter.
Computation 10 00038 g006
Figure 7. Heat, α optimal = 3.35 × 10 3 . Comparison of different DIAS regularizations and Tikhonov regularization. Top row: relative errors in the inverse solutions for different regularization methods for two values of regularization parameter. Middle row: inverse solutions with DIAS regularizations with one-dimensional active subspace and Tikhonov regularization for two values of regularization parameter. Bottom row: inverse solutions with DIAS regularizations with the optimal dimension active subspace ( r optimal = 22 ) and Tikhonov regularization for two values of regularization parameter.
Figure 7. Heat, α optimal = 3.35 × 10 3 . Comparison of different DIAS regularizations and Tikhonov regularization. Top row: relative errors in the inverse solutions for different regularization methods for two values of regularization parameter. Middle row: inverse solutions with DIAS regularizations with one-dimensional active subspace and Tikhonov regularization for two values of regularization parameter. Bottom row: inverse solutions with DIAS regularizations with the optimal dimension active subspace ( r optimal = 22 ) and Tikhonov regularization for two values of regularization parameter.
Computation 10 00038 g007
Figure 8. Phillips, α optimal = 0.175 . Comparison of different DIAS regularizations and Tikhonov regularization. Top row: relative errors in the inverse solutions for different regularization methods for two values of regularization parameter. Middle row: inverse solutions with DIAS regularizations with one-dimensional active subspace and Tikhonov regularization for two values of regularization parameter. Bottom row: inverse solutions with DIAS regularizations with the optimal dimension active subspace ( r optimal = 8 ) and Tikhonov regularization for two values of regularization parameter.
Figure 8. Phillips, α optimal = 0.175 . Comparison of different DIAS regularizations and Tikhonov regularization. Top row: relative errors in the inverse solutions for different regularization methods for two values of regularization parameter. Middle row: inverse solutions with DIAS regularizations with one-dimensional active subspace and Tikhonov regularization for two values of regularization parameter. Bottom row: inverse solutions with DIAS regularizations with the optimal dimension active subspace ( r optimal = 8 ) and Tikhonov regularization for two values of regularization parameter.
Computation 10 00038 g008
Figure 9. Foxgood, α optimal = 1.7 × 10 2 . Comparison of different DIAS regularizations and Tikhonov regularization. Top row: relative errors in the inverse solutions for different regularization methods for two values of regularization parameter. Middle row: inverse solutions with DIAS regularizations with one-dimensional active subspace and Tikhonov regularization for two values of regularization parameter. Bottom row: inverse solutions with DIAS regularizations with the optimal dimension active subspace ( r optimal = 3 ) and Tikhonov regularization for two values of regularization parameter.
Figure 9. Foxgood, α optimal = 1.7 × 10 2 . Comparison of different DIAS regularizations and Tikhonov regularization. Top row: relative errors in the inverse solutions for different regularization methods for two values of regularization parameter. Middle row: inverse solutions with DIAS regularizations with one-dimensional active subspace and Tikhonov regularization for two values of regularization parameter. Bottom row: inverse solutions with DIAS regularizations with the optimal dimension active subspace ( r optimal = 3 ) and Tikhonov regularization for two values of regularization parameter.
Computation 10 00038 g009
Figure 10. The spectrum of C (Uncentered Active Subspace) and C ˜ (Centered Active Subspace), noise level 1%.
Figure 10. The spectrum of C (Uncentered Active Subspace) and C ˜ (Centered Active Subspace), noise level 1%.
Computation 10 00038 g010
Figure 11. Eigenvectors of uncentered AS w i and centered AS v i , i = 1 , 2 , 6 , 9 and their corresponding eigenvalues, noise level 1%.
Figure 11. Eigenvectors of uncentered AS w i and centered AS v i , i = 1 , 2 , 6 , 9 and their corresponding eigenvalues, noise level 1%.
Computation 10 00038 g011
Figure 12. Solutions for X-ray tomography using Tikhonov regularization with regularization parameters α = 10 2 , 1.3 × 10 5 , 10 8 and the original image (right-most), noise level 1%.
Figure 12. Solutions for X-ray tomography using Tikhonov regularization with regularization parameters α = 10 2 , 1.3 × 10 5 , 10 8 and the original image (right-most), noise level 1%.
Computation 10 00038 g012
Figure 13. DIAS solutions for X-ray tomography problem with underregularization parameter α = 10 2 , noise level 1%.
Figure 13. DIAS solutions for X-ray tomography problem with underregularization parameter α = 10 2 , noise level 1%.
Computation 10 00038 g013
Figure 14. DIAS solutions for X-ray tomography problem with optimal regularization parameter α optimal = 1.3 × 10 5 , noise level 1%.
Figure 14. DIAS solutions for X-ray tomography problem with optimal regularization parameter α optimal = 1.3 × 10 5 , noise level 1%.
Computation 10 00038 g014
Figure 15. DIAS solutions for X-ray tomography problem with overregularization parameter α = 10 8 , noise level 1%.
Figure 15. DIAS solutions for X-ray tomography problem with overregularization parameter α = 10 8 , noise level 1%.
Computation 10 00038 g015
Table 1. Comparison of features between different approaches (✓ = possess, ✗ = not possess).
Table 1. Comparison of features between different approaches (✓ = possess, ✗ = not possess).
TSVDTikhonovDI [21]DIAS
Solving ill-posed inverse problem
Parameter robustness/Avoiding regularizing data-informed modes
Ordering data-informed modes based on observational data
Readily aplicable to nonlinear inverse problems
Taking the uncertain nature of the inverse problem into account
Table 2. cos w i , v j between first ten eigenvectors in W 1 and V 1 .
Table 2. cos w i , v j between first ten eigenvectors in W 1 and V 1 .
v 1 v 2 v 3 v 4 v 5 v 6 v 7 v 8 v 9 v 10
w 1 −0.010.970−0.0200.220−0.010−0.02
w 2 10.0100000000
w 3 0010000000
w 4 00.02010−0.010000
w 5 0000100000
w 6 0−0.2200.0100.970−0.010−0.01
w 7 000000−1000
w 8 00.010000.010100
w 9 0000000010
w 10 00.020000.010001
Table 3. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 10 2 , noise level 1%.
Table 3. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 10 2 , noise level 1%.
MethodNumber of Dimensions of Active Subspace
r = 1r = 500r = 2000r = 5000r = 10,000r = 12,000
DIAS-A79.4749.3136.1327.3922.7435.36
TSVD84.5856.5040.1729.1022.9835.36
DIAS-F389.15389.16389.15389.16389.13389.01
DI389.15389.16389.15389.17389.13389.01
Tikhonov389.15
Table 4. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for optimal regularization case with α optimal = 1.3 × 10 5 , noise level 1%.
Table 4. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for optimal regularization case with α optimal = 1.3 × 10 5 , noise level 1%.
MethodNumber of Dimensions of Active Subspace
r = 1r = 500r = 2000r = 5000r = 10,000r = 12,000
DIAS-A79.4749.3136.1327.3922.7435.36
TSVD84.5856.5040.1729.1022.9835.36
DIAS-F21.3521.2221.0821.0522.1735.39
DI21.3621.3321.2421.1922.2135.39
Tikhonov21.36
Table 5. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for overregularization case with α = 10 8 , noise level 1%.
Table 5. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for overregularization case with α = 10 8 , noise level 1%.
MethodNumber of Dimensions of Active Subspace
r = 1r = 500r = 2000r = 5000r = 10,000r = 12,000
DIAS-A79.4749.3136.1327.3922.7435.36
TSVD84.5856.5040.1729.1022.9835.36
DIAS-F74.2948.6535.9127.3222.7335.36
DI78.0855.639.8929.0122.9735.36
Tikhonov80.77
Table 6. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 10 2 , noise level 3%.
Table 6. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 10 2 , noise level 3%.
MethodNumber of Dimensions of Active Subspace
r = 1r = 500r = 2000r = 5000r = 10,000r = 12,000
DIAS-A79.4749.4237.3233.9548.4697.10
TSVD84.5856.5741.0535.1048.5397.10
DIAS-F371.83371.82371.81371.80371.81372.87
DI371.83371.83371.83371.81371.81372.88
Tikhonov371.83
Table 7. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 1.35 × 10 5 , noise level 3%.
Table 7. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 1.35 × 10 5 , noise level 3%.
MethodNumber of Dimensions of Active Subspace
r = 1r = 500r = 2000r = 5000r = 10,000r = 12,000
DIAS-A79.4749.4237.3233.9548.4697.10
TSVD84.5856.5741.0535.1048.5397.10
DIAS-F33.9932.2430.8932.7748.4297.10
DI34.1733.5632.1133.3048.4697.10
Tikhonov34.17
Table 8. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 10 8 , noise level 3%.
Table 8. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 10 8 , noise level 3%.
MethodNumber of Dimensions of Active Subspace
r = 1r = 500r = 2000r = 5000r = 10,000r = 12,000
DIAS-A79.4749.4237.3233.9548.4697.10
TSVD84.5856.5741.0535.1048.5297.10
DIAS-F78.7749.3437.3033.9548.4697.10
DI83.6856.4741.0235.0948.5297.10
Tikhonov95.19
Table 9. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 10 2 , noise level 5%.
Table 9. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 10 2 , noise level 5%.
MethodNumber of Dimensions of Active Subspace
r = 1r = 500r = 2000r = 5000r = 10,000r = 12,000
DIAS-A79.4749.5439.0343.4977.18161.65
TSVD84.5856.6942.6244.3977.22161.62
DIAS-F326.77326.77326.78326.77326.69332.71
DI326.77326.77326.78326.78326.69332.71
Tikhonov326.77
Table 10. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 1.35 × 10 5 , noise level 5%.
Table 10. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 1.35 × 10 5 , noise level 5%.
MethodNumber of Dimensions of Active Subspace
r = 1r = 500r = 2000r = 5000r = 10,000r = 12,000
DIAS-A79.4749.5439.0343.4977.18161.65
TSVD84.5856.6942.6244.3977.22161.62
DIAS-F43.5538.7335.9343.0877.19161.65
DI44.0541.7338.0743.7377.22161.62
Tikhonov44.06
Table 11. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 10 8 , noise level 5%.
Table 11. Relative error (%) of inverse solutions by DIAS and Tikhonov regularization for underregularization case with α = 10 8 , noise level 5%.
MethodNumber of Dimensions of Active Subspace
r = 1r = 500r = 2000r = 5000r = 10,000r = 12,000
DIAS-A79.4749.5439.0343.4977.18161.65
TSVD84.5856.6942.6244.3977.22161.61
DIAS-F79.2149.5239.0243.4977.18161.64
DI84.2556.6642.6144.3877.22161.61
Tikhonov98.06
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Nguyen, H.; Wittmer, J.; Bui-Thanh, T. DIAS: A Data-Informed Active Subspace Regularization Framework for Inverse Problems. Computation 2022, 10, 38. https://0-doi-org.brum.beds.ac.uk/10.3390/computation10030038

AMA Style

Nguyen H, Wittmer J, Bui-Thanh T. DIAS: A Data-Informed Active Subspace Regularization Framework for Inverse Problems. Computation. 2022; 10(3):38. https://0-doi-org.brum.beds.ac.uk/10.3390/computation10030038

Chicago/Turabian Style

Nguyen, Hai, Jonathan Wittmer, and Tan Bui-Thanh. 2022. "DIAS: A Data-Informed Active Subspace Regularization Framework for Inverse Problems" Computation 10, no. 3: 38. https://0-doi-org.brum.beds.ac.uk/10.3390/computation10030038

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop