- About
- Editorial board
- Articles
- Special issues
- Highlight articles
- Manuscript tracking
- Book reviews
- Subscribe to alerts
- Peer review
- For authors
- For reviewers
- EGU publications
- Imprint
- Data protection

Journal cover
Journal topic
**Nonlinear Processes in Geophysics**
An interactive open-access journal of the European Geosciences Union

Journal topic

- About
- Editorial board
- Articles
- Special issues
- Highlight articles
- Manuscript tracking
- Book reviews
- Subscribe to alerts
- Peer review
- For authors
- For reviewers
- EGU publications
- Imprint
- Data protection

- About
- Editorial board
- Articles
- Special issues
- Highlight articles
- Manuscript tracking
- Book reviews
- Subscribe to alerts
- Peer review
- For authors
- For reviewers
- EGU publications
- Imprint
- Data protection

**Research article**
26 Nov 2019

**Research article** | 26 Nov 2019

A fast approximation for 1-D inversion of transient electromagnetic data by using a back propagation neural network and improved particle swarm optimization

- The State Key Laboratory of Transmission Equipment & System Safety and Electrical New Technology, Chongqing University, Chongqing, 400044, China

- The State Key Laboratory of Transmission Equipment & System Safety and Electrical New Technology, Chongqing University, Chongqing, 400044, China

**Correspondence**: Huaiqing Zhang (zhanghuaiqing@cqu.edu.cn)

**Correspondence**: Huaiqing Zhang (zhanghuaiqing@cqu.edu.cn)

Abstract

Back to toptop
As one of the most active nonlinear inversion methods in transient electromagnetic (TEM) inversion, the back propagation (BP) neural network has high efficiency because the complicated forward model calculation is unnecessary in iteration. The global optimization ability of the particle swarm optimization (PSO) is adopted for amending the BP's sensitivity to its initial parameters, which avoids it falling into a local optimum. A chaotic-oscillation inertia weight PSO (COPSO) is proposed for accelerating convergence. The COPSO-BP algorithm performance is validated by two typical testing functions, two geoelectric models inversions and a field example. The results show that the COPSO-BP method is more accurate, stable and needs relatively less training time. The proposed algorithm has a higher fitting degree for the data inversion, and it is feasible to use it in geophysical inverse applications.

Download & links

How to cite

Back to top
top
How to cite.

Li, R., Zhang, H., Yu, N., Li, R., and Zhuang, Q.: A fast approximation for 1-D inversion of transient electromagnetic data by using a back propagation neural network and improved particle swarm optimization, Nonlin. Processes Geophys., 26, 445–456, https://doi.org/10.5194/npg-26-445-2019, 2019.

1 Introduction

Back to toptop
The transient electromagnetic (TEM) method applies the secondary receiving voltage induced by the rapid switching of pulse current, and it then deduces the geoelectrical parameters consisting of the resistivities and thicknesses of the layers. The later is a typical TEM inversion issue with nonlinear features. The linear inversion method was simple and widely used through linearization processes, yet it is extremely dependent on the selection of initial parameters and results in poor inversion accuracy. Hence, the nonlinear inversion methods have attracted more geophysicists' attention in recent years.

The artificial neural network (ANN) is one of the most active nonlinear inversion methods, and it has very high computation efficiency because the complicated forward model calculation is unnecessary in this iteration. All the geoelectrical parameters and the forward model relations are implied in the weight and threshold parameters of ANN. And it is different from the nonlinear Monte Carlo method with a global space search solution (He et al., 2018; Jha et al., 2008; Pekşen et al., 2014; Sharma, 2012; Tran and Hiltunen, 2012). Srinivas et al. (2012) compared the inversion performance of back propagation (BP), the radial basis function (RBF) and the generalized regression neural network (GRNN) in vertical electrical sounding data, then established a 1-D inversion model with BP and finally realized the parameter inversion. Maiti et al. (2012) proposed a Bayesian neural-network training method in 1-D electrical sounding. Jiang et al. (2018) improved the training method for the kernel principal-component wavelet neural network and achieved the resistivity imaging. Jiang et al. (2016a) produced a learning algorithm based on information criterion (IC) and particle swarm optimization for RBF network which improves the global search ability. Johnson and Aizebeokhai (2017) utilized neural-network method to invert multi-layer georesistivity sounding. Jiang et al. (2016b) presented a pruning Bayesian neural-network (PBNN) method for resistivity imaging and solved the instability and local minimization problems. Raj et al. (2014) solved nonlinear apparent resistivity inversion problems with ANN. The ANN has been widely applied in electric-prospecting data interpretation for its powerful fitting ability. However, the neural-network method is sensitive to its initial parameter settings and falls easily into a local minimum. Many improved methods were proposed for balancing the convergence rate and inversion quality. Zhang and Liu (2011) proposed ant colony optimization for ANN, and they applied high-density resistivity and acquired smaller inversion errors and higher determinant coefficients. Dai et al. (2014) suggested a differential evolution (DE) for BP which enhanced the global search ability. Rosas-Carbajal et al. (2014) introduced the genetic algorithm for ANN.

The particle swarm optimization (PSO) has a simple structure, fast convergence rate, high accuracy and global optimization ability. Fernández et al. (2010) successfully introduced the PSO in a 1-D resistivity inversion. Godio and Santilano (2018) applied it in a geophysical inversion and deduced a depth resistivity earth model. Due to the PSO's global searching performance, the BP's initial weights and thresholds can be trained by the PSO, and the BP's global optimization ability can be improved. Compared to the standard PSO (SPSO), a chaotic-oscillation inertia weight PSO (COPSO) can accelerate the convergence rate in the early stage, and it was proposed naturally (Shi et al., 2009).

The paper structure is as follows: the principles of the PSO algorithm with different inertia weights schemes, the BP neural network and the proposed COPSO-BP algorithm are given in Sect. 2. Then, the COPSO-BP algorithm performance is validated by two typical testing functions in Sect. 3. And in a later section, inversion simulations of three-layer and five-layer geoelectric models are carried out; the hidden-layer neuron numbers determining method is put forward; and algorithm performance is compared.

2 Principle of COPSO-BP algorithms

Back to toptop
For an *n*-dimensional optimization problem, it is supposed that the position (resistivity
and thickness for layered-model parameter inversion) and velocity (update
speed) of the *i*th particle (global search group number) at time *t* are
*x*_{i}= (*x*_{i1}, *x*_{i2}, …, *x*_{iN}) and
*v*_{i}= (*v*_{i1}, *v*_{i2}, …, *v*_{iN}), respectively.
Then, at time *t*+1, they can be calculated by the iterations as

$$\begin{array}{}\text{(1)}& {\displaystyle}{v}_{\mathrm{id}}^{t+\mathrm{1}}=\mathit{\omega}\cdot {v}_{\mathrm{id}}^{t}+{c}_{\mathrm{1}}{r}_{\mathrm{1}}({p}_{\mathrm{id}}^{t}-{x}_{\mathrm{id}}^{t})+{c}_{\mathrm{2}}{r}_{\mathrm{2}}({p}_{\mathrm{gd}}^{t}-{x}_{\mathrm{id}}^{t}),\text{(2)}& {\displaystyle}{x}_{\mathrm{id}}^{t+\mathrm{1}}={x}_{\mathrm{id}}^{t}+{v}_{\mathrm{id}}^{t+\mathrm{1}},\end{array}$$

where *r*_{1} and *r*_{2} are random values evenly distributed in the interval
(0,1), *c*_{1} and *c*_{2} are learning factors (usually equal to 2), and
*p*_{id} and *p*_{gd} are the individual and global maximum values.

The inertia weight parameter *ω* affects the algorithm performance
seriously. A fixed weight always was used in the early time, and then
various dynamic weights were proposed. Shi et al. (2010) have summarized
several methods as

$$\begin{array}{}\text{(3)}& {\displaystyle}{\mathit{\omega}}_{\mathrm{1}}\left(t\right)={\mathit{\omega}}_{\mathrm{s}}-({\mathit{\omega}}_{\mathrm{s}}-{\mathit{\omega}}_{\mathrm{e}})t/{T}_{max},\text{(4)}& {\displaystyle}{\mathit{\omega}}_{\mathrm{2}}\left(t\right)={\mathit{\omega}}_{\mathrm{s}}-({\mathit{\omega}}_{\mathrm{s}}-{\mathit{\omega}}_{\mathrm{e}})(t/{T}_{max}{)}^{\mathrm{2}},\text{(5)}& {\displaystyle}{\mathit{\omega}}_{\mathrm{3}}\left(t\right)={\mathit{\omega}}_{\mathrm{s}}-({\mathit{\omega}}_{\mathrm{s}}-{\mathit{\omega}}_{\mathrm{e}})\left[\mathrm{2}t/{T}_{max}-(t/{T}_{max}{)}^{\mathrm{2}}\right],\end{array}$$

where *ω*_{s} and *ω*_{e} are the start and end weight. The
*t* and *T*_{max} are the current and maximum iterations. The above weights
are smooth and monotonically decreasing. In this paper, we proposed a
decreasing oscillation weight scheme which was based on the chaotic logistic
equation. Its specific calculation formula is

$$\begin{array}{}\text{(6)}& {\displaystyle}{x}_{t+\mathrm{1}}=\mathit{\mu}{x}_{t}\left(\mathrm{1}-{x}_{t}\right)\phantom{\rule{1em}{0ex}}t=\mathrm{0},\mathrm{1},\mathrm{2},\mathrm{\dots},n,\text{(7)}& {\displaystyle}{\mathit{\omega}}_{\mathrm{c}}\left(t\right)={\mathit{\omega}}_{\mathrm{e}}+\left({\mathit{\omega}}_{\mathrm{s}}-{\mathit{\omega}}_{\mathrm{e}}\right)\left({\mathrm{0.99}}^{t}\cdot {x}_{t}\right),\end{array}$$

where *μ* is the control parameter. A complete chaos state is established
for *x*∈(0,1) and *μ*=4, and an inertia weight is then obtained from
Eq. (7). Numerical experiments were carried out correspondingly and showed
that the initial value of *x*_{0} has little effect on the inertia weight
*ω*. The inertia weight comparison was shown in Fig. 1, where *x*_{0}=0.234 and *μ*=4 for chaotic oscillation.

The BP neural network has a multi-layer feed-forward structure, and a typical three-layer network is shown in Fig. 2 (Li et al., 2009).

For Fig. 2, *x*_{1},*x*_{2}, …, *x*_{n} are the input values,
*y*_{1}, *y*_{2}, …, *y*_{m} are the predicted outputs, and *w*_{ij} and *w*_{jk} are the network
weights. The threshold parameter *α* is defined in the hidden layer with
its output,

$$\begin{array}{}\text{(8)}& {H}_{j}=f\left(\sum _{i=\mathrm{1}}^{n}{w}_{ij}{x}_{i}-{\mathit{\alpha}}_{j}\right)\phantom{\rule{1em}{0ex}}j=\mathrm{1},\mathrm{2},\phantom{\rule{0.25em}{0ex}}\mathrm{\dots},l,\end{array}$$

where *l* is the hidden-layer node numbers, *f* is the activation function with
different expressions, and the most widely used is a sigmoid-type function.
The predicted output for the *k*th unit is calculated by

$$\begin{array}{}\text{(9)}& {O}_{k}=\sum _{j=\mathrm{1}}^{l}{H}_{j}{w}_{jk}-{b}_{k},\end{array}$$

and parameter *b* is the output threshold. Then the prediction error can be
determined based on the predicted output *O*_{k} and the expected output
*T*_{k}, which is ${e}_{k}=({T}_{k}-{O}_{k}){O}_{k}(\mathrm{1}-{O}_{k})$. The updated formula
for weights and thresholds is the following:

$$\begin{array}{}\text{(10)}& \left\{\begin{array}{c}{w}_{ij}={w}_{ij}+\mathit{\eta}{H}_{j}(\mathrm{1}-{H}_{j}){x}_{i}{\sum}_{k=\mathrm{1}}^{m}{w}_{jk}{e}_{k}\\ {w}_{jk}={w}_{jk}+\mathit{\eta}{H}_{j}{e}_{k}\\ {\mathit{\alpha}}_{j}={\mathit{\alpha}}_{j}+\mathit{\eta}{H}_{j}(\mathrm{1}-{H}_{j}){\sum}_{k=\mathrm{1}}^{m}{w}_{jk}{e}_{k}\\ {b}_{k}={b}_{k}+{e}_{k}\end{array},\right.\end{array}$$

where *i* is 1,2, …, *n*, j is 1,2, …, *l*,
*k* is 1,2, …, *m* and *η* is the learning rate.

The initial parameters are chosen randomly, which affects the convergence rate, learning efficiency and perhaps falling into a local minimum. The chaotic-oscillation PSO (COPSO) has a much better global optimization capability; therefore, the COPSO algorithm is proposed to optimize the initial weight and threshold of the BP. The COPSO-BP pseudo-codes are briefly described in Algorithm 1.

The formula for calculating the *i*th particle fitness is defined as

$$\begin{array}{}\text{(11)}& {f}_{i}={\displaystyle \frac{\mathrm{1}}{S}}\sum _{s=\mathrm{1}}^{S}\sum _{j=\mathrm{1}}^{m}{\left({Y}_{sj}-{\widehat{Y}}_{sj}\right)}^{\mathrm{2}},\end{array}$$

where *S* is the number of training set samples, *m* is the output neurons number,
*Y*_{sj} is the *j*th true output of the *s*th sample, and ${\widehat{Y}}_{sj}$ is the
corresponding predicted output.

3 Algorithm testing

Back to toptop
In order to investigate the COPSO-BP performance and reliability, Rosenbrock and Bohachevsky testing functions were adopted, which are typical non-convex functions and mainly evaluate the performance of unconstrained algorithms. However, due to the random nature of the function, it is not easy to solve and has a global minimum function value of zero.

$$\begin{array}{}\text{(12)}& \begin{array}{rl}{f}_{\mathrm{1}}\left(x\right)& =\mathrm{100}\times {\left({x}_{\mathrm{1}}^{\mathrm{2}}-{x}_{\mathrm{2}}\right)}^{\mathrm{2}}+{\left(\mathrm{1}-{x}_{\mathrm{1}}\right)}^{\mathrm{2}},{x}_{i}\\ & \in \left[-\mathrm{10},\mathrm{10}\right],i=\mathrm{1},\mathrm{2}\end{array}\end{array}$$

$$\begin{array}{}\text{(13)}& \begin{array}{rl}{f}_{\mathrm{2}}\left(x\right)& ={x}_{\mathrm{1}}^{\mathrm{2}}+{x}_{\mathrm{2}}^{\mathrm{3}}-{x}_{\mathrm{1}}{x}_{\mathrm{2}}{x}_{\mathrm{3}}+{x}_{\mathrm{3}}-\mathrm{sin}\left({x}_{\mathrm{2}}^{\mathrm{2}}\right)\\ & -\mathrm{cos}\left({x}_{\mathrm{1}}{x}_{\mathrm{3}}^{\mathrm{2}}\right),{x}_{i}\in \left[-\mathrm{2}\mathit{\pi},\mathrm{2}\mathit{\pi}\right],i=\mathrm{1},\mathrm{2},\mathrm{3}\end{array}\end{array}$$

Using the standard PSO-BP (SPSO-BP) with linear decreasing inertia weight in
Eq. (3), the COPSO-BP tests were carried out. The three-layer BP with an
*n*-*s*-1 structure is constructed with different hidden nodes. The PSO parameters
are the population size *M*=60, learning factors ${c}_{\mathrm{1}}={c}_{\mathrm{2}}=\mathrm{2.0}$, maximum iteration *T*_{max}=30, inertia weight *ω*_{s}=0.9, *ω*_{e}=0.4, *x*_{0}=0.234 and *μ*=4 for
chaotic parameters, and the search dimension $D=n\times s+s\times \mathrm{1}+s+\mathrm{1}$, which includes all
the neuron weights and thresholds. For BP network, 150 training samples and
50 testing samples were randomly produced within the variable range. The
training error is defined as

$$\begin{array}{}\text{(14)}& E={\displaystyle \frac{\mathrm{1}}{S}}\sum _{s}^{S}{\left({T}_{s}-{O}_{s}\right)}^{\mathrm{2}},\end{array}$$

where *S* is the training samples number and *T*_{s} and *O*_{s} are the expected and
predicted outputs for training sample *s*, respectively. The network structures
with minimum training errors for the Rosenbrock and Bohachevsky functions are 2-7-1 and 3-6-1,
respectively. The simulation is performed 20 times for each testing function
with the SPSO-BP and COPSO-BP algorithms. The numerical result was shown in
Table 1. One of the evolutionary training error curves (randomly selected once in 20
times) was shown in Fig. 3, and the fitting curves of the COPSO-BP
algorithm were shown in Fig. 4.

It can be seen in Table 1 that although both the SPSO-BP and COPSO-BP algorithms can acquire relatively high accuracy for testing functions, the COPSO-BP algorithm is slightly better than the SPSO-BP algorithm. However, the COPSO-BP algorithm has a better convergence rate and optimization efficiency in the early stage in Fig. 3. Therefore, the SPSO-BP and COPSO-BP algorithms have a strong learning ability, good stability and generalization ability, which will be suitable for TEM inversion.

4 Layered-model and parameter analysis

Back to toptop
According to the derivation of Kaufman and Keller (1983), the frequency response of the central loop source for the layered model takes the following Hankel transform:

$$\begin{array}{}\text{(15)}& {H}_{z}(\mathit{\rho},\mathit{\omega})=Ia\underset{\mathrm{0}}{\overset{\mathrm{\infty}}{\int}}{\displaystyle \frac{{m}^{\mathrm{2}}}{m+{m}_{\mathrm{1}}/{R}_{\mathrm{1}}^{\ast}}}{J}_{\mathrm{1}}\left(m\mathit{\rho}\right)\mathrm{d}m,\end{array}$$

where *a* is the radius of the transmitting coil, *I* is the excitation current,
*ρ* is the center distance between the transmitting coil and the
receiving coil, *J*_{1}(*m**ρ*) is the 1st-order Bessel function, *m* is
the integral variable, *m*_{1} is equal to (*m*^{2}-${k}_{\mathrm{1}}^{\mathrm{2}}{)}^{\mathrm{1}/\mathrm{2}}$, *k*_{1} is the
conduction current, *σ*_{1} is the conductivity, *k*_{1} is equal to −*i**ω**μ**σ*_{1}, and ${R}_{\mathrm{1}}^{\ast}$ is the first-layer apparent
resistivity conversion function, which can be obtained by the following
recurrence formula:

$$\begin{array}{}\text{(16)}& \left\{\begin{array}{l}{R}_{n}^{\ast}=\mathrm{1}\\ {R}_{j}^{\ast}=\frac{{m}_{j}{R}_{j+\mathrm{1}}^{\ast}+{m}_{j+\mathrm{1}}\mathrm{th}\left({m}_{j}{h}_{j}\right)}{{m}_{j+\mathrm{1}}+{m}_{j}{R}_{j+\mathrm{1}}^{\ast}\mathrm{th}\left({m}_{j}{h}_{j}\right)}\end{array}.\right.\end{array}$$

There is no analytical solution for the time-domain response for the layered
model; it can only be solved by numerical calculation. The Hankel transform
in formula (15) is calculated by an improved digital filtering algorithm
with 47 points, the *J*_{1} filter coefficient, and then the time response can be
obtained using the Gaver–Stehfest transform as follows:

$$\begin{array}{}\text{(17)}& {H}_{z}(\mathit{\rho},t)={\displaystyle \frac{\mathrm{ln}\mathrm{2}}{t}}\sum _{n=\mathrm{1}}^{\mathrm{N}}{K}_{n}{H}_{z}(\mathit{\rho},{s}_{n}),\end{array}$$

where *s*_{n} is equal to $(\mathrm{ln}\mathrm{2}/t)\times n$, *K*_{n} is the coefficient and N is determined by the
computer bits, generally N =12.

The ramp excitation current of TEM is

$$\begin{array}{}\text{(18)}& I\left(t\right)=\left\{\begin{array}{c}\mathrm{0},t<\mathrm{0}\\ t/{T}_{\mathrm{1}},\mathrm{0}\le t<{T}_{\mathrm{1}}\\ \mathrm{1},{T}_{\mathrm{1}}<t\end{array},\right.\end{array}$$

where *T*_{1} is the turn-off time, and the
Laplace transform is

$$\begin{array}{}\text{(19)}& I\left(s\right)={\displaystyle \frac{\mathrm{1}}{{T}_{\mathrm{1}}{s}^{\mathrm{2}}}}-{\displaystyle \frac{\mathrm{1}}{{T}_{\mathrm{1}}{s}^{\mathrm{2}}}}{e}^{-{T}_{\mathrm{1}}s}={\displaystyle \frac{\mathrm{1}}{{T}_{\mathrm{1}}{s}^{\mathrm{2}}}}\left(\mathrm{1}-{e}^{-{T}_{\mathrm{1}}s}\right).\end{array}$$

Therefore, for a specific layered model, the apparent
resistivity conversion function ${R}_{\mathrm{1}}^{\ast}$ is
firstly calculated by the recurrence in Eq. (16) based on
geoelectric structure parameters. And then the frequency
response at fixed point *H*_{z}(*ω*) is
calculated by a Hankel transform as in Eq. (15). For ramp
excitation, the Laplace transform of *H*_{z}(*s*)
should multiplied by *I*(*s*). Finally, the
time response *H*_{z}(*t*) is obtained by a
Gaver–Stehfest transform in Eq. (17). So the
*H*_{z}(*t*) is obtained by a Gaver–Stehfest
transform, a Hankel transform and a recurrence calculation,
and it is somewhat computationally consuming.

However, the vertical magnetic field
*H*_{z}(*t*) is the actual observed signal in the
transient electromagnetic method in engineering applications.
It is the inversion input, and the outputs are geoelectric
structure parameters. A method which can avoid the
complicated forward model calculation is of great importance
in algorithm efficiency.

For BP structure, the output nodes are determined by
the number of inversion geoelectrical parameters, the input
nodes are determined by the samples number of
*H*_{z}(*t*), and the hidden nodes vary according
to approximation performance. As a three-layer or five-layer
geoelectric model, its geoelectrical parameters are five (three
resistivity and two thickness parameters) or nine (five
resistivity and four thickness parameters); the output nodes
are five or nine, correspondingly. The characteristic samplings
of *H*_{z}(*t*) are chosen as 10 or 20, which
are determined by the model's complexity, with more
layers meaning more sampling points needed. The 10 samplings
were selected in this paper, hence with 10 input nodes.
The hidden-layer neuron number is directly related
to the weight and threshold parameter amount
and greatly affects the BP performance. An appropriate
hidden-node number is necessary, and a determination
coefficient *R*^{2} is defined for evaluation as

$$\begin{array}{}\text{(20)}& {R}^{\mathrm{2}}={\displaystyle \frac{{\left(n\sum _{i=\mathrm{1}}^{n}{Y}_{i}{\widehat{Y}}_{i}-\sum _{i=\mathrm{1}}^{n}{Y}_{i}\sum _{i=\mathrm{1}}^{n}{\widehat{Y}}_{i}\right)}^{\mathrm{2}}}{\left(n\sum _{i=\mathrm{1}}^{n}{\widehat{Y}}_{i}^{\mathrm{2}}-{\left(\sum _{i=\mathrm{1}}^{n}{\widehat{Y}}_{i}\right)}^{\mathrm{2}}\right)\left(n\sum _{i=\mathrm{1}}^{n}{Y}_{i}^{\mathrm{2}}-{\left(\sum _{i=\mathrm{1}}^{n}{Y}_{i}\right)}^{\mathrm{2}}\right)}},\end{array}$$

where *Y*_{i} is the true value, ${\widehat{Y}}_{i}$
is the predicted value for the *i*th training
data and *n* is the training data number. A larger
determination coefficient means better approximation
performance. The simulations on hidden-nodes effects were
carried out for three-layer and five-layer geoelectric
models. The BP structure is 10-*s*-5 and
10-*s*-9, and its transfer, training and learning
functions are the “log sigmoidal”,
“Levenberg–Marquardt” and “gradient descent
momentum”, respectively. The average, minimum and
maximum values of *R*^{2} were obtained after
running 20 times for each simulation. The *R*^{2}
curves were shown in Fig. 5.

It can be seen that the optimal neural-network
structures were 10-2-5 and 10-5-9 for three- and five-layer
models based on the maximum *R*^{2} values. Then, the
PSO-BP algorithms with different inertia weight were
implemented and compared for the three-layer model. The BP
structure was chosen as 10-2-5, and four types of inertia
weight as in Eqs. (3)–(7) in the PSO were compared in Table 2.

The simulation was implemented on a Core (TM) i5-7500 processor with 8 GB of memory. It is obviously found in Table 2 that the COPSO algorithm has a much faster convergence rate and a lower iteration number and is time consuming.

A three-layer and five-layer geoelectric models were investigated, and the PSO parameter values are the same as those of the “Algorithm testing” section in this paper. In order to simulate actual TEM applications, the ramp turn-off is taken into account. Considering the probability distribution characteristic of the above algorithms, the average of 20 simulation results was chosen. The BP, SPSO-BP and COPSO-BP algorithms and a nonlinear programming genetic algorithm (NPGA) (Li et al., 2017) were compared.

The central loop TEM parameters were set as follows: the transmitting coil
radius was set to *a*=100 m; the ramp emission current was set to 100 A; and the turn-off time was set to
1 µs. In the geoelectric model, the resistivity is *ρ*_{1}=100, *ρ*_{2}=10 and *ρ*_{3}=100 Ω m and the thickness is *h*_{1}=100 and *h*_{2}=200 m.

The BP training samples, which are a series of *H*_{z}(*t*) for different
geoelectrical parameters, were generated by the TEM forward model. The
resistivity ranges were *ρ*_{1}∈(50, 150), *ρ*_{2}∈(5, 15) and *ρ*_{3}∈(50, 150) and the thickness ranges were *h*_{1}∈(50, 150) and *h*_{2}∈(100, 300); 1000 random groups were chosen. The
resistivity and thickness distributions of *ρ*_{1} and *h*_{1} were
shown in Fig. 6. The relative error is defined as

$$\begin{array}{}\text{(21)}& {\mathrm{Err}}_{\mathit{\_}\mathrm{rel}}=\left|\left.{\displaystyle \frac{{T}_{\mathit{\_}\mathrm{cal}}^{\ast}-{O}_{\mathit{\_}\mathrm{ref}}^{\ast}}{{O}_{\mathit{\_}\mathrm{ref}}^{\ast}}}\right|,\right.\end{array}$$

where ${T}_{\mathit{\_}\mathrm{cal}}^{\ast}$ and ${O}_{\mathit{\_}\mathrm{ref}}^{\ast}$ are the calculated and reference values for the geoelectric models.

The inversion results were shown in Table 3 and Figs. 7–8. The BP type algorithms were superior to the NPGA inversion in Table 3. Moreover, the inversion accuracy, convergence rate and optimization ability of the COPSO-BP algorithm were better than others.

Additional results showed that the solution range of *ρ*_{1} and
*h*_{1} in the 20 simulations for the above algorithms were *ρ*_{1}∈(97.980, 103.102) and *h*_{1}∈(96.962, 102.480) for BP, *ρ*_{1}∈(98.954, 101.137) and *h*_{1}∈(96.955, 101.829) for SPSO-BP, *ρ*_{1}∈(99.382, 100.989) and *h*_{1}∈(97.877, 101.044) for COPSO-BP,
respectively. Therefore, the COPSO-BP can acquire a higher accuracy and is
more stable.

A five-layer KHK-type geoelectric model was adopted, and its resistivities were
*ρ*_{1}=100, *ρ*_{2}=300, *ρ*_{3}=50,
*ρ*_{4}=200 and *ρ*_{5}=30 Ω m, and its thickness were *h*_{1}=100, *h*_{2}=200, *h*_{3}=300 and *h*_{4}=500 m.

The training samples with parameter ranges were *ρ*_{1}∈(50, 150),
*ρ*_{2}∈(150,450), *ρ*_{3}∈(25, 75), *ρ*_{4}∈(100, 300) and *ρ*_{5}∈(15, 45) for resistivity and *h*_{1}∈(50, 150), *h*_{2}∈(100, 300), *h*_{3}∈(150, 450) and *h*_{4}∈(250, 750) for thickness. The 1000 group training samples were generated
within the above ranges. The inversion results were shown in Table 4 and
Figs. 9–10. It can be seen that the COPSO-BP algorithm has
better global optimization performance.

Three kinds of BP methods, the traditional BP, SPSO-BP and COPSO-BP algorithms, were compared in Table 5. Hence, the training times of COPSO-BP were obviously less than those of SPSO-BP and were almost equal to BP; it can obtain better precision especially for its global optimization performance.

The inversion of COPSO-BP and NPGA were compared in Fig. 11. The fitting ability of COPSO-BP was much better than NPGA.

In order to verify the algorithm robustness, 5 % (26 dB) and 10 % (20 dB)
Gaussian random noise was added in TEM data for the three-layer geoelectric
model. Three kinds of inversions were implemented respectively. The results
and a comparison were shown in Table 6. The *H*_{z}(*t*) and data with 5 % noise
were shown in Fig. 12.

As can be seen in Table 2, after applying 5 % and 10 % Gaussian noise, the COPSO-BP inversion has a higher robust ability. The accuracy was obviously improved based on the total relative-error data.

In order to test the effectiveness of the method, a transient electromagnetic vertical magnetic field (Hz) with 10 measuring points at 380 to 1280 m of the no. 1 line from a mining area in Anhui Province were selected. After the data processing, the inversion was performed using the three-layer neural-network model in the previous section, and the results of the BP and COPSO-BP inversion were compared. Figure 13 shows the comparison between the surveyed data and the inversion data at 380 m of the no. 2 line in the mining area. Figure 14 displays the pseudo-sections of the 10 sets of inversion data combined with the geological data interpolation smoothing. It can be seen from Fig. 14 that the first layer has low resistivity (100–200 Ω m), which is inferred to be the second-layer (T2g22) gray dolomite of the Middle Triassic old Malague section, with a thickness of about 200 m. The second layer has the second-highest resistivity (300–400 Ω m), which is surmised to be the first-layer (T2g21) dolomite of the Middle Triassic old Malague section, with a thickness of about 400 m. The third layer has high resistivity (600–800 Ω m), which is speculated to be the sixth-layer (T2g16) limestone dolomite of the Middle Triassic old group. The results are basically consistent with the geological conditions of the mining area, indicating the feasibility and effectiveness of the neural-network method. And the results of COPSO-BP inversion are better than those of the BP, for which the inversion position is more accurate, the shape and spacing are clearer, and the resistivity of each layer is more consistent with the those of the actual geological model.

5 Discussion

Back to toptop
The inversion is performed for three-layer (H-type) and five-layer (KHK-type) geoelectric models in this paper. The results show that the BP neural network is better than the NPGA algorithm because the BP method does not need to use the forward algorithm repeatedly, and its calculation time is short, which is different from the nonlinear heuristic method based on a global space search solution.

The main advantage of the BP is that it can interpret the transient electromagnetic sounding results quickly after training the network. Furthermore, the BP algorithm could automatically obtain the “reasonable rules” between input and output data by learning, and it can adaptively store the learning content in the network weight, for which the BP neural network has a high self-learning and self-adaptation ability. In addition, the superior simulation results of the test function indicate that the BP algorithm can approximate any nonlinear continuous function with arbitrary precision, which means it has a strong nonlinear mapping ability; the inversion results of the layered geoelectric model with uncorrelated noise data prove that the BP algorithm has strong robustness, which means it has the ability to apply learning results to new knowledge. However, the BP neural-network weight is gradually adjusted by the direction of local improvement, which causes the algorithm to fall into a local extremum, and the weight converges to a local minimum that leads to the network training failure. Moreover, the BP is very sensitive to the initial network weight, and the initialization network with different weight values tends to converge at different local minimums, so it obtains different results each time. In addition, the BP algorithm is essentially a gradient descent method, which leads to a slow convergence rate.

From the results of the layered-model and parametric analysis part, it can be seen that the single BP algorithm has a higher error value than SPSO-BP because the BP method is sensitive to initial weight and easy to fall into local minimum values; thus a heuristic global search particle swarm optimization algorithm with a simple structure, rapid convergence and high precision is applied to optimize the weight and threshold of the BP neural network, which improves the global optimization performance of the algorithm. Furthermore, the PSO algorithm adjusts the inertia weight adaptively based on the chaotic-oscillation curve that is similar to the annealing process in the simulated annealing algorithm (SA), which jumps out the local extremum faster in the early stage and accelerates the convergence and reduces the training times. Therefore, compared with the SPSO-BP and BP algorithms, the inversion results of COPSO-BP are closer to the theoretical data with smaller error fluctuations, stronger anti-noise controls, a better generalization of performance and higher stability, which it is effective in solving geophysical inverse problems.

From the simulation experiment, it is not clear how the weight organization affects the BP neural-network weight learning process. It is necessary to conduct a more systematic study on this problem to improve our understanding of how the BP neural network handles training data.

6 Conclusions

Back to toptop
The nonlinear COPSO-BP method was proposed for TEM inversion. The BP's initial weight and threshold parameters were trained by the COPSO algorithm, which makes it easy to not fall into a local optimum. The chaotic-oscillation inertia weight for PSO was proposed so as to improve the PSO's global optimization ability and fast convergence in the early stage. The layered geoelectric model inversion showed that the COPSO-BP method is more accurate and stable and requires relatively less training time.

Code availability

Back to toptop
Code availability.

The PSOBP code was developed by Huaiqing Zhang and Ruiyou Li of Chongqing University in China, who are able to be contacted at the telephone number 13752954568 or e-mail address zhanghuaiqing@cqu.edu.cn. A computer with MATLAB R2016a is required is to run this code. The programming language of this code is C$++$, and its size is 10 KB. The source code is available at https://github.com/liruiyou/PSOBP (last access: 22 February 2019).

Author contributions

Back to toptop
Author contributions.

HZ conceived this paper. RL and HZ developed the main algorithmic idea and the mathematical part. RL and NY carried out the simulation and data analysis. QZ completed the writing and interpretation of this paper. All authors contributed to the writing of the paper and approved the final paper.

Competing interests

Back to toptop
Competing interests.

The authors declare that they have no conflict of interest.

Acknowledgements

Back to toptop
Acknowledgements.

This work was partly supported by the National Natural Science Foundation of China (grant no. 51377174).

Financial support

Back to toptop
Financial support.

This research has been supported by the National Natural Science Foundation of China (grant no. 51377174).

Review statement

Back to toptop
Review statement.

This paper was edited by Norbert Marwan and reviewed by two anonymous referees.

References

Back to toptop
Dai, Q., Jiang, F., and Dong, L.: Nonlinear inversion for electrical resistivity tomography based on chaotic DE-BP algorithm, J. Cent. South. Univ., 21, 2018–2025, https://doi.org/10.1007/s11771-014-2151-9, 2014.

Fernández Martínez, J. L., García Gonzalo, E., Fernández Álvarez, J. P., Kuzma, H. A., and Menéndez Pérez, C. O.: PSO: A powerful algorithm to solve geophysical inverse problems: Application to a 1D-DC resistivity case, J. Appl. Geophys., 71, 13–25, https://doi.org/10.1016/j.jappgeo.2010.02.001, 2010.

Godio, A. and Santilano, A.: On the optimization of electromagnetic geophysical data: Application of the PSO algorithm, J. Appl. Geophys., 148, 163–174, https://doi.org/10.1016/j.jappgeo.2017.11.016, 2018.

Jha, M. K., Kumar, S., and Chowdhury, A.: Vertical electrical sounding survey and resistivity inversion using genetic algorithm optimization technique, J. Hydrol., 359, 71–87, https://doi.org/10.1016/j.jhydrol.2008.06.018, 2008.

Jiang, F., Dai, Q., and Dong, L.: An ICPSO-RBFNN nonlinear inversion for electrical resistivity imaging, J. Cent. South. Univ., 23, 2129–2138, https://doi.org/10.1007/s11771-016-3269-8, 2016a.

Jiang, F., Dai, Q., and Dong, L.: Nonlinear inversion of electrical resistivity imaging using pruning Bayesian neural networks, J. Appl. Geophys., 13, 267–278, https://doi.org/10.1007/s11770-016-0561-1, 2016b.

Jiang, F., Dong, L., and Dai, Q.: Electrical resistivity imaging inversion: An ISFLA trained kernel principal component wavelet neural network approach, Neural. Networks., 104, 114–123, https://doi.org/10.1016/j.neunet.2018.04.012, 2018.

Johnson, O. L. and Aizebeokhai, A. P.: Application of Artificial Neural Network for the Inversion of Electrical Resistivity Data, Journal of Informatics and Mathematical Sciences, 9, 297–316, 2017.

Kaufman, A. A. and Keller, G. V.: Frequency and Transient Sounding, Elsevier Methods in Geochemistry & Geophysics No 16, 620 pp., 1983.

Li, F. P., Yang, H. Y., and Liu, X. H.: Nonlinear programming genetic algorithm in transient electromagnetic inversion, Geophysical and Geochemical Exploration, 41, 347–353, 2017.

Li, Y. Y., Chen, B. C., Zhao, Y. G., Yun, C., Ma, X. B., and Kong, X. R.: Nonlinear inversion for electrical resistivity tomography, Chinese J. Geophys., 52, 758–764, https://doi.org/10.1016/S1003-6326(09)60084-4, 2009.

Maiti, S., Erram, V. C., Gupta, G., and Tiwari, R. K.: ANN based inversion of DC resistivity data for groundwater exploration in hard rock terrain of western Maharashtra (India), J. Hydrol., 464, 294–308, https://doi.org/10.1016/j.jhydrol.2012.07.020, 2012.

Pekşen, E., Yas, T., and Kıyak, A.: 1-D DC Resistivity Modeling and Interpretation in Anisotropic Media Using Particle Swarm Optimization, Pure Appl. Geophys., 171, 2371–2389, https://doi.org/10.1007/s00024-014-0802-2, 2014.

Raj, A. S., Srinivas, Y., and Oliver, D. H.: A novel and generalized approach in the inversion of geoelectrical resistivity data using Artificial Neural Networks (ANN), J. Earth Syst. Sc., 123, 395–411, https://doi.org/10.1007/s12040-014-0402-7, 2014.

Rosas-Carbajal, M., Linde, N., Kalscheuer, T., and Vrugt, J. A.: Two-dimensional probabilistic inversion of plane-wave electromagnetic data: methodology, model constraints and joint inversion with electrical resistivity data, Geophys. J. Int., 196, 1508–1524, https://doi.org/10.1093/gji/ggt482, 2014.

Sharma, S. P.: VFSARES – a very fast simulated annealing FORTRAN program for interpretation of 1-D DC resistivity sounding data from various electrode arrays, Comput. Geosci., 42, 177–188, https://doi.org/10.1016/j.cageo.2011.08.029, 2012.

Shi, F., Wang, X. C., and Yun, L.: Matlab neural network case study, The Beijing University of Aeronautics & Astronautics Press, Beijing, 2010.

Shi, X. M., Xiao, M., Fan, J. K., Yang, G. S., and Zhang, X. H.: The damped PSO algorithm and its application for magnetotelluric sounding data inversion, Chinese J. Geophys., 52, 1114–1120, https://doi.org/10.3969/j.issn.0001-5733.2009.04.029, 2009.

Srinivas, Y., Raj, A. S., Oliver, D. H., Muthuraj, D., and Chandrasekar, N.: A robust behavior of Feed Forward Back propagation algorithm of Artificial Neural Networks in the application of vertical electrical sounding data inversion, Geosci. Front., 3, 729–736, https://doi.org/10.1016/j.gsf.2012.02.003, 2012.

Tran, K. T. and Hiltunen, D. R.: Two-Dimensional Inversion of Full Waveforms Using Simulated Annealing, J. Geotech. Geoenviron. Eng., 138, 1075–1090, https://doi.org/10.1061/(ASCE)GT.1943-5606.0000685, 2012.

Wang, H., Liu M. L., Xi, Z. Z., Peng, X. L., and He, H.: Magnetotelluric inversion based on BP neural network optimized by genetic algorithm, Chinese J. Geophys., 61, 1563–1575 https://doi.org/10.6038/cjg2018L0064, 2018.

Zhang, L. Y. and Liu, H. F.: The application of ABP method in high-density resistivity method inversion, Chinese J. Geophys., 54, 64–71, https://doi.org/10.1002/cjg2.1587, 2011.

Short summary

The chaotic-oscillation inertia weight back propagation (COPSO-BP) algorithm is proposed for transient electromagnetic inversion. The BP's initial weight and threshold parameters were trained by COPSO, overcoming the BP falling into a local optimum. Inversion of the layered geoelectric model showed that the COPSO-BP method is accurate and stable and needs less training time. It can be used in 1-D direct current sounding, 1-D magnetotelluric sounding, seismic-wave impedance and source detection.

The chaotic-oscillation inertia weight back propagation (COPSO-BP) algorithm is proposed for...

Nonlinear Processes in Geophysics

An interactive open-access journal of the European Geosciences Union