# Deep Learning Joint Inversion of Electrical Data for Ahead-Prospecting in Tunneling

Authors: Peng Jiang; Benchao Liu; Chuanwu Wang; Lei Chen; Yuting Tang
Journal: Advances in Civil Engineering
 (2023)

---

## Abstract

Water inrush has become one of the bottlenecks restricting tunnel construction. Among various advanced forecasting techniques, the direct current method is more cost-effective and sensitive to water-bearing structures. It has been widely used in exploring water inrush disasters in practical engineering. Although traditional resistivity linear inversion methods are reasonably practical, they usually suffer from volume effects and cannot accurately locate the location and morphology of water-bearing bodies. Therefore, nonlinear techniques such as deep learning have recently become popular to directly approximate the inversion function by learning the mapping of apparent resistivity data to the geoelectric model. This work presents a novel deep learning-based electrical approach that combines resistivity and polarizability to estimate water-bearing location and morphology. Specifically, we design an encoder-decoder network. A shared encoder extracts features from the input data, two encoders output resistivity, and polarizability models, respectively, and fine-tuned collinear regularization for both outputs reduces solutions’ multiplicity. Compared with traditional linear inversion methods and independent parameter inversion, our proposed joint inversion method shows superiority in locating and delineating anomalous bodies.


---

## Body

## 1. Introduction

Continued growth in the world’s population and economy has prompted the construction of tunnels worldwide to provide better energy for transportation, water, and electricity. However, many abnormal geological structures, such as fault zones, faults, aquifers, and karsts, may be located in front of the tunnel face during the excavation process. If the hazard source detection is not clear, geological disasters may occur due to excavation disturbances, resulting in delays in construction schedules, damage to instruments, economic losses, and heavy casualties. Ahead-prospecting has become a mandatory procedure in tunnel excavation, providing essential prior information to ensure safe, economical, and efficient tunnel excavation [1]. Over the past 50 years, many advanced exploration methods have been developed. Among them, geophysical exploration is more time-saving and cost-effective. This study focuses on electrical methods that inverse the observed electrical properties (resistivity and polarizability) to help estimate a water-bearing structure and water volume ahead tunnel face, avoiding water inrush hazards. Traditional linear inversion methods are sensitive to the initial model and easily trapped into local optimization. Therefore, the final inversion results are far from precise. Many works follow this way to solve electrical resistivity inversion for hydrogeological exploration [2, 3], environmental [4–6] and engineering investigations [7–10], and mineral prospecting [11, 12]. Although resistivity inversion has been widely used, traditional linear inversion algorithms for electrical property inversion still face some problems because the task is nonlinear and ill-posed [13, 14]. Rijo [15] proposed the Marquette fast inversion method for simple interpretation of 3D resistivity and polarizability data in the study of polarizability methods. Oldenburg and Li [16] proposed a 3D polarizability inversion algorithm using the linear equation of the excitation response, i.e., in the case of inversion to obtain the resistivity distribution, and the second step of inversion is performed to obtain the polarizability. However, this stepwise inversion method leads to good or lousy polarizability results depending on the resistivity inversion results. Li et al. [17] proposed a synchronous resistivity and polarizability inversion algorithm based on cross-gradient constraints to solve the stepwise inversion problem and improve the inversion accuracy in response to the problem that the polarizability results depend on the resistivity inversion results that exist in stepwise inversion. Geophysical inversion methods can infer geological structures using various physical parameters such as velocity, density, resistivity, and polarization rate. Because different geophysical exploration methods contain complementary information, joint inversion can improve the accuracy of geological structure imaging under complex geological conditions. Meju [18] proposed the cross-gradient operator. They implemented a joint inversion of seismic travel time data and resistivity data based on cross-gradient, which requires different physical parameters to vary in the same method or opposite directions.In recent years, nonlinear inversion methods represented by artificial neural networks have attracted widespread attention. They have gradually improved and evolved into deep neural networks (DNNs), which have significant nonlinear mapping capabilities and have been used to solve inversion problems [19–22]. Deep learning methods to capture global feature information and further introduce adversarial learning to extract invariant feature information, global features, and shared feature information can be fully utilized [23]. Yi et al. [24] achieved 85% accuracy in a benign/malignant classification task of breast tumours using a convolutional neural network approach, which was 14% more accurate than the manual feature-based approach at the time. Jin et al. [25] solved the ill-posed inverse problem in medical imaging by removing artefacts from the direct inversion results using CNNs. This brings a new solution to the ill-posed nonlinear problem of resistivity inversion. Araya-Polo et al. [26] used velocity models generated from raw seismic data as input to a convolutional neural network to achieve inverse imaging of subsurface geological wave velocities by laminar imaging, laying the foundation for the application of deep learning in the geophysical field. Yang and Ma [27] investigated a deep fully convolutional neural network (FCN)-based velocity model inversion method for effective velocity inversion of the salt mound model using multigun prestack seismic data as input and tested it on SEG public data. Wu and Li [28] constructed a convolutional neural network with an encoder-decoder structure to learn the correspondence between seismic data and subsurface. The correspondence of a velocity structure significantly improves the accuracy of velocity inversion of laminar and fault models. In addition, Li et al. [29] proposed an end-to-end seismic inversion network (SeisInvNet) to achieve the direct reconstruction of velocity models from seismic data by SeisInvNet. The network solved the problem of weak data spatial relativity through operations such as adjacent gunpoint data enhancement, global feature generation, and multichannel data decoding and achieved accurate prediction of laminar models with good generalization. Liu et al. [30] successfully applied a deep neural network-based seismic data in a tunnel observation environment based on the SeisInvNet inversion method to achieve geological inversion imaging in front of the tunnel face.Deep learning techniques were initially used in seismic detection, but many geophysicists are now attempting to apply them to address geophysical inversion problems. Resistivity inversion has a particular property different from normal image mapping; that is, the resistivity data have vertical variation characteristics. Therefore, traditional convolutional neural networks (CNNs) do not work well and can create ambiguities when building maps directly. To solve this problem, Liu designed a resistivity inversion network (ERSInvNet) with U-Net as the backbone, added the vertical position information of the resistivity data to the input, and added depth weighting constraints to the loss function, forming a high-quality inversion imaging method of resistivity based on deep learning. Instead of supplementing the features [31], Liu addresses this problem by modifying the convolution kernels to adapt them to the applied vertical position, further improving the inversion imaging quality of complex geological models [32]. Rijo [15] proposed the Marquette fast inversion method for simple interpretation of 3D resistivity and polarizability data in the study of polarizability methods. Oldenburg and Li [16] proposed a 3D polarizability inversion algorithm using the linear equation of the excitation response, i.e., in the case of inversion to obtain the resistivity distribution, the second step of inversion is performed to obtain the polarizability. However, this stepwise inversion method leads to excellent or lousy polarizability results depending on the resistivity inversion results.Geophysical inversion methods can infer geological structures using various physical parameters such as velocity, density, resistivity, and polarization rate. Because different geophysical exploration methods contain complementary information, joint inversion can improve the accuracy of geological structure imaging under complex geological conditions. Usually, for the same subsurface medium, the physical parameters of different geophysical methods have consistent spatial variation, and based on this feature, scholars have carried out research on joint inversion methods based on structural similarity, which define a metric based on structural similarity, and then continuously reduce the structural gap between the inversion results of different methods through repeated iterations, to enhance the structural similarity of the model. Jie and Morgan [33] and Haber and Oldenburg [34] constructed a structural similarity metric by using Laplace operators. This method reconstructs seismic velocity and resistivity simultaneously under the cross constraint of standardized structural curvature. An objective function is defined to quantify the structural difference between the two models, and then the objective function is minimized to solve the joint inversion under the condition of satisfying data constraints. Meju [18] proposed the cross-gradient operator. They implemented the joint inversion of seismic time-lapse data and resistivity data based on cross-gradient, which requires different physical parameters to vary in the same method or opposite directions. Moorkamp et al. [35] applied the joint cross-gradient inversion method in the in-well detection environment and obtained good results. When calculating the cross-gradient term, the discretization methods of different grids need to be considered. In 2020, Yan et al. [36] carried out a study on the 3D joint inversion of gravity, magnetic, and geomagnetic methods, divided the cross-gradient function grids into four categories, and performed discretization operations on the four types of grids, respectively; however, this discretization method has a large workload.Although the linear inversion method of an electrical method usually has the advantage of fast iteration speed, there are problems such as relying on the initial model and easily falling into local optimization. The idea of joint inversion has improved the inversion results, but it has not fundamentally solved the shortcomings of traditional linear inversion methods when dealing with nonlinear problems. In contrast, single-parameter deep learning algorithms are more sensitive to data perturbations and suffer from overfitting. In this work, in addition to the resistivity inversion, we also performed a joint resistivity and polarizability inversion for better accuracy and less solution multiplicity. Correspondingly, we design an encoder-decoder network with a shared encoder to extract features from the inputs of two channels (one channel for resistivity inversion results and one for polarizability inversion results) and the two encoders output resistivity and polarizability models, respectively, and a detailed collinear regularization is performed on both outputs to reduce solution multiplicity. In addition, to increase the computational rate, a two-channel parameter parallel computation method is implemented with the same training time as a single parameter. In order to further verify the stability of multiparameter inversion, the experimental noise tests are carried out in this study, and different tunnel resistivity models have good robustness and adaptability to noisy data.


## 2. Method

Moreover, the decoder takes that the prerequisite for applying the excitation polarization method is the difference in the excitation effect between the target geological body and the surrounding rock [37]. If it is assumed that in the absence of excitation effect, the resistivity value in front of the tunnel is ρ, the apparent resistivity forward response dr of the model can be expressed as follows:(1)dr=Fρ,where F denotes the forward operator. When there is an induced polarization effect in front of the tunnel, the polarizability is η, and the equivalent apparent resistivity d and apparent polarizability dp are as follows:(2)d=Fρ1−ηdp=d−drd=Fρ/1−η−FρFρ/1−η.As shown in Figure1, our task is to invert the resistivity model mr and polarizability model mp given dr and dp.Figure 1
Resistivity and polarizability inversion. (a) Resistivity data, (b) resistivity model, (c) polarizability data, and (d) polarizability model.The parameter-independent (resistivity or polarizability) linear inversion method uses iterative optimization to update the initial model to approximate the target model. Given, the electrical datad observed at step t and the model m, the model (mt+1) at step t+1 obtained by conventional linear inversion methods can be given by(3)dt⟵Fmt,mt+1⟵mt−ψ∂Ldt,d∂m,where L is the loss function, usually defined as the mismatch between synthetic data dt and observed data d. The synthetic data dt is obtained by forwarding modelling mt, and we usually set the initial model m0 as the unified background. ψ is usually an L-BFGS optimizer. A joint parameter inversion method is proposed to solve the ill-posed problem of electrical property inversion and improve the inversion accuracy. The main difference from independent parameter inversion is the loss function design; joint inversion has an additional cross-gradient term [38], defined as(4)Lr=12drt−dr‖2+λ2∇mr×∇mp‖2,Lp=12dpt−dp‖2+λ2∇mp×∇mr‖2,where ∇mr and ∇mp represent gradients related to mr and mp, respectively. λ is the weight for cross-gradient term.


### 2.1. Database Establishment and Neural Network Architecture

In this work, about 10,000 sets of numerical simulation data bureaus were designed. The complexity of the data set is mainly reflected in two aspects. There are three types of resistivity anomalies, including faults, single caves, and two caves. The background resistivity is set to 1000 Ωm, and the resistivity range of anomalies is 30–200 Ωm. We randomly divide the data set into a training set, a validation set, and a test set according to a ratio of 10 : 1 : 1.To achieve joint tunnel resistivity/polarizability inversion, our gradient-optimized deep neural network should also receive both resistivity and polarizability types of data and output both target inversion results. A shared encoder deep neural network based on multitask learning [39–41] is used to accomplish the task of fusing the two types of data. Our design uses five fully connected modules to form the encoder, extract data features, upscale the data, and use five conv-up blocks to form each decoder. Each conv-down block has two 3D convolutional layers and one max pooling layer. We use ReLU as an activation function, as shown in Figure 2.Figure 2
Neural network architecture.The decoder extracts the high-level semantic information features of the data and obtains global information by compressing the input data size and expanding the number of data channels. In this study, four conv-down blocks are used to form the encoder. Each conv-down block contains two convolution operations and one maximum pooling operation, where each convolution operation contains a 3D convolution operation, a batch normalization operation, and a ReLU activation function. Moreover, the decoder will upsample the global information obtained by the encoder to finally get the prediction result with the same size as the original input data. The lower local information obtained from the encoding layer is introduced into the decoding process through short-circuit connection. This method can aggregate global and local information to ensure that the decoding result does not lose detailed information. The decoder consists of four conv-up blocks, each of which first upsamples the features through a 3D transposed convolutional layer and then performs two additional common 3D convolutional operations through the convolutional layer. The specific network design is shown in Table1.Table 1
Setting for each convolutional layer module.
EncoderTypeFilterStrideDecoderTypeFilterStrideconv_down1Conv3 × 3 × 31conv_up1Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31conv_down2Max pool2 × 2 × 22conv_up2Up-conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31conv_down3Max pool2 × 2 × 22conv_up3Up-conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31conv_down4Max pool2 × 2 × 22conv_up4Up-conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Up-conv3 × 3 × 31Conv3 × 3 × 31


## 2.1. Database Establishment and Neural Network Architecture

In this work, about 10,000 sets of numerical simulation data bureaus were designed. The complexity of the data set is mainly reflected in two aspects. There are three types of resistivity anomalies, including faults, single caves, and two caves. The background resistivity is set to 1000 Ωm, and the resistivity range of anomalies is 30–200 Ωm. We randomly divide the data set into a training set, a validation set, and a test set according to a ratio of 10 : 1 : 1.To achieve joint tunnel resistivity/polarizability inversion, our gradient-optimized deep neural network should also receive both resistivity and polarizability types of data and output both target inversion results. A shared encoder deep neural network based on multitask learning [39–41] is used to accomplish the task of fusing the two types of data. Our design uses five fully connected modules to form the encoder, extract data features, upscale the data, and use five conv-up blocks to form each decoder. Each conv-down block has two 3D convolutional layers and one max pooling layer. We use ReLU as an activation function, as shown in Figure 2.Figure 2
Neural network architecture.The decoder extracts the high-level semantic information features of the data and obtains global information by compressing the input data size and expanding the number of data channels. In this study, four conv-down blocks are used to form the encoder. Each conv-down block contains two convolution operations and one maximum pooling operation, where each convolution operation contains a 3D convolution operation, a batch normalization operation, and a ReLU activation function. Moreover, the decoder will upsample the global information obtained by the encoder to finally get the prediction result with the same size as the original input data. The lower local information obtained from the encoding layer is introduced into the decoding process through short-circuit connection. This method can aggregate global and local information to ensure that the decoding result does not lose detailed information. The decoder consists of four conv-up blocks, each of which first upsamples the features through a 3D transposed convolutional layer and then performs two additional common 3D convolutional operations through the convolutional layer. The specific network design is shown in Table1.Table 1
Setting for each convolutional layer module.
EncoderTypeFilterStrideDecoderTypeFilterStrideconv_down1Conv3 × 3 × 31conv_up1Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31conv_down2Max pool2 × 2 × 22conv_up2Up-conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31conv_down3Max pool2 × 2 × 22conv_up3Up-conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31conv_down4Max pool2 × 2 × 22conv_up4Up-conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Conv3 × 3 × 31Up-conv3 × 3 × 31Conv3 × 3 × 31


## 3. Inversion Results

To demonstrate the performance of our proposed method, comparative experiments with traditional inversion methods, deep learning independent parameter inversion methods, and deep learning joint inversion methods are carried out. By minimizing the error in the observed data and the physical model parameters with smooth excess, the least-squares inversion objective function based on the smooth constraint [42] is constructed as follows:(5)Φ=ϕd+λϕm=e2+λ∂m∂x2+∂m∂y2+∂m∂z2=∆d−A∆mT∆d−A∆m+λC∆mTC∆m,where C is the smoothness matrix and λ is the Lagrangian daily number. From equation (5), it can be seen that the least-squares inversion objective function based on the smoothness constraint consists of two components, ϕd denoting the variance between the actual observed data and the theoretical orthogonal observed data, ϕm denoting the degree of difference between the adjacent grid resistivities, and λ denoting the weights of ϕd and ϕm. When the value of the objective function Φ is minimal, the model parameters under the smooth constraint make the minimum fit of the actual observed data to the theoretical observed data.Figure3 is the inversion result of the tunnel filling with water-bearing faults, Figure 3(a) is the geoelectric model of the tunnel, and y is the tunneling direction. Figure 3(b) shows the traditional linear inversion results, using the least-squares inversion method. Figure 3(c) shows the deep learning inversion results based on a single attribute parameter of tunnel resistivity. Figure 3(d) shows the joint inversion results of deep learning based on the two parameters of tunnel resistivity and polarizability. This study only shows the joint inversion results of tunnel resistivity. It can be seen from Figure 3 that the traditional linear inversion results can invert the existence of a low-resistance anomaly in front of the tunnel. Still, the shape and outline of the anomaly cannot be determined, and the difference in resistance values is significant. The inversion results of a single attribute of deep learning can accurately locate the fault location, but the whole is negligible. Compared with the inversion results of deep learning resistivity, the combined inversion results of deep learning can more accurately invert the location and volume shape of the water-bearing fault in front of the tunnel with the resistance value.Figure 3
Inversion results of faults. (a) Resistivity model, (b) traditional inversion results, (c) deep learning inversion results, and (d) deep learning joint inversion results.
(a)(b)(c)(d)To further study the universality of convolutional networks, we designed single caves. Compared with filling fault, the inversion of karst cave is more difficult because of its small volume. It can be seen from Figure4 that the general location of the karst cave in front of the tunnel can be found in the traditional linear inversion results. However, due to the volume effect, the volume and resistance of the inversion results are large, and the boundary description is not obvious, which is easy to affect the accurate identification of the disaster-causing structure in front of the tunnel. The deep learning inversion results are better than the traditional linear inversion results in overall performance, similar to the fault inversion results. The volume and shape of the inversion results of single resistivity parameters of deep learning are small, and the joint inversion results of deep learning can solve this problem.Figure 4
Inversion results of single caves. (a) Resistivity model, (b) traditional inversion results, (c) deep learning inversion results, and (d) deep learning joint inversion results.
(a)(b)(c)(d)In Figure5, we show the loss curves of the inversion results on the training and validation sets. It can be seen from the downward trend of the two curves that there is no overfitting during the training process. After 400 epoch trainings, the error dropped from 0.09 to less than 0.04, and the curve gradually converged.Figure 5
Loss curves of convolutional neural network on training and validation sets. (a) Train loss and (b) valid loss.
(a)(b)In this study, we quantitatively evaluate the inverse effect of large-scale arithmetic cases in the validation set and test set by using the weighted mean square error (WMSE) and weighted correlation coefficient (WR) [43]. The specific formulas for the weighted mean square error (WMSE) and weighted correlation coefficient (WR) are as follows:(6)WMSE=1N∑n=1Nωnm^n−mnTωnm^n−mn,WR=1N∑n=1N∑n=1Nωnm^n−m^¯nTωnm^n−m¯nωnm^n−m^¯n2ωnm^n−m¯n2.A quantitative comparison of WMSE and WR for the joint inversion method and the resistivity single-parameter inversion method for the test and validation sets is given in Table2. As seen from the table, the inversion results of the joint deep learning method have smaller WMSE values and larger WR values and show the same pattern in both the test and validation sets, which further demonstrates the superiority of the proposed method in this study.Table 2
Quantitative comparison of inversion results on test and validation sets.
TestValidationWMSE↓WR↑WMSE↓WR↑Deep learning joint inversion results0.0320900.9833120.0122100.98934Deep learning inversion results0.2102560.9042150.1803200.91230↑ indicates that the larger the value achieved, the better the performance is, while ↓ indicates the smaller the index value, the better the inversion result.


### 3.1. Noise Disturbance Test

During the numerical simulation, the background resistivity is relatively uniform for both the tunnel water-bearing cavern model and the filled fault model. However, in actual engineering probing, it is impossible to measure the apparent resistivity of a uniform stratigraphic background due to noise factors. Therefore, noise tests were performed on the input data. A Gaussian white (WGN) noise of 5 dBw was added to the synthetic data to investigate how the noise affects the output during training. The tunnel resistivity model is used as an example to show the tunnel resistivity inversion results under the influence of noise. The inversion results of the resistivity model, the deep learning method with individual parameter inversion, and the joint inversion method with 5 dBw Gaussian white noise are shown in Figure6 from left to right in the figure.Figure 6
Resistivity model and inversion results with 5 dBw noisy synthetic data. (a) Resistivity model, (b) deep learning inversion results, and (c) deep learning joint inversion results.Analyzing the inversion results of the two methods, we can find that the joint deep learning inversion method is superior to the single-parameter deep learning inversion method in terms of boundary inscription and resistance value degree for both faults and caves. In particular, for the geological model of two water-bearing caves, nearby, the joint inversion method can still identify the two water-bearing caves and trap the location distribution more accurately under noise interference.


## 3.1. Noise Disturbance Test

During the numerical simulation, the background resistivity is relatively uniform for both the tunnel water-bearing cavern model and the filled fault model. However, in actual engineering probing, it is impossible to measure the apparent resistivity of a uniform stratigraphic background due to noise factors. Therefore, noise tests were performed on the input data. A Gaussian white (WGN) noise of 5 dBw was added to the synthetic data to investigate how the noise affects the output during training. The tunnel resistivity model is used as an example to show the tunnel resistivity inversion results under the influence of noise. The inversion results of the resistivity model, the deep learning method with individual parameter inversion, and the joint inversion method with 5 dBw Gaussian white noise are shown in Figure6 from left to right in the figure.Figure 6
Resistivity model and inversion results with 5 dBw noisy synthetic data. (a) Resistivity model, (b) deep learning inversion results, and (c) deep learning joint inversion results.Analyzing the inversion results of the two methods, we can find that the joint deep learning inversion method is superior to the single-parameter deep learning inversion method in terms of boundary inscription and resistance value degree for both faults and caves. In particular, for the geological model of two water-bearing caves, nearby, the joint inversion method can still identify the two water-bearing caves and trap the location distribution more accurately under noise interference.


## 4. Conclusions

The nonlinear problem of tunnel resistivity inversion has always been a research hotspot. Deep learning is an essential means to achieve complex mapping between data because it can summarize laws based on massive data. The inversion process of inferring a geological model by analyzing the characteristics of apparent resistivity data can be regarded as a mapping between images. To further improve the robustness and generalization of deep learning algorithms and increase data constraint information, it is necessary to integrate joint inversion into deep learning methods.To verify the universality of this method, a fault and karst cave model is established based on the standard disaster-causing water body in front of the tunnel, and a comparative test of different inversion methods is carried out. It can be seen from the inversion results that the deep learning inversion results are far superior to the traditional linear inversion results, and the imaging quality of the deep learning joint inversion is also better than the inversion quality of a single resistivity parameter. In the future, we will consider integrating the physical laws of the tunnel electric field into deep learning to improve the generalization of deep learning inversion results and apply them in tunnel engineering practice.



---

*Source: 5639207-2023-03-15.xml*