Ground condition assessment based on Bayesian framework

Minimizing uncertainties is an important issue among the significant discussions pertaining to project design and planning. Usually, the uncertainties in subsurface projects arise from the unknown ground conditions that may cause the designer to fail to consider all potential issues prone to occur during the construction procedure. This paper presents a Bayesian framework to predict the ground conditions before excavation. The state transition probability is calculated by the difference between the simulation results and the observation data so that the model parameters gradually tend to a posterior distribution. Through numerical simulation, the algorithm is verified, and error analysis is carried out. Within the structure, the resistivity is mainly controlled by the structural characteristics, so the formation lithology has great uncertainty in the geological structure area. The inversion results of fracture degree and water content are accurate.


Introduction
The geological and geotechnical conditions are vital for engineering projects in the subsurface. At the preliminary investigation stage, boreholes are often used to detect the physical properties of underground media. However, the spatial coverage of boreholes is limited, and there are great difficulties and risks in the construction of boreholes in areas with complex terrain. As an indirect exploration method with large coverage and deep exploration depth, Geophysics can effectively make up for the above shortcomings [1,2]. The geophysical method is based on the difference of physical properties such as gravity [3], resistivity [4], and damping of the rock [5] to refer to its properties (lithology, fracture, water content, etc.). However, because the observed data are the physical properties of the rock mass, the obtained geological profile is inevitably full of multi-solution and uncertainty. Spatial uncertainty in geology and geotechnical engineering has become an increasingly popular topic of study in recent years [6]. Minimizing uncertainties is important for the prediction of factors that contribute to the success of a project [7].
Tunneling, similar, but more than other geotechnical endeavors, is characterized by the influence of uncertainty. Many algorithms have been proposed to determine subsurface properties, including lithology and geotechnical characteristics. Hui proposed a Monte Carlo simulation-based framework to quantify the uncertainties in the Q-system of rock mass classification in underground construction [8]. Buoyed proposed a Sequential Indicator Co-simulation for incorporating a geologist's interpreted cross-section to quantify spatial uncertainty of lithology [9]. Mahmoodzadeh [10] presents an updating procedure to refine predictions of ground conditions and time-cost scatter-gram during tunnel construction. The distinguishing characteristic of a Markov process is a single-step memory: history apart from the most recent event is neglected informing predictions [11]. A significant advantage of using single-step memory instead of multiple-step memory is that probability calculations are considerably simpler, and a full probability distribution can often be found. The Markov process is suitable for tunnel geology prediction models because it involves not only uncertain ground parameters along the tunnel line but also their locations [12,13]. Therefore, the concept of a stochastic process can be used.
This paper presents a probabilistic inversion method to determine the lithology, fracture and water content of underground space from the borehole log and high-density resistivity method (HDRE) data. We use herein a two-dimensional finite-difference time-domain (FDTD) model to simulate the apparent resistivity for a 2D numerical model. The posterior distributions of lithology, fracture, and water content are determined by fitting the FDTD model against the observed HDRE data and borehole information using MCMC simulation with the DREAM(ZS) method [14].

Methodology
In this section, we describe the methodology used herein to determine the distribution of lithology, fracture, and water content in the underground space using stochastic inversion. We first describe the framework of Bayesian inversion methodology., then discuss the DREAM(ZS) algorithm, followed by a description of the FDTD forward simulator.

Bayesian inversion framework
The relationship between the parameters and observed data could be described with the following equation: (1) Where  d is a vector with observed data (apparent resistivity and borehole log in this paper), ( ) f  describes (simulates) the physical relation between the model parameters,( m ) and  d . The u  represents the input data and e signifies a vector of system error. The measured data consists of the borehole log and measured apparent resistivity. The model parameters in this study define the twodimensional distribution of resistivity of the underground space, and the model input is the location of the borehole and electrode arrangement. The tilde operator~is used to denote measured entities.
The model parameters can generally not be derived from closed-form analytic solutions, and inverse methods are required to determine their values from the measured data. In the Bayesian framework, model parameters are represented by probability density functions whose posterior distribution ( ) p m  d can be derived from the data using Bayes' law [15]: denotes the likelihood function, and ( ) p d  is a normalizing constant (probability density of d  ) that ensures that the posterior parameter distribution integrates to unity. We can discard probability density of d  or ( ) p d  from Eq.(2) as all our inferences about the model parameters can be made from the unnormalized density The prior probability density, ( ) p m determined by our knowledge about the model parameters before getting the measured data d  . In most investigations, the spatial distribution of the apparent resistivity of the medium under consideration is unknown, and a uniform prior distribution is used.
The main problem now resides in the formulation of the likelihood function, ( ) L m d  used to summarize, in probabilistic terms, the level of agreement between the observed data and the FDTD The SSE term measures the difference between the observed and simulation result and thus evaluates the performance of the model for given parameter values. The smaller the value of the SSE, the larger the value of the likelihood function, and thus the better the FDTD model fits the observed data.
Once the prior distribution and likelihood function has been defined, what is left in Bayesian analysis is to summarize the posterior distribution. Confidence intervals construed from a classical first-order approximation around the values m derived from deterministic inversion methods can then only provide an approximate estimate of the posterior distribution. We, therefore, resort to MCMC simulation with the DREAM(ZS) algorithm to generate a sample of the posterior distribution.

Markov chain Monte Carlo simulation: DREAM(ZS)
The basis of MCMC simulation is a Markov chain that generates a random walk through the search space and successively visits solutions with stable frequencies stemming from a stationary distribution. To explore the target distribution, an MCMC algorithm generates trial moves from the current state of the Markov chain mt-1 to a new state mp. The earliest and most general MCMC approach is the random walk Metropolis (RWM) algorithm. Assuming that a random walk has already sampled points {m0,….,mt-1} this algorithm proceeds in the following three steps. First, a candidate point mp is sampled from a proposal distribution q(·) that is symmetric, q(mt-1,mp) = q(mp,mt-1) and may depend on the present location mt-1. Next, the candidate point is either accepted or rejected using the Metropolis acceptance probability. 1 1 (m ) (m m ) min 1, (m ) Where p(·) denotes the density of the target distribution, if the proposal is accepted, the chain moves to mt = mp. Otherwise, the chain remains at its current location, mt = mt-1.
The DREAM(zs) algorithm is an MCMC algorithm that has its roots within the DREAM and DE-MC algorithm, which automatically tunes scale and orientation. In DREAM(ZS), K different Markov chains are run simultaneously in parallel, therefore leads to high sampling efficiencies and rapid convergence.

HDRE forward simulation
As shown in Fig.1, the resistivity is obtained by two supply electrodes and two measurement electrodes. Arrange the four electrodes A, B, M, and N in an arbitrary arrangement, where A and B are power supply electrodes. The current starts from point a flows through the underground medium and returns to pole B. at the same time, observe the potential difference between M and N. at this time, the potential difference between M and N is Where I am the excited current intensity, the apparent resistivity(ρs) is the comprehensive value of the resistivity of the geological body on the current path.  Figure.1 The geologic model The electric current in the medium is governed by the well-known potential equation that describes current characteristics and potential field distribution under stable power. In the HDRE surveys, the supply electrode and measuring electrode are measured point by point according to a certain arrangement. When the supply electrode supplies a steady current to the ground, the potential underground field can be expressed by the following equation： The analytical solution of Eq.6 is very difficult for most practical situations, and we, therefore, resort to 2D numerical simulation using the FDTD simulator. This simulator discretizes the partial derivatives of Eq.6 in space using central differencing. The numerical solution is then obtained by solving in a leapfrog manner the resulting finite-difference equations, and this provides a simulated potential field for geometry, physical parameters, boundary conditions, time step, excitation, and electrode positions in the HDRE survey. The values of U at each position are subsequently derived from the FDTD model output.

Simulation and analysis
In this section, we establish a finite difference, numerical model. The MCMC simulation with DREAM(ZS) was used to estimate the posterior distribution of the physical parameters underground.

Markov chain Monte Carlo simulation: DREAM(ZS)
The numerical model is shown in Fig.2, which includes four different lithologies: interbedding of limestone and shale(ρ1=1317 Ω·m), shale(ρ2=743 Ω·m), limestone(ρ3=2145 Ω·m), and interbedding of sandy shale and limestone(ρ4=1000 Ω·m). The fracture and water content of background and geological structure are shown in Tab.1. Consequently, the resistivity of each pixel could be calculated by [16]: Where the  and w S are porosity and water content respectively, geological boreholes are arranged at 100m, 800m, and 1500m, respectively.  Figure.2 The geologic model , the electrode distance is 10m, and the detection depth is 200m. The forward result of apparent resistivity is shown in the Fig.3.

Result and analysis
The simulation results are used as the observation data to inverse the lithology, fracture degree, and water content of each 5m 5m  pixel. Fig.4 depicts the convergence diagnosis curve. Most parameters meet the condition of R<1.2 after 6e5 times of sampling process. These parameters converge to reliable posterior distributions. In order to obtain more posterior samples and make the form of posterior distribution more obvious, the Markov stochastic process is continued after stabilization.  Figure.4 Parameter convergence curve The inversion results of lithologic distribution are shown in Fig.5. The inversion results can basically accurately reflect the lithologic distribution area, but there are still many solutions in the lithologic interface and where the structure exists. This is because the comprehensive resistivity value in these areas is mainly caused by geological structure, so there is great uncertainty in lithologic inversion within these areas. The results show that the fracture distribution is similar to the numerical simulation model, and the crushing area can be seen at the location of the karst cave and water-bearing structure. However, there are a few abnormal areas between different lithologic interfaces. Fig.6(b) depicts the distribution of water content. The position of the water-bearing structure shows the abnormal area of great water content, the position of fault and cavity shows the abnormal area of low water content. Furthermore, the two edges of the area turn into opposite characteristics, resulting in local great water content on both sides of the cavity. There is some interference on the boundary of the model.  Figure.6 The inversion result By analyzing the inversion results of a single pixel, Fig.7(a) and Fig.7(b) show the lithology inversion results inside and outside the geological structure, respectively. It can be seen that the probability of lithology has more accurate inversion results outside the geological structure, and the lithology inversion results inside the geological structure are affected by other factors and have greater uncertainty. Fig.7(c) and Fig.7(d) are the inversion results of water content and fracture within waterbearing structures, respectively. It can be seen that they basically conform to the rule of normal distribution, and the expected value is close to the real value.   Figure.7 Posterior distributions

Conclusions
This paper introduces the different elements of a Bayesian inversion method for detection in twodimensional underground space. This framework uses the main building blocks of the twodimensional FDTD simulator and the DREAM(ZS) algorithm to reconstruct the lithology, fracture, and water content from HDRE inversion. The FDTD simulator solves numerically current equations and simulates the apparent resistivity iteratively for a given geometry, lithology, fracture, and other input data. The DREAM(ZS) algorithm is used to estimate the posterior distribution. The following conclusions can be obtained: inversion results can well identify the interface between different lithology, and within the geological structure, the uncertainty is inevitable; The inversion of fracture degree and water content is basically consistent with the forward model, and there is a little interference at the boundary of the model.