E Zio

N Pedroni

Chapter 1 Reliability estimation by advanced Monte Carlo simulation

Monte Carlo Simulation (MCS) offers a powerful means for evaluating the reliability of a system, due to the modeling flexibility that it offers indifferently of the type and dimension of the problem. The method is based on the repeated sampling of realizations of system configurations, which however are seldom of failure so that a large number of realizations must be simulated in order to achieve an acceptable accuracy in the estimated failure probability, with costly large computing times. For this reason, techniques of efficient sampling of system failure realizations are of interest, in order to reduce the computational effort.

In this paper, the recently developed Subset Simulation (SS) and Line Sampling (LS) techniques are considered for improving the MCS efficiency in the estimation of system failure probability. The SS method is founded on the idea that a small failure probability can be expressed as a product of larger conditional probabilities of some intermediate events: with a proper choice of the intermediate events, the conditional probabilities can be made sufficiently large to allow accurate estimation with a small number of samples. The LS method employs lines instead of random points in order to probe the failure domain of interest. An "important direction" is determined, which points towards the failure domain of interest; the high-dimensional reliability problem is then reduced to a number of conditional one-dimensional problems which are solved along the "important direction".

The two methods are applied on two structural reliability models of literature, i.e. the cracked plate model and the Paris-Erdogan model for thermal fatigue crack growth. The efficiency of the proposed techniques is evaluated in comparison to other stochastic simulation methods of literature, i.e., standard MCS, Importance Sampling (IS), Dimensionality Reduction (DR) and Orthogonal Axis (OA).

Introduction

In the performance-based design and operation of modern engineered systems, the accurate assessment of reliability is of paramount importance, particularly for civil, nuclear, aerospace and chemical systems and plants which are safety-critical and must be designed and operated within a risk-informed approach [START_REF] Thunnissen | Uncertainty quantification in estimating critical spacecraft component temperature[END_REF][START_REF] Patalano | Risk informed design changes in a passive decay heat removal system[END_REF].

The reliability assessment requires the realistic modelling of the structural/mechanical components of the system and the characterization of their material constitutive behaviour, loading conditions and mechanisms of deterioration and failure that are anticipated to occur during the working life of the system [START_REF] Schueller | Benchmark study on reliability estimation in higher dimension of structural systems -An overview[END_REF].

In practice, not all the characteristics of the system under analysis can be fully captured in the model. This is due to: i) the intrinsically random nature of several of the phenomena occurring during the system life; ii) the incomplete knowledge about some of these phenomena. Thus, uncertainty is always present in the hypotheses underpinning the model (model uncertainty) and in the values of its parameters (parameter uncertainty); this leads to uncertainty in the model output, which must be quantified for a realistic assessment of the system [START_REF] Nutt | Evaluations of nuclear safety from the outputs of computer codes in the presence of uncertainties[END_REF].

In mathematical terms, the probability of system failure can be expressed as a multi-dimensional integral of the form , otherwise. The failure domain F is commonly defined by a so-called Performance Function (PF) or Limit State Function (LSF) ( )

x x g
which is lower than or equal to zero if F ∈ x and greater than zero, otherwise.

In practical cases, the multi-dimensional integral (1.1) can not be easily evaluated by analytical methods nor by numerical schemes. On the other hand, Monte Carlo Simulation (MCS) offers an effective means for estimating the integral, because the method does not suffer from the complexity and dimension of the domain of integration, albeit it implies the nontrivial task of sampling from the multidimensional PDF. The MCS solution to (1.1) entails that a large number of samples of the values of the uncertain parameters x be drawn from ) ( x q and that these be used to compute an unbiased and consistent estimate of the system failure probability as the fraction of the number of samples that lead to failure. However, a large number of samples (inversely proportional to the failure probability) is necessary to achieve an acceptable estimation accuracy: in terms of the integral in (1.1) this can be seen as due to the high dimensionality n of the problem and the large dimension of the relative sample space compared to the failure region of interest [START_REF] Schueller | On the treatment of uncertainties in structural mechanics and analysis[END_REF]. This calls for new simulation techniques for performing robust estimations with a limited number of input samples (and associated low computational time).

In this respect, effective approaches are offered by Subset Simulation (SS) [START_REF] Au | Estimation of small failure probabilities in high dimensions by subset simulation[END_REF]Au and Beck 2003b) and Line Sampling (LS) (Koutsourelakis et al. 2004;[START_REF] Pradlwarter | Realistic and efficient reliability estimation for aerospace structures[END_REF].

In the SS method, the failure probability is expressed as a product of conditional failure probabilities of some chosen intermediate events, whose evaluation is obtained by simulation of more frequent events. The evaluation of small failure probabilities in the original probability space is thus tackled by a sequence of simulations of more frequent events in the conditional probability spaces. The necessary conditional samples are generated through successive Markov Chain Monte Carlo (MCMC) simulations [START_REF] Metropolis | Equations of state calculations by fast computing machines[END_REF]Hastings 1970;[START_REF] Fishman | Monte Carlo: concepts, algorithms, and applications[END_REF], gradually populating the intermediate conditional regions until the final target failure region is reached.

In the LS method, lines, instead of random points, are used to probe the failure domain of the high-dimensional problem under analysis [START_REF] Pradlwarter | Realistic and efficient reliability estimation for aerospace structures[END_REF]). An "important direction" is optimally determined to point towards the failure domain of interest and a number of conditional, one-dimensional problems are solved along such direction, in place of the high-dimensional problem [START_REF] Pradlwarter | Realistic and efficient reliability estimation for aerospace structures[END_REF]. The approach has been shown to perform always better than standard MCS; furthermore, if the boundaries of the failure domain of interest are not too rough (i.e., almost linear) and the "important direction" is almost perpendicular to them, the variance of the failure probability estimator could be ideally reduced to zero [START_REF] Koutsourelakis | Reliability of structures in high dimensions, Part I: algorithms and application[END_REF].

In this Chapter, SS and LS schemes are developed for application to two structural reliability models of literature, i.e., the cracked plate model [START_REF] Ardillon | Measures de sensibilitè dans les approaches probabilistes[END_REF] and the Paris-Erdogan thermal fatigue crack growth model [START_REF] Paris | A rational analytic theory of fatigue[END_REF]. The problem is rather challenging as it entails estimating failure probabilities of the order of 10 -7 . The effectiveness of SS and LS is compared to that of other simulation methods, e.g. the Importance Sampling (IS), Dimensionality Reduction (DR) and Orthogonal Axis (OA) methods [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]. In the IS method, the PDF ) ( x q in (1.1) is replaced with an Importance Sampling Distribution (ISD) arbitrarily chosen so as to generate samples that lead to failure more frequently (Au and Beck 2003a); in the DR method, the failure event is reexpressed in such a way as to highlight one important variable (say, xj) and the failure probability is then computed as the expected value of the Cumulative Distribution Function (CDF) of xj conditional on the remaining (n -1) variables; fi-nally, in the OA method, a sort of importance sampling is performed around the most likely point in the failure domain [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF].

The remainder of the Chapter is organized as follows. In Section 1.2, a general presentation of the SS and LS schemes implemented for this study is given. In Section 1.3, the IS, DR and OA methods taken as terms of comparison are briefly summarized. The results of the application of SS and LS to the cracked plate and thermal fatigue crack growth models are reported in Sections 1.4 and 1.5, respectively. Based on the results obtained, a critical discussion of the simulation techniques adopted and compared in this work is offered in the last Section. For completeness of the contents of the paper, detailed descriptions of the Markov Chain Monte Carlo (MCMC) simulation method used for the development of the SS and LS algorithms are provided in Appendices 1 and 2, respectively.

Simulation methods implemented in this study

The Subset Simulation method

Subset Simulation (SS) is an adaptive stochastic simulation method originally developed for efficiently computing small failure probabilities in structural reliability analysis [START_REF] Au | Estimation of small failure probabilities in high dimensions by subset simulation[END_REF]. The underlying idea is to express the (small) failure probability as a product of (larger) probabilities conditional on some intermediate events. This allows converting a rare event simulation into a sequence of simulations of more frequent events. During simulation, the conditional samples are generated by means of a Markov chain designed so that the limiting stationary distribution is the target conditional distribution of some adaptively chosen event; by so doing, the conditional samples gradually populate the successive intermediate regions up to the final target (rare) failure region (Au and Beck 2003b).

The basic principles

For a given target failure event F of interest, let

F F F F m = ⊃ ⊃ ⊃ ... 2 1 be a se- quence of intermediate events, so that i k i k F F 1 = ∩ = , k = 1, 2, …, m.
By sequentially conditioning on the event Fi, the failure probability P(F) can be written as

∏ - = + = = 1 1 1 1 ) | ( ) ( ) ( ) ( m i i i m F F P F P F P F P (1.2)
Notice that even if P(F) is small, the conditional probabilities involved in (1.2) can be made sufficiently large by appropriately choosing m and the intermediate events {Fi, i = 1, 2, …, m -1}.

The original idea of SS is to estimate the failure probability P(F) by estimating P(F1) and

}. 1 ..., , 2 , 1 : ) | ( { 1 - = + m i F F P i i
Considering for example P(F) ≈ 10 -5 and choosing m = 4 intermediate events such that P(F1) and

1 . 0 } 4 , 3 , 2 , 1 : ) | ( { 1 ≈ = + i F F P i i
, the conditional probabilities can be evaluated efficiently by simulation of the relatively frequent intermediate events [START_REF] Au | Estimation of small failure probabilities in high dimensions by subset simulation[END_REF].

Standard MCS can be used to estimate P(F1). On the contrary, computing the conditional probabilities in (1.2) by MCS entails the non-trivial task of sampling from the conditional distributions of x given that it lies in Fi, i = 1, 2, ..., m -1, i.e. from

) ( / ) ( ) ( ) | ( F P I q F q i F i x x x = .
In this regard, Markov Chain Monte Carlo (MCMC) simulation provides a powerful method for generating samples condi-tional on the intermediate regions Fi, i = 1, 2, ..., m -1 [START_REF] Au | Estimation of small failure probabilities in high dimensions by subset simulation[END_REF]Au and Beck 2003b). For completeness of the paper, the related algorithm is presented in Appendix 1.

The algorithm

In the actual SS implementation, with no loss of generality it is assumed that the failure event of interest can be defined in terms of the value of a critical system response variable Y being lower than a specified threshold level y, i.e., F = {Y < y}.

The sequence of intermediate events

} ..., , 2 , 1 : { m i F i =
can then be correspondingly defined as

{ } i i y Y F < = , m i ..., , 2 , 1 =
where y1 > y2 > ... > yi > ... > ym = y > 0 is a decreasing sequence of intermediate threshold values [START_REF] Au | Estimation of small failure probabilities in high dimensions by subset simulation[END_REF]Au and Beck 2003b).

The choice of the sequence { }

m i y i ..., , 2 , 1 : = affects the values of the condi- tional probabilities } 1 ..., , 2 , 1 : ) | ( { 1 - = + m i F F P i i in (1.
2) and hence the efficiency of the SS procedure. In particular, choosing the sequence { }

m i y i ..., , 2 , 1 : = a pri- ori makes it difficult to control the values of the conditional probabilities } 1 ..., , 2 , 1 : ) | ( { 1 - = + m i F F P i i
. For this reason, in this work, the intermediate threshold values are chosen adaptively in such a way that the estimated conditional probabilities are equal to a fixed value p0 [START_REF] Au | Estimation of small failure probabilities in high dimensions by subset simulation[END_REF]Au and Beck 2003b).

Schematically, the SS algorithm proceeds as follows (Figure 1.1):

1. Sample N vectors { } N k k ..., , 2 , 1 : 0 = x
by standard MCS, i.e., from the original probability density function q(•). The subscript '0' denotes the fact that these samples correspond to 'Conditional Level 0'; 2. Set i = 0;

Compute the values of the response variable

{ } N k Y k i ..., , 2 , 1 : ) ( = x
; 4. Choose the intermediate threshold value yi+1 as the (1 -p0)N th value in the de-

creasing list of values { } N k Y k i ..., , 2 , 1 : ) ( = x
(computed at step 3. above) to define Fi+1 = {Y < yi+1}. By so doing, the sample estimate of P(Fi+1|Fi) = P(Y < yi+1|Y < yi) is equal to p0 (note that it has been implicitly assumed that p0N is an integer value); 5. If yi+1 ym, proceed to 10. below; 6. Viceversa, i.e. if yi+1 > ym, with the choice of yi+1 performed at step 4. above,

identify the p0N samples { } N p u u i 0 ..., , 2 , 1 : = x among { } N k k i ..., ,
+ + ∈ = i k i F N k x
, at 'Conditional level i + 1'; 8. Set i ← i + 1; 9. Return to step 3. above; 10.Stop the algorithm.

For clarity sake, a step-by-step illustration of the procedure for Conditional levels 0 and 1 is provided in Figure 1.2 by way of example. Notice that the procedure is such that the response values { } m i y i ..., , 2 , 1 : = at the specified probability levels

0 1 ) ( p F P = , 2 0 1 1 2 2 ) ( ) | ( ) ( p F P F F p F P = = , …, m m p F P 0 ) ( = are
estimated, rather than the event probabilities

) ( 1 F P , ) | ( 1 2 F F P , …, ) | ( 1 - m m F F P
, which are a priori fixed at p0. In this view, SS is a method for generating samples whose response values correspond to specified probability levels, rather than for estimating probabilities of specified failure events. As a result, it produces information about ) ( y Y P < versus y at all the simulated values of Y rather than at a single value of y. This feature is important because the whole trend of

) ( y Y P <
versus y provides much more information than a point estimate [START_REF] Au | Reliability-based design sensitivity by efficient simulation[END_REF].

The Line Sampling method

Line Sampling (LS) was also originally developed for the reliability analysis of complex structural systems with small failure probabilities [START_REF] Koutsourelakis | Reliability of structures in high dimensions, Part I: algorithms and application[END_REF]. The underlying idea is to employ lines instead of random points in order to probe the failure domain of the high-dimensional system under analysis [START_REF] Pradlwarter | Realistic and efficient reliability estimation for aerospace structures[END_REF].

In extreme synthesis, the problem of computing the multidimensional failure probability integral (1.1) in the original "physical" space is transformed into the so-called "standard normal space", where each random variable is represented by an independent central unit Gaussian distribution. In this space, a unit vector α (hereafter also called "important unit vector" or "important direction") is determined, pointing towards the failure domain F of interest (for illustration purposes, two plausible important unit vectors, α 1 and α 2 , pointing towards two different failure domains, F 1 and F 2 , are visually represented in Figure 1.3, left and right, respectively, in a two-dimensional uncertain parameter space). The problem of computing the high-dimensional failure probability integral (1.1) is then reduced to a number of conditional one-dimensional problems, which are solved along the "important direction" α in the standard normal space. The conditional onedimensional failure probabilities (associated to the conditional one-dimensional problems) are readily computed by using the standard normal cumulative distribution function [START_REF] Pradlwarter | Realistic and efficient reliability estimation for aerospace structures[END_REF]. . For problems where the dimension n is not so small, the parameter vector x can be transformed into the vector

n ℜ ∈ θ
, where each element of the vector θ j , j = 1, 2, …, n, is associated with a central unit Gaussian standard distribution [START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF] ( )

x θ θ x T = (1.4) ( ) θ x x T θ = (1.5)
Transformations (1.4) and (1.5) are in general nonlinear and are obtained by applying Rosenblatt's or Nataf's transformations, respectively [START_REF] Rosenblatt | Remarks on multivariate transformations[END_REF][START_REF] Nataf | Determination des distribution dont les marges sont donnees[END_REF][START_REF] Huang | A robust design method using variable transformation and Gauss-Hermite integration[END_REF]. They are linear only if the random vector x is jointly Gaussian distributed. By transformation (1.4), also the Performance Function (PF) or Limit State Function (LSF)

( ) ⋅ x g
defined in the physical space (Section 1) can be transformed into ( ) ⋅ θ g in the standard normal space:

However, the design points, and their neighborhood, do not always represent the most important regions of the failure domain, especially in high-dimensional spaces [START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF]). Moreover, the computational cost associated with the calculation of the design point can be quite high, in particular if long-running numerical codes are required to simulate the response of the system to its uncertain input parameters [START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF]), as it is frequently the case in structural reliability.

In [START_REF] Pradlwarter | Realistic and efficient reliability estimation for aerospace structures[END_REF], the direction of α is taken as the normalized gradient of the performance function in the standard normal space. Since the unit vector

{ } n j α α α α ..., , ..., , , 2 1 = α
points towards the failure domain F, it can be used to draw information about the relative importance of the random parameters

{ } n j j ..., , 2 , 1 : = θ
with respect to the failure probability P(F): the more relevant a random variable in determining the failure of the system, the larger the corre-sponding component of the unit vector α will be [START_REF] Pradlwarter | Realistic and efficient reliability estimation for aerospace structures[END_REF]. Such quantitative information is obtained from the gradient of the performance function

() θ θ g
in the standard normal space, ( ) [START_REF] Pradlwarter | Realistic and efficient reliability estimation for aerospace structures[END_REF].

θ θ g ∇ : () () () () () T n j g g g g g ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ ∂ ∂ ∂ ∂ ∂ ∂ ∂ θ θ g g ∇ ∇ = (
On the other hand, when the performance function is defined on a highdimensional space, i.e. when many parameters of the system under analysis are random, the computation of the gradient ( ) Finally, the important unit vector α can also be computed as the normalized "center of mass" of the failure domain F of interest [START_REF] Koutsourelakis | Reliability of structures in high dimensions, Part I: algorithms and application[END_REF]

θ θ g ∇ in (1.
∑ = ⋅ = s N u u u s N 1 2 1 θ θ α (Fig- ure 1.4
). This direction is not optimal, but it is clear that it provides a good approximation of the important regions of the failure domain (at least as the sample size Ns is large). On the other hand, it should be noticed that the procedure implies Ns additional system analyses by the deterministic model simulating the system, which substantially increase the computational cost associated to the simulation method. In the implementation of LS for this work, the method based on the normalized "center of mass" of the failure domain F has been employed, because it relies on a "map" approximating the failure domain F under analysis (given by the failure samples generated through a Markov chain) and thus it provides in principle the most realistic and reliable estimate for the LS important direction α.

For completeness, a thorough description of the Line Sampling algorithm and its practical implementation issues is given in Appendix 2 at the end of the paper.

Simulation methods considered for comparison

The performances of Subset Simulation (Section 1.2.1) and Line Sampling (Section 1.2.2) will be compared to those of the Importance Sampling (IS) (Section 1.3.1), Dimensionality Reduction (DR) (Section 1.3.2) and Orthogonal Axis (OA) (Section 1.3.3) methods; the comparison will be made with respect to the results reported in [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF] for the two literature case studies considered, of the cracked plate and thermal fatigue crack growth models.

The Importance Sampling method

The concept underlying the Importance Sampling (IS) method is to replace the original PDF ) ( x q with an Importance Sampling Distribution (ISD)

) ( ~x q
arbitrarily chosen by the analyst so as to generate a large number of samples in the "important region" of the sample space, i.e. the failure region F (Au and Beck 2003a;[START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF].

The IS algorithm proceeds as follows [START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF]):

1. Identify a proper Importance Sampling Density (ISD), ( ) ⋅ q ~, in order to increase the probability of occurrence of the failure samples.

2. Express the failure probability P(F) in (1.1) as a function of the ISD ( )

⋅ q ~: ( ) ( ) ( ) ()() () () ()() () ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ = ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ = = ∫ ∫ x x x x x x x x x x x q q I E d q q q I d q I F P F q F F ~~ (1.8) 3. Draw NT independent and identically distributed (i.i.d.) samples { } T k N k ..., , 2 , 1 : = x from the ISD ( ) ⋅ q ~; if a good choice for the ISD () ⋅ q ~ has been made, the samples { } T k N k ..., , 2 , 1 : = x
should be concentrated in the failure region F of interest.

Compute an estimate ( )

F P ˆ
for the failure probability P(F) in (1.1) by resorting to the last expression in (1.8):

() ( ) ( ) ( ) ∑ = = T N k k k k F T q q I N F P 1 1 ˆx x x (1.9) 5. The variance [ ] ) ( ˆF P V
of the estimator ( )

F P ˆ in (1.9) is given by () [ ] ( ) ( ) () ()() () () () ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - = ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ = ∫ 2 2 2 ~1 1 ˆF P d q q q I N q q I V N F P V F T F q T x x x x x x x x (1.10)
It is straightforward to verify that the quantity (1.10) becomes zero when

() () ( ) ( ) () F P q I q q F opt x x x x = = ~ (1.11)
This represents the optimal choice for the importance sampling density which is practically unfeasible since it requires the a priori knowledge of P(F). Several techniques have been developed in order to approximate the optimal sampling density (1.11) or to at least find one giving small variance of the estimator (1.9). Recent examples include the use of engineering judgment [START_REF] Pagani | The impact of uncertainties on the performance of passive systems[END_REF]), design points [START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF]) and kernel density estimators (Au and Beck 2003a).

The Dimensionality Reduction method

Objective of the Dimensionality Reduction (DR) method is to reduce the variance associated to the failure probability estimates by exploiting the property of conditional expectation [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]. In extreme synthesis, the failure event

( ) 0 ≤ x x g
is re-expressed in such a way as to highlight one of the n uncertain input variables of x (say, xj); then, the failure probability estimate is computed as the expected value of the CDF of xj conditional on the remaining (n -1) input variables. By so doing, the zero values contained in the standard MCS estimator (i.e., IF(x) = 0, if x ∈ F) are removed: this allows to i) reach any level of probability (even very small) and ii) reduce the variance of the failure probability estimator [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF].

The DR algorithm proceeds as follows [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]:

1. Write the failure event ( ) ( ) 0 ..., , ..., , , 2 1 ≤ = n j x x x x x x g g x
in such a way as to highlight one of the n uncertain input variables (e.g., xj):

( )

j x j h x - ≤ x , j = 1, 2, …, n (1.12)
where hx(·) is a function defined on 1 -ℜ n which takes values on the set of all (measurable) subsets of ℜ and j -

x is a vector containing all the uncertain input variables except xj, i.e., ( )

n j j j x x x x x ..., , , ..., , , 1 1 2 1 + - -= x ;
2. Write the failure probability P(F) as follows:

( ) ( ) [ ] ( ) [ ] ( ) [ ] { } j x x j x j x h F E h x P g P F P j j j - - - - = ≤ = ≤ = x x x x x | 0 (1.13)
where

() ⋅ -j j x F x | is the Cumulative Distribution Function of xj conditional on j - x , i.e., (
)

n j j j x x x x x ..., , , ..., , , 1 1 2 1 + - -= x ; 3. Draw NT samples { } T k j N k ..., , 2 , 1 : = - x
, where ( )

k n k j k j k k k j x x x x x ..., , , ..., , , 1 1 2 1 + - -= x ,
from the (n -1)-dimensional marginal probability density function

( ) j m q - x , i.e., ( ) ( )= = + - - n j j m j m x x x x x q q ..., , , ..., , , 1 1 2 1 x ( ) ∫ j x j n j dx x x x x q
..., , ..., , , 2 1 ; 4. Using the last expression in (1.13), compute an unbiased and consistent estimate () F P ˆ

for the failure probability P(F) as follows:

() ( ) [ ] ∑ = - - = T j j N k k j x x T h F N F P 1 | 1 ˆx x (1.14)
It is worth noting that in (1.14) the failure probability estimate is computed as the expected value of the cumulative distribution function

( ) ⋅ -j j x F x | of xj condi-
tional on the remaining (n -1) input variables. Since this quantity takes values between 0 and 1, the zero values contained in the standard MCS estimator (i.e., IF(x) = 0, if x ∈ F) are removed: this allows to i) reach any level of failure probability (even very small) and ii) reduce the variance of the failure probability estimator. However, such method can not always be applied: first, the performance function gx(•) must be known analytically; second, it must have the property that one of the uncertain input variables can be separated from the others to allow re-writing the failure condition ( )

0 ≤ x x g
in the form of (1.12) [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF].

Finally, notice that DR can be considered a very special case of LS (Section 1.2.2) where the performance function gx(•) is analytically known and the important direction α coincides with the "direction" of the variable xj, i.e., ( )

0 , 0 ..., , ..., , 0 , 0 j x = α .

The Orthogonal Axis method

The Orthogonal Axis (OA) method combines the Fist Order Reliability Method (FORM) approximation (Der [START_REF] Kiureghian | The geometry of random vibrations and solutions by FORM and SORM[END_REF] and Monte Carlo Simulation (MCS) in a sort of importance sampling around the "design point" of the problem (see Section 1.2.2.2).

The OA algorithm proceeds as follows [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]:

1. Transform { } n n j x x x x ℜ ∈ = ..., , ..., , , 2 1 x
, i.e., the vector of uncertain parameters defined in the original physical space

n ℜ ∈ x , into the vector n ℜ ∈ θ
, where each element of the vector θj, j = 1, 2, …, n, is associated with a central unit Gaussian standard distribution [START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF]) (see Section 1.2.2.1). Thus, the joint probability density function of θ can simply be written as

( ) ( ) ∏ = = n j j n 1 θ φ ϕ θ (1.15)
where ( ) ( ) ( ) [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF].

It is worth noting that the idea underlying the OA method is essentially the same as that of LS (Section 1.2.2). However, in OA the "important direction" is forced to coincide with that of the design point of the problem; moreover, OA employs a rotation of the coordinate system which can be difficult to define in very high-dimensional problems.

Application 1: the cracked plate model

The cracked plate model is a classical example in Fracture Mechanics and its relative simplicity allows a detailed and complete study of different simulation techniques. A thorough description of this model can be found in [START_REF] Ardillon | Measures de sensibilitè dans les approaches probabilistes[END_REF].

The mechanical model

A metal plate of infinite length with a defect of initial length equal to a [m] is considered. The plate is supposed to be subject to a uniform normal loading (i.e., stress) s [MPa]. The intensity factor K [MPa m ], determined by the uniform loading in the neighborhood of the defect is defined as follows:

a Fs K π ∞ = (1.19)
where F is the shape factor of the defect. The plate is supposed to break (i.e., fail) when the intensity factor K in (1.19) becomes greater than or equal to a critical value Kc, i.e.:

c K a Fs K ≥ = ∞ π (1.20)
The variables of the mechanical model are summarized in Table 1.1.

The structural reliability model

From the point of view of a structural reliability analysis, the cracked plate mechanical model of Section 1.4.1 is analyzed within a probabilistic framework in which the variables Kc, a, F and s are uncertain (for simplicity of illustration with respect to the notation of the previous Sections, the four variables are hereafter named x1, x2, x3 and x4, respectively).

Referring to (1.20), the performance function gx( x ) of the system is ( ) ( )

2 4 3 1 4 3 2 1 , , , x x x x x x x x g g x x π - = = x (1.21)
The failure region F is then expressed as

( ) { } ( ) { } 2 4 3 1 4 3 2 1 : , , , 0 : x x x x x x x x g F x π ≤ = ≤ = x x (1.22)
Finally, the probability of system failure P(F) is written as follows:

( ) ( ) ( ) [ ] (
)

2 4 3 1 0 x x x x P g P F P F P x π ≤ = ≤ = ∈ = x x
.

(1.23)

Case studies

Four case studies, namely Case 0 (Reference case), 1, 2 and 3, are considered with respect to the structural reliability model of the previous Section 1.4.2. Each case study is characterized by different PDFs for the uncertain variables x1, x2, x3 and x4 and by different failure probabilities P(F): these features are summarized in Table 1.2. Notice that in Cases 0, 1 and 2 the random variables are independent and normally distributed, whereas in Case 3 they are independent and lognormally distributed. Moreover, it is worth noting that the exact (i.e., analytically computed) failure probabilities P(F) approximately range from 10 -3 to 10 -7 , allowing a deep exploration of the capabilities of the simulation algorithms considered and a meaningful comparison between them [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF].

x2 (a)

N(5·10 -3 , 10 -3 ) N(5·10 -3 , 10 -3 ) N(5·10 -3 , 10 -3 )

LG(5·10 -3 , 10 -3 ) x3 (F) LG(600, 60)

P(F)

1.165·10 -3 4.500·10 -7 4.400·10 -7 3.067·10 -4

Results

In this Section, the results of the application of SS and LS for the reliability analysis of the cracked plate model of Section 1.4.1 are illustrated with reference to Case studies 0, 1, 2 and 3 described in the previous Section 1.4.3.

For fair comparison, all methods have been run with a total of NT = 50000 samples in all four cases. The efficiency of the simulation methods under analysis is evaluated in terms of four quantities: the failure probability estimate () F P ˆ, the sample standard deviation σˆ of the failure probability estimate ( ) F P ˆ, the coefficient of variation (c.o.v.) δ of ( )

F P ˆ
(defined as the ratio of the sample standard deviation σˆ to the estimate ( ) , where tcomp is the computational time required by the simulation method). The closer is the estimate ( )

F P ˆ
to the exact (i.e., analytically computed) failure probability P(F), the more accurate is the simulation method.

The sample standard deviation σˆ and the c.o.v. δ of ( )

F P ˆ
are used to quantify the variability of the failure probability estimator; in particular, the lower are the values of σˆ and δ , the lower is the variability of the corresponding failure probability estimator and thus the higher is the efficiency of the simulation method adopted. Finally, the FOM is introduced to take into account the computational time required by the method. The value of the FOM increases as the sample variance 2 σ of the failure probability estimate ( )

F P ˆ
and the computational time tcomp required by the method decrease; thus, in this case the higher is the value of the index, the higher is the efficiency of the method [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF].

The different simulation methods are also compared with respect to two direct performance indicators relative to standard MCS. First, the ratio of the sample standard deviation MC σˆ obtained by Standard MCS to that obtained by the simula- tion method under analysis meth σˆ is computed. This ratio only quantifies the improvement in the precision of the estimate achieved by using a given simulation method instead of standard MCS. Then, the ratio of the FOM of the simulation method in object, namely FOMmeth, to that of standard MCS, namely FOMMC, is considered to quantify the overall improvement in efficiency achieved by a given simulation method with respect to standard MCS, since it takes into account also the computational time required. Obviously, the higher are the values of these two indices for a given method, the higher is the efficiency of that method [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]. It can be seen that SS performs consistently better than standard MCS and its performance significantly grows as the failure probability to be estimated decreases: for instance, in Case 0 (Reference), where P(F) ~ 10 -3 , the FOM of SS, namely FOMSS, is only four times larger than that of Standard MCS, namely FOMMC; whereas in Case 1, where P(F) ~ 10 -7 , the ratio FOMSS/FOMMC is about 557. On the other hand, LS outperforms SS with respect to both (equation (1.5') in Appendix 2) require much more than NT system analyses by the model, this is significantly overweighed by the accelerated convergence rate that can be attained by the LS method with respect to SS.

Comparison with other stochastic simulation methods

The results obtained by SS and LS are compared to those obtained by the Importance Sampling (IS), Dimensionality Reduction (DR), Orthogonal Axis (OA) methods and by a combination of IS and DR (Section 1.3) [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]. For DR, the variable x1 is explicited.

The values of the performance indicators

meth MC σ σ / ˆ
and FOMmeth/FOMMC obtained by the four methods in Cases 0, 1, 2 and 3 are summarized in Table 1.4. ) are more than one order of magnitude larger than those reported in [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF] for IS (630, 376 and 643, 1.5•10 5 for Cases 1 and 2, respectively) and DR (856, 7.3•10 5 and 242, 242 for Cases 1 and 2, respectively). Moreover, it is worth noting that in the reference studies by [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF] a significant number of simulations has been run to properly tune the parameters of the ISDs for the IS method (in particular, 8, 6, 6 and 8 simulations have been performed for Cases 0, 1, 2 and 3, respectively), with a significant increase in the associated computational effort.

LS is found to perform slightly worse than OA in all the case studies considered: actually, the values of both LS MC σ σ ˆ and FOMLS/FOMMC are slightly lower than those reported in [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF] for OA. However, it should be considered that in these studies the OA method has been applied to a simplified version of the problem described in Sections 1.4.1 and 1.4.2; actually, only three uncertain variables (i.e., x1, x2 and x4) have been considered by keeping variable x3 (i.e., F) fixed to its mean value (i.e., 0.99): this certainly reduces the variability of the model output and contributes to the reduction of the variability of the associated failure probability estimator.

Further, LS performs consistently better than the combination of IS and DR in the task of estimating failure probabilities around 10 -3 ÷ 10 -4 (for instance, in Case 0

ˆ

= 8077 and LS MC σ σ ˆ = 6746. However, it has to be noticed again that in the reference studies by [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]) a significant number of simulations has been run to properly tune the parameters of the ISDs for the IS method (in particular, 4, 8, 8 and 10 simulations have been performed in Cases 0, 1, 2 and 3, respectively).

Finally, it is worth noting that in these cases SS performs worse than the other methods proposed.

Application 2: thermal fatigue crack growth model

The thermal fatigue crack growth model considered in this study is based on the deterministic Paris-Erdogan model which describes the propagation of a manufacturing defect due to thermal fatigue [START_REF] Paris | A rational analytic theory of fatigue[END_REF].

The mechanical model

The evolution of the size a of a defect satisfies the following equation:

() ( ) m c K R f C dN da Δ ⋅ ⋅ = (1.24)
where Nc is the number of fatigue cycles, C and m are parameters depending on the properties of the material, f(R) is a correction factor which is a function of the material resistance R and K is the variation of the intensity factor, defined as

( ) a a Y s K π ⋅ ⋅ Δ = Δ (1.25)
In (1.25), s is the variation of the uniform loading (stress) applied to the system and Y(a) is the shape factor of the defect. Let Si = si be the variation of the uniform normal stress at cycle i = 1, 2, …, Nc. The integration of equation (1.24) gives

() ( ) () ( ) ∫ ∑ = ⋅ ⋅ = c N c a a N i m i m S R f C a a Y da 0 1 π (1.26)
where a0 and c N a are the initial and final size of the defect, respectively. In (1.26) the following approximation can be adopted

() ( ) ( ) ( ) ( ) ∑ = ⋅ ⋅ ⋅ - ≈ ⋅ c N i m c m i S R f N T T S R f 1 0 (1.27)
where T and T0 are the initial and final times of the thermal fatigue treatment (of Nc cycles).

The system is considered failed when the size aNc of the defect at the end of the Nc cycles exceeds a critical dimension i.e.: Using (1.27), a safety margin ( )

T M
can then be defined as follows:

() () ( ) ( )( ) ( ) m c a a m S R f N T T C a a Y da T M c ⋅ ⋅ ⋅ - ⋅ - ⋅ = ∫ 0 0 π (1.31)
The failure criterion can then be expressed in terms of the safety margin (1.31):

( ) 0 ≤ T M (1.32)
The variables of the model are summarized in Table 1.5.

The structural reliability model

For the purpose of a structural reliability analysis, the thermal fatigue crack growth model is framed within a probabilistic representation of the uncertainties affecting the nine variables a0, ac, T0, T, C, m, f(R), Nc and S (hereafter named x1,x2,x3,x4,x5,x6,x7,x8 and x9,respectively).

From (1.32), the probability of system failure P(F) is written as

() ( ) [ ] () ( ) ( )( ) ( ) ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ ≤ ⋅ ⋅ ⋅ - ⋅ - ⋅ = ≤ = ∫ 0 0 0 0 m c a a m S R f N T T C a a Y da P T M P F P c π (1.33) or () ( ) [ ] () ( ) ( ) ( ) ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ ≤ ⋅ ⋅ ⋅ - ⋅ - ⋅ = ≤ = ∫ 0 0 6 2 1 6 9 7 8 3 4 5 x x x x x x x x x x a a Y da P T M P F P π (1.34)
It is worth noting the highly nonlinear nature of expressions (1.33) and (1.34) which increases the complexity of the problem.

Case studies

Two different case studies, namely Case 1 and Case 2, are built with reference to the structural reliability model of the previous Section 1.5.2. The characteristics of the PDFs of the uncertain variables of Table 1.5 are summarized in Table 1.6; the values of the exact (i.e., analytically computed) failure probabilities, P(F), for both Cases 1 and 2 are also reported in the last row of Table 1.6.

Table 1.6. Probability distributions and parameters (i.e., means and standard deviations) of the uncertain variables x1, x2, …, x9 of the thermal fatigue crack growth model of Section 1.5.2 for the two case studies considered (i.e., Cases 1 and 2); the last row reports the values of the corresponding exact (i.e., analytically computed) failure probabilities, P(F) [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]. Exp=exponential distribution; LG=Lognormal distribution; N=Normal distribution

Case 1

Case 2

x1 (a0) Exp(0.61·10 -3 ) Exp(0.81·10 -3 ) x2 (ac) N(21.4·10 -3 , 0.214·10 -3 ) N(21.4·10 -3 , 0.214·10 -3 ) x3 (T0) 0 0 x4 (T) 40 40 x5 (C)
LG(6.5·10 -13 , 5.75·10 -13 )

LG(1.00·10 -12 , 5.75·10 -13 )

x6 (m) 3.4 3.4 x7 (f(R)) 2 2 x8 (Nc) N(20, 2) N(20, 2) x9 (S)
LG(300, 30)

LG(200, 20)

P(F)

3.3380•10 -4 1.780·10 -5

Results

In this Section, the results of the application of SS and LS for the reliability analysis of the thermal fatigue crack growth model of Sections 1.5.1 and 1.5.2 are illustrated with reference to Cases 1 and 2 (Table 1.6 of Section 1.5.3).

Again for fair comparison all simulation methods have been run with the same total number of samples (NT = 40000) in both Cases 1 and 2. The efficiency of the methods has been evaluated in terms of the same indices and performance indicators defined in Section 1.4.4. Also in this application, the LS methodology is found to outperform SS in both Cases 1 and 2: for example, in Case 2, where the failure probability P(F) to be estimated is around 10 -5 , the ratio FOMLS/FOMMC is about 1500 times larger than the ratio FOMSS/FOMMC.

Comparison with other stochastic simulation methods

As done for the previous application of Section 1.4, the results obtained by SS and LS have been compared to those obtained by other literature methods, in particular the Importance Sampling (IS) and a combination of Importance Sampling and Dimensionality Reduction (Section 1.3) which have turned out to give the best results in the case studies considered [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]. Notice that the Orthogonal Axis (OA) method has not been implemented for this application in the reference study [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]: this is due to the high dimensionality of the problem which makes the definition of a proper rotation matrix very difficult (step 3. in In this application, LS is found to outperform both IS and the combination of IS and DR: for example, in Case 2, the ratio FOMLS/FOMMC is 65 and 35 times larger than FOMIS/FOMMC and FOMIS+DR/FOMMC, respectively. This confirms the capability of the LS method to efficiently probe complex high-dimensional domains of integration.

Summary and critical discussion of the techniques

One of the major obstacles in applying simulation methods for the reliability analysis of engineered systems and structures is the challenge posed by the estimation of small failure probabilities: the simulation of the rare events of failure occurrence implies a significant computational burden [START_REF] Schueller | On the treatment of uncertainties in structural mechanics and analysis[END_REF].

In order to overcome the rare-event problem, the Importance Sampling (IS) method has been introduced (Au and Beck 2003a;[START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF]). This technique amounts to replacing the original PDF of the uncertain random variables with an Importance Sampling Distribution (ISD) chosen so as to generate samples that lead to failure more frequently (Au and Beck 2003). IS has the capability to considerably reduce the variance compared with standard MCS, provided that the ISD is chosen similar to the theoretical optimal one (equation (1.11) of Section 1.3.1). However, generally substantial insights on the system stochastic behaviour and extensive modelling work is needed to identify a "good" ISD, e.g. by identifying "design points" [START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF], setting up complex kernel density estimators (Au and Beck 2003a) or simply by tuning the parameters of the ISD based on expert judgment and trial-and-error [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF][START_REF] Pagani | The impact of uncertainties on the performance of passive systems[END_REF]. Overall, this greatly increases the effort associated to the simulation for accurate failure probability estimation. Furthermore, there is always the risk that an inappropriate choice of the ISD may lead to worse estimates compared to Standard MCS [START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF].

Subset Simulation (SS) offers a clever way out of this problem by breaking the small failure probability evaluation task into a sequence of estimations of larger conditional probabilities. During the simulation, more frequent samples conditional to intermediate regions are generated from properly designed Markov chains. The method has been proven much more effective than standard MCS in the very high-dimensional spaces characteristic of structural reliability problems in which the failure regions are just tiny bits [START_REF] Au | Estimation of small failure probabilities in high dimensions by subset simulation[END_REF].

The strength of Subset Simulation lies in the generality of its formulation and the straightforward algorithmic scheme. In contrast to some of the alternative methods (e.g., Line Sampling and Orthogonal Axis), it is not restricted to standard normal spaces and can provide equally good results irrespectively of the joint distribution of the uncertain variables as long as one can draw samples from it. Furthermore, a single run of the SS algorithm leads to the calculation of the probabilities associated with all the conditional events considered: if for example, the probability of exceeding a critical level by a system response statistic of a stochastic system (e.g., the mean or a percentile of the displacement, stress, temperature, etc) is sought, then by appropriate parametrization of the intermediate conditional events, a single run can provide the probabilities of exceedance associated with a wide range of values of the response statistic of interest irrespectively of their magnitude [START_REF] Au | Reliability-based design sensitivity by efficient simulation[END_REF].

On the other hand, a word of caution is in order with respect to the fact that the conditional samples generated during the Markov Chain Monte Carlo (MCMC) simulation are correlated by construction. Since it is demonstrated that a high correlation among conditional samples increases the variance of the SS estimates, a good choice/tuning of the SS parameters (i.e., the conditional probability p0 and the proposal PDFs for MCMC simulation) is required to avoid it (Au and Beck 2003b). Finally, another drawback of the SS method is the need to express the failure event F in terms of a real valued parameter crossing a given threshold (i.e., F = {Y < y}). This parametrization is natural for the cases of practical interest in structural reliability and otherwise specific for other system reliability problems [START_REF] Zio | Reliability analysis of discrete multi-state systems by means of subset simulation[END_REF].

An alternative way to perform robust estimations of small failure probabilities without the extensive modelling effort required by IS is offered by Line Sampling (LS). The LS method employs lines instead of random points in order to probe the high-dimensional failure domain of interest. An "important direction" is optimally determined to point towards the failure domain of interest and a number of conditional, one-dimensional problems are solved along such direction, in place of the original high-dimensional problem [START_REF] Pradlwarter | Realistic and efficient reliability estimation for aerospace structures[END_REF]. In case the boundaries of the failure domain of interest are not too rough (i.e., approximately linear) and the "important direction" is almost perpendicular to them, only few simulations suffice to arrive at a failure probability with acceptable confidence. The determination of the important direction requires additional evaluations of the system performance which increases the computational cost (Section 1.2.2.2). Further, for each random sample (i.e., system configuration) drawn, two or three evaluations of the system performance are necessary to estimate the conditional one-dimensional failure probability estimates by linear or quadratic interpolation (equation (1.5') in Appendix 2). In case the "important direction" is not the optimal one, the variance of the estimator will increase. Of particular advantage of Line Sampling is its robustness: in the worst possible case where the "important direction" is selected orthogonal to the (ideal) optimal direction, line sampling performs at least as well as standard Monte Carlo simulation [START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF].

Finally, the Dimensionality Reduction (DR) method and the Orthogonal Axis (OA) method employ simulation concepts similar to those of LS, but with important limitations [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]. In the DR method, the failure event of interest is re-expressed in such a way as to highlight one (say, xj) of the input random variables, recognized as more important; then, the failure probability estimate is computed as the expected value of the CDF of xj conditional on the remaining (n -1) input variables. By so doing, the zero values contained in the standard MCS estimator (i.e., IF(x) = 0, if x ∈ F) are removed: this allows to i) reach any level of probability (even very small) and ii) reduce the variance of the failure probability estimator [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]. Notice that DR can be considered a very special case of LS where the important direction α coincides with the "direction" of the variable xj, i.e., ( )

0 , 0 ..., , ..., , 0 , 0 j x = α
. However, such method can not always be applied: first, the performance function of the system must be analytically known (which is never the case for realistic systems simulated by detailed computer codes); second, the performance function must have the characteristic that one of the variables can be separated from the others [START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF].

Finally, the Orthogonal Axis (OA) method performs a sort of importance sampling around the design point of the problem in the standard normal space. Thus, if the design point is actually representative of the most important regions of the failure domain, the OA leads to an impressive reduction in the variance of the failure probability estimator. However, it is worth noting that the design points and their neighbors do not always represent the most important regions of the failure domain, especially in high-dimensional problems. Moreover, the computational cost associated with the identification of the design points may be quite relevant which adversely affect the efficiency of the method [START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF]. Finally, the implementation of the OA method requires the definition of a rotation matrix in order to modify the coordinate system, which can be very difficult for high-dimensional problems.

A synthetic comparison of the stochastic simulation methods considered in this work is given in Table 1.9 (the second column, namely "Decisions", refers to parameters, distributions and other characteristics of the methods that have to be chosen or determined by the analyst in order to perform the simulation). -high computational cost (in presence of long-running models for determining system response and small failure probabilities)

SS

-express a small probability as a product of larger conditional probabilities -no assumptions about regularity of the limit state function (robustness)

-if limit state function is almost linear, few simulations suffice to achieve acceptable estimation accuracies -no necessity to estimate important direction α with excessive accuracy -even in the worst possible case (α orthogonal to optimal direction) the performance is at least comparable to standard MCS -determination of important direction α requires additional evaluation of system performance (with increase in the computational cost)

-for each sample drawn, two or three evaluations of system performance are necessary to estimate failure probability (with increase in the computational cost)

-essential restriction to standard normal space (Rosenblatt's or Nataf's transformations are required) [START_REF] Rosenblatt | Remarks on multivariate transformations[END_REF][START_REF] Nataf | Determination des distribution dont les marges sont donnees[END_REF] . -any probability level can be reached (also the very small ones of rare events)

-analytical expression for the system performance function is required -performance function must have the characteristics that one of the variables can be separated out from the others (1.6')

OA

The variance of the estimator (1.6') is In general, a relatively low number NT of simulations has to be carried out to obtain a sufficiently accurate estimate. A single evaluation would suffice for the ideal case in which the limit state function is linear and a Line Sampling direction α perpendicular to it has been identified [START_REF] Koutsourelakis | Reliability of structures in high dimensions, Part I: algorithms and application[END_REF]. [START_REF] Pradlwarter | Realistic and efficient reliability estimation for aerospace structures[END_REF]

Fig. 1 . 1 .

 11 Fig. 1.1. Sketch of the SS algorithm

Fig. 1 . 3 .

 13 Fig. 1.3. Examples of possible important unit vectors α 1 (left) and α 2 (right) pointing towards the corresponding failure domains F 1 (left) and F 2 (right) in a two-dimensional uncertain parameter space

  performance analyses are required to determine accurately the gradient at a given point of the performance function merical differentiation, e.g. the secant method(Ahammed and Melchers 2006;[START_REF] Fu | Monte Carlo sampling methods using Markov chains and their applications[END_REF].

  . A point 0 θ is taken in the failure domain F. This can be done by traditional Monte Carlo sampling or by engineering judgment when possible. Subsequently, 0 θ is used as the initial point of a Markov chain which lies entirely in the failure domain F. For that purpose a MCMC Metropolis-Hastings algorithm is employed to generate a sequence of Ns points { } …, Ns, are then averaged in order to obtain the LS important unit vector as

Fig. 1 . 4 .

 14 Fig.1.4. Line Sampling important unit vector α taken as the normalized "center of mass" of the failure domain F in the standard normal space. The "center of mass" of F is computed as an average of Ns failure points generated by means of a Markov chain starting from an initial failure point θ 0[START_REF] Koutsourelakis | Reliability of structures in high dimensions, Part I: algorithms and application[END_REF]

  the "design point" * θ of the problem (see Section 1.2.2.2); 3. Rotate the coordinate system (i.e., by means of a proper rotation matrix R ) so that the new coordinate n θ is in the direction of the axis defined by the design point * 1, 2, …, NT, are evaluated with an iterative algorithm which searches for the roots of the equation

F

  P ˆ) and the Figure Of Merit (FOM) of the method (defined as ( )

  in all the Cases considered. For instance, in Case 2, where the failure probability P(F) to be estimated is very small, i.e., P(F) = 4.4•10 -7 , the ratio is 11750 times larger than the ratio FOMSS/FOMMC. Notice that for the LS method even though the determination of the sampling important direction α (Section 1.2.2.2) and the calculations of the conditional one-dimensional failure probability estimates

ˆ

  = 426). In addition, LS performs comparably to the combination of IS and DR in the estimation of failure probabilities around 10 -7 : actually, in Case 1

-

  high-dimensional problem in the physical space into one-dimensional problems in the standard normal space -project the problem onto a line α pointing at the important regions of the failure domain -use line α almost perpendicular to the failure domain to reduce the variance of the estimates -one failure point to start the Markov chain for the determination of α

  -identification of the design point -rotation of system coordinates -solve one-dimensional problems along direction of design point / -if the design point is representative of the most important regions of the failure domain, then the variance is significantly reduced -design point frequently not representative of the most important regions of the failure domain (high-dimensional problems)-high computational cost associated to design point the independent conditional "one-dimensional" failure probability estimates

Fig. 1

 1 Fig. 1.1'. The Line Sampling procedure[START_REF] Pradlwarter | Realistic and efficient reliability estimation for aerospace structures[END_REF]

Table 1 .

 1 1. Names, descriptions and units of measure of the variables of the cracked plate model

Variables of the cracked plate model Name Description Unit of measure

	Kc	Critical stress intensity factor	MPa m
	a	Initial length of the defect	m
	F	Shape factor of the defect	/
	s	Uniform normal loading (stress) to which the plate is subject MPa

Table 1 .2. Probability

 1

		Case 0	Case 1	Case 2	Case 3
	x1 (K)	N(149.3, 22.2)	N(149.3, 22.2)	N(160, 18)	LG(149.3, 22.2)

distributions and parameters (i.e., means and standard deviations) of the uncertain variables x1, x2, x3, and x4 of the cracked plate model of Section 1.4.2 for the four case studies considered (i.e., Cases 0, 1, 2 and 3); the last row reports the values of the corresponding exact (i.e., analytically computed) failure probabilities, P(F)

[START_REF] Gille | Evaluation of failure probabilities in structural reliability with Monte Carlo methods[END_REF][START_REF] Gille | Probabilistic numerical methods used in the applications of the structural reliability domain[END_REF]

. N = Normal distribution; LG = Lognormal distribution

Table 1

 1 Nsys of system response analyses (i.e., model evaluations) is also reported. Notice that for both SS and LS the actual number Nsys of system analyses does not coincide with the total number NT of random samples drawn (i.e., NT = 50000). In particular, in the SS method, the presence of repeated conditional samples in each Markov chain (used to gradually populate the intermediate event regions) allows a reduction in the number of model evaluations required: actually, one evaluation is enough for all identical samples (see Appendix 1). In the LS method, instead, the actual number Nsys of system analyses is given by Nsys = Ns + 2•NT: in particular, Ns = 2000 analyses are performed to generate the Markov chain used to compute the important unit vector α as the normalized "center of mass" of the failure domain F (Section 1.2.2.2); the 2•NT analyses are carried out to compute the NT conditional one-dimensional probability estimates Results of the application of standard MCS, SS and LS to the reliability analysis of Cases 0 (Reference), 1, 2 and 3 of the cracked plate model of Section 1.4.2; the values of the performance indicators used to compare the effectiveness of the methods (i.e.,

	.3 reports the values of ( ) F P ˆ, σˆ, δ , FOM,	σ	MC σ ˆ	meth	and FOM-
	meth/FOMMC obtained by Standard MCS, SS and LS in Cases 0, 1, 2 and 3 (Section

Table 1.3.

.v., δ N sys FOM meth MC

								σ ˆ	σ /	FOMmeth/FOMMC
	Standard MCS	4.400•10 -7	3.000•10 -6	6.667	50000	2.222•10 6	1				1
	SS	4.679•10 -7	6.890•10 -8	1.473•10 -1	49888	4.222•10 9	43.54		1.9•10 3
	LS	4.381•10 -7	4.447•10 -10	1.015•10 -3	102000 4.959•10 13 6746.7	2.2•10 7
	Case 3
		P ˆ	() F	σ ˆ	c.o.v., δ	N sys	FOM	σ ˆ	MC	σ /	meth	FOMmeth/FOMMC
	Standard MCS	3.000•10 -4	7.745•10 -5	2.582•10 -1	50000	3.334•10 3 1				1
	SS	3.183•10 -4	2.450•10 -5	7.697•10 -2	49907	3.339•10 4 3.16			10.01
	LS	3.068•10 -4	1.817•10 -7	5.923•10 -4	102000 3.028•10 8 426.16		9.1•10 4

Table 1 . 4 .

 14 Values of the performance indicatorsComparing Table1.3 and Table 1.4, it can be seen that LS performs significantly better than IS and DR in all the case studies considered: in particular, in Cases 1 and 2 the values of the performance indicators

	σ	MC σ ˆ	meth	and FOMmeth/FOMMC obtained by IS,

Table 1 .

 1 5. Names, descriptions and units of measure of the variables of the thermal fatigue crack growth model

	Variables of the thermal fatigue crack growth model
	Name	Description	Unit of measure
	a0	Initial size of the defect	[m]
	ac	Critical size of the defect	[m]
	T0	Initial time	[years]
	T	Final time	[years]
	C	Parameter of the material	/
	m	Parameter of the material	/
	f(R)	Correction factor	/
	Nc	Number of cycles per year	/
	S	Stress per cycle	[MPa]

Table 1 .

 1 7 reports the values of ( ) /FOMMC obtained by Standard MCS, SS and LS in the Cases 1 and 2 of Section 1.5.3; the actual number Nsys of system response analyses (i.e., model evaluations) is also reported. Results of the application of standard MCS, SS and LS to the reliability analysis of Cases 1 and 2 of the thermal fatigue crack growth model of Section 1.5.2; the values of the performance indicators used to compare the effectiveness of the methods (i.e.,

						F P ˆ, σˆ, δ , FOM,	σ	MC σ ˆ	meth	and FOM-
	Table 1.7. meth MC σ σ ˆ	and FOM-
	meth/FOMMC) are highlighted in bold
	Case 1
		P ˆ	() F	σ ˆ	c.o.v., δ	N sys	FOM			σ ˆ	MC	σ /	meth	FOMmeth/FOMMC
	Standard MCS	2.500•10 -4	7.905•10 -5	3.162•10 -1	40000 4.001•10 3 1	1
	SS	3.006•10 -4	3.214•10 -5	1.069•10 -1	40019 2.419•10 4 2.46	6.05
	LS	3.768•10 -4	4.610•10 -7	1.223•10 -3	82000 5.737•10 7 171.46	1.434•10 4
	Case 2
		P ˆ	() F	σ ˆ	c.o.v., δ	N sys	FOM		σ ˆ	MC	σ /	meth	FOMmeth/FOMMC
	Standard MCS 1.780•10 -5	2.269•10 -5	1.102	40000 4.860•10 4		1		1
	SS	1.130•10 -5	1.653•10 -6	1.462•10 -1	39183 9.341•10 6		13.73	192.36
	LS	1.810•10 -5	2.945•10 -8	1.627•10 -3	81999 1.188•10 13 770.02	2.892•10 5

meth

  Section 1.3.3). and IS and DR for Cases 1 and 2 of the thermal fatigue crack growth model of Sections 1.5.1 and 1.5.2 are summarized in Table1.8.

	The values of the performance indicators	σ ˆ	meth / MC σ	and FOMmeth/FOMMC ob-
	tained by IS Table 1.8. Values of the performance indicators	σ ˆ	MC	σ /	meth	and FOMmeth/FOMMC obtained by IS
	and IS + DR when applied for the reliability analysis of Cases 1 and 2 of the thermal fatigue
	crack growth model of Section 1.5.2 (Gille 1998 and 1999)
	Case 1
		σ ˆ	MC	σ /	meth		FOMmeth/FOMMC
	IS	16.9				424.36
	IS + DR 65.4				864.36
	Case 2
		σ ˆ	MC	σ /	meth		FOMmeth/FOMMC
	IS	41.1				4.396•10 3
	IS + DR 172.4			8.317•10 3

Table 1 . 9 .

 19 Synthetic comparison of the stochastic simulation methods considered in this work (Part I)

	Method	Simulation concepts	Decisions	Advantages	Drawbacks
	Standard	-repeat random sampling of	/	-samples the full range of each input variable
	MCS	possible system configurations		-consistent performance in spite of complexity
				and dimension of the problem
				-accuracy easily assessed
				-no need for simplifying assumptions nor
				surrogate models
				-no complex elaborations of the original model
				-identification of nonlinearities, thresholds and
				discontinuities
				-simplicity

Table 1 .

 1 9. Synthetic comparison of the stochastic simulation methods considered in this work(Continued)

	Method	Simulation concepts	Decisions	Advantages	Drawbacks
	IS	-repeated random sampling of	-construction/choice of	-if the ISD is similar to optimal one: significant	-many system behavior insights and
		possible system configurations	the ISD	increase in estimation accuracy (or, conversely,	and much modeling work needed for
		-sample from Importance		reduction in sample size for given accuracy)	identification of good ISD
		Sampling Density (ISD) to			-inappropriate ISD leads to worse
		generate more samples in the region			estimates compared to Standard MCS
		of interest (e.g., low probability of
		occurrence)
	DR	-express failure event in such a way	-random variable to be	-remove zero values included in the Standard
		as to highlight one random variable	separated from others	MCS estimator (reduced variance)
		-estimate failure probability as expected
		value of the CDF of the chosen
		variable conditional on the remaining
		(n -1) variables

Appendix 1. Markov Chain Monte Carlo (MCMC) Simulation

Markov Chain Monte Carlo (MCMC) simulation comprises a number of powerful simulation techniques for generating samples according to any given probability distribution [START_REF] Metropolis | Equations of state calculations by fast computing machines[END_REF].

In the context of the reliability assessment of interest in the present work, MCMC simulation provides an efficient way for generating samples from the multidimensional conditional PDF ) | ( F q x

. The distribution of the samples thereby generated tends to the multidimensional conditional PDF ) | ( F q x as the length of the Markov chain increases. In the particular case of the initial sample 1 x being distributed exactly as the multidimensional conditional PDF ) | ( F q x , then so are the subsequent samples and the Markov chain is always stationary [START_REF] Au | Estimation of small failure probabilities in high dimensions by subset simulation[END_REF].

In the following it is assumed without loss of generality that the components of x are independent, that is,

, where ) ( j j x q denotes the onedimensional PDF of j

x (Au and Beck 2001).

To illustrate the MCMC simulation algorithm with reference to a generic failure region Fi, let

Sample a precandidate value

a. Compute the acceptance ratio:

x of the j th element of 1 ~+ u x as follows: x is taken as the next sample

x , depending on whether the candidate 1 ~+ u x lies in the failure region Fi or not.

Finally, notice that in this work, the one-dimensional proposal PDF * j p , j = 1, 2, …, n, is chosen as a symmetrical uniform distribution centered at the current sample xj, j = 1, 2, …, n, with width 2lj, where lj is the maximum step length, i.e. the maximum allowable distance that the next sample can depart from the current one. The choice of lj is such that the standard deviation of * j p is equal to that of qj, j = 1, 2, …, n.

Appendix 2. The Line Sampling algorithm

The LS algorithm proceeds as follows [START_REF] Pradlwarter | Realistic and efficient reliability estimation for aerospace structures[END_REF]:

. Any of the methods summarized in Section 1.2.2.2 can be employed to this purpose.

Notice that the computation of α implies additional system analyses, which substantially increase the computational cost associated to the simulation method (Section 1.2.2.2). 2. From the original multidimensional joint probability density function () is the scalar product between α and k θ , k = 1, 2, ..., NT.

Finally, it is worth noting that since the standard Gaussian space is isotropic, both the scalar c k and the vector ⊥ , k θ are also standard normally distributed [START_REF] Pradlwarter | Application of line sampling simulation method to reliability benchmark problems[END_REF]). b. Compute the value k c as the intersection between the limit state function


P Baraldi

E Zio

G Gola

D Roverso

M Hoffmann

Robust nuclear signal reconstruction by a novel ensemble model aggregation procedure

Monitoring of sensor operation is important for detecting anomalies and reconstructing the correct values of the signals measured. This can be done, for example, with the aid of auto-associative regression models. However, in practical applications difficulties arise because of the need of handling large numbers of signals. To overcome these difficulties, ensembles of reconstruction models can be used. Each model in the ensemble handles a small group of signals and the outcomes of all models are eventually combined to provide the final outcome. In this work, three different methods for aggregating the model outcomes are investigated and a novel procedure is proposed for obtaining robust ensembleaggregated outputs. Two applications are considered concerning the reconstruction of 920 simulated signals of the Swedish Forsmark-3 Boiling Water Reactor (BWR) and 215 signals measured at the Finnish Pressurized Water Reactor (PWR) situated in Loviisa.

Introduction

Sensors are placed at various locations in a production plant to monitor its state and consequently operate its control and protection systems. Early detection of sensor malfunctions and reconstruction of the signals measured is then necessary [START_REF] Hoffmann | On-line Monitoring for Calibration Reduction[END_REF][START_REF] Hoffmann | Signal Grouping Algorithm for an Improved On-line Calibration Monitoring System[END_REF].

In real plants there are thousands of sensors whose signals cannot be effectively handled by a single reconstruction model [START_REF] Hoffmann | Signal Grouping Algorithm for an Improved On-line Calibration Monitoring System[END_REF][START_REF] Roverso | Solutions for plant-wide on-line calibration monitoring[END_REF][START_REF] Zio | Genetic Algoritms for Grouping of Signals for System Monitoring and Diagnostics[END_REF][START_REF] Baraldi | Genetic algorithms for signal grouping in sensor validation: a comparison of the filter and wrapper approaches[END_REF]. A possible way to overcome this problem is by subdividing the set of signals into small overlapping groups, building a corresponding number of reconstruction models and aggregating their outcomes to provide the ensemble-aggregated output [START_REF] Perrone | When networks disagree: ensemble methods for hybrid neural networks[END_REF][START_REF] Krogh | Neural network ensembles, cross-validation and active learning[END_REF][START_REF] Sharkey | On combining artificial neural nets[END_REF][START_REF] Gola | Signal Grouping for Sensor Validation: a Multi-Objective Genetic Algorithm Approach[END_REF][START_REF] Baraldi | Reconstruction of faulty signals by an ensemble of principal component analysis models optimized by a multi-objective genetic algorithm[END_REF][START_REF] Gola | Reconstructing signals for sensor validation by a GAoptimized ensemble of PCA models[END_REF]. This latter task is crucial for the good performance of the ensemble.

In this work, three methods of aggregating the outcomes of the models in an ensemble are first investigated: Simple Mean (SM), Median (MD) and Trimmed Mean (TM) [START_REF] Polikar | Ensemble based systems in decision making[END_REF]. SM amounts to using all the outcomes of the individual models in the ensemble; MD considers only the central value in the distribution of the models outcomes; TM discards the outcomes in the tails of the distribution. A novel procedure based on the combination of the MD and TM aggregation approaches is then proposed.

The groups at the basis of the ensemble are generated by randomly selecting their signals [START_REF] Polikar | Ensemble based systems in decision making[END_REF][START_REF] Bryll | Attribute bagging: improving accuracy of classifiers ensembles by using random feature subsets[END_REF][START_REF] Baraldi | A randomized model ensemble approach for reconstructing signals from faulty sensors[END_REF]. This allows injecting high diversity in the group structure; for each group, a regression model based on Principal Component Analysis (PCA) [START_REF] Jolliffe | Principal Component Analysis[END_REF][START_REF] Diamantaras | Principal component neural networks: theory and applications[END_REF][START_REF] Scholkopf | Kernel principal component analysis[END_REF][START_REF] Moore | Principal component analysis in linear systems: Controllability, observability, and model reduction[END_REF] is developed.

The paper is organized as follows. In Section 2, the three methods for aggregating the models outcomes are described. Two applications are illustrated in Section 3: the first one concerns the reconstruction of a data set of 920 signals obtained with the HAMBO code [START_REF] Karlsson | The HAMBO BWR Simulator of HAMMLAB[END_REF] which simulates the Forsmark-3 Boiling Water Reactor (BWR) located in Sweden; the second addresses the reconstruction of 215 signals measured at the Finnish Pressurized Water Reactor (PWR) situated in Loviisa.

The novel procedure is applied and tested on both applications. Conclusions on the advantages and limitations of the aggregation methods are drawn in the last Section.

Methods of aggregation of the outcomes of the models in the ensemble

Given the set of The K diverse signal groups generated are used as bases for developing a corresponding number of PCA reconstruction models. To do this, the data set X of N signal patterns available is partitioned into a training set X TRN (made of TRN N patterns) and a test set X TST (made of TST N patterns). The former is used to train the individual models, whereas the latter is used to verify the ensemble performance in the signal reconstruction task.

Each signal i is present in a number i K of groups and thus a corresponding number of individual PCA models provide its reconstruction. Three different methods are here considered to aggregate the outcomes of these individual models: Simple Mean (SM), Median (MD) and Trimmed Mean (TM).

SM amounts to considering all the i K available model outcomes in a simple average [12, 14, 20, 21]; for the generic pattern 1, 2,..., TST tN  :

,

1 1 ˆ( ) ( ) 1, 2,..., i K E SM k ii i k f t f t i n K    (1)
where , ˆ() E SM i ft is the ensemble aggregate of the reconstructions ˆ() k i ft of pattern t of signal i provided by the individual models 1, 2,..., i kK  containing signal i .

The MD aggregation approach amounts to considering for the generic pattern t of signal i only the single outcome ˆ() C k i ft lying in the centre of the distribution of the i K model outcomes for that pattern, i.e.:

,

ˆ( ) ( ) 1, 2,..., C k E MD ii f t f t i n  (2)
where C k denotes the index of the model whose outcome is central with respect to the reconstructed values of the i K models including signal i .

The TM approach amounts to discarding a fraction TM  of the i K outcomes of signal i at the tails of the distribution and then simply averaging the  

12 TM i TM i KK   remaining outcomes: , 1 1 ˆ( ) ( ) 1, 2,..., TM i K E TM k ii TM i k f t f t i n K    (3)
In a sense, the TM approach represents a compromise between the SM and MD methods, for it allows discarding the tails of the distribution while still considering multiple outcomes in the ensemble reconstruction of the signal. In fact, 11 0, 1 2

TM i K         
, the lower limit corresponding to considering all the available outcomes (SM), the upper limit representing the case of using the single central outcome (MD).

Finally, to evaluate the performance of the ensemble aggregates, first the absolute signal reconstruction error is computed1 :

, (

TST N E SM MD TM E SM MD TM i i i TST t f t f t N     (4) , , ) , ( , , ) 1 1 ( ) ( )
Then, the ensemble performance index is retained as the average of the absolute signal reconstruction errors (Eq. 4):

, , , ) , ( , , ) 1 1
n E SM MD TM E SM MD TM i i n     ( (5)

Applications

The applications concern the reconstruction of a set of 920 simulated signals of the Swedish Forsmark-3 BWR and a set of 215 signals measured at the Finnish PWR located in Loviisa. Table 1 reports the main characteristics and parameters adopted in the two case studies.

The PCA models have been constructed with the code http://lib.stat.cmu.edu/multi/pca, adapted to perform the signal reconstruction task of interest here.

The robustness of the ensemble has been tested on the reconstruction of signals when in presence of sensor failures, e.g. random noises or offsets. Within the proposed ensemble approach, a faulty sensor sends a faulty signal in input to the PCA models which contain that signal; each PCA model reconstructs the value of the signal by auto-associating the information of the non-faulty signals in the model; finally, the ensemble aggregate of the individual outcomes of the models is obtained.

To verify the performance of the ensemble, disturbs are introduced in the test set of both case studies to challenge the corresponding ensemble in the reconstruction of the true, undisturbed signal. More

 

for MD outperforming TM or viceversa) and the corresponding average error gain in using the best performing one. For the Forsmark-3 case study (Table 3a), approximately half of the signals are better reconstructed by MD and half by TM; nevertheless, on disturbed signals, the average error reduction achieved by MD is larger, i.e. the improvements of using MD are more relevant than those obtained with TM, even though on a slightly smaller number of signals; on the other hand, in the Loviisa case study (Table 3b), results show that on disturbed signals MD performs better on a large number of signals (more than two thirds) with a considerably higher average error reduction than TM. On the basis of these insights, a novel procedure is proposed for combining MD and TM with the aim of exploiting the advantages of both methods. The idea is to define a reconstruction scheme finalized to discern automatically which method is the most effective for reconstructing a signal. The procedure (sketched in Figure 1) is based on the ensemble reconstruction errors obtained on the training set by the MD and TM methods. The conjecture is that if the training patterns of a signal are better reconstructed by one method (MD or TM), then the same method will better reconstruct also the signal's test patterns. The procedure is tailored with respect to the disturbed training signals for application to the reconstruction of the test signals (be it undisturbed or disturbed).

For each signal i , the ensemble reconstruction errors

Conclusions

This work addresses the problem of reconstructing the correct signal values measured by faulty sensors in nuclear power plants. The task is rather complex due to the large number of signals involved. A feasible approach to treat the high dimensionality of the problem is to resort to an ensemble of models for signal reconstruction.



  by the MD and TM methods, respectively, are computed using the disturbed training set: then, the method by which reconstructing the test patterns of signal i is MD if

Figure 1 .

 1 Figure 1. Sketch of the novel procedure for combining the MD and TM aggregation methods

  To construct the ensemble, in this work the set of signals is first subdivided into small, overlapping groups by random feature selection. Then, one PCA-based auto-associative reconstruction model is developed based on the signals of each group. The outcomes of the models have finally been aggregated by three methods to obtain the ensemble reconstruction output: Simple Mean (averaging all the available model outcomes), Median (taking only the single outcome in the centre of the distribution of outcomes) and Trimmed Mean (discarding the model outcomes lying in the tails of the distribution). These methods have been applied to two high-dimensional problems regarding the reconstruction of 920 simulated signals of a nuclear boiling water reactor and 215 signals measured at a pressurized water reactor. The advantages of discarding the outlying outcomes in the ensemble aggregation have been proved in both cases. A novel procedure has then been developed for combining the Median and Trimmed Mean approaches in a way to fully exploit their benefits. The procedure has been applied to both case studies showing an increased robustness in reconstructing signals affected by disturbs.

Table 1 .

 1 Main characteristics and parameters of the two case studiesTable2reports the results for the two case studies in terms of the values of the ensemble reconstruction performance indexes (Eq. 5) on the test patterns obtained by the three aggregation methods on undisturbed ( ,

	precisely, the signals of a test pattern

TST

X .

3 Once m and i K are set, the number of groups to generate by random feature selection is obtained from the relation i mK nK  . For a more detailed explanation of the procedure, the interested reader may refer to

[START_REF] Baraldi | A randomized model ensemble approach for reconstructing signals from faulty sensors[END_REF]

.

Table 2 .

 2 Ensemble reconstruction performance indexes obtained with SM, MD and TM aggregations on undisturbed and disturbed signals for the Forsmark-3 and Loviisa case studies.To further delve into the MD and TM approaches, their reconstruction capabilities on each signal are analyzed. For the Forsmark-3 and Loviisa applications, respectively, Tables3a and 3bsummarize numerically the comparison on the training and test sets between MD and TM in terms of the number of signals for which one method outperforms the other (i.e.

	,, E MD E TM
	ii

Table 3b .

 3b Loviisa case study: comparison of MD and TM aggregation methods on undisturbed and disturbed training and test signals

	Aggregation	Number of signals for which the aggregation method is better / Average error reduction (	3 10  	)
	method	Training		Test
		Undisturbed	Disturbed	Undisturbed	Disturbed
	MD	356 / 0.14	451 / 2.92	474 / 1.58	420 / 2.97
	TM	564 / 0.13	469 / 1.60	446 / 1.34	500 / 1.92
	Table 3a. Forsmark-3 case study: comparison of MD and TM aggregation methods on undisturbed and disturbed training
			and test signals
	Aggregation	Number of signals for which the aggregation method is better / Average error reduction (	3 10  	)
	method	Training		Test
		Undisturbed	Disturbed	Undisturbed	Disturbed
	MD	78 / 0.79	155 / 5.68	113 / 4.55	163 / 6.12
	TM	137 / 0.73	60 / 2.88	102 / 4.89	52 / 2.91

Table 4 .

 4 Ensemble reconstruction errors obtained with MD, TM and MX approaches on undisturbed and disturbed signals for the two case studies

	2 10  	Forsmark-3 Boiling Water Reactor (BWR) Loviisa Pressurized Water Reactor (PWR)
			MD	TM	MX	MD	TM	MX
		, EU	1.878	1.894	1.890	0.956	0.957	0.958
		, ED	7.607	7.638	7.298	5.169	5.261	5.023

In the application that follows, each signal of the validation set has been previously normalized in the range [0.2, 1], for convenience.