Open AccessArticle

Global Reliability Sensitivity Analysis Based on Maximum Entropy and 2-Layer Polynomial Chaos Expansion

Jianyu Zhao

¹,

Shengkui Zeng

^1,2

Jianbin Guo

^1,2,* and

Shaohua Du

School of Reliability and Systems Engineering, Beihang University, Beijing 100191, China

Science and Technology on Reliability and Environmental Engineering Laboratory, Beihang University, Beijing 100191, China

CRRC ZIC Research Institute of Electrical Technology & Material Engineering, Zhuzhou 412001, China

Author to whom correspondence should be addressed.

Entropy 2018, 20(3), 202; https://doi.org/10.3390/e20030202

Submission received: 5 February 2018 / Revised: 13 March 2018 / Accepted: 14 March 2018 / Published: 16 March 2018

(This article belongs to the Special Issue Maximum Entropy and Its Application II)

Download

Browse Figures

Versions Notes

Abstract

To optimize contributions of uncertain input variables on the statistical parameter of given model, e.g., reliability, global reliability sensitivity analysis (GRSA) provides an appropriate tool to quantify the effects. However, it may be difficult to calculate global reliability sensitivity indices compared with the traditional global sensitivity indices of model output, because statistical parameters are more difficult to obtain, Monte Carlo simulation (MCS)-related methods seem to be the only ways for GRSA but they are usually computationally demanding. This paper presents a new non-MCS calculation to evaluate global reliability sensitivity indices. This method proposes: (i) a 2-layer polynomial chaos expansion (PCE) framework to solve the global reliability sensitivity indices; and (ii) an efficient method to build a surrogate model of the statistical parameter using the maximum entropy (ME) method with the moments provided by PCE. This method has a dramatically reduced computational cost compared with traditional approaches. Two examples are introduced to demonstrate the efficiency and accuracy of the proposed method. It also suggests that the important ranking of model output and associated failure probability may be different, which could help improve the understanding of the given model in further optimization design.

Keywords:

global reliability sensitivity analysis; polynomial chaos expansion; Sobol’s indices; the maximum entropy method

1. Introduction

Sensitivity analysis (SA) is one of the most important methodologies dealing with optimization in engineering practice. It aims to quantify importance in the output of a model with regards to some input parameters or design variables [1,2,3,4]. With the results of SA, important variables of given model are chosen for further study when less important variables may be treated as constant. This could help reduce the dimension of inputs (number of input variables) and improve the model understanding.

Generally, SA methodologies can be classified into local sensitivity analysis (LSA) and global sensitivity analysis (GSA) [1,5]. The former only concentrate on the influence of single input parameter on model output by gradient estimating in a small neighborhood around a normal point. It is often used to provide deterministic direction to be close to the optimal object [6]. Meanwhile, the latter one tries to quantify the output uncertainty when all input parameters are varying over entire value domain. And it does not depend on the choice of a nominal point [5,7]. Thus, GSA is suitable for uncertainty quantification (UQ) and reliability based design optimization (RBDO) [8].

Numerous methodologies are developed for GSA, among which variance-based methods are widely used [1,5]. Variance-based methods take the well-known variance decomposition formula V(Y) = V(E(Y|X)) + E(V(Y|X) as their foundation. This formula has variance separated by using the variance of conditional expectation and associated residual variance. Thus, the contribution or importance of X on Y could be further studied. Many researchers take Sobol’s indices because they could provide a clear, accurate and robust indication of the selected uncertain inputs both quantitatively and qualitatively [9]. Although Sobol’s indices are traditionally applied for non-statistical model output, some researchers have extended Sobol’s indices to identify uncertainty effects of inputs on reliability. Morio only studied the effects of some variable distribution parameters, such as variable means, but it did not consider the variable uncertainties themselves [10]. Borgonovo proposed the moment-independent global sensitivity indices which are derived from conditional probability density function (PDF) [11,12]. Lu accepted this idea and further developed the variance-based importance measurement for global reliability sensitivity analysis (GRSA) with associated global reliability sensitivity indices [13,14]. Cui and Li proved the equivalence relation between these global reliability sensitivity indices and Sobol’s indices when taking the indicator function as the assumed model output.

Unlike traditional global sensitivity indices calculation, the global reliability sensitivity indices aim at calculating importance of input variables on reliability that is a statistical parameter. These indices are specifically concerned about the variance of conditional expectation of indicator function according to variance decomposition formula. However, the indicator function maps every single model output to only two real numbers 0 or 1 discretely, so the analytical method and surrogate method, which are popular with uncertainty analysis, have not been developed for GRSA because it is difficult to obtain an analytical formula of the indicator function. In this way, the Crude Monte Carlo Simulation (MCS) is always employed to calculate such measurement [15,16], which leads that the computation cost is unacceptable when the given model is complex. Although Li et al. used the importance sampling technologies to improve sampling efficiency [13], quite a few samples are inevitable.

Fortunately, the conditional expectation of indicator function is equal to the conditional failure probability (CFP) that could be regarded as a continuous random variable. Thus, the variance of conditional expectation of indicator function could be replaced by the variance of CFP. In this way, if CFP could be assumed as the model output, there might be an unknown expression between CFP and corresponding conditional variable. So, the surrogate method could be used to approximate the CFP. Then, the statistical properties of CFP, such as mean and variance, may be calculated by surrogate model efficiently.

As for the issue of the statistical properties of a surrogate model, polynomial chaos expansions (PCE) technology provides a powerful tool for this application [17,18]. The PCE is a surrogate modelling technology with orthogonal polynomial basis functions, and it has statistical moments calculation simplified. In other words, variance information of surrogate model could be readily obtained by post-processing PCE’s coefficients analytically [19,20]. Therefore, it is suitable to make a surrogate model for CFP and calculate its variance. However, a challenge for PCE method in such application is that quite a few specific probability values are required. In addition, since the CFP is dependent on a set of conditional variables, each global reliability sensitivity index corresponding to different conditional variables requires a specific surrogate model. This may be time consuming when computing several indices. To circumvent this problem, the maximal entropy (ME) method is introduced since it could compute probabilities precisely and analytically only with prior moment information of model output [21,22,23,24]. Considering the orthogonal characteristics, the PCE could not only provide the first two order moments directly, but also provide higher order moments by PCE multiplication algorithm without additional sampling. In this way, the combination of PCE and ME could make GRSA efficiently. To the author’s knowledge, this idea has not been developed for global reliability sensitivity indices calculation.

Therefore, a new 2-layer PCE algorithm coupling PCE and ME is developed. That is, the first layer PCE is to generate a surrogate model by traditional algorithm. The second layer PCE is to construct several surrogate models for CFPs with the help of ME, which are referred as reliability polynomial chaos expansions (R-PCEs). The moment information that ME requires is provided by the first layer PCE when the conditional variables are given and fixed. Then, the variance of each CFP comes out by analytically post processing R-PCE coefficients, followed by global reliability sensitivity indices.

The presentation of the work is structured as follows: Section 2 provides a brief introduction of GSA and GRSA respectively. Section 3 presents the proposed 2-layer PCE framework and detailed calculation procedure. A numerical example and an engineering example are presented to test the performance of the proposed method in Section 4, while Section 5 summarizes the main contributions of this paper. The symbols used in this paper are listed in Table A1 in Appendix A.

2. Related Work

2.1. Global Sensitivity Analysis and Sobol’s Indices

In this section, the traditional variance based GSA is briefly reviewed [16,19]. Let us consider a model y = g(x) with n-dimensional independent identical distribution random input x = (x₁, …, x_n) ∈Ωⁿ and a scalar output y. Denote f_x(x) as the joint probability density function (PDF), then we have

f_{x} (x) = \prod_{i = 1}^{n} f_{x} (x_{i})

. In addition, we assume that

x \in Ω^{n} \mapsto g (x) \in L^{2} (Ω^{n}, f_{x} (x))

The Sobol’s decomposition of g(x) can be expressed as:

g (x) = g_{0} + \sum_{i = 1}^{n} g_{i} (x_{i}) + \sum_{1 \leq i < j \leq n} g_{i j} (x_{i}, x_{j}) + \dots + g_{1, 2, \dots, n} (x_{1}, \dots, x_{n}),

(1)

where:

\begin{array}{l} g_{0} = \int_{Ω^{n}} g (x) f_{x} (x) d x, \\ g_{i} (x_{i}) = \int_{Ω^{n - 1}} g (x) f_{x} (x_{~ [i]}) d x_{~ [i]} - g_{0}, \\ g_{i j} (x_{i}, x_{j}) = \int_{Ω^{n - 2}} g (x) f_{x} (x_{~ [i, j]}) d x_{~ [i, j]} - g_{i} (x_{i}) - g_{j} (x_{j}) - g_{0}, \end{array}

(2)

and so on x_~[i] in Equation (2) means the vector without the element of i, that is, x_~[i] = (x₁, …, x_i−1, x_i+1, …, x_n).

Furthermore, the expansion items above are orthogonal to each other, that is:

\int_{Ω^{n}} g_{1, 2, \dots, i_{s}} (x_{1}, x_{2}, \dots, x_{i_{s}}) g_{1, 2, \dots, j_{t}} (x_{1}, x_{2}, \dots, x_{j_{t}}) f_{x} (x) d x = 0 \forall (x_{1}, x_{2}, \dots, x_{i_{s}}) \neq (x_{1}, x_{2}, \dots, x_{j_{t}}) .

(3)

And the variance of output y = g(x) is defined as:

V (g (x)) = \int_{Ω^{n}} {(g (x) - g_{0})}^{2} f_{x} (x) d x .

(4)

Sobol’s proposed a variance-based representation to measure the importance of input variable. With Equations (1), (3) and (4), the Sobol’s indices are given as:

S_{i_{1}, \dots, i_{s}} = \frac{V (g_{i_{1}, \dots, i_{s}} (x_{i_{1}}, \dots, x_{i_{s}}))}{V (g (x))} = \frac{\int_{Ω^{s}} g_{i_{1}, \dots, i_{s}}^{2} (x_{i_{1}}, \dots, x_{i_{s}}) f_{x} (x_{i_{1}}, \dots, x_{i_{s}}) d x_{i_{1}} \dots x_{i_{s}}}{\int_{Ω^{n}} {(g (x) - g_{0})}^{2} f_{x} (x) d x},

(5)

where

V (g_{i_{1}, \dots, i_{s}} (x_{i_{1}}, \dots, x_{i_{s}}))

is called as the partial variance, which is distinguished from V(g(x)) above.

Generally, S_i is the first-order sensitivity index, S_ij is the second-order sensitivity index, etc. There are 2ⁿ − 1 Sobol’s sensitivity indices, and they are subjected to the following formula:

\sum_{i = 1}^{n} S_{i} + \sum_{1 \leq i < j \leq n} S_{i j} + \dots S_{1, 2, \dots, n} = 1 .

(6)

Moreover, since S_i represents the single contribution to the output uncertainty only due to an input x_i, it is also called as the main effect. And the total effect index

S_{T_{i}}

has been defined to measure the total contribution of an input x_i. It could be calculated as follows:

S_{T_{i}} = \sum_{L_{i}} S_{i_{1}, \dots, i_{s}} = 1 - S_{~ [i]}, \forall L_{i} = {(i_{1}, \dots, i_{s}) : \exists k, 1 \leq k \leq s, i_{k} = i},

(7)

where S_~[i] is the sum of all

S_{i_{1}, \dots, i_{s}}

without subscribe index i.

Traditionally, the Sobol’s indices could be solved by MCS [16], which may be computationally demanding. Sudret proposed a PCE-based method to mitigate this problem with much less computational costs [18].

PCE, which was originally introduced by Wiener [25], employs the Hermite polynomials in the random space to approximate the Gaussian stochastic processes [19,26]. In [20], Xiu and Karniadakis developed the PCE under Wiener-Askey scheme that could be applied to non-Gaussian scenarios. It can uniformly approximate any random process with finite second-order moments. Specifically, the PCE approximation of model y = g(x) is:

\begin{array}{l} y \approx g_{P C E} (ξ) = & c_{0} + \sum_{i = 1}^{n} \sum_{α \in φ_{i}} c_{α} ψ_{α} (ξ_{i}) + \sum_{1 \leq i_{1} < i_{2} \leq n} \sum_{α \in φ_{i_{1} i_{2}}} c_{α} ψ_{α} (ξ_{i_{1}}, ξ_{i_{2}}) + \dots + \\ \sum_{1 \leq i_{1} < \dots < i_{s} \leq n} \sum_{α \in φ_{i_{1} \dots i_{s}}} c_{α} ψ_{α} (ξ_{i_{1}}, \dots, ξ_{i_{s}}) + \dots + \sum_{α \in φ_{1, 2, \dots, n}} c_{α} ψ_{α} (ξ_{1}, \dots, ξ_{n}), \end{array}

(8)

where ξ = (ξ₁, ξ₂, …, ξ_n) is independent standardized orthogonal random variables associated with input x = (x₁, …, x_n), the subscript α is a tuple defined as α = (α₁, …, α_n), and

φ_{i_{1}, \dots, i_{s}}

is defined as a realization of α tuple in which only the indices {i₁, …, i_s} are nonzero:

φ_{i_{1}, \dots, i_{s}} = {α : \begin{matrix} α_{k} > 0 \forall k = 1, 2, \dots, n . k \in (i_{1}, \dots, i_{s}) \\ α_{k} = 0 \forall k = 1, 2, \dots, n . k \notin (i_{1}, \dots, i_{s}) \end{matrix}} .

(9)

Moreover,

ϕ_{k}

is the kth one-dimensional orthogonal polynomials from the Askey scheme.

ψ_{α} (ξ_{i_{1}}, \dots, ξ_{i_{s}})

is the multi-dimensional PCE basis, which is constructed by tensor products of the corresponding one-dimensional polynomials, namely,

ψ_{α} (ξ_{i_{1}}, \dots, ξ_{i_{s}}) = \prod_{k = 1}^{i_{s}} ϕ_{α_{k}} (ξ_{k})

. Besides, the highest degree of orthogonal polynomials in Equation (8) is a predefined PCE degree p. Therefore, we have |α| ≤ p. It should be noted that p could be determined by considering the balance of accuracy and efficiency. Specifically, PCEs with small degrees p and p + 1 could be obtained respectively, which is followed by comparing the difference (usually the first two order moments). If the difference is smaller than a pre-defined threshold, the PCE with p + 1 degree is adopted. Otherwise, the PCE with p + 2 degree is used to repeat this procedure.

Generally, Equation (8) can be rewritten as:

y \approx g_{P C E} (ξ) = \sum_{j = 0}^{N - 1} c_{j} ψ_{j} (ξ), f o r ξ = (ξ_{1}, \dots, ξ_{n}),

(10)

where coefficient c_j and expansion base ψ_j(ξ) are corresponding to Equation (8) sequentially,

N = (\begin{matrix} n + p \\ p \end{matrix}) = \frac{(n + p)!}{n! p!}

is the number of PCE items. Thus, the main effect index S_i and the total effect index

S_{T_{i}}

could be solved directly with PCE coefficients as:

S_{i} = \frac{\sum_{α \in φ_{i}} c_{α}^{2} E [ψ_{α}^{2}]}{\sum_{j = 1}^{N - 1} c_{j}^{2} E (ψ_{j}^{2} (ξ))},

(11)

S_{T_{i}} = \sum_{(i_{1}, \dots, i, \dots, i_{s}) \in φ_{i_{1}, \dots, i_{s}}} S_{i_{1}, \dots, i, \dots i_{s}} .

(12)

For more details see Appendix B.

2.2. Global Reliability Sensitivity Analysis

If the given model g(x) above is the limit-state function, the input x falls in the failure region if g(x) returns a negative value. Denote P_f as the failure probability, then we have:

\begin{array}{l} P_{f} & = P (g (x) < 0) \\ = \int \dots \int_{g (x <)} f_{x} (x_{1}, x_{2}, \dots, x_{n}) d x_{1} d x_{2} \dots d x_{n} \\ = \int \dots \int_{Ω^{n}} I_{F} f_{x} (x_{1}, x_{2}, \dots, x_{n}) d x_{1} d x_{2} \dots d x_{n} \\ = E [I_{F}], \end{array}

(13)

where

I_{F} = {\begin{matrix} 1, & g (x) \leq 0, \\ 0, & g (x) > 0 . \end{matrix}

is the indicator function of the failure domain, and E(I_F) is the mean of indicator function.

Similarly, the CFP is defined as follows:

P_{f | x_{I}} = E [I_{F} | x_{I}],

(14)

where x_I is the vector of conditional input variable with subscript I = {i₁, …, i_m} (1 ≤ i₁ ≤…≤ i_m ≤ n),

I_{F} | x_{I} = {\begin{matrix} 1, & g (x | x_{I}) \leq 0, \\ 0, & g (x | x_{I}) > 0 . \end{matrix}

is the conditional indicator function of the failure domain of the conditional limit-state function g(x|x_I) and it can be obtained by assuming x_I at a fixed realization value while other inputs are free.

Cui and Lü put forward the moment-based importance measure

δ_{I}^{P}

of the basic variable on the distribution as [14]:

δ_{I}^{P} = \int_{- \infty}^{+ \infty} {(P_{f} - P_{f | x_{I}})}^{2} f_{x} (x_{I}) d x_{I} = E {[P_{f} - P_{f | x_{I}}]}^{2},

(15)

where

δ_{I}^{P}

eflects the effect of conditional variable x_I on failure probability when x_I takes its value according to its PDF. Furthermore, Cui proved that:

δ_{I}^{P} = V [E (I_{F} | x_{I})] .

(16)

In this way,

δ_{I}^{P}

is equal to the Sobol’s index when I_F is taken as the model output. Therefore, Cui proposed that the main effect index and the total effect index of x_i can be defined as:

S_{i}^{(R)} = \frac{V (E (I_{F} | x_{i}))}{V (I_{F})},

(17)

S_{T_{i}}^{(R)} = 1 - \frac{V (E (I_{F} | x_{~ [i]}))}{V (I_{F})},

(18)

where V(I_F), V(E(I_F|x_i)) and V(E(I_F|x_~[i])) are expressed as follows:

V (I_{F}) = E (I_{F}^{2}) - E {(I_{F})}^{2} = E (I_{F}) - E {(I_{F})}^{2} = P_{f} - P_{f}^{2},

(19)

V (E (I_{F} | x_{i})) = E (E {(I_{F} | x_{i})}^{2}) - E {(E (I_{F} | x_{i}))}^{2} = E (E {(I_{F} | x_{i})}^{2}) - P_{f}^{2},

(20)

V (E (I_{F} | x_{~ [i]})) = E (E {(I_{F} | x_{~ [i]})}^{2}) - E {(E (I_{F} | x_{~ [i]}))}^{2} = E (E {(I_{F} | x_{~ [i]})}^{2}) - P_{f}^{2} .

(21)

The MCS method and associated improved Importance Sampling (IS) method for solving

S_{i}^{(R)}

and

S_{T_{i}}^{(R)}

are presented in [27], where the IS method uses the advanced first-order second-moment method [28] to help construct importance sampling density functions within MCS framework. However, they are both time-consuming due to the sampling principle, and the computational burden becomes worse when given model are complex. In this way, a new algorithm is proposed to mitigate this problem.

3. R-PCE Framework for GRSA

3.1. PCE-Based GRSA Method

Take the limit-state function in Section 2.2 into consideration. As mentioned above, the model output for global reliability sensitivity indices calculation is I_F, and it is difficult to compute I_F related value except for MCS. This is due to the fact that I_F only maps every single model output to two discrete real numbers 0 or 1. Either analytical methods or surrogate methods are efficient when they are applied to continuous model output rather than the discrete case. However, E(I_F|x_i) and E(I_F|x_~[i]) could be replaced by two CFP values P_f(x_i) and P_f(x_~[i]) respectively, which could be treated as continuous random variables. In this way, the key point of global reliability sensitivity index calculation is to compute variance values V(P_f(x_i)) and V(P_f(x_~[i])). Equations (17) and (18) are replaced by Equations (22) and (23):

S_{i}^{(R)} = \frac{V (E (I_{F} | x_{i}))}{V (I_{F})} = \frac{E (P_{f} {(x_{i})}^{2}) - P_{f}^{2}}{P_{f} - P_{f}^{2}},

(22)

S_{T_{i}}^{(R)} = 1 - \frac{V (E (I_{F} | x_{~ [i]}))}{V (I_{F})} = 1 - \frac{E (P_{f} {(x_{~ [i]})}^{2}) - P_{f}^{2}}{P_{f} - P_{f}^{2}} .

(23)

Since the PCE technology provides an efficient surrogate modeling tool to determinate the mean as well as the variance of the model output, a PCE-based method is developed for this implement. However, several probability values required during surrogate modeling procedure may cause computational demanding. So, the analytical ME method is integrated with PCE to solve this problem. Thus, the surrogate model of CFP could be developed, which is referred as R-PCE in this work.

The PCE-based GRSA method involves a 2-layer PCE surrogate model construction process. The first layer PCE is to generate traditional surrogate model approximating the limit-state function. It is used to obtain mean and variance when given conditional variables are fixed in specific values. And the second layer is to construct R-PCEs for different CFPs. To obtain probability response during the surrogate modeling of R-PCE, the ME method is taken to analytically compute these probability values with the first layer PCE. Besides, since each CFP is dependent on a different conditional variable, each global reliability sensitivity index has its own R-PCE formula. The algorithm framework is shown in Figure 1, and the whole procedure consists of four phases.

The first phase is to prepare inputs or data for GRSA.

In the second phase, the PCE technology is an efficient and accurate tool to approximate the given model with finite samples. Thus, the original model is replaced by the 1st-layer PCE in the following procedures. Besides, the first and second moment of the model output are obtained directly without traditional Mote Carlo sampling.

In the third phase, the PCE technology is developed to make the surrogate model for CFP. Unlike the surrogate model in the second phase, quite a few failure probabilities of given model require to be calculated. Unlike the model output which may be obtained by simulating directly, it is usually difficult to obtain these statistical values. Fortunately, this problem could be overcome by the ME method. According to ME, the probability distribution that best represents the current state of knowledge, which is usually presented in terms of the first n-order moments, is the one having the largest entropy [24]. Under this principle, the probability distribution of the model output could be analytically obtained with moment information based on the 1st-layer PCE. The details of this procedure are presented as follows.

First of all, for each global reliability sensitivity index, a set of conditional variable samples is generated by the probabilistic collocation method (PCM). Second, due to the regularity of orthogonal polynomials, the conditional variable is removed from the 1st-layer PCE. Also, the specific formula of the 1st-layer PCE with other variables could remain the form of linear combination of orthogonal polynomials, so does the statistical characteristic. For simplicity, we call it conditional PCE (C-PCE). Thus, the associated moment evaluations are also available. Third, the ME method is taken to solve failure probabilities with regards to specific conditional variable. Fourth, a set of conditional variable samples and associated conditional failure probabilities are prepared for following 2nd-layer PCE (R-PCE) construction, and this procedure is similar to PCE in the second phase. In this way, R-PCE corresponding to specific conditional variable is generated.

Finally, with the coefficients of R-PCE, we could obtain variance of each CFP and calculate global reliability sensitivity index. During the procedure, the number of sampling of both the 1st-layer PCE and the 2nd-layer PCE of proposed method is limited according to PCE sample principle. This makes proposed method efficient. The following section will provide details of how to construct the proposed 2-layer PCE.

3.2. Conditional Polynomial Chaos Expansion and Moments Calculation

Without loss of generality, PCE formula of given limit-state function is the same as Equation (8). If the value of ξ_I is given as the conditional variable, Equation (8) could be simplified as:

\begin{array}{l} g_{P C E} (ξ | ξ_{I}) = & {c^{'}}_{0} + \sum_{i = 1}^{n - | I |} \sum_{α \in φ_{i \notin I}} {c^{'}}_{α} ψ_{α} (ξ_{i}) \\ + \sum_{1 \leq i_{1} < \dots < i_{s} \leq n - | I |} \sum_{α \in φ_{{i_{1} \dots i_{s}} \notin I}} {c^{'}}_{α} ψ_{α} (ξ_{i_{1}}, \dots, ξ_{i_{s}}) + \dots \\ + \sum_{α \in φ_{{1, 2, \dots, n - | I |} \notin I}} {c^{'}}_{α} ψ_{α} (ξ_{1}, \dots, ξ_{n - | I |}), \end{array}

(24)

where |I| is the number of elements in subscript I, {i₁, …, i_s}∉I represents that

\forall i_{k} \in {i_{1}, \dots, i_{s}}, i_{k} \notin I

. The coefficients are subject to:

{c^{'}}_{0} = c_{0} + \sum_{\forall I^{'} \in I, I^{'} \notin \emptyset, α \in φ_{I^{'}}} c_{α} ψ_{α} (ξ_{I^{'}}) = c_{0} + \sum_{\forall I^{'} = {i_{1}, \dots, i_{k}} \in I, I^{'} \notin \emptyset, α \in φ_{I^{'}}} c_{α} ψ_{α} (ξ_{i_{1}}, \dots, ξ_{i_{k}})

(25)

and:

{c^{'}}_{α} ψ_{α} (ξ_{i_{1}}, \dots, ξ_{i_{s}}) = \sum_{\forall I^{'} \in I, I^{'} \notin \emptyset, α \in φ_{i_{1} \dots i_{s}, I^{'}}} c_{α} ψ_{α} (ξ_{i_{1}}, \dots, ξ_{i_{s}}, ξ_{I^{'}}) .

(26)

The total number of coefficients is

N^{'} = (\begin{matrix} n + p - | I | \\ p \end{matrix}) = \frac{(n + p - | I |)!}{(n - | I |)! p!}

. Equation (24) is called as C-PCE. Thus, the mean and variance could be obtained by the following equations:

E (g_{P C E} (ξ | ξ_{I})) = {c^{'}}_{0},

(27)

\begin{array}{l} V (g_{P C E} (ξ | ξ_{I})) = & \sum_{i = 1}^{n - | I |} \sum_{α \in φ_{i \notin I}} {({c^{'}}_{α})}^{2} E (ψ_{α}^{2} (ξ_{i})) + \dots \\ + \sum_{1 \leq i_{1} < \dots < i_{s} \leq n - | I |} \sum_{α \in φ_{{i_{1} \dots i_{s}} \notin I}} {({c^{'}}_{α})}^{2} E (ψ_{α}^{2} (ξ_{i_{1}}, \dots, ξ_{i_{s}})) + \dots \\ + \sum_{α \in φ_{{1, 2, \dots, n - | I |} \notin I}} {({c^{'}}_{α})}^{2} E (ψ_{α}^{2} (ξ_{1}, \dots, ξ_{n - | I |})) . \end{array}

(28)

To simplify the expression, C-PCE could be rewritten as:

g_{P C E} (ξ | ξ_{I}) = \sum_{j = 0}^{N^{'} - 1} {c^{'}}_{j} ψ_{α} (ξ | ξ_{I}),

(29)

where

{c^{'}}_{j}

and ψ_α(ξ|ξ_I) are corresponding to Equation (28) sequentially, notation “ξ|ξ_I” in ψ_α(ξ|ξ_I) means that expansion base does not include conditional variable ξ_I.

The mean and variance of C-PCE are given as Equations (30) and (31) respectively:

E (g_{P C E} (ξ | ξ_{I})) = {c^{'}}_{0},

(30)

V (g_{P C E} (ξ | ξ_{I})) = \sum_{j = 1}^{N^{'} - 1} {c^{'}}_{j}^{2} E (ψ_{j}^{2} (ξ | ξ_{I})) .

(31)

Specifically, for main effect indices, since we have I = {i}, Equation (29) can be further rewritten as:

g_{P C E} (ξ_{~ [i]}) = \sum_{j = 0}^{N^{'} - 1} {c^{'}}_{j} ψ_{α} (ξ_{~ [i]}) .

(32)

and for total effect indices, since we have I = {1, …, n}/{i}, Equation (29) can be further rewritten as:

g_{P C E} (ξ_{i}) = \sum_{j = 0}^{N^{'} - 1} {c^{'}}_{j} ϕ_{α} (ξ_{i}) .

(33)

The C-PCE could provide the first two order moments accurately. However, there are some difficulties for high order moment calculation. Consider the fact that the C-PCE consists of orthogonal bases, the products of two C-PCEs could also be expanded as a linear combination of orthogonal bases [29]. Therefore, the orthogonal characteristic of PCE multiplication remains. The generalization form of orthogonal bases are presented in [29].

Take the C-PCE with Hermite polynomial bases for example, suppose two C-PCEs

u = \sum_{| α | \leq p_{α}} {c^{'}}_{α} ψ_{α} (ξ | ξ_{I})

and

v = \sum_{| β | \leq p_{β}} {c^{'}}_{β} ψ_{β} (ξ | ξ_{I})

. If

E ({| u v |}^{2}) < \infty

, then the product of u and v has the PCE formula as:

u v = \sum_{| θ | \leq p_{α} + p_{β}} \sum_{0 \leq β \leq θ} \sum_{\begin{array}{l} | θ - β + r | \leq p_{α}, \\ | β + r | \leq p_{β} \end{array}} C (θ, β, r) u_{θ - β + r} v_{β + r} ψ_{θ} (ξ),

(34)

where:

C (θ, β, r) = {[(\begin{matrix} θ - β + r \\ r \end{matrix}) (\begin{matrix} β + r \\ r \end{matrix}) (\begin{matrix} θ \\ θ - β \end{matrix})]}^{\frac{1}{2}} .

(35)

The subscripts α, β, r, and θ in Equations (34) and (35) are tuples associated with C-PCE terms, e.g., α = (α₁, α₂, …, α_n). We say β ≤ θ if β_i ≤ θ_i for all i = 1, 2, …, n. The operation of these subscripts, such as + or −, is also defined as component-wise. Especially, the factorial of tuples is defined like α! = ∏_i α_i!.

Then, the mean of uv could be obtained as:

E (u v) = \sum_{| r | \leq \min (p_{α}, p_{β})} C (θ = 0, β = 0, r) u_{r} v_{r} = \sum_{| r | \leq \min (p_{α}, p_{β})} u_{r} v_{r} .

(36)

Based upon the preparation above, we could calculate the high-order moments of PCE by replacing u and v with two PCEs respectively. That is, u and v in Equation (34) could be replaced by

g_{P C E} (ξ | ξ_{I})

g_{P C E}^{2} (ξ | ξ_{I})

g_{P C E}^{3} (ξ | ξ_{I})

, …, respectively as Equation (37):

\begin{array}{l} E (y^{2} (ξ | ξ_{I})) \approx E (g_{P C E}^{2} (ξ | ξ_{I})) = E (g_{P C E} (ξ | ξ_{I}) \cdot g_{P C E} (ξ | ξ_{I})) \\ E (y^{3} (ξ | ξ_{I})) \approx E (g_{P C E}^{3} (ξ | ξ_{I})) = E (y_{p}^{2} (ξ | ξ_{I}) \cdot g_{P C E} (ξ | ξ_{I})) \\ ⋮ \\ E (y^{k} (ξ | ξ_{I})) \approx E (g_{P C E}^{k} (ξ | ξ_{I})) = E (g_{P C E}^{⌊ \frac{k}{2} ⌋} (ξ | ξ_{I}) \cdot g_{P C E}^{⌈ \frac{k}{2} ⌉} (ξ | ξ_{I})), \end{array}

(37)

where ⌊•⌋ is rounded down and ⌈•⌉ is rounded up. Thus, the high order moments could be computed analytically.

3.3. The Maximum Entropy Method

Let f_y(y) be the unknown PDF of the output of limit-state function and assume H(f_y(y)) to be the entropy of f_y(y). Then we have:

H (f_{y} (y)) = - \int_{y \in Ω_{y}} f_{y} (y) \ln (f_{y} (y)) d y,

(38)

where Ω_y is the domain of y. According to the principle of ME, the failure probability distribution which best represents the current information, such as the moments, is the one with the largest entropy [24]. Thus, the ME method for estimating f_y(y) can be formulated as:

\begin{array}{l} \max_{f_{y} (\cdot)} H (f_{y} (y)) \\ s . t . \int_{y \in Ω_{y}} f_{y} (y) d y = 1 \\ \int_{y \in Ω_{y}} y^{k} f_{y} (y) d y = μ_{k}, k = 1, 2, \dots, N_{M E}, \end{array}

(39)

where N_ME is the number of moments, and μ_k is the kth moment, i.e., μ₁ = E(y), μ₂ = E(y²), etc.

It can be proved that the optimal solution for the PDF formula is given as:

f_{y} (y) = \exp (- \sum_{i = 0}^{N_{M E}} λ_{i} y^{i}) .

(40)

Thus, λ_i in Equation (40) should satisfy the following equation:

\begin{array}{l} \int_{y \in Ω_{y}} e^{- \sum_{i = 0}^{N} λ_{i}} d y = 1 \\ \int_{y \in Ω_{y}} y^{k} e^{- \sum_{i = 0}^{N} λ_{i}} d y = μ_{k}; k = 1, 2, \dots, N . \end{array}

(41)

Specifically, if the first two moments of the output are available, we have:

λ_{0}^{*} = \frac{E {(y)}^{2}}{2 V (y)} + \frac{1}{2} \ln (2 π V (y)), λ_{1}^{*} = - \frac{E (y)}{V (y)}, λ_{2}^{*} = \frac{1}{2 V (y)},

(42)

where

V (y) = E (y^{2}) - E {(y)}^{2}

and the failure probability can be calculated as:

P_{f} = \frac{1}{2} \frac{\sqrt{π} \exp (\frac{1}{4} \frac{{(λ_{1}^{*})}^{2}}{λ_{2}^{*}}) (e r f (\frac{1}{2} \frac{λ_{1}^{*}}{\sqrt{λ_{2}^{*}}}) + 1) \exp (- λ_{0}^{*})}{\sqrt{λ_{2}^{*}}} .

(43)

In the proposed method, if taking Equations (27) and (28) into Equations (42) and (43), the failure probability of C-PCE could be obtained directly without additional sampling.

For more general cases, the first four moments are usually used for PDF calculation because they are corresponding to mean, variance, skewness, kurtosis respectively, which are clear in physics meaning [30]. The 3rd moment and the 4th moment could be obtained efficiently by Equation (37) instead of traditional sampling method. Besides, since the ME method may not have the closed form solution, some optimization techniques, such as Newton method or its improved versions [31,32], could be utilized to solve the failure probability.

3.4. Global Reliability Sensitivity Indices Calculation

Let

{ξ_{I}^{(i)}}_{i = 1}^{M_{R}}

denote a set of conditional variable samples, where M_R is the number of samples. As described above, when the PCE for limit-state function is given, if taking each conditional variable sample into Equation (24), it returns the mean and variance of C-PCE. Then, the CFP is calculated by ME method. In this way, we could get a set of surrogate model relationship between

{{[\begin{matrix} ψ_{0} (ξ_{I}^{(i)}) & \dots & ψ_{N} (ξ_{I}^{(i)}) \end{matrix}]}^{T}}_{i = 1}^{M_{R}}

and

{P_{f} (ξ_{I}^{(i)})}_{i = 1}^{M_{R}}

. Similar to the traditional PCE procedure as Equation (A7) in Appendix B, we could solve the coefficients of R-PCE. Then, the result of R-PCE is:

P_{f | X_{I}} = E [I_{F} | x_{I}] = P_{f} (x_{I}) = P_{f} (ξ_{I}) = \sum_{j = 0}^{N_{R} - 1} c_{R, j} ψ_{j} (ξ_{I}),

(44)

where c_R,j is jth coefficient of R-PCE, and the subscript “R” means the coefficient is used for GRSA. ξ_I is independent standardized random variables corresponding to original variables x_I, N_R is the total number of coefficients, which is determined by:

N_{R} = \frac{(| I | + p_{R})!}{| I |! p_{R}!},

(45)

where p_R is the predefined degree of R-PCE. Then, the CFP is:

P_{f} = E (P_{f} (ξ_{I})) = c_{R, 0} .

(46)

And the variance of CFP is:

V (P_{f} (ξ_{I})) = \sum_{j = 1}^{N_{R} - 1} c_{R, j}^{2} E (ψ_{j}^{2} (ξ_{I})) .

(47)

Substituting Equations (46) and (47) into Equations (22) and (23), we could calculate the main effect and the total effect of global reliability sensitivity indices.

3.5. Implementation

For a given limit-state function g = g(x), we assume here that the input random variables are independent. With preparation mentioned above, the complete algorithm description is introduced as follows:

Step 1:: Represent distributions of the input into the associated standardized random variable. If the distribution types are not unified, we could map different distributions into the normal distribution. The following computation will be based upon standardized random variables.
Step 2:: Compute PCE for limit-state function by Equation (A7) in Appendix B, and it refers as g_PCE(ξ).
Step 3:: Based upon PCE above, compute R-PCE. First of all, generate a set of conditional variable samples. Then, integrate g_PCE(ξ) and ME method to obtain the CFP value P_f(ξ|ξ_I), where I = {i} for main effect indices and I = {1, 2, ..., n}/{i} for total effect indices. Specifically, for each sample set ${ξ_{I}^{(j)}}$ , substitute it into Equations (24)–(26) to obtain C-PCE $g_{P C E}^{(j)} (ξ | ξ_{I})$ . With the moments given by C-PCE, we could calculate the $P_{f} (ξ_{I}^{(j)})$ with ME method analytically. Third, the R-PCE comes out as Equation (44).
Step 4:: Compute the global reliability sensitivity indices by Equation (22), Equation (23) and Equations (46) and (47).

The step-by-step procedure of algorithm for the main effect indices and the total effect indices by the proposed 2-layer PCE method are presented in Algorithms A1 and A2 in Appendix C respectively. The total evaluations of the limit-state function is

\frac{2 (n + p)!}{n! p!} + 1

, where n and p are input variable dimension and the PCE degree respectively. If taking times of PCE evaluation in R-PCE construction process into consideration, the computation costs are

\frac{2 (n + p)!}{n! p!} + 2 n p_{R} + \frac{2 (n + p_{R} - 1)! n}{(n - 1)! p_{R}!} + 3 + n

, where p_R is the degree of R-PCE. This computation burden is usually 1–2 times less than IS procedure and 2–3 times less than MCS procedure, which could be seen in the following test examples. And the efficiency of the proposed method could be further improved by adopting sparse algorithms [6,17,33]. However, such work is not the scope of this paper.

4. Test Examples

Implementation of the proposed method for global reliability sensitivity indices calculation is illustrated by two test examples. The first example is a numerical example, while the second one is a structural analysis of a bridge crane. The results obtained by the proposed method are compared with solutions obtained by the full-scale MCS and the IS method from the perspective of accuracy and efficiency respectively.

4.1. A Numerical Example

Consider a limit-state function:

g (x) = 10 - 18 x_{1} + x_{2}^{2} + x_{2} + x_{3}^{2} + 5 x_{3},

(48)

where the input variables follow the normal distribution or the lognormal distribution, i.e., lnx₁~

N

(0.5,0.1), lnx₂~

N

(0.6,0.1), and lnx₃~

N

(0.5,0.1). The proposed method is employed to compute the main effect indices

S_{i}^{(R)}

and total effect indices

S_{T_{i}}^{(R)}

, when the PCE degree for the limit-state function is 2, and the R-PCE degree for CFP is 7. The results are compared with the results obtained by the full-scale MCS method and the IS method in Table 1. Moreover, the total times of limit-state function evaluation by different methods are also shown in the last column. It could be seen that the relative error of both methods are less than 3%, and about half of the results calculated by the proposed method is better than those obtained by the IS method. Moreover, the computational cost of proposed method is reduced sufficiently.

Also, the main effect indices and the total effect indices of GRSA are schematically shown in Figure 2. We could see that both the main effect indices

S_{i}^{(R)}

and the total effect indices

S_{T_{i}}^{(R)}

produce the same importance ranking: x₁ > x₃ > x₂, and there are strong interaction effects among these input variables since the total effects are much larger than the main effects. This suggests that the risk may be reduced by decreasing the uncertainty of x₁.

4.2. An Engineering Example

The second example is a bridge crane with respect to its rigidity design in Figure 3. The crane is made by universal beam with movable sheaves and hook. When the weight is in the middle of the beam, the maximal deflection

ω_{\max}

of such a crane beam could be computed as:

ω_{\max} (q, E, I, L, F_{p}) = \frac{L^{3}}{E I} (\frac{F_{p}}{48} + \frac{5 q L}{384}),

(49)

where q is the magnitude of a normally distributed load that represents the dead-weight of crane beam, E and I denote Young’s modulus and inertia moment of universal beam respectively, L is the span length, and F_p is maximal weight including listed load, sheaves and hook. Two distribution types are taken for this example, namely, the normal distribution and the lognormal distribution. Besides, two different variances are also considered to test this example, which are labeled as case1 and case 2 respectively. The detail distribution parameters are given in Table 2.

According to design requirement, the admissible ratio of deflection to span length is lower than 0.0018. Thus, the limit-state function is given as:

g (q, E, I, L, F_{p}) = 0.0018 - \frac{ω_{\max}}{L} = 0.0018 - \frac{L^{2}}{E I} (\frac{F_{p}}{48} + \frac{5 q L}{384}) .

(50)

The results obtained by the full-scale MCS method, the IS method and the proposed method are given in Table 3 and Table 4. The full-scale MSC results with large sample size are seen as the accurate results. In the calculation procedure of the proposed method, the PCE degree for the limit-state function and the PCE degree for CFP are 5 and 7 respectively.

It is shown that both the IS method and the proposed method could approximate

S_{i}^{(R)}

and

S_{T_{i}}^{(R)}

with much less computational burden, and the latter one presents better efficiency than the former one. In case 1, the average relative error of global reliability sensitivity indices obtained by proposed method is about 0.0196, while that obtained by the IS method is about 0.0283. In case 2, the average relative error of global reliability sensitivity indices obtained by proposed method is about 0.0206, while that obtained by the IS method is about 0.0295. The computational costs of our method are about two orders of magnitude less than the IS method, and four orders of magnitude less than the MCS method. Thus, the proposed method may be better than the other two methods.

Moreover, the global reliability sensitivity indices are also schematically shown in Figure 4a and Figure 5a. In case 1, the importance ranking based upon global reliability sensitivity indices

S_{i}^{(R)}

and

S_{T_{i}}^{(R)}

are given as L > I > F_p > E > q. It indicates that the reliability could be increase mostly by decreasing the uncertainty of the span length. Besides, since the total effect indices of input variables are much higher than their main effect indices, there are strong interaction effects among these variables. In case 2, the importance ranking based upon global reliability sensitivity indices

S_{i}^{(R)}

and

S_{T_{i}}^{(R)}

is given as E > L > I > F_p > q. Obviously, the influence of Young’s modulus on failure probability is higher than other random variables, which is different from case 1. Moreover, strong interaction effects also exist in this case. However, the interaction effects make more contributions than those in case 1 even when input variances are smaller compared with those of case 1. Therefore, the importance ranking should be paid more attention during further optimization design.

The global sensitivity indices of maximal deflection of cases 1 and 2 are also computed in Table 5 by traditional GSA method mentioned in Section 2.1, and they are plotted in Figure 4b and Figure 5b respectively. It is shown that the results of global sensitivity indices and global reliability sensitivity indices may be different. First of all, although the importance ranking based on GSA of case 1 is the same as that based on GRSA, the importance ranking obtained by GSA of case 2 is L > E > I > Fp > q, which is different from that obtained by GRSA. Secondly, the total effect indices of both cases are similar to the main effect indices, which means that there are few interaction effects among input variables. This is also different from the results of GRSA.

Comparing the results of GRSA with those of GSA, it suggests that the importance of input variables on the failure probability is not equivalent to that on the model output. So, both global reliability sensitivity indices and global sensitivity indices may be considered in engineering practice when choosing importance variables. And the proposed method could provide a useful tool for such purpose.

5. Conclusions

Global reliability sensitivity analysis aims at measuring the importance statistical parameter with regards to the input variable. The variance of expectation of conditional indicator function in global reliability sensitivity indices is difficult to calculate. This paper develops a PCE-based GRSA method to solve the problem. The proposed method takes conditional failure probability to replace the expectation of conditional indicator function, and then utilizes a 2-layer PCE construction technology to obtain the variance of conditional failure probability. The efficiency of proposed method is illustrated by a numerical example and an engineering model example. Both main effect indices and total effect indices present good performance in accuracy and significantly less computational burden in limit-state function evaluations compared with the traditional full-scale MCS method and IS method. Through global reliability sensitivity analysis, we could determine which uncertainties of inputs variables should be controlled and pave the way for further design improvement or optimization.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (NSFC: 61304218) and Beijing Natural Science Foundation (BNSF: 3153027). The authors would like to extend their sincere thanks to anonymous reviewers for their valuable comments.

Author Contributions

Jianbin Guo conceived the framework and structured the whole paper; Jianyu Zhao performed the experiments, wrote and revised the paper; Shengkui Zeng checked the results and revised the paper; Shaohua Du helped write programs.

Conflicts of Interest

The authors declare no conflict of interest.

Appendix A

Table A1. Notation.

$n$	the dimension of a variable
$x$	$x = (x_{1}, x_{2}, \dots, x_{n})$ , the input of given model
$ξ$	$ξ = (ξ_{1}, ξ_{2}, \dots, ξ_{n})$ , the independent standardized orthogonal random variable
$y$	the output of given model
$x_{~ [i]} (ξ_{~ [i]})$	the vector without the ith element
$f_{x} (x_{i})$	the PDF of x_i
$f_{x} (x)$	$f_{x} (x) = \prod_{i = 1}^{n} f_{x} (x_{i})$ , the joint PDF of $x = (x_{1}, x_{2}, \dots, x_{n})$
$f_{y} (y)$	the PDF of y
$E (\cdot)$	expectation
$V (\cdot)$	variance
$S_{T_{i}}$	the main effect of GSA
$S_{i}$	the total effect of GSA
$S_{i}^{(R)}$	the main effect of GRSA
$S_{T_{i}}^{(R)}$	the total effect of GRSA
$p$	the degree of PCE
$p_{R}$	the degree of R-PCE
$c_{j}$	the coefficient of jth PCE item
${c^{'}}_{j}$	the coefficient of jth CPCE item
$c_{R, j}$	the coefficient of jth R-PCE item
$N$	the number of PCE items
$N^{'}$	the number of CPCE items
$N_{R}$	the number of R-PCE items
$ϕ_{k}$	the kth one-dimensional orthogonal polynomial
$ψ_{α} (ξ)$	$ψ_{α} (ξ) = \prod_{k = 1}^{n} ϕ_{α_{k}} (ξ_{k})$ where $\sum_{k = 1}^{j} α_{k} \leq p$ and $α = (α_{1}, α_{2}, \dots, α_{n})$ , the base of PCE
$I_{F}$	the indicator function
$P_{f}$	the failure probability
$λ_{i}$	the unknown parameter of ME formula
$N_{M E}$	the number of unknown parameter of ME formula
$μ_{k}$	the kth moment of $y$

Appendix B

Consider the model y = g(x) in Section 2.1 again. Then, let

{ξ_{i}}_{i = 0}^{\infty}

be independent standardized orthogonal random variables associated with

{x_{i}}_{i = 1}^{n}

. For example, if x_i~

N

(2,0.5), then we have x_i = 0.5ξ_i + 2 where ξ_i is subject to the standardized normal distribution

N

(0,1).

Then, the PCE of the limit-state function is:

y = g_{P C E} (ξ) = c_{0} ψ_{0} + \sum_{i_{1} = 1}^{n} c_{i_{1}} ψ_{1} (ξ_{i_{1}}) + \sum_{i_{1} = 1}^{n} \sum_{i_{2} = 1}^{i_{1}} c_{i_{1} i_{2}} ψ_{2} (ξ_{i_{1}}, ξ_{i_{2}}) + \sum_{i_{1} = 1}^{n} \sum_{i_{2} = 1}^{i_{1}} \sum_{i_{3} = 1}^{i_{2}} c_{i_{1} i_{2} i_{3}} ψ_{3} (ξ_{i_{1}}, ξ_{i_{2}}, ξ_{i_{3}}) + \dots,

(A1)

where

{c_{j}}_{j = 0}^{\infty}

is the coefficient, and expansion base

ψ_{s} (ξ_{i_{1}}, \dots, ξ_{i_{s}})

denotes the sth multi-dimensional hyper-geometric polynomials in terms of the multi-dimensional standardized random variable

ξ = {(ξ_{1}, \dots, ξ_{n})}^{T}

. It is defined as tensor products of the corresponding one-dimensional orthogonal polynomials

{ϕ_{k}}_{k = 0}^{\infty}

ψ_{s} (ξ_{1}, ξ_{2}, \dots, ξ_{s}) = \prod_{k = 1}^{s} ϕ_{α_{k}} (ξ_{k}),

(A2)

where ϕ_i is an one-dimensional orthogonal basis with orthogonality relation

〈 ϕ_{i}, ϕ_{j} 〉 = 〈 ϕ_{i}^{2} 〉 δ_{i j},

(A3)

where δ_ij is the Kronecker delta and

〈 \cdot, \cdot 〉

denotes the ensemble average which is the inner product in Hilbert space. And α_k is the vector of index, which is subjected to a nonnegative integer.

The type of polynomials depends on the distributed type of input variables, e.g., the orthogonal polynomials corresponding to normal distribution are Hermite polynomials. For more details about interpretation of polynomial terms see [32]. It should be noted that, if the distribution types of input variables are not the same, Isukapalli [26] provided several distribution transformation approaches that map different distributions as functions of normal random variables. With this pre-processing, the problem of inconsistent distribution could be overcome.

In engineering practice, the PCE is truncated to finite terms. So, considering an n-dimensional orthogonal polynomial with degree not exceeding p, Equation (A1) can be rewritten as another form with limited input uncertainties, that is,

\begin{array}{l} y \approx g_{P C E} (ξ) = & c_{0} + \sum_{i = 1}^{n} \sum_{α \in φ_{i}} c_{α} ψ_{α} (ξ_{i}) + \sum_{1 \leq i_{1} < i_{2} \leq n} \sum_{α \in φ_{i_{1} i_{2}}} c_{α} ψ_{α} (ξ_{i_{1}}, ξ_{i_{2}}) + \dots \\ + \sum_{1 \leq i_{1} < \dots < i_{s} \leq n} \sum_{α \in φ_{i_{1} \dots i_{s}}} c_{α} ψ_{α} (ξ_{i_{1}}, \dots, ξ_{i_{s}}) + \dots + \sum_{α \in φ_{1, 2, \dots, n}} c_{α} ψ_{α} (ξ_{1}, \dots, ξ_{n}), \end{array}

(A4)

where the subscript α is a tuple defined as α = (α₁, …, α_n), and

φ_{i_{1}, \dots, i_{s}}

is defined as a realization of α tuple so that only the indices {i₁, …, i_s} are nonzero:

φ_{i_{1}, \dots, i_{s}} = {α : \begin{matrix} α_{k} > 0 \forall k = 1, 2, \dots, n . k \in (i_{1}, \dots, i_{s}) \\ α_{k} = 0 \forall k = 1, 2, \dots, n . k \notin (i_{1}, \dots, i_{s}) \end{matrix}} .

(A5)

Correspondingly, the expansion base

ψ_{α} (ξ) = \prod_{k = 1}^{n} ϕ_{α_{k}} (ξ_{k})

, where

\sum_{k = 1}^{j} α_{k} \leq p

. Denote N as the total number of polynomials, then we have

N = (\begin{matrix} n + p \\ p \end{matrix}) = \frac{(n + p)!}{n! p!} .

Generally, Equation (A4) can be further simplified as

y \approx g_{P C E} (ξ) = \sum_{j = 0}^{N - 1} c_{j} ψ_{j} (ξ), f o r ξ = (ξ_{1}, \dots, ξ_{n}),

(A6)

where coefficient c_j and expansion base ψ_j(ξ) are corresponding to Equation (A4) sequentially.

Let

{ξ^{i}}_{i = 1}^{M}

denote a set of the random variable samples, and it could be determined by the probabilistic collocation method (PCM). Again Let

{g (ξ^{i})}_{i = 1}^{M}

denote the corresponding set of model output, where M is the number of samples. Denote

c = {(c_{0}, \dots, c_{N - 1})}^{T}

, thus an approximation

\hat{c}

is given by:

\hat{c} = \arg \min_{c} \sum_{i = 1}^{M} {(g (ξ^{i}) - \sum_{j = 0}^{N} c_{j} ψ_{j} (ξ^{i}))}^{2} = {(ψ^{T} ψ)}^{- 1} ψ^{T} y (ξ^{1 : M}),

(A7)

where

ψ = {{[\begin{matrix} ψ_{0} (ξ^{(i)}) & \dots & ψ_{N} (ξ^{(i)}) \end{matrix}]}^{T}}_{i = 1}^{M}

denotes the coefficient matrix and

y (ξ^{1 : M}) = {(g (ξ^{1}), \dots, g (ξ^{M}))}^{T}

. Since

ψ^{T} ψ

may be ill-conditioned, M is suggested to be selected as M = 2(N + 1) [18].

Due to the orthogonality of the basis, the mean value and the variance of y in Equation (A6) can be calculated as:

\begin{array}{l} E (y) = c_{0}, \\ V (y) = \sum_{j = 1}^{N - 1} c_{j}^{2} E (ψ_{j}^{2} (ξ)) . \end{array}

(A8)

It could be proved that the PCE-based main effect index S_i and total effect index

S_{T_{i}}

are derived as [6]:

S_{i} = \frac{\sum_{α \in φ_{i}} c_{α}^{2} E [ψ_{α}^{2}]}{V (y (ξ))}

(A9)

and:

S_{T_{i}} = \sum_{(i_{1}, \dots, i, \dots, i_{s}) \in φ_{i_{1}, \dots, i_{s}}} S_{i_{1}, \dots, i, \dots i_{s}} .

(A10)

Therefore, the PCE could obtain the expectation and variance of the model output efficiently. Moreover, compared with the traditional MCS method which requires a large number of samples [16], the PCE method could calculate Sobol’s indices with quit a few computation costs. In this way, PCE is widely employed for GSA application.

Appendix C

Algorithm A1: Algorithm for main effect indices using proposed method.

Algorithm A2: Algorithm for Total Effect Indices Using Proposed Method.

References

Borgonovo, E.; Plischke, E. Sensitivity analysis: A review of recent advances. Eur. J. Oper. Res. 2016, 248, 869–887. [Google Scholar] [CrossRef]
Li, L.; Lu, Z. A new algorithm for importance analysis of the inputs with distribution parameter uncertainty. Int. J. Syst. Sci. 2016, 47, 3065–3077. [Google Scholar] [CrossRef]
Song, X.; Zhang, J.; Zhan, C.; Xuan, Y.; Ye, M.; Xu, C. Global sensitivity analysis in hydrological modeling: Review of concepts, methods, theoretical framework, and applications. J. Hydrol. 2015, 523, 739–757. [Google Scholar] [CrossRef]
Jensen, H.A.; Mayorga, F.; Valdebenito, M.A. Reliability sensitivity estimation of nonlinear structural systems under stochastic excitation: A simulation-based approach. Comput. Methods Appl. Mech. Eng. 2015, 289, 1–23. [Google Scholar] [CrossRef]
Saltelli, A.; Ratto, M.; Andres, T.; Campolongo, F.; Cariboni, J.; Gatelli, D.; Saisana, M.; Tarantola, S. Global Sensitivity Analysis: The Primer; John Wiley & Sons: Chichester, UK, 2008. [Google Scholar]
Blatman, G.; Sudret, B. Efficient computation of global sensitivity indices using sparse polynomial chaos expansions. Reliab. Eng. Syst. Saf. 2010, 95, 1216–1229. [Google Scholar] [CrossRef]
Borgonovo, E. Measuring Uncertainty Importance: Investigation and Comparison of Alternative Approaches. Risk Anal. 2006, 26, 1349–1361. [Google Scholar] [CrossRef] [PubMed]
Park, H.; Chung, J.; Behdinan, K.; Lee, J. Multidisciplinary wing design optimization considering global sensitivity and uncertainty of approximation models. J. Mech. Sci. Technol. 2014, 28, 2231–2242. [Google Scholar] [CrossRef]
Reuter, U.; Liebscher, M.; Müllerschön, H. Global sensitivity analysis in structural optimization. In Proceedings of the 7th European LS-DYNA Conference, Salzburg, Austria, 14–15 May 2009. [Google Scholar]
Morio, J. Influence of input PDF parameters of a model on a failure probability estimation. Simul. Model. Pract. Theory 2011, 19, 2244–2255. [Google Scholar] [CrossRef]
Borgonovo, E.; Castaings, W.; Tarantola, S. Model emulation and moment-independent sensitivity analysis: An application to environmental modelling. Environ. Model. Softw. 2012, 34, 105–115. [Google Scholar] [CrossRef]
Borgonovo, E.; Castaings, W.; Tarantola, S. Moment independent importance measures: New results and analytical test cases. Risk Anal. 2011, 31, 404–428. [Google Scholar] [CrossRef] [PubMed]
Li, L.; Lu, Z.; Feng, J.; Wang, B. Moment-independent importance measure of basic variable and its state dependent parameter solution. Struct. Saf. 2012, 38, 40–47. [Google Scholar] [CrossRef]
Cui, L.; Lü, Z.; Zhao, X. Moment-independent importance measure of basic random variable and its probability density evolution solution. Sci. China Technol. Sci. 2010, 53, 1138–1145. [Google Scholar] [CrossRef]
Wei, P.; Lu, Z.; Yuan, X. Monte Carlo simulation for moment-independent sensitivity analysis. Reliab. Eng. Syst. Saf. 2013, 110, 60–67. [Google Scholar] [CrossRef]
Sobol, I.M. Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates. Math. Comput. Simul. 2001, 55, 271–280. [Google Scholar] [CrossRef]
Blatman, G.; Sudret, B. Adaptive sparse polynomial chaos expansion based on least angle regression. J. Comput. Phys. 2011, 230, 2345–2367. [Google Scholar] [CrossRef]
Sudret, B. Global sensitivity analysis using polynomial chaos expansions. Reliab. Eng. Syst. Saf. 2008, 93, 964–979. [Google Scholar] [CrossRef]
Crestaux, T.; le Maıtre, O.; Martinez, J. Polynomial chaos expansion for sensitivity analysis. Reliab. Eng. Syst. Saf. 2009, 94, 1161–1172. [Google Scholar] [CrossRef]
Xiu, D.; Karniadakis, G.E. Modeling uncertainty in steady state diffusion problems via generalized polynomial chaos. Comput. Methods Appl. Mech. Eng. 2002, 191, 4927–4948. [Google Scholar] [CrossRef]
Chakraborty, S.; Chowdhury, R. Sequential experimental design based generalised ANOVA. J. Comput. Phys. 2016, 317, 15–32. [Google Scholar] [CrossRef]
Lasota, R.; Stocki, R.; Tauzowski, P.; Szolc, T. Polynomial chaos expansion method in estimating probability distribution of rotor-shaft dynamic responses. Bull. Pol. Acad. Sci. Tech. Sci. 2015, 63, 413–422. [Google Scholar] [CrossRef]
Chakraborty, S.; Chowdhury, R. A semi-analytical framework for structural reliability analysis. Comput. Methods Appl. Mech. Eng. 2015, 289, 475–497. [Google Scholar] [CrossRef]
Xi, Z.; Hu, C.; Youn, B.D. A comparative study of probability estimation methods for reliability analysis. Struct. Multidiscip. Optim. 2012, 45, 33–52. [Google Scholar] [CrossRef]
Wiener, N. The homogeneous chaos. Am. J. Math. 1938, 60, 897–936. [Google Scholar] [CrossRef]
Isukapalli, S.S. Uncertainty Analysis of Transport-Transformation Models. Ph.D. Thesis, Rutgers University, New Brunswick, NJ, USA, 1999. [Google Scholar]
Wei, P.; Lu, Z.; Hao, W.; Feng, J.; Wang, B. Efficient sampling methods for global reliability sensitivity analysis. Comput. Phys. Commun. 2012, 183, 1728–1743. [Google Scholar] [CrossRef]
Melchers, R.E. Structural Reliability Analysis and Prediction; John Wiley & Son Ltd.: Chichester, UK, 1999. [Google Scholar]
Savin, E.; Faverjon, B. Higher-order moments of generalized polynomial chaos expansions for intrusive and non-intrusive uncertainty quantification. In Proceedings of the 19th AIAA Non-Deterministic Approaches Conference, Grapevine, TX, USA, 9–13 January 2017. [Google Scholar]
Zellner, A.; Highfield, R.A. Calculation of maximum entropy distributions and approximation of marginalposterior distributions. J. Econom. 1988, 37, 195–209. [Google Scholar] [CrossRef]
Abramov, R.V. An improved algorithm for the multidimensional moment-constrained maximum entropy problem. J. Comput. Phys. 2007, 226, 621–644. [Google Scholar] [CrossRef]
Bandyopadhyay, K.; Bhattacharya, A.K.; Biswas, P.; Drabold, D.A. Maximum entropy and the problem of moments: A stable algorithm. Phys. Rev. E Stat. Nonlinear Soft Matter Phys. 2005, 71, 1–4. [Google Scholar] [CrossRef] [PubMed]
Blatman, G.; Sudret, B. An adaptive algorithm to build up sparse polynomial chaos expansions for stochastic finite element analysis. Probab. Eng. Mech. 2010, 25, 183–197. [Google Scholar] [CrossRef]

Figure 1. PCE based GRSA method framework.

Figure 2. The importance ranking based on global reliability sensitivity indices.

Figure 3. The beam structure of example 2.

Figure 4. The importance ranking of case 1. (a) The importance ranking based on GRSA; (b) The importance ranking based on GSA.

Figure 5. The importance ranking of case 2. (a) The importance ranking based on GRSA; (b) The importance ranking based on GSA.

Table 1. The global reliability sensitivity indices of example 1.

Method		x₁	x₂	x₃	Times of Function Evaluations
Full-scale MCS	$S_{i}^{(R)}$	0.7120	0.0005274	0.005265	1.5 × 10⁷
Full-scale MCS	$S_{T_{i}}^{(R)}$	0.9946	0.1130	0.2602	1.5 × 10⁷
IS	$S_{i}^{(R)}$	0.7206	0.0005147	0.005327	1.2 × 10⁵
IS	$S_{T_{i}}^{(R)}$	0.9945	0.1124	0.2583	1.2 × 10⁵
Proposed method	$S_{i}^{(R)}$	0.7181	0.0005426	0.005139	182
Proposed method	$S_{T_{i}}^{(R)}$	0.9932	0.1119	0.2625	182

Table 2. The distribution parameters of the input variables of example 2.

Input Variables		q/(N/M)	E/(kN/m²)	I/m⁴	L/m	F_p
Distribution type		Normal	Normal	Normal	Lognormal	Lognormal
Mean value		662.5	2 × 10¹¹	2.172 × 10⁻⁴	10	3 × 10⁴
variances	Case 1	0.1	0.08	0.1	0.06	0.1
variances	Case 2	0.05	0.06	0.05	0.03	0.04

Table 3. The global reliability sensitivity indices of case 1.

Method		q	E	I	L	F_p	Times of Function Evaluations
Full-scale MCS	$S_{i}^{(R)}$	0.001283	0.05398	0.1029	0.1446	0.05583	1.5 × 10⁷
Full-scale MCS	$S_{T_{i}}^{(R)}$	0.07086	0.4701	0.5967	0.6891	0.4895	1.5 × 10⁷
IS	$S_{i}^{(R)}$	0.001267	0.05103	0.09755	0.1397	0.05129	1.2 × 10⁵
IS	$S_{T_{i}}^{(R)}$	0.06961	0.4727	0.6026	0.6985	0.4905	1.2 × 10⁵
Proposed method	$S_{i}^{(R)}$	0.001118	0.05404	0.1031	0.1445	0.05591	3882
Proposed method	$S_{T_{i}}^{(R)}$	0.07251	0.4714	0.5879	0.6865	0.4983	3882

Table 4. The global reliability sensitivity indices of case 2.

Method		q	E	I	L	F_p	Times of Function Evaluations
Full-scale MCS	$S_{i}^{(R)}$	0.000486	0.03123	0.01150	0.02078	0.002947	1.5 × 10⁷
Full-scale MCS	$S_{T_{i}}^{(R)}$	0.09433	0.8405	0.7389	0.8240	0.5270	1.5 × 10⁷
IS	$S_{i}^{(R)}$	0.000502	0.03111	0.01194	0.02260	0.002806	1.2 × 10⁵
IS	$S_{T_{i}}^{(R)}$	0.09845	0.8509	0.7356	0.8333	0.5336	1.2 × 10⁵
Proposed method	$S_{i}^{(R)}$	0.000532	0.03064	0.01183	0.02047	0.003107	3882
Proposed method	$S_{T_{i}}^{(R)}$	0.09512	0.8497	0.7298	0.8199	0.5277	3882

Table 5. The global sensitivity indices of example 2.

		q	E	I	L	F_p
Case 1	$S_{i}$	0.003538	0.15818	0.253083	0.387629	0.183083
Case 1	$S_{T_{i}}$	0.003718	0.16375	0.260962	0.397182	0.188962
Case 2	$S_{i}$	0.003182	0.317395	0.218915	0.349925	0.106503
Case 2	$S_{T_{i}}$	0.003227	0.319898	0.220892	0.352517	0.107553

© 2018 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

Zhao, J.; Zeng, S.; Guo, J.; Du, S. Global Reliability Sensitivity Analysis Based on Maximum Entropy and 2-Layer Polynomial Chaos Expansion. Entropy 2018, 20, 202. https://doi.org/10.3390/e20030202

AMA Style

Zhao J, Zeng S, Guo J, Du S. Global Reliability Sensitivity Analysis Based on Maximum Entropy and 2-Layer Polynomial Chaos Expansion. Entropy. 2018; 20(3):202. https://doi.org/10.3390/e20030202

Chicago/Turabian Style

Zhao, Jianyu, Shengkui Zeng, Jianbin Guo, and Shaohua Du. 2018. "Global Reliability Sensitivity Analysis Based on Maximum Entropy and 2-Layer Polynomial Chaos Expansion" Entropy 20, no. 3: 202. https://doi.org/10.3390/e20030202

APA Style

Zhao, J., Zeng, S., Guo, J., & Du, S. (2018). Global Reliability Sensitivity Analysis Based on Maximum Entropy and 2-Layer Polynomial Chaos Expansion. Entropy, 20(3), 202. https://doi.org/10.3390/e20030202

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Menu

Global Reliability Sensitivity Analysis Based on Maximum Entropy and 2-Layer Polynomial Chaos Expansion

Abstract

1. Introduction

2. Related Work

2.1. Global Sensitivity Analysis and Sobol’s Indices

2.2. Global Reliability Sensitivity Analysis

3. R-PCE Framework for GRSA

3.1. PCE-Based GRSA Method

3.2. Conditional Polynomial Chaos Expansion and Moments Calculation

3.3. The Maximum Entropy Method

3.4. Global Reliability Sensitivity Indices Calculation

3.5. Implementation

4. Test Examples

4.1. A Numerical Example

4.2. An Engineering Example

5. Conclusions

Acknowledgments

Author Contributions

Conflicts of Interest

Appendix A

Appendix B

Appendix C

References

Share and Cite

Article Metrics

Article Access Statistics

Further Information

Guidelines

MDPI Initiatives

Follow MDPI