A Differential Evaluation Markov Chain Monte Carlo algorithm for Bayesian Model
Updating
M. Sherri a, I. Boulkaibet b, T. Marwala b, M. I. Friswell c,
a

Department of Mechanical Engineering Science, University of Johannesburg, PO Box 524, Auckland
Park 2006, South Africa.
b
Institute of Intelligent Systems, University of Johannesburg, PO Box 524, Auckland Park 2006, South
Africa.
c
College of Engineering, Swansea University, Bay Campus, Swansea SA1 8EN, United Kingdom.

Abstract:
The use of the Bayesian tools in system identification and model updating paradigms has been increased in the last
ten years. Usually, the Bayesian techniques can be implemented to incorporate the uncertainties associated with
measurements as well as the prediction made by the finite element model (FEM) into the FEM updating procedure. In this
case, the posterior distribution function describes the uncertainty in the FE model prediction and the experimental data.
Due to the complexity of the modeled systems, the analytical solution for the posterior distribution function may not exist.
This leads to the use of numerical methods, such as Markov Chain Monte Carlo techniques, to obtain approximate
solutions for the posterior distribution function. In this paper, a Differential Evaluation Markov Chain Monte Carlo (DEMC) method is used to approximate the posterior function and update FEMs. The main idea of the DE-MC approach is to
combine the Differential Evolution, which is an effective global optimization algorithm over real parameter space, with
Markov Chain Monte Carlo (MCMC) techniques to generate samples from the posterior distribution function. In this
paper, the DE-MC method is discussed in detail while the performance and the accuracy of this algorithm are investigated
by updating two structural examples.
Keywords: Bayesian model updating; Markov Chain Monte Carlo; differential evolution; finite element model; posterior
distribution function.
1.

Introduction

During the last thirty years, the application of the finite element method (FEM) [1-3] has exponentially increased
where this numerical technique has become one of the most popular engineering tools in systems modelling and prediction.
In the domain of structural dynamics, the FEM tools are widely applied to model complex systems where this technique
can produce results with high accuracy, especially when the modelled system is simple. However, the results attained by
the FEM can be relatively inaccurate and the mismatches between the FEM results and the results attained from
experimental studies are relatively significant. This is due to the errors associated with the modelling process as well as
the complexity of modelled structure, which may reduce the accuracy of the modelling process. Consequently, the model
obtained by an FEM needs to be updated to reduce the errors between the experimental and modelled outputs. The
procedure of minimizing the differences between the numerical results and the measured data is known as the finite
element model updating (FEMU) [4, 5], where the FEMU methods can be divided into two main classes. In the first class,
which is the direct methods, we equate the experimental data directly to the FEM outputs resulting in a procedure that
constrains the updating to the FE system matrices (mass, stiffness) only. This kind of approach may produce non-realistic
results where the resulting updating parameters may not have physical meaning. In the second class, which is also known
as the iterative (or indirect) approaches, the FEM outputs are not directly equated to the experimental data, but instead, an
objective function is introduced and iteratively minimised to reduce the errors between the analytical and experiential
results. Thus, we vary the system matrices and the model output during the minimisation process, and realistic results are
often expected at the end of the updating process.

Generally, several sources of uncertainty are associated with the modelling process, such as the mathematical
simplifications made during the modelling, where this kind of uncertainty may affect the accuracy of the modelling
process. Moreover, the noise that contaminates the experimental results may also have a significant impact on the updating
process. To deal with such uncertainty problems, another class of methods called the uncertainty quantification methods
accomplishes the updating process. The most common uncertainty quantification method is known as the Bayesian
approach in which the unknown parameters and their uncertainty are identified by defining each unknown parameter with
a probability density distribution (PDF). Recently, the use of the Bayesian methodology has massively increased in the
domain of system identification and uncertainty quantification. In this approach, the uncertainties associated with the
modelled structure are expressed in terms of probability distributions where the unknown parameters are defined as a
random vector with a multi-variable probability density function, and the resulting function is known as the posterior PDF.
Solving the posterior PDF helps in identifying the unknown parameters and their uncertainties. Unfortunately, the posterior
PDF cannot be solved in an analytical way for sufficiently complex problems which is the case for the FEMU problems
since the search space is usually nonlinear and high dimensional. In this case, sampling techniques are employed to identify
these uncertain parameters. The most recognised sampling methods are these related to Markov chain Monte Carlo
(MCMC) methods.
Generally, the MCMC methods are very useful tools that can efficiently cope with large search spaces and generate
samples from complex distributions. These methods draw samples with an element of randomness while being guided by
the values of the posterior distribution function. Then, the drawn samples are accepted or rejected according to the
Metropolis criterion. Unfortunately, the updated models, with relatively large complexities, may have multiple optimal
(or near optimal) simple MCMC algorithms cannot easily identify solution, and this. In this paper, another version of the
MCMC algorithms, known as the Differential Evolution Markov Chain (DE-MC) [6, 7] algorithm, is used to update FEMs
of structural systems. The DE-MC algorithm combines the abilities of the differential evaluation algorithm [8, 9], which
is one of the genetic algorithms for global optimization, with the Metropolis-Hastings algorithm. In this algorithm, multiple
chains are run in parallel, and the exploration and exploitation of the search space in the current chain are achieved by the
difference of two randomly selected chains, multiplied by the value of the difference with a preselected factor and then
the result is added to the value of the current chain. The value of the current chain is then accepted or rejected according
to the Metropolis criterion. In this paper, the efficiency, reliability and the limitations of the DE-MC algorithm are
investigated when the Bayesian approach is applied for FEMU. This paper is organized as follows: in the next section, the
Bayesian formulations are introduced. Section 3 describes the DE-MC algorithm while section 4 presents the results when
a simple mass-spring structure is updated. Section 5 presents the updating results of an unsymmetrical H-shaped Structure.
The paper is concluded in section 6.
2.

Bayesian formulations

In this paper, the Bayesian approach is adopted to compute the posterior distribution function in order to update the
FEMs. The posterior function can be represented by Bayes rule [10-14]:
𝑃(𝜽|𝒟, ℳ) 𝑃(𝒟|𝜽, ℳ) 𝑃(𝜽|ℳ)

(2.1)

where ℳ describes the model class for the target system where each model class ℳ is defined by certain updating
parameters 𝜽 ∈ 𝚯 ⊂ ℛ 𝑑 . The experimental data 𝒟 of the structural system, which is represented by the natural frequencies
𝑓𝑖𝑚 and mode shapes 𝝓𝑚
𝑖 , are used to improve the FEM results. 𝑃(𝜽|ℳ) is the prior probability distribution function
(PDF) that represents the initial knowledge of the uncertain parameters given a specific model ℳ, and in the absence of
the measured data 𝒟. The function 𝑃(𝒟|𝜽, ℳ) is known as the likelihood function and represents the difference between
the experimental data and the FEM results. Finally, the probability distribution function 𝑃(𝜽|𝒟, ℳ) is the posterior
function of the unknown parameters given a model class ℳ and the measured data 𝒟. The model class ℳ is used only
when several classes are investigated for both model updating and model selection. In this paper, only one model class is
considered, and therefore, the term ℳ is omitted in order to simplify the Bayesian formulations.
In this paper, the likelihood function is given by:

2

1

𝑃(𝒟|𝜽) =
(

2𝜋 𝑁𝑚
)
𝛽𝑐

⁄2

exp (−
𝑚
𝑚
∏𝑁
𝑖=1 𝑓𝑖

𝑁𝑚 𝑓 𝑚 − 𝑓
𝛽𝑐
𝑖
𝑖
∑ (
) )
2
𝑓𝑖𝑚
𝑖

(2.2)

where 𝑁𝑚 is the number of measured modes, 𝛽𝑐 is a constant, 𝑓𝑖𝑚 and 𝑓𝑖 are the 𝑖th analytical and measured natural
frequencies. The initial knowledge of the updating parameters 𝜽, which is defined by a prior PDF, is given by the following
Gaussian distribution:
𝑃(𝜽) =

𝑄𝛼

1
(2𝜋)𝑄⁄2 ∏𝑄
𝑖=1

1
√𝛼𝑖

exp (− ∑

𝑖

𝑖

2

2

‖𝜃 𝑖 − 𝜃0𝑖 ‖ ) =

1
(2𝜋)𝑄⁄2 ∏𝑄
𝑖=1

1
√𝛼𝑖

1
exp (− (𝜽 − 𝜽0 )𝑇 Σ−1 (𝜽 − 𝜽0 ))
2

(2.3)

where 𝑄 is the number of the uncertain parameters, 𝜽0 represents the mean value of the updating parameters, 𝛼𝑖 , 𝑖 =
1, … , 𝑄 are the coefficients of the updating parameters and the Euclidean norm is given by the notation: ‖∗‖.
After substituting Eqs. (2.2) and (2.3) into the Bayesian inference defined by Eq. (2.1), the posterior 𝑃(𝜽|𝒟) of the
unknown parameters 𝜽 given the experimental data 𝒟 is characterized by:

2

𝑃(𝜽|𝒟)

𝑁𝑚 𝑓 𝑚 − 𝑓
𝑄𝛼
1
𝛽𝑐
2
𝑖
𝑖
𝑖
exp (−
∑ (
‖𝜃 𝑖 − 𝜃0𝑖 ‖ )
𝑚 ) −∑
𝑍𝑠 (𝛼, 𝛽𝑐 )
2
𝑓𝑖
𝑖
𝑖 2

(2.4)

where
𝑁𝑚
𝑄
2𝜋 𝑁𝑚⁄2
1
𝑍𝑠 (𝛼, 𝛽𝑐 ) = ( )
∏ 𝑓𝑖𝑚 (2𝜋)𝑄/2 ∏
𝛽𝑐
𝑖=1
𝑖=1 √𝛼𝑖

(2.5)

Generally, the complexity of the posterior PDF, which depends on the modal parameters of the analytical model, is
related to the complexity of the analytical model, and for certain relatively complex structural models the analytical results
for the posterior distribution are difficult to obtain due to the high dimensionality of the search space. In this case, sampling
techniques [5, 10, 11, 13, 14] are the only practical approaches in order to approximate the posterior PDF. The main idea
of sampling techniques is to generate a 𝑁𝑠 sequence of vectors {𝜽1 , 𝜽2 , … , 𝜽𝑁𝑠 } and use these samples to approximate the
future response of the unknown parameters at different time instances. The most recognized sampling techniques are
Markov Chain Monte Carlo (MCMC) methods [5, 13-18]. In this paper, the combination of one of the basic MCMC
algorithms, known as the Metropolis-Hasting algorithm, with one of the genetic algorithms, known as differential
evolution (DE), is used to generate samples from the posterior PDF in order to update structural models.
3.

The Differential Evolution Markov Chain Monte Carlo (DE-MC) method

In this paper, the DE and MCMC methods, which are extremely popular methods in several scientific domains, are
combined to improve the convergence of the sampling procedure. In this approach, multiple chains are run in parallel in
order to improve the accuracy of the updating parameters, while these chains learn from each other instead of running all
the chains independently. This may improve the efficiency of the searching procedure and avoid sampling in the vicinity
of a local minimum. The new chains are then accepted or rejected according to the Metropolis-Hastings criterion.
The Metropolis-Hastings (M-H) [18, 19, 20] algorithm is one of the common MCMC methods that can be used to
draw samples from multivariate probability distributions. To sample from the posterior PDF 𝑃(𝜃|𝐷), where 𝜽 =
{𝜃1 , 𝜃2 , … , 𝜃𝑑 } is a 𝑑-dimensional parameters vector, a proposal density distribution 𝑞(𝜽|𝜽𝑡−1 ) is used to generate a

proposed random vector 𝜽∗ given the value at the previous accepted vector 𝜽𝑡−1 at the iteration 𝑡 − 1 of the algorithm.
Next, the Metropolis criterion is used to accept or reject the proposed sample 𝜽∗ as follows:
𝛼(𝜽∗ |𝜽𝑡−1 ) = min {1,

𝑃(𝜽∗ |𝐷) 𝑞(𝜽𝑡−1 |𝜽∗ )
}
𝑃(𝜽𝑡−1 |𝐷)𝑞(𝜽∗ |𝜽𝑡−1 )

(3.1)

On the other hand, the Differential Evolution (DE) [8] is a very effective genetic algorithm in solving various realworld global optimization problems. As one of the genetic algorithms, the DE algorithm begins by randomly initialising
the population within certain search area, and then these initial values are evolved over the generations in order to find the
global minimum. This can be achieved using genetic operators such as: mutation, selection, and crossover.
By integrating the Metropolis-Hastings criterion within the search abilities of the DE algorithm, the resulted MCMC
method can be more efficient in determining where other chains can be employed to create the new candidates for the
current chain. In the DE-MC algorithm, the new value of the chain is obtained by a simple mutation operation where the
difference between two randomly selected chains (different from the current chain) are added to the current chain. Thus,
the proposal for each chain depends on a weighted combination of other chains which can be easily defined as [6, 7]:
𝜽∗ = 𝜽𝑖 + 𝛾(𝜽𝑎 − 𝜽𝑏 ) + 𝜺

(3.2)

where 𝜽∗ represents the new proposed vector, 𝜽𝑖 is the current state of the 𝑖-th chain, 𝜽𝑎 and 𝜽𝑏 are randomly selected
chains, 𝛾 is a tuning factor that always take a positive value and can be set to vary between [0.4, 1]. Note that the vectors:
𝜽𝑖 ≠ 𝜽𝑎 ≠ 𝜽𝑏 . Finally, the noise 𝜺, which is defined as a Gaussian distribution 𝜺~𝑁𝑝 (𝟎, 𝝈2 ) with a very small variance
vector 𝝈2 , is added to the proposed vector to avoid degeneracy problems. The factor 𝛾 can be seen as the magnitude that
controls the jumping distribution. The main idea of the DE-MC algorithm can be illustrated in Figure 3.1b.

(a): Metropolis-Hastings

(b): DE-MC

Figure 3.1: Proposed vector generation in the M-H and DE-MC methods

Figure 3.1 explains the way to generate proposed vectors for the M-H method (Figure 3.1a) and for the DE-MC
method (Figure 3.1b). As illustrated, the difference vector between the two randomly selected chains 𝜽𝑎 and 𝜽𝑏 represents
the direction of the new proposed vector, where this difference is multiplied by the factor 𝛾 to define the moving distance.
The moving distance is then added to the current chain 𝜽𝑖 to create the proposed vector. Note that, the DE-MC method
has only one tuning factor 𝛾 in comparing to other versions of evolutionary MCMC methods. Finally, the new proposal
𝜽∗ of the 𝑖-th chain is accepted or rejected according to the Metropolis criterion which is given as:

𝑟 = min {1,

𝑃(𝜽∗ |𝐷)
}
𝑃(𝜽𝒊 |𝐷)

(3.3)

The steps to update FEMs using the DE-MC algorithm are summarized as follows:
1234-

Initialize the population 𝜽𝑖,𝜊 , 𝑖 ∈ {1,2, … , 𝑁}.
Set the tuning factor 𝛾. In this paper, 𝛾 = 2.38/√2𝑑 and 𝑑 is the dimension of the updating parameters.
Calculate the Posterior PDF for all chains.
For all chains 𝑖 ∈ {1,2, … , 𝑁}:
4.1 Sample uniformly two random vectors 𝜽𝑎 ,𝜽𝑏 where 𝜽𝑎 ≠ 𝜽𝑏 ≠ 𝜽𝑖 .
4.2 Sample the random value 𝜺 with small variance 𝜺~𝑁𝑝 (𝟎, 𝝈2 ).
4.3 Calculate the proposed vector 𝜽∗ = 𝜽𝑖 + 𝛾(𝜽𝑎 − 𝜽𝑏 ) + 𝜺.
4.4 Calculate the Posterior PDF for the vector 𝜽∗ .
4.5 Calculate the Metropolis ratio 𝑟 = min {1,

𝑃(𝜽∗ |𝐷)
𝑃(𝜽𝒊 |𝐷)

}

4.6 Accept the proposed vector 𝜽𝑖 ← 𝜽∗ with probability 𝑚𝑖𝑛 (1, 𝑟), otherwise 𝜽𝑖 is unchanged.
5- Repeat the steps 4.1 to 4.6 until the number of samples required is achieved.
In next two sections, the DE-MC performance is highlighted when two structural examples are updated.
4.

Application 1: Simple Mass-Spring system

In this section, a five degrees of freedom mass-spring linear system, as presented in Figure 4.1, is updated using the DEMC algorithm.

Figure 4.1: The five degrees of freedom mass-spring system

The system contains 5 masses connected to each other using 10 springs (see Figure 4.1). The deterministic values of
the masses are: 𝑚1 = 2.7 kg, 𝑚2 = 1.7 kg, 𝑚3 = 6.1 kg, 𝑚4 = 5.3 kg and 𝑚5 = 2.9 kg. The stiffness of the springs are:
𝑘3 = 3200 N/m, 𝑘5 = 1840 N/m, 𝑘7 = 2200 N/m, 𝑘9 = 2800 N/m and 𝑘10 = 2000 N/m. The spring stiffnesses
𝑘1 , 𝑘2 , 𝑘4 , 𝑘6 , and 𝑘8 are considered as the uncertain parameters where the updating vector is: 𝜽 = {𝜃1 , 𝜃2 , 𝜃3 , 𝜃4 , 𝜃5 } =
{𝑘1 , 𝑘2 , 𝑘4 , 𝑘6 , 𝑘8 }.

Since the DE-MC method is used for the updating procedure, the population used by the algorithm is selected to be
𝑁 = 10. The updating vectors are bounded by 𝜽𝑚𝑎𝑥 and 𝜽𝑚𝑖𝑛 which are set to {4800, 2600,2670, 3400, 2750} and
{3200, 1800,1600, 1800, 2050}, respectively. The tuning factor is set to 𝛾 = 2.38/√2𝑑 while 𝑑 = 5, the initial vector
of 𝜽 is set to 𝜽0 = {4600, 2580, 1680, 3100, 2350} and the number of generations (number of samples) is set to 𝑁𝑠 =
10000. The obtained samples are illustrated in Figure 4.2 while the updating parameters, as well as the initial and updated
natural frequencies, are shown in Tables 4.1 and 4.2, respectively.
2600

2350
2300

2500

2250

2400
2200

2300

3

2

2150
2100
2050

2200
2100

2000

2000
1950

1900

1900
1850
3400

3600

3800

4000

4200

1800
1800

4400

1

2000

2

2200

2400

2600

3100
3000

2500

2900
2800

2400

5

4

2700
2600

2300

2500

2200

2400
2300

2100

2200
2100
1800

2000

2200

2400

2600

3

2000
2000

2500

4

3000

3500

Figure 4.2: The scatter plots of the samples using the DE-MC algorithm

Figure 4.2 shows the scatter plots of the uncertain parameters using the DE-MC algorithm. The confidence ellipses
(error ellipse) of the samples are also shown in the same figures (in red) where these ellipses visualize the regions that
contain 95% of the obtained samples. As expected, the figure shows that the DE-MC algorithm has found the high
probability area after only few iterations. Table 4.1 contains the initial values, the nominal values and the updated values
of the uncertain parameters. The coefficient of variation (c.o.v) values, which are estimated by dividing the standard
deviations 𝜎𝑖 by the updated vectors 𝜽𝑖 (or 𝝁𝑖 ), are also presented in Table 4.1 and used to measure the errors in the
updating. It is clear that the obtained values of the c.o.v when the DE-MC algorithm is used to update the structure are
small and less than 2.5% which means that the DE-MC algorithm performed well and was able to identify the areas with
high probability. This also can be verified from the same table where the updating parameters are close the nominal values.
Table 4.1: The updating parameters using DE-MC technique
Unknown parameters (N/m)
Nominal values

Error (%)

𝜃1

4600

4010

14.71

4004.4

0.14

1.03

𝜃2

2580

2210

16.74

2197.6

0.56

1.71

𝜃3

1680

2130

21.13

2109.4

0.97

2.07

𝜃4

3100

2595

19.46

2600.9

0.23

2.32

𝜃4

2350

2398

02.00

2410.4

0.52

1.77

DE-MC (𝜇𝑖 )

Error (%)

𝜎𝑖

Initial

𝜇𝑖

c.o.v %)

Table 4.2 contains the initial, nominal and updated natural frequencies. Furthermore, the absolute errors, which are
estimated by

|𝑓𝑖𝑚−𝑓𝑖 |

, the total average error (𝑇𝐴𝐸), which is computed by 𝑇𝐴𝐸 =

𝑓𝑖𝑚

1
𝑁𝑚

𝑚
∑𝑁
𝑖=1

|𝑓𝑖𝑚 −𝑓𝑖 |
𝑓𝑖𝑚

, 𝑁𝑚 = 5, and the c.o.v

values are also displayed. Obviously, the updated frequencies obtained by the DE-MC are better than the initial
frequencies, and almost equal to the nominal frequencies.
Table 4.2: The updated natural frequencies and the errors obtained using the DE-MC
Modes

Nominal
Frequency

Initial Frequency

Error

Frequency DE-MC

c.o.v

Error

(Hz)

(%)

(Hz)

(%)

(%)

(Hz)
1

3.507

3.577

1.97

3.507

0.118

0.00

2

5.149

5.371

4.30

5.149

0.126

0.00

3

7.083

7.239

2.21

7.082

0.119

0.02

4

8.892

9.030

1.56

8.894

0.140

0.03

5

9.426

9.412

0.16

9.426

0.117

0.00

TAE

_______

_______

1.98

_______

_____

0.012

The total average error of the FEM output was reduced from 1.98% to 0.012%. On the other hand, the values of the
c.o.v for all updated frequencies are smaller than 0.15% which indicates that the DE-MC technique efficiently updated the
structural system. Figure 4.3 shows the evaluation of the total average error at each iteration. The 𝑇𝐴𝐸 in Figure 4.3 is
1
obtained as follows: first, the mean value of the samples at each iteration is computed as 𝜃̂ = 𝐸(𝜃) =
̃ ∑𝑖𝑗=1 𝜃 𝑖 where 𝑖
𝑁𝑠

is the current iteration. Next, the mean value is used to compute the analytical frequencies of the FEM, and then the total
average error is calculated as: 𝑇𝐴𝐸(𝑖) =

1
𝑁𝑚

𝑚
∑𝑁
𝑗=1

|𝑓𝑗𝑚−𝑓𝑗 |
𝑓𝑗𝑚

. As a result, it is clear that the DE-MC algorithm converges

efficiently after the first 2000 iterations.

Total average error

10

10

10

0

-1

-2

0

1000

2000

3000

4000

5000

6000

7000

8000

9000

10000

Iteration

Figure 4.3: The evaluation of the TAE using the DE-MC method
Figure 4.4 illustrates the correlation between the updating parameters where all parameters are correlated (the values
are different from zero). Moreover, the majority of these parameters are weakly correlated (small values <0.3) except the
pairs (𝜃1 , 𝜃2 ) and (𝜃4 , 𝜃5 ) which are highly correlated (values >0.7).

1.5

Correlation

1
0.5
0
-0.5
-1
1
2
5

3

4

4

3
5

2
1

i

i

Figure 4.4: The correlation between the updating parameters.
In the next section, the DE-MC method is used to update an unsymmetrical H-shaped aluminum structure with real
experimental data.
5.

Application 2: The FEMU of the unsymmetrical H-Shaped structural system.

In this section, the performance of the DE-MC algorithm is examined by updatintg an unsymmetrical H-shaped
aluminum structure with real measured data. The FEM model of the H-shaped structure is presented in figure (5.1) where
the structure is divided into 12 elements, and each element is modelled as an Euler-Bernoulli beam. The location displayed
by a double arrow at the middle beam indicates the position of excitation which is produced by an electromagnetic shaker.
An accelerometer was used to measure the set of frequency-response functions. The initial analytical natural frequencies
are 53.9, 117.3, 208.4, 254.0 and 445.0 Hz. In this example, the moments of inertia 𝐼𝑥𝑥 and the cross-sectional areas 𝐴𝑥𝑥
of the left, middle and right subsections of the H-shaped beam are selected to be updated in order to improve the analytical
natural frequencies. Thus, the updating parameters are: 𝜽 = {𝐼𝑥1 , 𝐼𝑥2 , 𝐼𝑥3 , 𝐴𝑥1 , 𝐴𝑥2, 𝐴𝑥3, }.

Figure 5.1: The Unsymmetrical H-Shaped Structure

The rest of the H-shaped structure parameters are given as follows: The Young’s modulus is set to 7.2 × 1010 N/m2
and the density is set to 2785 kg/m3. The updating parameters 𝜽 are bounded by maximum and minimum vectors given
by: [4.73 × 10−8 , 4.73 × 10−8 , 4.73 × 10−8 , 5.16 × 10−4 , 5.16 × 10−4 , 5.16 × 10−4 ] and [0.73 × 10−8 , 0.73 ×
10−8 , 0.73 × 10−8 , 1.16 × 10−4 , 1.16 × 10−4 , 1.16 × 10−4 ], respectively. These boundaries are used to ensure that the
updating parameters are physically realistic. The number of samples is set to 𝑁𝑠 = 5000, the factor 𝛽𝑐 of the likelihood

function is set equal to 10, the coefficients 𝛼𝑖 of the prior PDF are set to

1
𝜎𝑖2

where 𝜎𝑖2 is the variance of the 𝑖th uncertain

parameters, and 𝜎 = [5 × 10−8 , 5 × 10−8 , 5 × 10−8 , 5 × 10−4 , 5 × 10−4 , 5 × 10−4 ].
Figure 5.2 illustrates the scatter plots of the updating parameters. The confidence ellipse that contains 95% of samples
are also included in the figure. The updating parameters were normalized by dividing the parameters by 𝒌 =
[10−8 , 10−8 , 10−8 , 10−4 , 10−4 , 10−4 ]. As expected, the DE-MC algorithm was able to find the area with high probably
after a few iterations. The rest of the updating parameters are shown in Table 5.1 as well as the initial values, the c.o.v
values and the updating parameters obtained by the M-H algorithm.

4

4.2

4

3.5

3.8

3.6

4 / k4

3 / k3

3
3.4

3.2

2.5
3

2.8

2

2.6

1.5
1.5

2

2.5

3

3.5

1.5

2

2.5

3

3.5

4

3 / k3

1 / k1

Figure 5.2: The scatter plots of the samples using the DE-MC algorithm

Table 5.1: The initial parameters, the c.o.v values and the updating parameters using the DE-MC and M-H algorithms
Initial

DE-MC

𝜎𝑖
𝜇𝑖

(%)

(𝜇𝑖 )

M-H

𝜎𝑖
𝜇𝑖

(%)

(𝜇𝑖 )

𝜃1

2.7265 × 10−8

2.8965 × 10−8

5.76

2.31 × 10−8

22.59

𝜃2

2.7265 × 10−8

2.9739 × 10−8

1.91

2.68 × 10−8

15.25

𝜃3

2.7265 × 10−8

1.7676 × 10−8

1.67

2.17 × 10−8

13.96

𝜃4

3.1556 × 10−4

3.8966 × 10−4

0.65

2.85 × 10−4

14.36

𝜃5

3.1556 × 10−4

2.1584 × 10−4

2.93

2.83 × 10−4

14.36

𝜃6

3.1556 × 10−4

2.9553 × 10−4

0.026

2.77 × 10−4

13.08

The results in Table 5.1 indicate that the updating parameters obtained by the DE-MC and M-H algorithms are
different from the initial values which means that the uncertain parameters have been successfully updated. Furthermore,
the c.o.v values of the updating parameters obtained by the DE-MC algorithm are relatively small (<2.5%) with verifies
that the algorithm was able to identify the areas with high probability in a reasonable amount of time; however, the c.o.v
obtained by the M-H algorithm are relatively high (>13.08%) which means that the M-H algorithm does not have the
efficiency of the DE-MC algorithm.
Table 5.2 illustrates the updating frequencies using the DE-MC and M-H algorithms, the errors and the c.o.v values.
As expected, the analytical frequencies obtained by the DE-MC algorithm are better than the initial frequencies as well as
the frequencies obtained by the M-H algorithm. The DE-MC method has improved all natural frequencies and reduced
the total average error (TAE) from 5.37% to 1.53%. Also, the c.o.v values obtained by the DE-MC method are relatively
small (<0.65%).
Table 5.2: Natural frequencies, c.o.v values and errors when DE-MC and M-H techniques are used for FEMU
Modes

Measured
Frequency
(Hz)

Initial
Frequency
(Hz)

Error
(%)

Frequency
DE-MC
(Hz)

c.o.v
(%)

Error
(%)

Frequency
M-H
(Hz)

c.o.v
(%)

Error
(%)

1

53.90

51.04

5.31

52.56

0.30

2.49

53.92

3.96

0.04

2

117.30

115.79

1.29

119.42

0.35

1.81

122.05

4.28

4.05

3

208.40

199.88

4.09

210.46

0.54

0.99

210.93

4.95

1.22

4

254.00

245.76

3.25

253.37

0.41

0.25

258.94

4.81

1.94

5

445.00

387.53

12.92

435.71

0.63

2.09

410.33

4.74

7.79

TAE

_______

_______

5.37

_______

_____

1.53

_____

_____

3.01

1

Correlation

0.5
0
-0.5
-1
1
2
6

3
5

4

4

5

3
6

i

2
1



i

Figure 5.3: The correlation between the updating parameters

Figure 5.3 shows the correlation of the updating parameters using the DE-MC parameters where the majority of these
parameters are weakly correlated except the pairs (𝜃2 , 𝜃5 ) and (𝜃4 , 𝜃5 ), where the correlation between these pairs are
relatively high (<0.7%).

The evaluation of total average error after each accepted (or rejected) sample is illustrated in Figure 5.4. The result
indicates that the DE-MC has a fast convergence rate and was able converge after 500 iterations.

10

Total average error

10

10

10

10

10

10

10

0.9

0.8

0.7

0.6

0.5

0.4

0.3

0.2

0

500

1000

1500

2000

2500

3000

3500

4000

4500

5000

Iteration

Figure 5.4: The evaluation of the TAE using the DE-MC method

6.

Conclusion

In this paper, the Differential Evolution Markov Chain (DE-MC) algorithm is used to approximate the Bayesian
formulations in order to perform a finite element model updating procedure. In the DE-MC method, multi-chains are run
in parallel which allows the chains to learn from each other in order to improve the sampling process where the jumping
step depends on the difference between randomly selected chains. This method is investigated by updating two structural
system: the first one is a five DOF mass-spring linear system and the second one is the unsymmetrical H-shaped aluminum
structure. In the first case, the total average error was reduced from 1.98% to 0.012%, while in the second case, the FEM
updating of the unsymmetrical H-shaped structure, the total average error was reduced from 5.37% to 1.53%. Also the
DE-MC algorithm appeared to have better results than the M-H algorithm when the unsymmetrical H-shaped structure is
updated. In further work, the DE-MC algorithm will be modified and improved to include several steps such as the cross
over and the exchange between the parallel chains. These changes may further improve the sampling procedure.

References
[1]
[2]
[3]
[4]
[5]
[6]
[7]
[8]
[9]
[10]

Zienkiewicz, O.C., & Taylor, R.L. The finite element method, Vol. 3, London: McGraw-hill, 1977.
Hutton, D., Fundamentals of finite element analysis, McGraw-Hill, 2004.
Cook, R.D., Concepts and applications of finite element analysis. John Wiley & Sons, 2007.
Friswell, M. I. and J.E. Mottershead, Finite element model updating in structural dynamics. Vol. 38, Springer
Science & Business Media, 2013.
Marwala, T., Finite element model updating using computational intelligence techniques: applications to
structural dynamics, Springer Science & Business Media, 2010.
Ter Braak, C.J., A Markov Chain Monte Carlo version of the genetic algorithm Differential Evolution: easy
Bayesian computing for real parameter spaces. Statistics and Computing, 16(3), pp. 239-249, 2006.
ter Braak, C.J. and J.A. Vrugt, Differential evolution Markov chain with snooker updater and fewer chains.
Statistics and Computing, 18(4), pp. 435-446, 2008.
Storn, R. and K. Price, Differential evolution–a simple and efficient heuristic for global optimization over
continuous spaces. Journal of global optimization, 11(4), pp. 341-359, 1997.
Price, K., R.M. Storn, and J.A. Lampinen, Differential evolution: a practical approach to global optimization,
Springer Science & Business Media, 2006.
Bishop, C.M., Pattern recognition and machine learning, springer, 2006.

[11]

[12]
[13]

[14]
[15]
[16]
[17]
[18]
[19]

[20]
[21]

Boulkaibet, I., Marwala, T., Mthembu, L., Friswell, M. I., & Adhikari, S. Sampling techniques in Bayesian finite
element model updating, in Topics in Model Validation and Uncertainty Quantification, Vol 4. Springer. pp. 7583, 2012.
Yuen, K.-V., Bayesian methods for structural dynamics and civil engineering. John Wiley & Sons, 2010.
Boulkaibet, I., Mthembu, L., Marwala, T., Friswell, M. I., & Adhikari, S. Finite element model updating using
the shadow hybrid Monte Carlo technique. Mechanical Systems and Signal Processing, Vol. 52, pp. 115-132,
2015.
Boulkaibet, I., Mthembu, L., Marwala, T., Friswell, M. I., & Adhikari, S. Finite element model updating using
Hamiltonian Monte Carlo techniques. Inverse Problems in Science and Engineering, 2017. 25(7): p. 1042-1070.
Marwala, T. and S. Sibisi, Finite Element Model Updating Using Bayesian Framework and Modal Properties.
Journal of Aircraft, Vol. 42(1), pp. 275-278, 2005.
Bishop, C.M., Neural networks for pattern recognition. Oxford university press, 1995.
Ching, J. & S.-S. Leu, Bayesian updating of reliability of civil infrastructure facilities based on condition-state
data and fault-tree model. Reliability Engineering & System Safety, Vol. 94(12), pp. 1962-1974, 2009.
Metropolis, N., Rosenbluth, A. W., Rosenbluth, M. N., Teller, A. H., & Teller, E., Equation of state calculations
by fast computing machines. The journal of chemical physics, Vol. 21(6), pp. 1087-1092, 1953.
Boulkaibet, I., Mthembu, L., Marwala, T., Friswell, M. I., & Adhikari, S., Finite element model updating using
an evolutionary Markov Chain Monte Carlo algorithm, in Dynamics of Civil Structures, Vol. 2., Springer. pp.
245-253, 2015.
Marwala, T., I. Boulkaibet, and S. Adhikari, Probabilistic Finite Element Model Updating Using Bayesian
Statistics: Applications to Aeronautical and Mechanical Engineering. John Wiley & Sons, 2016.
Boulkaibet, I., Finite element model updating using Markov Chain Monte Carlo techniques, PhD thesis,
University of Johannesburg, 2014.

