Conditional Karhunen-Lo\`eve expansion for uncertainty quantification and active learning in partial differential equation models

TL;DR

This paper introduces a conditional Karhunen-Loève expansion approach for quantifying and reducing uncertainty in stochastic PDEs with partially known coefficients, employing active learning to optimize observation locations.

Contribution

It compares two conditioning approaches for KL expansions and demonstrates their effectiveness in uncertainty reduction and active learning in PDE models.

Findings

Approach 1 yields more accurate conditional coefficients and solutions.

Active learning based on Method 2 achieves greater variance reduction.

Conditioning reduces the dimensionality of the KL representation.

Abstract

We use a conditional Karhunen-Lo\`eve (KL) model to quantify and reduce uncertainty in a stochastic partial differential equation (SPDE) problem with partially-known space-dependent coefficient, $Y(x)$. We assume that a small number of $Y(x)$ measurements are available and model $Y(x)$ with a KL expansion. We achieve reduction in uncertainty by conditioning the KL expansion coefficients on measurements. We consider two approaches for conditioning the KL expansion: In Approach 1, we condition the KL model first and then truncate it. In Approach 2, we first truncate the KL expansion and then condition it. We employ the conditional KL expansion together with Monte Carlo and sparse grid collocation methods to compute the moments of the solution of the SPDE problem. Uncertainty of the problem is further reduced by adaptively selecting additional observation locations using two active learning methods. Method 1 minimizes the variance of the PDE coefficient, while Method 2 minimizes the variance of the solution of the PDE. We demonstrate that conditioning leads to dimension reduction of the KL representation of $Y(x)$. For a linear diffusion SPDE with uncertain log-normal coefficient, we show that Approach 1 provides a more accurate approximation of the conditional log-normal coefficient and solution of the SPDE than Approach 2 for the same number of random dimensions in a conditional KL expansion. Furthermore, Approach 2 provides a good estimate for the number of terms of the truncated KL expansion of the conditional field of Approach 1. Finally, we demonstrate that active learning based on Method 2 is more efficient for uncertainty reduction in the SPDE's states (i.e., it leads to a larger reduction of the variance) than active learning using Method 2.

Equations

$$0<k_l\le k(\mathbf{x},\omega )\le k_u<\infty \text{a.e.}\text{in}D\times \Omega .$$

$$\nabla \cdot \left[k(\mathbf{x},\omega )\nabla u(\mathbf{x},\omega )\right]$$

$$u(\mathbf{x},\omega )$$

$$\nabla u(\mathbf{x},\omega )\cdot n$$

$$g\sim \mathcal{G}\mathcal{P}(0,C_g(\mathbf{x},{\mathbf{x}}^{\prime }\mid \mathbf{\theta })),$$

$$C_g(\mathbf{x},{\mathbf{x}}^{\prime }\mid \mathbf{\theta })={\sigma }^2\exp \left[-\frac{(x_1-x_1^{\prime }{)}^2}{l_1^2}-\frac{(x_2-x_2^{\prime }{)}^2}{l_2^2}\right],$$

$$\mathbf{y}=\mathbf{g}+\mathbf{ϵ},\mathbf{ϵ}\sim \mathcal{N}(0,{\sigma }_{ϵ}^2\mathbf{I}),\mathbb{E}[\mathbf{g}{\mathbf{ϵ}}^{\top }]=0.$$

$$(\mathbf{\theta },{\sigma }_{ϵ})\equiv \underset{\mathbf{\theta },{\sigma }_{ϵ}}{\mathrm{argmax}}L(\mathbf{\theta },{\sigma }_{ϵ},X,\mathbf{y}),$$

$$L(\mathbf{\theta },{\sigma }_{ϵ},X,\mathbf{y})=-\frac{1}{2}{y}^{\top }|C_s(\mathbf{\theta })+{\sigma }_{ϵ}^{2}I{|}^{-1}y-\frac{1}{2}\ln |C_s(\mathbf{\theta })+{\sigma }_{ϵ}^{2}I|-\frac{N_s}{2}\ln 2\pi ,$$

$${\mu }_g^c(\mathbf{x})$$

$$C_g^c(\mathbf{x},{\mathbf{x}}^{\prime })$$

$$g^c(\mathbf{x},\omega )={\mu }_g^c(\mathbf{x})+\sum _{i=1}^{d^c}\sqrt{{\lambda }_i^c}{\varphi }_i^c(\mathbf{x}){\xi }_i,$$

$${\int }_D{C}_g^c(\mathbf{x},{\mathbf{x}}^{\prime }){\varphi }^c({\mathbf{x}}^{\prime })\mathrm{d}{\mathbf{x}}^{\prime }={\lambda }^c{\varphi }^c(\mathbf{x}),$$

$$g(\mathbf{x},\omega )=\sum _{i=1}^{\infty }\sqrt{{\lambda }_i}{\varphi }_i(\mathbf{x}){\xi }_i(\omega ),$$

$$g(\mathbf{x},\omega )={\Phi }^{\top }(\mathbf{x}){\Lambda }^{1/2}\mathbf{\xi }(\omega ),$$

$$C_g(\mathbf{x},{\mathbf{x}}^{\prime })=\sum _{i=1}^{\infty }{\lambda }_i{\varphi }_i(\mathbf{x}){\varphi }_i({\mathbf{x}}^{\prime }),$$

$${\int }_D{C}_g(\mathbf{x},{\mathbf{x}}^{\prime })\varphi ({\mathbf{x}}^{\prime })\mathrm{d}{\mathbf{x}}^{\prime }=\lambda \varphi (\mathbf{x}).$$

$$C_g(\mathbf{x},{\mathbf{x}}^{\prime })={\Phi }^{\top }(\mathbf{x})\Lambda \Phi (\mathbf{x}).$$

$$\mathbf{y}={\Phi }^{\top }(X){\Lambda }^{1/2}\mathbf{\xi }+\mathbf{ϵ}.$$

$$\left[\begin{array}{c}\mathbf{y}\\ \mathbf{\xi }\end{array}\right]\sim \mathcal{N}\left(\left[\begin{array}{c}0\\ 0\end{array}\right],\left[\begin{array}{cc}{C}_s+{\sigma }_{ϵ}^{2}I& {\Phi }^{\top }(X){\Lambda }^{1/2}\\ {\Lambda }^{1/2}{\Phi }^{\top }(X)& I\end{array}\right]\right),$$

$$\mathbf{\mu }$$

$$M$$

$$\tilde{g}(\mathbf{x},\omega )={\Phi }^{\top }(\mathbf{x}){\Lambda }^{1/2}\tilde{\mathbf{\xi }}(\omega ).$$

$$g^d(\mathbf{x},\omega )$$

$$=({\Phi }^d{)}^{\top }(\mathbf{x})({\Lambda }^d{)}^{1/2}{\mathbf{\xi }}^d(\omega ),$$

$$C_g^d(\mathbf{x},{\mathbf{x}}^{\prime })=({\Phi }^d{)}^{\top }(\mathbf{x}){\Lambda }^d{\Phi }^d({\mathbf{x}}^{\prime }),$$

$${\tilde{g}}^d(\mathbf{x},\omega )=({\Phi }^d{)}^{\top }(\mathbf{x})({\Lambda }^d{)}^{1/2}{\tilde{\mathbf{\xi }}}^d,$$

$${\mathbf{\mu }}^d$$

$$M^d$$

$$\mathrm{rk}(M^d)=r\equiv d-N_s,$$

$$D=\left[\begin{array}{cc}{D}^r& 0\\ 0& 0\end{array}\right],$$

$${\tilde{g}}^d(\mathbf{x},\omega )={\hat{g}}^d(\mathbf{x})+({\Phi }^d{)}^{\top }(\mathbf{x})({\Lambda }^d{)}^{1/2}QD{Q}^{-1}{\mathbf{\zeta }}^d(\omega ),$$

$$\begin{gathered} {\tilde{g}}^d(\mathbf{x},\omega )-{\hat{g}}^d(\mathbf{x})= \\ \left[\begin{array}{cc}({\Phi }_r^d{)}^{\top }(\mathbf{x})& ({\Phi }^d{)}_{r^{\prime }}^{\top }(\mathbf{x})\end{array}\right]{\left[\begin{array}{cc}{\Lambda }_r^d& 0\\ 0& {\Lambda }_{r^{\prime }}^d\end{array}\right]}^{1/2}\left[\begin{array}{cc}{Q}_r& Q_{r{r}^{\prime }}\\ Q_{r^{\prime }r}& Q_{r^{\prime }}\end{array}\right]\left[\begin{array}{c}{D}_r{\mathbf{\eta }}_r\\ 0\end{array}\right]. \end{gathered}$$

$${\tilde{g}}^d(\mathbf{x},\omega )-\hat{g}(\mathbf{x})$$

Full text

We employ the conditional KL expansion together with Monte Carlo and sparse grid collocation methods to compute the moments of the solution of the SPDE problem. Uncertainty of the problem is further reduced by adaptively selecting additional observation locations using two active learning methods. Method 1 minimizes the variance of the PDE coefficient, while Method 2 minimizes the variance of the solution of the PDE.

We demonstrate that conditioning leads to dimension reduction of the KL representation of $Y(x)$. For a linear diffusion SPDE with uncertain log-normal coefficient, we show that Approach 1 provides a more accurate approximation of the conditional log-normal coefficient and solution of the SPDE than Approach 2 for the same number of random dimensions in a conditional KL expansion. Furthermore, Approach 2 provides a good estimate for the number of terms of the truncated KL expansion of the conditional field of Approach 1. Finally, we demonstrate that active learning based on Method 2 is more efficient for uncertainty reduction in the SPDE’s states (i.e., it leads to a larger reduction of the variance) than active learning using Method 2.

††journal: Journal of Computational Physics

1 Introduction

Uncertainty quantification in partial differential equations (PDE) problems with partially known parameters (e.g., coefficients and source terms) is often performed by modeling these partially known quantities as random variables with appropriate probability distribution functions. Spectral methods such as Polynomial chaos (PC)-based stochastic Galerkin [1, 2] and stochastic collocation [3, 4] are commonly used for solving PDEs with random parameters. In spectral methods, random fields are represented in terms of their Karhunen-Loève (KL) expansions. While infinite KL expansions are necessary to exactly represent the two-point statistics of a random field, numerical treatment of PDEs requires truncating KL expansions. These KL expansions are truncated based on the decay of their eigenvalues, so that the random fields can be reconstructed with sufficient accuracy using the retained KL terms.

The computational cost of spectral methods exponentially increases with the dimensionality of the stochastic problem (i.e., the number of terms in the truncated KL expansion) [5, 6, 7, 8, 9]. Various approaches have been recently proposed to address this issue: by finding the solution in reduced dimensional spaces using basis adaption [10], domain decomposition methods [11, 12], sliced inverse regression [13, 14], and sparsity enhancing together with the active subspace method [15], among others. In this work, we propose reducing the computational cost of spectral methods by conditioning the KL expansion on available data. We also demonstrate that conditioning on data reduces uncertainty of predictions, i.e., reduces the variance of quantities of interest of the governing stochastic PDEs.

Most of existing UQ methods, including all work referenced above, is based on stochastic models of unknown fields with constant variance and stationary covariance functions, that is, covariance functions of the form $C(x,y)=C(|x-y|)$. The unconditional statistics of fields is estimated from data using the so-called ergodicity assumption, where fields are treated as realizations of a random process, and the spatial statistics of the fields is assumed to be the same as the ensemble statistics of the generating random process. However, there is no reason to assume that the variance at locations where measurements are available should be the same as the variance in locations with no measurements.

Gaussian process (GP) regression, also known as kriging, has been used in geostatistics and hydrogeological modeling since its introduction in the sixties [16] to represent partially observed properties of materials (e.g., permeability of geological porous media) as random field conditioned on observations [17, 18]. A characteristic feature of GP models of random fields is that their conditional statistics, i.e., statistics conditioned on data, are not stationary. For example, the conditional variance of a GP model is a function of space (e.g., in the absence of measurement errors, it is zero at the observation locations), and the conditional covariance function $C^{c}(x,y)$ depends explicitly on both $x$ and $y$. There are very few studies of PDEs with non-stationary random fields conditioned on data, including the conditional Moment Method [19, 20] limited to parameters with small variances, and a few papers on conditional PC methods [21, 22, 23].

Two main approaches have been proposed for conditional spectral methods: (i) first truncating the KL expansion of the random parameters and then conditioning the resulting truncated expansion on data [22], and (ii) first conditioning the infinite KL expansion and then truncating it [23]. Here, we demonstrate that the truncating and conditioning operations do not commute as the two approaches produce different results. A detailed analysis of spectral methods for PDEs with non-stationary random inputs is clearly lacking. In this work, we study the application of conditional KL models to quantify and reduce uncertainty in physical systems modeled with stochastic PDEs (SPDEs). We compare the solution of the SPDE in terms of its conditional mean and variance obtained using both the conditioning first and truncating first approaches, and discuss the merits of both constructions. Our results show that the approach of conditioning first and then truncating (with $N$ terms) approach is more accurate than the approach of truncating first (with $N$ terms) and then conditioning. We also show that for the truncating-first approach the final dimension of the conditional model is reduced from $N$ to $N-N_{s}$, where $N_{s}$ is the number of measurements.

Furthermore, we adopt active learning [24, 25] to identify additional observation locations in order to efficiently reduce predictive uncertainty in terms of the variance of the solution of the SPDE. We consider two criteria for identifying observation locations. In the first criterion, we choose the location that minimizes a norm of the variance of the conditional KL expansion of the model parameters. In the second criterion, we propose a novel, GP regression-based approximation to the conditional variance of the solution of the SPDE, and choose locations that minimize this approximation. We demonstrate that the second strategy leads to higher reduction of predictive uncertainty for the same number of additional measurements.

This manuscript is organized as follows: In Section 2 we formulate a steady-state stochastic diffusion equation with random coefficient and the GP model for the random coefficient. In Section 3, we describe two approaches for constructing finite-dimensional conditional KL models of the random coefficient. In Section 4, we present active learning criteria for further reducing uncertainty in conditional KL models. Numerical examples are given in Section 5 and conclusions are presented in Section 6.

2 Governing Equations

We study a two-dimensional steady-state diffusion equation with a random diffusion coefficient $k(\mathbf{x},\omega):D\times\Omega\to\mathbb{R}$, $D\subset\mathbb{R}^{2}$ that is bounded and strictly positive:

[TABLE]

We seek the stochastic solution $u(\mathbf{x},\omega):D\times\Omega\to\mathbb{R}$ to the problem

[TABLE]

where the boundary conditions $u_{\Gamma}$ and $u_{\Gamma^{\prime}}$ are deterministic and known, and $n$ denotes the outward-pointing unit vector normal to $\partial D$. In many practical applications, the full probabilistic characterization of the coefficient $k(\mathbf{x},\omega)$ is not known, but measurements of $k$ are available at a few spatial locations. In this work, we assume that the distribution of $k$ is known and is log-normal [26, 17], i.e., $k(\mathbf{x},\omega)\equiv\exp[g(\mathbf{x},\omega)]$, where $g(\mathbf{x},\omega)\equiv\log k(\mathbf{x},\omega)$ is a Gaussian random field. We construct the GP prior or unconditional (i.e., not conditioned on measurements) model for $g$ employing the following two-step approach common in geostatistics [18]: we first select a parameterized GP prior covariance kernel $C_{g}(\mathbf{x},\mathbf{x}^{\prime}\mid\bm{\theta}):D\times D\to\mathbb{R}$, with hyperparameters $\bm{\theta}$; and next, we compute an estimate $\hat{\bm{\theta}}$ for the hyperparameters from available measurements of $\log k$ by type-II maximum likelihood estimation [18, 27]. In the first step, we assume that $g$ is wide-sense stationary with zero mean and squared exponential covariance function

[TABLE]

where $\bm{\theta}=\{\sigma,l_{1},l_{2}\}$ is the set of hyperparameters of the covariance kernel: $\sigma$ is the standard deviation, and $l_{1}$ and $l_{2}$ are the correlation lengths along the $x_{1}$ and $x_{2}$ spatial coordinates, respectively. To estimate $\bm{\theta}$, we assume that measurements are of the form $y_{i}=g(\mathbf{x}^{*}_{i})+\epsilon_{i}$, $i\in[1,N_{s}]$, where $N_{s}$ is the number of observations, $\mathbf{x}^{*}_{i}$ is the measurement location for the $i$th measurement, and the $\epsilon_{i}\sim\mathcal{N}(0,\sigma_{\epsilon})$ are iid measurement errors with standard deviation $\sigma_{\epsilon}$ independent of $g(\mathbf{x},\omega)$. The observations and observation locations are arranged into the vector $\mathbf{y}=(y_{1},\dots,y_{N_{s}})^{\top}$ and the matrix $X=(\mathbf{x}^{*}_{1},\dots,\mathbf{x}^{*}_{N_{s}})$, respectively. Similarly, the values of $g$ at the observation locations are arranged into the vector $\mathbf{g}=(g(\mathbf{x}^{*}_{1}),\dots,(\mathbf{x}^{*}_{N_{s}}))^{\top}$, so that $\mathbf{y}$ and $\mathbf{g}$ are related by

[TABLE]

As stated above, in this work we employ type II maximum likelihood estimation [18, 27] to estimate the hyperparameters of the GP prior and of the likelihood of the observations:

[TABLE]

where $L(\bm{\theta},\sigma_{\epsilon},X,\mathbf{y})$ is the log-marginal likelihood function

[TABLE]

and $C_{s}(\bm{\theta})$ is the observation covariance matrix with $ij$th entry given by $C_{g}(\mathbf{x}^{*}_{i},\mathbf{x}^{*}_{j}\mid\bm{\theta})$.

To compute the solution of the SPDE (2) conditioned on the observations $\mathbf{y}$, we adopt the stochastic collocation approach, which requires a finite-dimension stochastic representation of the random field $g(\mathbf{x},\omega)$. For constructing this representation, we will consider two strategies that rely on GP regression and the KL expansion of random fields:

In the first approach [23], the conditional random field $g^{c}(\mathbf{x},\omega)$ is obtained using GP regression and then discretized by calculating its KL expansion in terms of standard Gaussian random variables. Then, the KL expansion is truncated to an appropriate number of random dimensions $d^{c}$. 2. 2.

In the second approach [22], the unconditioned random field $g(\mathbf{x},\omega)$ is first discretized in terms of its KL expansion and unconditioned standard Gaussian random variables $\bm{\xi}$. Then, conditional Gaussian random variables $\tilde{\bm{\xi}}$, conditioned on the observations $\mathbf{y}$, are obtained by projection.

We describe these two approaches in details in Section 3.

3 Conditional KL models

3.1 Approach 1: truncated KL expansion of the conditioned GP field

In the first approach, the conditional random field $g^{c}(\mathbf{x},\omega)$ is approximated with a KL expansion written in terms of standard Gaussian random variables. We assume that the hyperparameters of the prior covariance function $C_{g}(\mathbf{x},\mathbf{x}^{\prime})$ have been estimated and are thus dropped from the notation. The mean and covariance of conditioned Gaussian random field $g^{c}(\mathbf{x},\omega)\sim\mathcal{GP}\left(\mu_{g}^{c},C_{g}^{c}(\mathbf{x},\mathbf{x}^{\prime})\right)$ are computed using GP regression [18, 27]:

[TABLE]

The conditional field is then expanded using a truncated KL expansion [1] as

[TABLE]

where $\{\xi_{i}\sim\mathcal{N}(0,1)\}^{d^{c}}_{i=1}$ are standard Gaussian random variables, and $\{\lambda_{i}^{c},\phi_{i}^{c}\}^{d^{c}}_{i=1}$ are the first $d^{c}$ pairs of eigenvalues and eigenfunctions stemming from the eigenvalue problem

[TABLE]

where $C^{c}_{g}$ is given by (9). The conditional solution of the SPDE, $u^{c}$, can then be computed using MC or sparse grid collocation methods by sampling $g^{c}(\mathbf{x},\omega)$ using (10), or by the PC method by constructing a spectral approximation of $u^{c}$ in terms of the $\xi_{i}$, $i\in[1,d^{c}]$.

3.2 Approach 2. Conditioning truncated KL expansion of the unconditioned field

In this approach, introduced in [22], the KL expansion of $g(\mathbf{x},\omega)$ is first truncated, and the resulting set of random variables are conditioned on the observations $\mathbf{y}$. Here, we demonstrate that this approach reduces the number of random dimensions of the representation of $g^{c}$ by the number of observations and propose a method for rewriting the representation of $g^{c}$ in terms of the reduced number of random variables. Effectively, this reduces dimensionality of the SPDE solution.

To present our approach we first summarize the conditional KL construction presented in [22] and next describe our proposed reduced-dimension conditional KL representation. We start with the KL expansion of the unconditional Gaussian random field $g(\mathbf{x},\omega)$:

[TABLE]

where the $\{\xi_{i}\sim\mathcal{N}(0,1)\}^{\infty}_{i=1}$ are iid standard Gaussian random variables. In matrix notation,

[TABLE]

where $\Phi(\mathbf{x})=(\phi_{1}(\mathbf{x}),\phi_{2}(\mathbf{x}),\dots)^{\top}$, $\Lambda=\operatorname{diag}(\lambda_{1},\lambda_{2},\dots)$, and $\bm{\xi}(\omega)$ is the infinite-dimensional random vector $\bm{\xi}(\omega)=(\xi_{1}(\omega),\xi_{2}(\omega),\dots)^{\top}$. By Mercer’s theorem, the covariance function of $g$ can be represented as the infinite sum

[TABLE]

where the eigenpairs $\{\lambda_{i},\phi_{i}\}^{\infty}_{i=1}$ are the solutions to the eigenvalue problem

[TABLE]

In matrix notation, the Mercer expansion (13) can be written as

[TABLE]

In [22], the KL expansion (11) is conditioned on the observations $\mathbf{y}$ by conditioning $\bm{\xi}$ on these observations as described in the following Eqs (15)–(22). Evaluating (12) at $X$ and substituting into (5) yields

[TABLE]

Therefore, the joint distribution of $\mathbf{y}$ and $\bm{\xi}$ is given by

[TABLE]

where we have employed the relation $\mathbb{E}[\mathbf{g}\bm{\epsilon}^{\top}]=0$ and introduced the notation $C_{s}\equiv\mathbb{E}[\mathbf{g}\mathbf{g}^{\top}]$.

It follows that the distribution of $\bm{\xi}$ conditioned on the measurements is $\bm{\xi}\mid X,\mathbf{y}\sim\mathcal{N}(\bm{\mu},M)$, where

[TABLE]

For simplicity, we denote the GP conditional on $(X,\mathbf{y})$ by $\tilde{g}$, and the conditional random vector by $\tilde{\bm{\xi}}=(\tilde{\xi}_{1},\tilde{\xi}_{2},\dots)^{\top}$. The conditional GP then reads

[TABLE]

We now apply the process of conditioning random variables to the truncated KL expansion of the unconditional random field. The KL expansion of the unconditional random field $g(\mathbf{x},\omega)$, truncated to $d$ terms, reads

[TABLE]

where $\Phi^{d}(\mathbf{x})=(\phi_{1}(\mathbf{x}),\dots,\phi_{d}(\mathbf{x}))^{\top}$, $\Lambda^{d}=\operatorname{diag}(\lambda_{1},\dots,\lambda_{d})$, and $\bm{\xi}^{d}(\omega)=(\xi_{1}(\omega),\dots,\xi_{d}(\omega))^{\top}$. This expansion corresponds to the truncated covariance

[TABLE]

obtained by substituting $\Phi^{d}(\mathbf{x})$ and $\Lambda^{d}$ for $\Phi$ and $\Lambda$ in Eq (13). The truncated representation (18) can be conditioned on the data $(X,\mathbf{y})$ by following the procedure outlined above, resulting in the conditional model

[TABLE]

where $\tilde{\bm{\xi}}^{d}\equiv\bm{\xi}^{d}\mid X,\mathbf{y}\sim\mathcal{N}(\bm{\mu}^{d},M^{d})$, with $\bm{\mu}^{d}$ and $M^{d}$ given by

[TABLE]

and where $C^{d}_{s}\equiv(\Phi^{d})^{\top}(X)\Lambda^{d}\Phi^{d}(X)$ is the truncated measurement covariance matrix.

Note that both (17) and (21) employ the same set of eigenpairs derived from the unconditioned covariance function $C_{g}$. Nevertheless, $\tilde{\bm{\xi}}^{d}$ and $(\tilde{\xi}_{1},\dots,\tilde{\xi}_{d})$ in (17) have different joint distribution, which implies that the conditional truncated model (21) is different than the model (17) after truncating to $d$ terms.

Due to the conditioning on $N_{s}$ measurements, the rank of $M^{d}$ is

[TABLE]

so that the model (21) is effectively of dimension $r$. In other words, conditioning the truncated model results in dimension reduction of the GP model for $g$, and reduction of the stochastic dimensionality of the SPDE problem.

To leverage this dimension reduction, we propose rewriting the model (21) in terms of $r$ random variables. We write the eigendecomposition of $M^{d}$ as $M^{d}=QDQ^{-1}$, where $D$ is the diagonal matrix of the form

[TABLE]

since $\operatorname{rk}(M^{d})=r$. Substituting the eigendecomposition of $M^{d}$ into (21), we obtain

[TABLE]

where $\hat{g}^{d}(\mathbf{x})\equiv(\Phi^{d})^{\top}(\mathbf{x})(\Lambda^{d})^{1/2}\bm{\mu}^{d}$ and $\bm{\zeta}^{d}\sim\mathcal{N}(0,I_{d})$.

Let $\bm{\eta}=Q^{-1}\bm{\zeta}^{d}$; then, substituting $\bm{\eta}$ into (26), we obtain

[TABLE]

Similarly, let $\tilde{\bm{\eta}}_{r}=D_{r}\bm{\eta}_{r}$; then, $\tilde{g}^{d}$ can be expanded in terms of the reduced set of random variables $\tilde{\bm{\eta}}_{r}$ as

[TABLE]

Here, the components of the random vector $\tilde{\bm{\eta}}_{r}$ are correlated with one another; nevertheless, $\tilde{\bm{\eta}}_{r}$ can be converted to a set of uncorrelated random variables using an orthogonalization method such as Gram-Schmidt process, so that $\bm{\zeta}^{r}=A\tilde{\bm{\eta}}_{r}$ and $\tilde{\bm{\eta}}_{r}=A^{-1}\bm{\zeta}^{r},$ where $\bm{\zeta}^{r}\sim\mathcal{N}(0,I_{r})$. Now, $\tilde{g}^{d}$ can be expanded in terms of the new set of random variables $\bm{\zeta}^{r}$ as

[TABLE]

where

[TABLE]

The $r$ components of $\bm{\zeta}^{r}$ are uncorrelated Gaussian random variables and hence are also independent. Therefore, the conditional KL expansion (28) can be used in combination with standard PC-based or stochastic collocation-based methods for solving stochastic PDEs.

4 Active learning for uncertainty reduction

The conditioning of $g(\mathbf{x},\omega)$ KL expansion presented in Section 3 leads to reduced uncertainty in $u(\mathbf{x},\omega)$. This uncertainty can be further reduced by collecting additional measurements of $g$. Giving limited ability to collect additional $g$ measurements, it is important to identify locations for additional measurements that optimally minimize the uncertainty in $u$. This process of optimal data acquisition is referred to as active learning [28, 24, 25].

In this section, we discuss the standard active learning method (Method 1) based on minimizing the variance of $g^{c}$  [24, 25] and propose an alternative approach (Method 2) based on minimizing the variance of $u^{c}$. In Method 1, the variance of $u^{c}$ is reduced by minimizing the variance of $g^{c}$, but, as we show below, it does not lead to the maximum reduction of the $u^{c}$ variance.

It is worth noting that active learning also reduces the cost of computing the conditional solution $u^{c}$, as solving SPDE problems with coefficients having smaller variance requires smaller number of MC simulations, lower order of polynomial chaos in the stochastic Galerkin method, and lower sparse grid level in the stochastic collocation method.

4.1 Method 1: minimization of the conditional variance of $g^{c}$

In the standard active learning method [24, 25], a new location $\mathbf{x}^{*}$ for sampling $g$ is selected to minimize the variance of the conditional model $g^{c}$ conditioned on the full set of observations, including the new observation. This can be done in closed form as, per Eq. (9), the conditional variance $C^{c}_{g}$ depends only on the observation locations and not the observed values. Therefore, the new observation location $\mathbf{x}^{*}$ is selected following the acquisition policy

[TABLE]

where $C^{c}_{g}(\cdot,\cdot\mid[X,\mathbf{x}^{\prime}])$ denotes the covariance of $g$ conditioned on the original set of observation locations, $X$, and a new location $\mathbf{x}^{\prime}$, and is computed using (9). In practice, the minimization problem (30) is approximately solved by choosing $\mathbf{x}^{*}$ as $\operatorname*{arg\,max}_{\mathbf{x}}C^{c}_{g}\left(\mathbf{x},\mathbf{x}\mid X\right)$ [24, 25].

Reducing the variance of $g^{c}$ also reduces the variance of $u^{c}$. Nevertheless, there is no reason to assume that the observation locations provided by the acquisition strategy (30) will lead to optimal variance reduction for $u^{c}$.

4.2 Method 2: By minimizing the conditional covariance of $u^{c}$

In Method 2, we chose location for measurements that minimize the conditional variance of $u$, that is,

[TABLE]

where $C^{c}_{u}(\cdot,\cdot\mid[X,\mathbf{x}^{\prime}])$ denotes the covariance of $u$ conditioned on the original set of observation locations, $X$, and a new location $\mathbf{x}^{\prime}$. Note that solving this minimization problem is significantly more challenging than solving (30), as it requires (i) constructing the conditional model $g^{c}$ conditioned on the new observation, and (ii) forward uncertainty propagation of $g^{c}$ through the SPDE problem to estimate $C^{c}_{u}$. Here, step (i) is the most challenging as it requires knowledge of $g(\mathbf{x}^{\prime})$ in order to construct the conditioned model $g^{c}$  [28]. Note that $g(\mathbf{x})$ is only partially known; otherwise, this problem would not be uncertain.

As an alternative, we propose modeling $g^{c}$ and $u^{c}$ as components of a bivariate Gaussian field $h^{c}$ in order to derive an approximation to $C^{c}_{u}(\mathbf{x},\mathbf{x}\mid[X,\mathbf{x}^{\prime}])$ that does not require knowledge of $g(\mathbf{x}^{\prime})$. Specifically, let $g^{c}$ and $u^{c}$ denote the fields conditioned on the original data set $(\mathbf{x},\mathbf{y})$, and let $h^{c}(\mathbf{x},\omega)=[g^{c}(\mathbf{x},\omega),u^{c}(\mathbf{x},\omega)]^{\top}$ be a bivariate Gaussian field. It is possible to compute the conditional covariance of $h^{c}$ conditioned on an additional measurement location $x^{\prime}$ by employing (9). After marginalizing the $u^{c}$ component, we obtain the approximation

[TABLE]

where $C^{c}_{ug}(\mathbf{x},\mathbf{x}^{\prime}\mid X)$ denotes the $u$–$g$ cross-covariance conditioned on $X$. Note that this approximation only requires the observation location $x^{\prime}$ and not the observation value $g(\mathbf{x}^{\prime})$.

In order to apply the approximation (32), it is necessary to compute the conditional covariances $C^{c}_{u}(\mathbf{x},\mathbf{x}^{\prime}\mid X)$ and $C^{c}_{ug}(\mathbf{x},\mathbf{x}^{\prime}\mid X)$, which, unlike $C^{c}_{g}(\mathbf{x},\mathbf{x}^{\prime}\mid X)$ given by (9), are not available in closed form; therefore, we compute sample approximations of these covariances. For this purpose, we draw $M$ realizations of $g^{c}(\mathbf{x},\omega)$ with conditioned mean and covariance given by (8) and (9), resulting in the ensemble of synthetic fields $\{g^{c}_{(i)}(\mathbf{x})\equiv g^{c}(\mathbf{x},\omega_{(i)})\}^{M}_{i=1}$. For each member of the ensemble we solve the corresponding deterministic PDE problem, resulting in the ensemble of solutions $\{u^{c}_{(i)}\}^{M}_{i=1}$. Employing both ensembles, we compute the sample covariances $\hat{C}^{c}_{g}(\cdot,\cdot)$, $\hat{C}^{c}_{u}(\cdot,\cdot)$ and $\hat{C}^{c}_{ug}(\cdot,\cdot)$. Substituting these sample covariances into (32), and the result into (31), leads to the acquisition policy

[TABLE]

5 Numerical Experiments

In this section, we apply the conditional KL modeling approaches presented in Section 3 and the active learning methods presented in Section 4 to solve the stochastic diffusion equation problem (2). We consider the following two-dimensional steady state diffusion equation with a random diffusion coefficient over the domain $D=[0,2]\times[0,1]$, subject to Dirichlet and Neumann boundary conditions:

[TABLE]

The reference $g=\log k$ field is constructed synthetically by drawing a realization of the GP process (3) and (4) with the parameters $\sigma_{g}=0.65,l_{x}=0.15$ and $l_{y}=0.2$, shown in Figure 1. From this reference field, we draw $40$ observations at random locations to be used for constructing the conditional GP model $g^{c}$ and computing conditional solution of Eq (34).

5.1 Conditional GP models

We construct finite-dimensional conditional GP models for the field $g$ based on the synthetic dataset using the two approaches introduced in Section 3 and employ these models to compute the mean and standard deviation of the conditional solution $u^{c}$ of the SPDE problem (34). We employ MC simulations to compute reference unconditional and conditional mean and standard deviation of $g$ and $u$. We use Approach 2 of Section 3.2, which provides a quantitative measure of dimension reduction due to conditioning, to estimate the dimensionality of the reduced conditional model $g^{c}$. We construct finite-dimensional models for $g^{c}$ of the estimated dimension using both Approaches 1 and 2 and propagate uncertainty through the SPDE problem (34) using the stochastic collocation method [3, 4].

5.1.1 Reference Monte Carlo unconditional and conditional solutions

To sample the reference unconditional $g$ field, we employ the unconditional KL expansion (11) truncated to $d=60$ terms, where $d$ is chosen such as to retain $99\%$ of the unconditional variance, i.e., $\sum_{i=1}^{d}\lambda_{i}\geq 0.99\operatorname{Tr}(C_{g})$, where $\operatorname{Tr}(C_{g})$ is the trace of $C_{g}$. The unconditional mean and covariance of $u$ are computed from MC with $15\,000$ samples. The resulting unconditional mean and standard deviation of $g$ and $u$ are presented in Figure 2. In Figure 3 we present a convergence study of the MC estimators of the unconditional mean and standard deviation of $u$, where the $L_{2}$ norm of the estimators is plotted against the number of MC samples. This figure demonstrates that $15\,000$ samples are sufficient to compute these estimators.

Similarly, to sample the reference conditional field, we compute the mean and covariance using Eqs (8) and (9), compute the KL expansion of the conditional covariance, and truncate the resulting conditional KL expansion to $d=53$ terms to retain $99\%$ of the conditional variance. The conditional mean and covariance of $g$ and $u$ are computed using MC with $15\,000$ samples. The results are presented in Figure 4. The convergence study of the MC estimators of the conditional mean and standard deviation of $u$ is shown in Figure 5 and shows that 15000 samples are sufficient to compute these estimators.

The comparison of Figures 2(b) and 4(b) show that the standard deviation of $g$ is reduced after conditioning on $g$ measurements. Similarly, the standard deviation of $u$, shown in Figures 2(d) and 4(d), also is reduced after conditioning on $g$ measurements.

5.1.2 Conditional solution computed using Approach 2: Conditioning truncated KL expansion of the unconditioned field

In this section we employ Approach 2, presented in Section 3.2, for constructing a finite-dimensional conditional KL model of the random coefficient of the SPDE problem (34). In Section 5.1.1, we determined that 60 dimensions are required to represent the unconditional $g$ field. By virtue of (24), it follows that conditioning these $d=60$ dimensions on $N_{s}=40$ observations reduces the dimensionality of the $g^{c}$ KL representation to $r=20$. As this approach for constructing a discretized conditional model disregards information provided by the higher eigenpairs of the unconditional expansion, we expect the resulting conditional moments (mean and standard deviation) of $g^{c}$ to differ from the reference moments computed in Section 5.1.1. The absolute point-wise errors in the mean and standard deviation of $g^{c}$ are shown in Figures 6(a) and 6(b), respectively.

Next, we employ the conditional KL expansion to estimate the mean and variance of $u^{c}$ using the sparse grid collocation method [3, 4] with $41$, $841$, and $11\,561$ quadrature points. Figure 5 shows the $L_{2}$ norm of the estimators. It can be seen that $841$ points are sufficient to obtain a convergent solution for both mean and standard deviation of $u^{c}$. Absolute point-wise errors of the mean and standard deviation of $u^{c}$ with respect to the reference conditional solution are shown in Figures 6(c) and 6(d), respectively.

5.1.3 Conditioned solution computed using Approach 1: Truncated KL expansion of the conditional field

Next, we apply Approach 1 for constructing a finite-dimensional conditional KL model for the random coefficient of the SPDE problem (34). In the previous section, it was determined that the conditional KL expansion obtained with Approach 2 has 20 dimensions. To compare the accuracy of the two conditional KL approaches, here we truncate the conditional KL expansion of $g^{c}$ to $20$ dimensions. Note that in Section 5.1.1 we demonstrated that $53$ dimensions are necessary to retain $99\%$ of the variance of the exact conditional field $g^{c}$. Therefore, by truncating the KL expansion to $20$ dimensions, we incur in the absolute point-wise errors in the conditional mean and standard deviation of $g^{c}$ shown in Figs. 7(a) and 7(b).

As in the previous section, we employ the conditional KL expansion to estimate the mean and variance of $u^{c}$ using the sparse grid collocation method [3, 4]. Figures 7(c) and 7(d) show the absolute point-wise error in the mean and standard deviation of $u^{c}$ with respect to the reference moments computed in Section 5.1.1.

Comparing Figures 6 and 7, it can be seen that, for a set number of stochastic dimensions of the $g^{c}$ KL expansion, Approach 1 provides a more accurate approximation of the moments of $g^{c}$ and $u^{c}$. On the other hand, through Eq. (24), Approach 2 provides a priori estimate of the dimensionality of the conditional KL model sufficient to obtain an accurate solution.

In order to study why Approach 1 provides a more accurate approximation of conditional moments, we compare the eigendecompositions of $g^{c}(x,\omega)$ provided by Approaches 1 and 2. Here, we note that Approach 2 does not provide an explicit eigendecomposition (cf. Eq. (28)). Therefore, we compute the corresponding eigendecomposition by first computing the covariance matrix induced by Approach 2, given by $(\tilde{\Psi}^{r})^{\top}(\tilde{\Psi}^{r})$, and then compute its eigendecomposition.

Figure 8 shows the eigenvalues of $g^{c}(x,\omega)$ resulting from Approaches 1 and 2, together with the eigenvalues of $g(x,\omega)$. It can be seen that the magnitude of the eigenvalues of $g^{c}(x,\omega)$ is smaller than those of $g(x,\omega)$, which follows from the fact that the variance of the conditional KL expansions is smaller than the variance of the unconditional KL expansion. It can also be seen that the conditional eigenvalues decay faster than the unconditional eigenvalues, especially for the first ten eigenvalues.

Figure 9 shows the first and second eigenfunctions of $g^{c}(x,\omega)$ resulting from Approaches 1 and 2, together with the eigenvalues of $g(x,\omega)$. Here, we note that, by construction, the eigenfunctions of Approach 2 are calculated from the first 20 eigenfunctions of $g(\mathbf{x},\omega)$. In contrast, in Approach 1, all the eigenfunctions of $g(\mathbf{x},\omega)$ contribute to the first 20 eigenfunctions of $g^{c}(\mathbf{x},\omega)$; therefore, the eigenfunctions of Approach 1 can resolve finer-scale features than the eigenfunctions of Approach 2. We attribute the superior accuracy of Approach 1 for approximating conditional moments to its superior capacity for resolving fine-scale features of $g^{c}(\mathbf{x},\omega)$.

5.2 Active learning

Here, we apply the two active learning data acquisition methods presented in Section 4 to identify additional measurement locations for the $g$ field. We use these additional observations together with the previously available observations to construct the conditional KL expansion of $g$ using Approach 1, as it was shown in Section 5.1.3 to be more accurate than Approach 2 for the considered application.

Here, we consider again the SPDE problem (34), and we aim to explore the behavior of the presented data acquisition methods for two choices of $\sigma_{g}$, namely $0.65$ and $1.3$. For $\sigma_{g}=0.65$, Figure 10 shows the $L^{2}$-norm of the standard deviation of $g^{c}$ and $u^{c}$ as a function of the number of additional measurements, $N_{am}$ identified with both active learning methods. In Method 2, the conditional covariances $\hat{C}^{c}_{u}$ and $\hat{C}^{c}_{ug}$ are computed using $200$ MC realizations. This small number of MC realizations is justified as empirical observation shows that very accurate estimates of the covariances are not necessary for obtaining good estimates of the new observation locations from the minimization problem (33).

Figure 10(a) shows that the decay rate of the norm of the $g^{c}$ standard deviation is approximately the same for both methods. As expected, the variance reduction in Method 1 is larger (by approximately $5\%$) than in Method 2, because by construction, Method 1 reduces the variance of $g^{c}$ optimally. On the other hand, Figure 10(b) shows that Method 2 reduces the $L^{2}$-norm of $\sigma^{c}_{u}$ more than Method 1 for most $N_{am}$, and by more than $15\%$ for $N_{am}>11$. For some $N_{am}<11$, Method 2 has larger norm of $\sigma^{c}_{u}$ than Method 1, which we attribute to the error in the $\sigma^{c}_{u}$ approximation (32).

We obtain qualitatively similar (but more pronounced) results for $\sigma_{g}=1.3$, as shown in Figure 11. Method 1 leads to a sharper decrease of the $\sigma^{c}_{g}$ norm and Method 2 results in a sharper decrease of the $\sigma^{c}_{u}$ norm than what we observe for $\sigma_{g}=0.65$. From Figures 10 and 11, we conclude that Method 2 is more efficient than Method 1 for reducing uncertainty in $u^{c}$, and that (32) provides a sufficiently accurate approximation for solving the minimization problem (31).

For the case $\sigma_{g}=0.65$, Figures 12 and 13 show the $\sigma^{c}_{g}$ and $\sigma^{c}_{u}$ fields for $N_{am}=1$, $5$, and $10$ obtained using both active learning methods. Similarly, $\sigma^{c}_{g}$ and $\sigma^{c}_{u}$ for $\sigma_{g}=1.3$ are shown in Figures 14 and 15, respectively. As expected, the locations of the additional observations obtained with Methods 1 and 2 are different. Our results show that if the main objective is to predict states, i.e., $u(\mathbf{x})$, rather than coefficients, i.e., $k(\mathbf{x})$, than Method 2 is more efficient than Method 1.

6 Conclusions

We presented two methods for constructing finite-dimensional conditional Karhunen-Loève (KL) expansions of partially known parameters in PDE problems. We demonstrated that conditioning on data reduces the dimensionality of KL expansions and, most importantly, reduces uncertainty in the solution of PDE problems. Finally, conditioning on data reduces the computational cost of solving stochastic PDEs. We also present a new active learning strategy for acquiring new observations of the data based on minimizing variance of the conditional PDE solution (referred to as Method 2 in the paper) and compared it with the standard active learning method based on minimizing the conditional variance of the partially known parameters (referred to in the paper as Method 1).

In the first approach for constructing finite-dimensional conditional GPs, presented in [23], the parameter field is conditioned first on data and then discretized by computing its KL expansion. In the second approach, presented in [22], the parameter field is discretized first by computing its KL expansion, and then the resulting KL expansion is conditioned on data. For the second approach, we demonstrated that conditioning leads to dimension reduction of the conditional representation, and we proposed a method for constructing a reduced representation in terms of the effective number of random dimensions. For a linear diffusion SPDE with uncertain log-normal coefficient, we show that Approach 1 provides a more accurate approximation of the conditional log-normal coefficient and solution of the SPDE than Approach 2 for the same number of random dimensions in a conditional KL expansion. Furthermore, Approach 2 provides a good estimate for the number of terms of the truncated KL expansion of the conditional field of Approach 1.

Finally, we demonstrate that the proposed active learning method (Method 2) is more efficient for reducing uncertainty in the solution of the SPDE under consideration (i.e., it leads to a larger reduction of the variance) than the standard active learning method (Method 1). The difference between two methods increases as the variance of the partially known coefficient increases.