4.2 Regularization in Electrical Imaging Inversion

Most tomographic problems in geophysics are solved with an excess number of model parameters and use regularization to create a mathematically stable solution (e.g., Constable et al., 1987). Due to these issues, this non-linear problem is solved using iterative inversion (Tripp et al., 1984; Daily and Owen, 1991). The solution to the ER inverse problem is a cross section or volume of electrical conductivity values and is usually based on non-linear least-squares minimization of a two-part objective function, F, which is given by Equation 9a.

[latex]\displaystyle F=\left\| C{_{D}}^{-0.5}\left ( d_{obs}-g(m) \right )\right\|^{2}+\varepsilon \left\| Dm\right\|^{2}[/latex] (9a)

This can also be written as Equation 9b.

[latex]\displaystyle F=\left ( d_{obs}-g(m) \right )^{\textrm{T}}C{_{D}}^{-1}\left ( d_{obs}-g(m) \right )+\varepsilon m^{\textrm{T}}D^{\textrm{T}}Dm[/latex] (9b)

where (for the case of ER inversion alone):

dobs = vector of electrical resistance (V/I) or apparent resistivity measurements, where the vector length is equal to the total number of quadripoles.
g( ) = forward model for electric potential (Equation 1 or 5), which produces the simulated estimates of the measurements.
m = vector of parameter estimates, (log electrical conductivity).
CD = covariance matrix of data uncertainty or errors (often where the diagonal of the matrix is defined by reciprocal or stacked measurements, and the off-diagonal values are zero), which defines how certain we are in our measurements. This matrix often includes some measure of modeling errors, as described in Section 4.3 below.
ε = regularization parameter that determines the importance given to the smooth appearance of the electrical conductivity field relative to the misfit between calculated and observed resistances. An overly small ε will minimize the residual error between measured and modeled resistances and may overfit the data, producing spurious heterogeneity in the solution. In contrast, an overly large ε will identify an overly smooth electrical conductivity field that may not fit the measured field data (resistances) well (see, for example, Tikhonov and Arsenin, 1977).
D = the model-weighting regularization matrix, which can either be defined by a discretized spatial-derivative operator or be based on the covariance of the model parameters (Tarantola, 1987; Gouveia and Scales, 1997; Kitanidis, 1997; Vasco et al., 1997; Day-Lewis et al., 2003). This matrix defines how each pixel is related to one another.

The first part of the objective function is the data misfit term, which minimizes the discrepancy between field resistance data (or resistance and IP data for IP inversion) and the computed resistances (or resistances and IP measurements for IP inversion) based on Equation 1 or Equation 5, within measurement errors. The second part is the regularization term, often called the model roughness if smoothing is sought, which typically minimizes the roughness (or maximizes the smoothing) of the electrical conductivity field and allows for a well-posed inverse problem. This term is required due to the overparameterization of the inverse problem, meaning that the information provided by the measurements cannot uniquely resolve each of the conductivity parameters. It is possible to be creative with the regularization term, depending on prior information available to develop a conceptual model (e.g., Caterina et al., 2014; Nguyen et al., 2016; Hermans et al., 2016).

In a non-linear inverse problem, model parameters are updated iteratively by repeated solution of a linearized system of equations for Δm at successive iterations. Such an approach results in the regularization changing throughout the iterative process. This process makes it difficult to map the effect of regularization throughout the inversion process, and consequently impairs quantitative inference from the images. The update appears as shown in Equation 10.

[latex]\displaystyle \left [ J^{\textrm{T}}C{_{D}}^{-1}J+\varepsilon D^{\textrm{T}}D \right ]\Delta m=J^{\textrm{T}}C{_{D}}^{-1}\left ( d_{obs}-g(m_{k-1}) \right )-\varepsilon D^{\textrm{T}}Dm_{k-1}[/latex] (10a)
[latex]\displaystyle m_{k}=m_{k-1}+\Delta m[/latex] (10b)


J = Jacobian matrix at iteration k, with elements Jij = ∂dsim,i/∂mj
mk = vector of parameter estimates after updating in iteration k
Δm = vector of parameter updates for iteration k

Although tomographic inversion with regularization is useful for imaging large-scale (low-spatial-frequency) structures, it yields poor results when attempting to infer quantitative values from the recovered images (e.g., Binley et al., 2002; Singha and Gorelick, 2005; Day-Lewis et al., 2007) due to uncertainty in the inversions for reasons discussed in the next section.

3-D acquisition and inversion are increasingly possible and appropriate, although many practitioners still use 2-D inversion. As noted earlier, commercially available software for 2-D inversion commonly invokes the 2.5-D assumption for computational efficiency, where 2-D ER cross sections are constructed by simulating 3-D current flow under the assumption of 2-D heterogeneity (Dey and Morrison, 1979; LaBrecque et al., 1996). The 2.5-D assumption amounts to assuming that all heterogeneity in the 2-D imaged cross section extends infinitely out the 2-D plane. Clearly, this assumption is violated in the presence of strong 3-D heterogeneity, which requires 3-D acquisition and inversion approaches to image accurately.


Electrical Imaging for Hydrogeology Copyright © 2022 by Kamini Singha, Timothy C. Johnson, Frederick D. Day Lewis and Lee D. Slater. All Rights Reserved.