K Konakli

B Sudret

Global sensitivity analysis using low-rank tensor approximations

Keywords:

In the context of global sensitivity analysis, the Sobol' indices constitute a powerful tool for assessing the relative significance of the uncertain input parameters of a model. We herein introduce a novel approach for evaluating these indices at low computational cost, by postprocessing the coefficients of polynomial meta-models belonging to the class of low-rank tensor approximations. Meta-models of this class can be particularly efficient in representing responses of high-dimensional models, because the number of unknowns in their general functional form grows only linearly with the input dimension. The proposed approach is validated in example applications, where the Sobol' indices derived from the meta-model coefficients are compared to reference indices, the latter obtained by exact analytical solutions or Monte-Carlo simulation with extremely large samples. Moreover, low-rank tensor approximations are confronted to the popular polynomial chaos expansion meta-models in case studies that involve analytical rank-one functions and finite-element models pertinent to structural mechanics and heat conduction. In the examined applications, indices based on the novel approach tend to converge faster to the reference solution with increasing size of the experimental design used to build the meta-model.

Introduction

Robust predictions via computer simulation necessitate accounting for the prevailing uncertainties in the parameters of the computational model. Uncertainty quantification provides the mathematically rigorous framework for propagating the uncertainties surrounding the model 1 input to a response quantity of interest. It comprises three fundamental steps [START_REF] Sudret | Uncertainty propagation and sensitivity analysis in mechanical modelscontributions to structural reliability and stochastic spectral methods[END_REF][START_REF] De Rocquigny | Modelling under risk and uncertainty -An introduction to statistical, phenomenological and computational methods[END_REF]: First, the model representing the physical system under consideration is defined; the model maps a given set of input parameters to a unique value of the response quantity of interest. The second step involves the probabilistic description of the input parameters by incorporating available data, expert judgment or a combination of both. In the third step, the uncertainty in the input is propagated upon the response quantity of interest through repeated evaluations of the model for appropriate combinations of the input parameters. In cases when the uncertainty in the response proves excessive or when it is of interest to reduce the dimensionality of the model, sensitivity analysis may be employed to rank the input parameters with respect to their significance for the response variability. Accordingly, important parameters may be described in further detail, whereas unimportant ones may be fixed to nominal values.

Methods of sensitivity analysis can be generally classified into local and global methods.

Local methods are limited to examining effects of variations of the input parameters in the vicinity of nominal values. Global methods provide more complete information by accounting for variations of the input parameters in their entire domain. Under the simplifying assumption of linear or nearly linear behavior of the model, global sensitivity measures can be computed by fitting a linear-regression model to a set of input samples and the respective responses (see e.g. [START_REF] Iooss | A review on global sensitivity analysis methods[END_REF]; [START_REF] Wei | Variable importance analysis: A comprehensive review[END_REF] for definitions of such measures). The same methods can be employed in cases with models that behave nonlinearly but monotonically, after applying a rank transformation on the available data. Variance-based methods represent a more powerful and versatile approach, also applicable to nonlinear and non-monotonic models. These methods, known as functional ANOVA (denoting ANalysis Of VAriance) techniques in statistics, rely upon the decomposition of the response variance as a sum of contributions from each input parameter or their combinations [START_REF] Efron | The Jacknife estimate of variance[END_REF]. The Sobol' indices, originally proposed in Sobol' (1993), constitute the most popular tool thereof. Although these indices have proven powerful in a wide range of applications, their definition is ambiguous in cases with dependent input variables, which has led to different extensions of the original framework [START_REF] Da Veiga | Local polynomial estimation for sensitivity analysis on models with correlated inputs[END_REF][START_REF] Li | Global sensitivity analysis for systems with independent and/or correlated inputs[END_REF][START_REF] Kucherenko | Estimation of global sensitivity indices for models with dependent variables[END_REF][START_REF] Mara | Variance-based sensitivity indices for models with dependent inputs[END_REF][START_REF] Zhang | A new framework of variance based global sensitivity analysis for models with correlated inputs[END_REF].

An alternative perspective is offered by the distribution-based indices, which are well-defined regardless of the dependence structure of the input [START_REF] Borgonovo | A new uncertainty importance measure[END_REF][START_REF] Liu | A new importance measure for sensitivity analysis[END_REF][START_REF] Borgonovo | Transformations and invariance in the sensitivity analysis of computer experiments[END_REF][START_REF] Zhou | Sparse grid integration based solutions for momentindependent importance measures[END_REF][START_REF] Greegar | Global response sensitivity analysis of uncertain structures[END_REF]. The key idea is to use an appropriate distance measure to evaluate the effect of suppressing the uncertainty in selected variables on the Probability Density Function (PDF) or the Cumulative Distribution Function (CDF) of the response. These indices are especially useful when consideration of the variance only is not deemed sufficient to characterize the response uncertainty. However, contrary to the Sobol' indices, they do not sum up to unity, which may hinder interpretation. For further information on global sensitivity analysis methods, the interested reader is referred to the review papers [START_REF] Saltelli | Global Sensitivity Analysis -The Primer[END_REF][START_REF] Iooss | A review on global sensitivity analysis methods[END_REF][START_REF] Wei | Variable importance analysis: A comprehensive review[END_REF][START_REF] Borgonovo | Sensitivity analysis: A review of recent advances[END_REF].

It should be mentioned that different classifications of sensitivity analysis techniques can be found in the literature. In cases when one needs to perform a fast exploration of the model behavior with respect to a possibly large number of uncertain input parameters, the so-called screening methods may be employed. These methods can provide a preliminary ranking of the importance of the various input parameters at low computational cost before more precise and costly methods are applied. The Cotter method [START_REF] Cotter | A screening design for factorial experiments with interactions[END_REF] and the Morris method [START_REF] Morris | Factorial sampling plans for preliminary computational experiments[END_REF] are widely used screening methods, with the latter covering the input space in a more exhaustive manner. The more recently proposed derivative-based global sensitivity measures can also be classified into this category, while they also provide upper bounds for the Sobol' indices [START_REF] Sobol | Derivative based global sensitivity measures and their link with global sensitivity indices[END_REF][START_REF] Lamboni | Derivative-based global sensitivity measures: general links with Sobol' indices and numerical tests[END_REF][START_REF] Sudret | Computing derivative-based global sensitivity measures using polynomial chaos expansions[END_REF].

The focus of the present paper is on sensitivity analysis by means of Sobol' indices. We limit our attention to cases with independent input and address the computation of these indices for high-dimensional expensive-to-evaluate models, which are increasingly used across engineering and sciences. Various methods have been investigated for computing the Sobol' indices based on Monte Carlo simulation [START_REF] Archer | Sensitivity measures, ANOVA-like techniques and the use of bootstrap[END_REF][START_REF] Sobol | Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates[END_REF][START_REF] Saltelli | Making best use of model evaluations to compute sensitivity indices[END_REF][START_REF] Sobol | Global sensitivity indices for nonlinear mathematical models[END_REF][START_REF] Saltelli | Variance based sensitivity analysis of model output. Design and estimator for the total sensitivity index[END_REF]; because of the large number of model evaluations required, these methods are not affordable for computationally costly models. To overcome this limitation, more efficient estimators have recently been proposed [START_REF] Sobol | Estimating the approximation error when fixing unessential factors in global sensitivity analysis[END_REF][START_REF] Janon | Asymptotic normality and efficiency of two Sobol index estimators[END_REF][START_REF] Kucherenko | Application of the control variate technique to estimation of total sensitivity indices[END_REF][START_REF] Owen | Variance components and generalized Sobol' indices[END_REF]. A different approach is to substitute a complex model by a meta-model, which has similar statistical properties while maintaining a simple functional form (see e.g. [START_REF] Oakley | Probabilistic sensitivity analysis of complex models: a Bayesian approach[END_REF]; [START_REF] Marrel | Calculations of Sobol indices for the Gaussian process metamodel[END_REF]; [START_REF] Storlie | Implementation and evaluation of nonparametric regression procedures for sensitivity analysis of computationally demanding models[END_REF]; [START_REF] Zuniga | Metamodelling with independent and dependent inputs[END_REF]; [START_REF] Zhang | An effective approximation for variance-based global sensitivity analysis[END_REF]; Le [START_REF] Gratiet | Metamodel-based sensitivity analysis: polynomial chaos expansions and Gaussian processes[END_REF] for global sensitivity analysis with various types of meta-models). [START_REF] Sudret | Global sensitivity analysis using polynomial chaos expansions[END_REF] proposed to compute the Sobol' indices by postprocessing the coefficients of Polynomial Chaos Expansion (PCE) meta-models. The key concept in PCE is to expand the model response onto a basis made of orthogonal multivariate polynomials in the input variables. The computational cost of the associated Sobol' indices essentially reduces to the cost of estimating the PCE coefficients, which can be curtailed by using sparse PCE (Blatman and Sudret, 2010b). The PCE-based approach for computing the Sobol' indices is employed by a growing number of researchers in various fields including hydrogeology [START_REF] Fajraoui | Use of global sensitivity analysis and polynomial chaos expansion for interpretation of nonreactive transport experiments in laboratory-scale porous media[END_REF][START_REF] Formaggia | Global sensitivity analysis through polynomial chaos expansion of a basin-scale geochemical compaction model[END_REF][START_REF] Deman | Using sparse polynomial chaos expansions for the global sensitivity analysis of groundwater lifetime expectancy in a multi-layered hydrogeological model[END_REF], geotechnics [START_REF] Al Bittar | Bearing capacity of strip footings on spatially random soils using sparse polynomial chaos expansion[END_REF], ocean engineering [START_REF] Alexanderian | Global sensitivity analysis in an ocean general circulation model: a sparse spectral projection approach[END_REF], biomedical engineering [START_REF] Huberts | Applicability of the polynomial chaos expansion method for personalization of a cardiovascular pulse wave propagation model[END_REF], hybrid dynamic simulation (Abbiati, Marelli, Bursi, Sudret, and Stojadinovic, Abbiati et al.) and electromagnetism [START_REF] Kersaudy | Stochastic analysis of scattered field by building facades using polynomial chaos[END_REF][START_REF] Yuzugullu | Global sensitivity analysis of a morphology based electromagnetic scattering model[END_REF]. Unfortunately, the PCE approach faces the so-called curse of dimensionality, meaning the exploding size of the candidate basis with increasing dimension.

The goal of this paper is to derive a novel approach for solving global sensitivity analysis problems in high dimensions. To this end, we make use of a recently emerged technique for building meta-models with polynomial functions based on Low-Rank Approximations (LRA) [START_REF] Nouy | Proper generalized decompositions and separated representations for the numerical solution of high dimensional stochastic problems[END_REF]Chevreuil et al., 2013;[START_REF] Doostan | Non-intrusive low-rank separated approximation of high-dimensional stochastic models[END_REF][START_REF] Hadigol | Partitioned treatment of uncertainty in coupled domain problems: A separated representation approach[END_REF][START_REF] Validi | Low-rank separated representation surrogates of high-dimensional stochastic functions: Application in Bayesian inference[END_REF]Konakli and Sudret, 2015b). LRA express the model response as a sum of a small number or rank-one tensors, which are products of univariate functions in each of the input parameters.

Because the number of unknown coefficients in LRA grows only linearly with the input dimension, this technique is particularly promising for dealing with cases of high dimensionality. We herein derive analytical expressions for the Sobol' sensitivity indices based on the general functional form of LRA with polynomial bases. As in the case of PCE, the computational cost of the LRA-based Sobol' indices reduces to the cost of estimating the coefficients of the meta-model.

Once a LRA meta-model of the response quantity of interest is available, the Sobol' indices can be calculated with elementary operations at nearly zero additional computational cost.

The paper is organized as follows: In Section 2, we review the basic concepts of Sobol' sensitivity analysis and define the corresponding sensitivity indices. In Section 3, we describe the mathematical setup of non-intrusive meta-modeling and define error measures that characterize the meta-model accuracy. After reviewing the computation of Sobol' indices using PCE metamodels in Section 4, we introduce the LRA-based approach in Section 5. In this, we first detail a greedy algorithm for the construction of LRA in a non-intrusive manner and then, use their general functional form to derive analytical expressions for the Sobol' indices. In Section 6, we demonstrate the accuracy of the proposed method by comparing the LRA-based indices to reference ones, with the latter representing the exact solution or Monte-Carlo estimates relying on a large sample of responses of the actual model. Furthermore, we examine the comparative performance of the LRA-and PCE-based approaches in example applications that involve analytical rank-one functions and finite-element models pertinent to structural mechanics and heat conduction. The main findings are summarized in Section 7.

Sobol' Sensitivity analysis

We consider a computational model M describing the behavior of a physical or engineered system of interest. Let X = {X 1 , . . . , X M } denote the M -dimensional input random vector of the model with prescribed joint PDF f X . Due to the input uncertainties embodied in X, the model response becomes random. By limiting our focus to a scalar response quantity Y , the computational model represents the map:

X ∈ D X ⊂ R M -→ Y = M(X) ∈ R, (1)
where D X denotes the support of X.

Sobol' sensitivity analysis aims at apportioning the uncertainty in Y , described by its variance, to contributions arising from the uncertainty in individual input variables and their interactions. As explained in the Introduction, the theoretical framework described in the sequel is confined to the case when the components of X are independent. Under this assumption, the joint PDF f X is the product of the marginal PDF f Xi of each input parameter.

Sobol' decomposition

Assuming that M is a square-integrable function with respect to the probability measure associated with f X , its Sobol' decomposition in summands of increasing dimension is given by Sobol'

(1993):

M(X) = M 0 + M i=1 M i (X i ) + 1≤i<j≤M M i,j (X i , X j ) + . . . + M 1,2, ... ,M (X) = M 0 + u⊂{1, ... ,M } u =∅ M u (X u ), (2)
where u = {i 1 , . . . , i s }, 1 ≤ s ≤ M , denotes a subset of {1, . . . , M } and X u = {X i1 , . . . , X is } is the subvector of X containing the variables of which the indices comprise u.

The uniqueness of the decomposition is ensured by choosing summands that satisfy the conditions:

M 0 = E [M(X)] (3)
and

E [M u (X u ) M v (X v )] = 0 ∀ u, v ⊂ {1, . . . , M }, u = v. (4)
Note that the above condition implies that all summands {M u , u = ∅}, in Eq. ( 2) have zero mean values. A recursive construction of summands satisfying the above conditions is obtained as:

M i (X i ) =E [M(X)|X i ] -M 0 M i,j (X i , X j ) =E [M(X)|X i , X j ] -M i (X i ) -M j (X j ) -M 0
(5) and so on. By introducing:

M u (X u ) = E [M(X)|X u ] , (6)
the above recursive relationship is written in the general form:

M u (X u ) = M u (X u ) - v⊂u v =u M v (X v ). (7)

Sobol' sensitivity indices

The uniqueness and orthogonality properties of the Sobol' decomposition allow the following decomposition of the variance D of M(X):

D = Var [M(X)] = M i=1 D i + 1≤i<j≤M D i,j + . . . + D 1,2, ... ,M = u⊂{1, ... ,M } u =∅ D u , (8)
where D u denotes the partial variance:

D u = Var [M u (X u )] = E (M u (X u )) 2 . ( 9
)
The Sobol' index S u , defined as:

S u = D u D = Var [M u (X u )] Var [M(X)] , (10)
represents the fraction of the total variance that is due to the interaction between the components of X u , i.e. S i1, ... ,is describes the influence from the interaction between variables {X i1 , . . . , X is }. By definition, the Sobol' indices satisfy:

M i=1 S i + 1≤i<j≤M S i,j + . . . + S 1,2, ... ,M = u⊂{1, ... ,M } u =∅ S u = 1. (11)
Accordingly, the first-order index for a single variable X i is defined as:

S i = D i D = Var [M i (X i )] Var [M(X)] (12)
and describes the influence of X i considered separately, also called main effect. It is noted that the first-order Sobol' indices are equivalent to the first-order indices obtained by the Fourier amplitude sensitivity test FAST method [START_REF] Cukier | Nonlinear sensitivity analysis of multiparameter model systems[END_REF][START_REF] Saltelli | A quantitative, model independent method for global sensitivity analysis of model output[END_REF]. The total Sobol' index, denoted by S T i , represents the total effect of X i , accounting for its main effect and all interactions with other input variables. It is derived from the sum of all partial indices S u that involve the variable X i :

S T i = u⊂{1, ... ,M } i∈u S u . (13)
First-order and total Sobol' indices are also defined for groups of variables [START_REF] Sobol | Sensitivity estimates for nonlinear mathematical models[END_REF][START_REF] Kucherenko | Estimation of global sensitivity indices for models with dependent variables[END_REF][START_REF] Owen | Variance components and generalized Sobol' indices[END_REF]. We herein respectively denote by S u and S T u the firstorder and total Sobl' indices of the subvector X u of X. The first-order index S u describes the influence of the elements of X u only, including their main effects and all interactions with other components of X u :

S u = v⊂u v =∅ S v . ( 14
)
The total index S T u describes the total effect of the components of X u , including their main effects and all interactions with other variables in X:

S T u = v⊂{1, ... ,M } v∩u =∅ S v . ( 15
)
The total index S T u can equivalently be obtained as:

S T u = 1 -S u , (16)
where u denotes the complementary set of u so that {X u , X u } = X. For u = {i}, Eq. ( 14)

and Eq. ( 15) respectively reduce to Eq. ( 12) and Eq. (13).

To elaborate on the above definitions of the Sobol' indices, let us consider a model with three input variables i.e. X = {X 1 , X 2 , X 3 }. Setting u = {1, 2} for instance, the first-order and total Sobol' indices of X u = {X 1 , X 2 } are respectively given by:

S 1,2 = S 1 + S 2 + S 1,2 (17)
and

S T 1,2 = S 1 + S 2 + S 1,2 + S 1,3 + S 2,3 + S 1,2,3 = 1 -S 3 . (18)
Eq. ( 17) follows from the specialization of Eq. ( 14) to the considered case. In Eq. ( 18), the first part is consistent with the definition in Eq. ( 15), whereas the second part follows from Eq. ( 16).

Note that the first-order index of the subvector X u can alternatively be obtained as:

S u = Var [E [M(X)|X u ]] Var [M(X)] = Var M u (X u ) Var [M(X)] . ( 19
)
The equivalence between the above equation and Eq. ( 14) can be verified easily by considering the recursive relationship in Eq. ( 7) and the orthogonality property in Eq. ( 4). For the total index of the subvector X u , Eq. ( 16) in conjunction with Eq. ( 19) leads to:

S T u = 1 - Var [E [M(X)|X u ]] Var [M(X)] = 1 - Var M u (X u ) Var [M(X)] . ( 20
)
The Sobol' indices can be computed by estimating the variance of the summands in Eq. ( 7)

using sampling-based formulas [START_REF] Sobol | Sensitivity estimates for nonlinear mathematical models[END_REF][START_REF] Saltelli | Making best use of model evaluations to compute sensitivity indices[END_REF]. However, this computation often requires a large number of model evaluations (say more than 10 4 per index) to obtain accurate estimates; it thus becomes cumbersome or even intractable in cases when a single evaluation of the model is time consuming. In Sections 4 and 5, we demonstrate that the Sobol' indices can be evaluated analytically in terms of the coefficients of LRA or PCE meta-models. Consequently, once a PCE or LRA representation of M(X) is available, the Sobol' indices indices can be obtained at nearly zero additional computational cost. In the case of PCE meta-models, the derivation is based on Eq. ( 14) and Eq. ( 15), following the original idea presented in [START_REF] Sudret | Global sensitivity analysis using polynomial chaos expansions[END_REF]. In the case of LRA meta-models, analytical expressions for the Sobol' indices are herein derived by relying on Eq. ( 19) and Eq. ( 20).

3 Non-intrusive meta-modeling and error estimation A meta-model M is an analytical function that mimics the behavior of M in Eq. (1). In nonintrusive approaches, which are of interest herein, the meta-model is developed using a set of realizations of the input vector E = {χ (1) , . . . , χ (N ) }, called Experimental Design (ED), and the corresponding responses of the original model Y = {M(χ (1) ), . . . , M(χ (N ) )}. Thus, nonintrusive meta-modeling treats the original model as a black box.

In the following, we describe measures of accuracy of the meta-model response Y = M(X).

To this end, we introduce the discrete

L 2 semi-norm of a function x ∈ D X -→ a(x) ∈ R,
given by:

a X = 1 n n i=1 a 2 (x i ) 1/2 , ( 21
)
where X = {x 1 , . . . , x n } ⊂ D X denotes a set of realizations of the input vector.

A good measure of the accuracy of the meta-model is the generalization error Err G , which represents the mean-square of the difference (Y -Y ). Err G can be estimated by:

Err G = M -M 2 X val , (22)
where X val = {x 1 , .

Err E = M -M 2 E , (23)
in which the subscript E indicates that the semi-norm is evaluated at the points of the ED. The relative empirical error, denoted by err E , is obtained by normalizing Err E with the empirical variance of Y = {M(χ (1) ), . . . , M(χ (N ) )}, the latter representing the set of model responses at the ED. Unfortunately, err E tends to underestimate the actual generalization error, which might be severe in cases of overfitting. Indeed, it can even decrease down to zero if the obtained meta-model interpolates the data at the ED, while it is not necessarily accurate at other points.

By using the information contained in the ED only, a fair approximation of the generalization error can be obtained with cross-validation techniques. In k-fold cross-validation, (i) the ED is randomly partitioned into k sets of approximately equal size, (ii) for i = 1, . . . , k, a metamodel is built considering all but the i-th partition and the excluded set is used to evaluate the respective generalization error, (iii) the generalization error of the meta-model built with the full ED is estimated as the average of the errors of the k meta-models obtained in (ii). Leave-one-out cross-validation corresponds to the case when k = N , i.e. one point of the ED being set apart in turn [START_REF] Allen | The prediction sum of squares as a criterion for selecting predictor variables[END_REF].

4 Sobol' sensitivity analysis using polynomial chaos expansions

Formulation and construction of polynomial chaos expansions

A meta-model of Y = M(X) in Eq. ( 1) belonging to the class of PCE has the form [START_REF] Xiu | The Wiener-Askey polynomial chaos for stochastic differential equations[END_REF]:

Y PCE = M PCE (X) = α∈A y α Ψ α (X), ( 24
)
where A denotes a set of multi-indices α = (α 1 , . . . , α M ), {Ψ α , α ∈ A} is a set of multivariate polynomials that are orthonormal with respect to f X and {y α , α ∈ A} is the set of polynomial coefficients. The orthonormality condition reads:

E [Ψ α (X)Ψ β (X)] = δ αβ , (25)
where δ αβ is the Kronecker delta, equal to one if α = β and zero otherwise.

The multivariate polynomials that comprise the PCE basis are obtained by tensorization of appropriate univariate polynomials:

Ψ α (X) = M i=1 P (i) αi (X i ). ( 26
)
In the above equation, P

αi (X i ) is a polynomial of degree α i in the i-th input variable belonging to a family of polynomials that are orthonormal with respect to f Xi , thus satisfying:

E P (i) j (X i )P (i) k (X i ) = δ jk . ( 27
)
For standard distributions, the associated family of orthonormal polynomials is well-known, e.g.

a standard normal variable is associated with the family of Hermite polynomials, whereas a uniform variable over [-1, 1] is associated with the family of Legendre polynomials. A general case can be treated through an isoprobabilistic transformation of the input random vector X to a basic random vector, e.g. a vector with independent standard normal components or independent uniform components over [-1, 1].

The set of multi-indices A in Eq. ( 24) is determined by an appropriate truncation scheme.

A typical scheme consists in selecting multivariate polynomials up to a total degree p t , i.e.

A = {α ∈ N M | α 1 ≤ p t }, with α 1 = M i=1 α i .
The corresponding number of terms in the truncated series is:

cardA = M + p t p t = (M + p t )! M !p t ! , (28)
which grows exponentially with M , thus giving rise to the so-called curse of dimensionality.

Because the terms that include interactions between input variables are usually less significant, Blatman and Sudret (2010a) proposed to use a hyperbolic truncation scheme, where the set of retained multi-indices is defined as A = {α ∈ N M | α q ≤ p t }, with the q-norm given by:

α q = M i=1 α i q 1/q , 0 < q < 1. ( 29
)
According to Eq. ( 29), the lower the value of q, the smaller is the number of interaction terms in the PCE basis. The case q = 1 corresponds to the common truncation scheme using a maximal total degree p t , whereas the case q = 0 corresponds to an additive model.

Once the basis has been specified, the set of coefficients y = {y α , α ∈ A} may be computed by minimizing the mean-square error of the approximation over the ED. More efficient solution schemes can be devised by considering the respective regularized problem:

y = arg min υ∈R cardA M - α∈A υ α Ψ α 2 E + λP(υ), (30)
in which P(υ) is an appropriate regularization functional of υ = {υ 1 , . . . , υ cardA }. If P(υ) is selected as the L 1 norm of υ, i.e. P(υ) = cardA i=1

|υ i |, insignificant terms may be disregarded from the set of predictors, leading to sparse solutions. [START_REF] Blatman | Adaptive sparse polynomial chaos expansion based on least angle regression[END_REF] proposed to use the hybrid Least Angle Regression (LAR) method for building sparse PCE. This method employs the LAR algorithm [START_REF] Efron | Least angle regression[END_REF] to select the best set of predictors and subsequently, estimates the coefficients using Ordinary Least Squares (OLS). Other techniques to derive sparse expansions can be found in e.g. [START_REF] Doostan | A non-adapted sparse approximation of pdes with stochastic inputs[END_REF]; [START_REF] Yang | Reweighted 1 minimization method for stochastic elliptic differential equations[END_REF].

A good measure of the PCE accuracy is the leave-one out error. As mentioned in Section 3, this corresponds to the cross-validation error for the case k = N . Using algebraic manipulations, the leave-one-out error can be computed based on a single PCE that is built with the full

ED. Let h(χ (i) ) denote the i-th diagonal term of matrix Ψ(Ψ T Ψ) -1 Ψ T , with Ψ = {Ψ ij = Ψ j (χ (i) ), i = 1, . . . , N ; j = 1, . . . , cardA}.
The leave-one-out error can then be computed as [START_REF] Blatman | Adaptive sparse polynomial chaos expansions for uncertainty propagation and sensitivity analysis[END_REF]:

Err LOO = M -M PCE 1 -h 2 E . (31)
The relative leave-one-out error, denoted by err LOO , is obtained by normalizing Err LOO with the empirical variance of Y = {M(χ (1) ), . . . , M(χ (N ) )}, the latter representing the set of model responses at the ED. Because err LOO may be too optimistic, the following corrected estimate may be used instead [START_REF] Chapelle | Model selection for small sample regression[END_REF]:

err * LOO = err LOO 1 - cardA N -1 1 + tr((Ψ T Ψ) -1 ) . (32)
This corrected leave-one-out error is a good compromise between fair error estimation and affordable computational cost.

Computation of Sobol' sensitivity indices using polynomial chaos expansions

Let us consider the PCE representation Y PCE = M PCE (X) in Eq. ( 24). It is straightforward to obtain the Sobol' decomposition of M PCE in an analytical form by observing that the summands in Eq. ( 2) can be written as:

M PCE u (X u ) = α∈Au y α Ψ α (X), A u = {α ∈ A | α k = 0 ⇔ k ∈ u} (33)
and

M PCE 0 = y 0 . (34)
Because of the orthonormality of the PCE basis, the conditions in Eq. ( 3) and Eq. ( 4) are satisfied, ensuring the uniqueness of the decomposition. Note that Eq. ( 3) in conjunction with

Eq. ( 34) lead to:

E M PCE (X) = y 0 . (35)
Furthermore, the orthonormality of the PCE basis allows to obtain the total variance D PCE and any partial variance D PCE u of M PCE (X) by a mere combination of the squares of the PCE coefficients. These quantities are respectively given by:

D PCE = Var M PCE (X) = α∈A\{0} y α 2 (36)
and

D PCE u = Var M PCE u (X u ) = α∈Au y α 2 . (37)
By utilizing Eq. ( 36) and Eq. ( 37), the first-order and total Sobol' indices of a subvector X u of X are respectively obtained as:

S PCE u = α∈ Au y α 2 /D PCE , A u = {α ∈ A\{0} | α i = 0 ∀ i / ∈ u} (38)
and

S T, PCE u = α∈ A T u y α 2 /D PCE , A T u = {α ∈ A | ∃ i ∈ u : α i > 0}. ( 39
)
To elaborate the above equations, let us consider a model with three input variables i.e. X = {X 1 , X 2 , X 3 } and focus on the case u = {1, 2}. Then, the set A u includes the multi-indices of the form (α 1 , 0, 0), (0, α 2 , 0) or (α 1 , α 2 , 0), where α i , i = 1, 2, 3, denotes a non-zero element.

The set A T u is a superset of A u , additionally including the multi-indices of the form (α 1 , 0, α 3 ), (0, α 2 , α 3 ) or (α 1 , α 2 , α 3 ). Note that the set A u in Eq. ( 33) includes only the multi-indices of the form (α 1 , α 2 , 0). By specializing Eq. ( 38) and Eq. ( 39) to the case u = {i}, the first-order and total Sobol' indices of a single variable X i are respectively obtained as:

S PCE i = α∈Ai y α 2 /D PCE , A i = {α ∈ A | α i > 0, α j =i = 0} (40) and S T, PCE i = α∈A T i y α 2 /D PCE , A T i = {α ∈ A | α i > 0}. ( 41
)
5 Sobol' sensitivity analysis using low-rank tensor approximations

Formulation and construction of low-rank tensor approximations

Let us consider again the mapping in Eq. ( 1). A rank-one function of the input vector X has the form:

w(X) = M i=1 v (i) (X i ), ( 42
)
where v (i) denotes a univariate function of X i . A representation of M as a sum of a finite number of rank-one functions constitutes a canonical decomposition with rank equal to the number of rank-one components. Naturally, of interest are representations where the exact response is approximated with sufficient accuracy by using a relatively small number of terms, leading to the name low-rank approximations.

A meta-model of Y = M(X) belonging to the class of LRA can be written as:

Y LRA = M LRA (X) = R l=1 b l M i=1 v (i) l (X i ) , ( 43
)
where R is the rank of the approximation, v

l is a univariate function of X i in the l-th rank-one component and {b l , l = 1, . . . , R} are scalars that can be viewed as normalizing constants.

In order to obtain a representation in terms of polynomial functions, we expand v (i) l onto a polynomial basis that is orthonormal with respect to f Xi , i.e. satisfies the condition in Eq. ( 27). This leads to the expression:

v (i) l (X i ) = pi k=0 z (i) k,l P (i) k (X i ), ( 44
)
where

P (i)
k is the k-th degree univariate polynomial in the i-th input variable, p i is the maximum degree of

P (i) k and z (i) k,l is the coefficient of P (i)
k in the l-th rank-one term. Appropriate families of univariate polynomials satisfying the orthonormality condition can be determined as discussed in Section 4.1. By substituting Eq. ( 44) into Eq. ( 43), we obtain:

Y LRA = M LRA (X) = R l=1 b l M i=1 pi k=0 z (i) k,l P (i) k (X i ) . ( 45
)
Disregarding the redundant parameterization arising from the normalizing constants, the number of unknowns in the above equation is R • M i=1 (p i + 1). Note that this number grows only linearly with the input dimension M . A representation of Y = M(X) in the form of LRA drastically reduces the number of unknowns compared to PCE. In order to emphasize this, we consider PCE with the candidate basis determined by the truncation scheme A = {α ∈ N M | α i ≤ p i , i = 1, . . . , M }, so that the expansion relies on the same polynomial functions as those used in Eq. ( 45). For the case when p i = p, i = 1, . . . , M , the resulting number of

unknowns is (p + 1) M in PCE versus R • M • (p + 1) in LRA.
Considering a typical engineering problem with M = 10 and low-degree polynomials with p = 3, the aforementioned formulas yield 1, 048, 576 unknowns in PCE versus 40R unknowns in LRA; for a low rank, say R ≤ 10, the number of unknowns in LRA does not exceed a mere 400. In other words, the LRA metamodel constitutes a compressed representation of a PCE meta-model, with the basis of the latter defined by constraining the maximum polynomial degree in each dimension. As will be seen in the following, the compressed formulation leads to a fundamentally different construction algorithm.

Non-intrusive algorithms proposed in the literature for building LRA in the form of Eq. ( 45) rely on Alternated Least-Squares (ALS) minimization for the computation of the polynomial coefficients (Chevreuil et al., 2013,?;[START_REF] Doostan | Non-intrusive low-rank separated approximation of high-dimensional stochastic models[END_REF][START_REF] Rai | Sparse Low Rank Approximation of Multivariate Functions -Applications in uncertainty quantification[END_REF]. The ALS approach consists in sequentially solving a least-squares minimization problem along each dimension 1, . . . , M separately, while "freezing" the coefficients in all remaining dimensions. We herein employ the greedy algorithm proposed in Chevreuil et al. (2013) and further investigated in Konakli and Sudret (2015a), which involves progressive increase of the rank by successively adding rank-one components.

Let Y LRA r = M LRA r (X) denote the rank-r approximation of Y = M(X): Y LRA r = M LRA r (X) = r l=1 b l w l (X), (46)
with:

w l (X) = M i=1 v (i) l (X i ) = M i=1 pi k=0 z (i) k,l P (i) k (X i ) . ( 47
)
The employed algorithm comprises a sequence of pairs of a correction step and an updating step, so that the r-th correction step yields the rank-one term w r and the r-th updating step yields the set of coefficients {b 1 , . . . , b r }. These steps are detailed in the sequel.

Correction step: Let R r (X) denote the residual after the completion of the r-th iteration:

R r (X) = M(X) -M LRA r (X). ( 48
)
The sequence is initiated by setting M LRA 0 (X) = 0 leading to R 0 (X) = M(X). Based on the available experimental design, E, the rank-one tensor w r is obtained as the solution of:

w r = arg min ω∈W R r-1 -ω 2 E , (49)
where W represents the space of rank-one tensors. By employing an ALS scheme, the minimization problem in Eq. ( 49) is substituted by a series of smaller minimization problems, each involving the coefficients along one dimension only:

z (j) r = arg min ζ∈R p j R r-1 -   i =j v (i) r   pj k=0 ζ k P (j) k 2 E , j = 1, . . . , M. (50)
The correction step is initiated by assigning arbitrary values to v

(i)
r , i = 1, . . . , M , e.g. unity values, and may involve several iterations over the set of dimensions {1, . . . , M }. The stopping criterion proposed in Konakli and Sudret (2015a) combines the number of iterations, denoted by I r , with the decrease in the relative empirical error in two successive iterations, denoted by ∆ err r , where the relative empirical error is given by: ii. For i = 1, . . . , M , repeat steps A-B (correction step):

err r = R r-1 -w r 2 E Var [Y] . (51)
A. Determine z

(i) r = {z (i) 1,r . . . z (i)
pi,r } (Eq. ( 50)).

B. Use the current values of z

(i) r = {z (i) 1,r . . . z (i) pi,r } to update v (i)
r (Eq. ( 44)).

iii. Use the current functional forms of v

(i)
r , i = 1, . . . , M , to update w r (Eq. ( 47)).

iv. Compute err r (Eq. ( 51)) and update ∆ err r .

(c) Determine b = {b 1 . . . b r } (Eq. ( 52), updating step).

(d) Evaluate M LRA r at the points of the ED (Eq. ( 46)).

(e) Evaluate R r at the points of the ED (Eq. ( 48)).

Note that the above algorithm relies on the solution of several small-size minimization problems: each iteration in a correction step involves M minimization problems of size {p i + 1, i = 1, . . . , M } (usually p i < 20), whereas the r-th updating step involves one minimization problem of size r (recall that small ranks are of interest in LRA). Because of the small number of unknowns, Eq. ( 50) and Eq. ( 52) can be efficiently solved with OLS, as shown in Konakli and Sudret (2015a). An alternative approach employed in Chevreuil et al. ( 2013) is to substitute these equations with respective regularized problems.

In a typical application, the optimal rank R is not known a priori. As noted earlier, the progressive construction of LRA results in a set of decompositions of increasing rank. Thus, one may set r = 1, . . . , r max in Step 2 of Algorithm 1, where r max is a maximum allowable candidate rank, and at the end, select the optimal-rank approximation using error-based criteria. In the present work, we select the optimal rank R ∈ {1, . . . , r max } based on the 3-fold cross-validation error, as proposed by Chevreuil et al. (2013) (see Section 3 for details on k-fold cross-validation).

Konakli and Sudret (2015a) investigated the accuracy of the approach in a number of applications and showed that it leads to optimal or nearly optimal rank selection in terms of the relative generalization error estimated with a large validation set.

Computation of Sobol' sensitivity indices using low-rank tensor approximations

In this section, we consider the LRA meta-model Y LRA = M LRA (X) in Eq. ( 45) and derive analytical expressions for the Sobol' indices in terms of the polynomial coefficients z

(i)
k,l and the normalizing constants b l . To this end, we rely on the definitions of the first-order and total

Sobol' indices for a subset X u of X given in Eq. ( 19) and Eq. ( 20), respectively. In the specific case u = {i}, these equations yield the first-order and total Sobol' indices of a single variable X i .

Employing the definition of variance and the property:

E M u (X u ) = E [E [M(X)|X u ]] = E [M(X)] , (53)
Eq. ( 19) is equivalently written as:

S u = E M u (X u ) 2 -(E [M(X)]) 2 Var [M(X)] . ( 54
)
We next derive analytical expressions for each of the quantities in the right-hand side of Eq. ( 54) considering M LRA in place of M. In these derivations, we make use of the equalities:

E v (i) l (X i ) = z (i) 0,l (55)
and

E v (i) l (X i )v (i) l (X i ) = pi k=0 z (i) k,l z (i) k,l , (56)
which follow from Eq. ( 44) in conjunction with the orthonormality condition in Eq. ( 27).

Employing the expression for the LRA response in Eq. ( 43), its mean value can be computed

as:

E M LRA (X) = R l=1 b l M i=1 E v (i) l (X i ) = R l=1 b l M i=1 z (i) 0,l . (57)
In the first part of the above equation, we have utilized the property

E [A • B] = E [A]•E [B], which
holds under the condition that A and B are independent; note that because the components of X are assumed independent, the quantities {v

(i)
l (X i ), i = 1, . . . , M } are also independent. The second part of Eq. ( 57) follows directly from Eq. ( 55).

The mean-square of the LRA response in Eq. ( 43) can be computed as:

E M LRA (X) 2 = R l=1 R l =1 b l b l M i=1 E v (i) l (X i )v (i) l (X i ) = R l=1 R l =1 b l b l M i=1 pi k=0 z (i) k,l z (i) k,l . (58)
The first part of the above equation relies on the independence of the quantities {v

(i) l (X i )v (i) l (X i ), i = 1, . . . , M }, which allows the use of the property E [A • B] = E [A] • E [B]
, whereas the second part follows directly from Eq. ( 56).

Substituting Eq. ( 57) and Eq. ( 58) into the definition of the variance, we have:

Var M LRA (X) =E M LRA (X) 2 -E M LRA (X) 2 = R l=1 R l =1 b l b l M i=1 pi k=0 z (i) k,l z (i) k,l - M i=1 z (i) 0,l z (i) 0,l . (59)
In the following, we derive an analytical expression for E M LRA u (X u )

2

, beginning with

the case u = {i}. Noting that v (i)
l (X i ) represents a constant in the quantity (M LRA (X)|X i ), substitution of Eq. ( 43) into the definition of M LRA i (X i ) (Eq. ( 6)) yields:

M LRA i (X i ) = R l=1 b l v (i) l (X i )   M j =i E v (j) l (X j )   = R l=1 b l   M j =i z (j) 0,l   v (i) l (X i ). ( 60
)
Using the above expression, we have:

E M LRA i (X i ) 2 = R l=1 R l =1 b l b l   M j =i z (j) 0,l z (j) 0,l   E v (i) l (X i )v (i) l (X i ) = R l=1 R l =1 b l b l   M j =i z (j) 0,l z (j) 0,l   pi k=0 z (i) k,l z (i) k,l . (61)
It is straightforward to extend Eq. ( 60) to M LRA u (X u ) considering that the components of X u represent constants in the quantity (M LRA (X)|X u ). We thereby obtain:

M LRA u (X u ) = R l=1 b l i∈u v (i) l (X i )   j / ∈u E v (j) l (X j )   = R l=1 b l   j / ∈u z (j) 0,l   i∈u v (i) l (X i ) (62)
and finally:

E M LRA u (X u ) 2 = R l=1 R l =1 b l b l   j / ∈u z (j) 0,l z (j) 0,l   i∈u E v (i) l (X i )v (i) l (X i ) = R l=1 R l =1 b l b l   j / ∈u z (j) 0,l z (j) 0,l   i∈u pi k=0 z (i) k,l z (i) k,l . (63)
Eq. ( 57), Eq. ( 59) and Eq. ( 63) provide all the elements required to compute the first-order index for any subvector X u according to Eq. ( 54). As explained above, Eq. ( 63) simplifies to Eq. ( 61) for the case of a single variable X i .

Respective analytical expressions for the total Sobol' indices can be derived using Eq. ( 20), which is equivalently written as:

S T u = 1 - E M u (X u ) 2 -(E [M(X)]) 2 Var [M(X)] . ( 64
)
By considering that u is the complementary set of u with respect to {1, . . . , M }, we can obtain an expression for E M LRA u (X u )

2 by simply interchanging the indices i and j in Eq. ( 63):

E M LRA u (X u ) 2 = R l=1 R l =1 b l b l i∈u z (i) 0,l z (i) 0,l   j / ∈u pj k=0 z (j) k,l z (j) k,l   . ( 65
)
For the special case u = {i}, Eq. ( 65) reduces to:

E M LRA i (X i ) 2 = R l=1 R l =1 b l b l z (i) 0,l z (i) 0,l   j =i pj k=0 z (j) k,l z (j) k,l   . (66)
Note that by appropriate combinations of first-order indices for single variables and groups of variables, we can obtain any higher-order index representing the effect from the interaction between a set of variables. For instance, the second-order index representing the effect from the interaction between X i and X j can be obtained as S i,j = S i,j -S i -S j . A general expression can be derived by dividing each part of Eq. ( 7) by Var [M(X)], leading to:

S u = S u - v⊂u v =u S v . ( 67
)

Example applications

In this section, we perform global sensitivity analysis for the responses of four models with different characteristics and dimensionality. The first two models correspond to analytical functions, The meta-models are built using two types of EDs, based on Sobol pseudo-random sequences [START_REF] Niederreiter | Random number generation and quasi-Monte Carlo methods[END_REF] and the so-called maximin Latin Hypercube Sampling (LHS). Each ED of the latter type represents the best among 5 random LHS-based designs, where the selection criterion is the maximum of the minimum distance between the points. The LRA meta-models are built by implementing Algorithm 1 in Section 5.1. A common polynomial degree p 1 = . . . = p M = p is considered with its optimal value selected as the one leading to the minimum 3-fold crossvalidation error, denoted by err CV3 (see Konakli and Sudret (2015a) for an investigation of the accuracy of the approach). The involved minimization problems are solved using the OLS method. In building the PCE meta-models, a candidate basis is first determined by employing a hyperbolic truncation scheme and then, a sparse expansion is obtained by evaluating the PCE coefficients with the hybrid LAR method. The optimal combination of the maximum total polynomial degree p t and the parameter q controlling the truncation, with q ∈ {0.25, 0.50, 0.75, 1.0}, is selected as the one leading to the minimum corrected leave-one-out error err * LOO (see Section 4.1 for details).

Analytical functions

Sobol function

In the first example, we consider the Sobol function:

Y = M i=1 |4X i -2| + c i 1 + c i , (68)
where X = {X 1 , . . . , X M } are independent random variables uniformly distributed over [0, 1] and c = {c 1 , . . . , c M } are non-negative constants. We examine the case with M = 20 and the constants given by [START_REF] Kersaudy | A new surrogate modeling technique combining kriging and polynomial chaos expansions-application to uncertainty analysis in computational dosimetry[END_REF]: 2,5,10,20,50,100,500,500,500,500,500,500,500,500,500,500,500, 500, 500} T . (69

c = {1,
)
This function is a commonly used benchmark for global sensitivity analysis with well-known analytical solutions for the Sobol' sensitivity indices. In particular, the index S i1, ... ,is is obtained as [START_REF] Sudret | Global sensitivity analysis using polynomial chaos expansions[END_REF]:

S i1, ... ,is = 1 D is i=i1 D i , (70)
with the partial and total variances respectively given by:

D i = 1 3(1 + c i ) 2 (71)
and

D = M i=1 (D i + 1) -1. (72)
Then, the first-order or total indices for any single variable or groups of variables can be computed using Eq. ( 12)-( 15). Because the input variables follow a uniform distribution, we build In the sequel, we compare LRA-and PCE-based sensitivity indices, obtained by postprocessing the coefficients of the meta-models according to Sections 4.2 and 5.2, with the corresponding analytical solutions for the actual model. To further assess the effect of the ED size, we examine the convergence of the LRA-and PCE-based indices of the two most important variables, i.e. X 1 and X 2 , while N varies from 100 to 2,000; as in the above investigations, the considered EDs are obtained with Sobol sequences.

X 1 X 2 X 3 X 4 X 5 S i N=200 Analytical LRA PCE 0.1 0.2 0.3 0.4 0.5 0.6 0.7 X 1 X 2 X 3 X 4 X 5 S i N=500 Analytical LRA PCE
.7 X 1 X 2 X 3 X 4 X 5 S T i N=200 Analytical LRA PCE 0.1 0.2 0.3 0.4 0.5 0.6 0.7 X 1 X 2 X 3 X 4 X 5 S T i N=500 Analytical LRA PCE
Figure 3 shows the differences between the first-order indices obtained from the meta-model We conclude this example by investigating the accuracy of the LRA-and PCE-based Sobol'

indices when the meta-models are built with maximin LHS designs. Figure 5 shows boxplots of the differences between the total Sobol' indices obtained with the meta-models and their exact values, considering 20 EDs of size N = 200 and N = 500, i.e. for the same ED sizes examined in Figure 2. Note that in all cases, the absolute median errors do not exceed 0.01. As expected, the error spread decreases with increasing ED size.

-0.04 -0.02 0 0.02 0.04

X 1 X 2 X 3 X 4 X 5 S T, M i -S T i N=200 LRA PCE -0.04 -0.02 0 0.02 0.04 X 1 X 2 X 3 X 4 X 5 S T, M i -S T i N=500 LRA PCE

Beam deflection

In this example, we consider a simply supported beam with a uniform rectangular cross-section subjected to a concentrated load at the mid-span. The response quantity of interest is the midspan deflection, which is obtained through basic structural mechanics as:

U = P L 3 4Ebh 3 , ( 73
)
where b and h respectively denote the width and height of the cross section, L is the length of the beam, E is the Young's modulus and P is the load. The aforementioned parameters are modeled as independent random variables with their distributions, mean and coefficient of variation (CoV) values listed in Table 3. The dimensionality of the problem is thus M = 5. LRA and PCE metamodels are built using Hermite polynomials after an isoprobabilistic transformation of the input variables into standard normal variables. In Table 4, we assess the accuracy of LRA and PCE meta-models in estimating the mean µ U and standard deviation σ U of the response by post-processing the meta-model coefficients.

To this end, we consider two EDs of size N = 30 and N = 50 obtained with Sobol sequences (see the Appendix for the parameters of the respective meta-models). Because U in Eq. ( 73) is a product of lognormal random variables, it easy to obtain the mean and standard deviation of the actual model response analytically. The relative errors ε of the LRA-and PCE-based estimates with respect to the analytical values are given in parentheses. N = 30 is sufficient to obtain highly accurate estimates of both µ U and σ U with the LRA approach; for both EDs, the PCE approach appears slightly inferior in the estimation of σ U . 6.2 Finite element models

Truss deflection

We consider the truss structure shown in Figure 9 (also studied in Sudret ( 2007)), which is subjected to six vertical loads. The mid-span deflection, denoted by u, represents the response quantity of interest and is computed with an in-house finite-element analysis code developed in Matlab environment. The random input includes M = 10 independent random variables:

the vertical loads, denoted by P 1 , . . . , P 6 , the cross-sectional area and Young's modulus of the horizontal bars, respectively denoted by A 1 and E 1 , and the cross-sectional area and Young's modulus of the vertical bars, respectively denoted by A 2 and E 2 . The distributions of the input random variables are listed in Table 6. LRA and PCE meta-models are built using Hermite polynomials after an isoprobabilistic transformation of the input variables into standard normal variables. computed with a MCS approach using n = 10 6 samples for each index. Figure 10 shows the total Sobol' indices in order of total importance. For N = 50, the indices obtained with the meta-models are overall in fair agreement with the reference ones. In particular, the dominant contributions of A 1 and E 1 are estimated with good accuracy with the LRA approach, but with slightly lesser accuracy with the PCE approach. For N = 200, both meta-modeling approaches yield practically perfect estimates of the reference indices. Because contributions from higherorder effects are negligible in this problem, the first-order indices are nearly equal to the total indices and are thus not shown. The values of both the total and the first-order indices for the cases examined in Figure 10 are listed in the Appendix. Note that in the reference solution, total indices of less significant variables may be marginally smaller than the corresponding first-order indices. This contradiction results from a small bias in the employed Monte-Carlo estimator, which is addressed in [START_REF] Owen | Variance components and generalized Sobol' indices[END_REF].

We next examine the errors (differences) of the LRA-and PCE-based indices with respect to influential variables, i.e. A 1 and E 1 , Figure 13 shows similar boxplots considering EDs of size varying from N = 30 to N = 500. For N ≤ 100, the LRA-based indices exhibit a smaller spread, while for larger N , the spread of the estimates is nearly zero for both types of meta-models.

-0.08 -0.04 0 0.04 0.08

A 1 E 1 P 3 P 4 P 2 P 5 E 2 A 2 P 6 P 1 S T, M i -S T, Ref i N=50 LRA PCE

Heat conduction with spatially varying diffusion coefficient

We consider stationary heat conduction over the two-dimensional domain D = (-0.5, 0.5) m × (-0.5, 0.5) m shown in Figure 14, with the temperature field T (z), z ∈ D, described by the partial differential equation:

-∇ • (κ(z)∇T (z)) = QI A (z) (74)
and boundary conditions: T = 0 on the top boundary; ∇T • n = 0 on the left, right and bottom boundaries, where n is the vector normal to the boundary. In Eq. ( 74), Q = 500 W/m 3 , A = (0.2, 0.3) m × (0.2, 0.3) m is a square domain within D (see Figure 14) and I A is the indicator function equal to 1 if z ∈ A and 0 otherwise. The diffusion coefficient, κ(z), is a lognormal random field defined as:

κ(z) = exp[a κ + b κ g(z)], (75)
where g(z) denotes a standard Gaussian random field with autocorrelation function:

ρ(z, z ) = exp (-z -z 2 / 2 ). ( 76
)
In Eq. ( 75), the parameters a κ and b κ are such that the mean and standard deviation of κ are µ κ = 1 W/C • m and σ κ = 0.3 W/C • m, respectively, while in Eq. ( 76), = 0.2 m.

The response quantity of interest is the average temperature T in the square domain B = (-0.3, -0.2) m × (-0.3, -0.2) m (see Figure 14):

T = 1 |B| z∈B T (z)dz. ( 77
)
Using the Expansion Optimal Linear Estimation (EOLE) method [START_REF] Li | Optimal discretization of random fields[END_REF], the random field g(z) is approximated by:

g(z) = M i=1 ξ i √ l i φ T i C zζ . (78)
In the above equation, {ξ 1 , . . . , ξ M } are independent standard normal variables; C zζ is a vector with elements

C (k) zζ = ρ(z, ζ k ), k = 1, . . . , n,
where {ζ 1 , . . . , ζ n } are the points of an appropriately defined grid in D; (l i , φ i ) are the eigenvalues and eigenvectors of the correlation matrix [START_REF] Sudret | Stochastic finite element methods and reliability: a state-of-the-art report[END_REF] that for a square-exponential autocorrelation function, the size of the element in the EOLE grid must be 1/2 -1/3 of . Accordingly, in the present numerical application, we use a square grid with element size 0.01 m, thus comprising n = 121 points. The number of terms in the EOLE series is determined according to the rule:

C ζζ with elements C (k,l) ζζ = ρ(ζ k , ζ l ), k, l = 1, . . . , n. It is recommended in
M i=1 l i / n i=1 l i ≥ 0.99, (79)
which herein leads to M = 53. The shapes of the first 20 basis functions {φ T i C zζ (z), i = 1, . . . , 20}, are shown in Figure 15. Because the random input in this problem comprises independent standard normal variables, LRA and PCE meta-models are built using Hermite polynomials.

For a given realization of {ξ 1 , . . . , ξ M }, the model response considered "exact" in the metamodeling application is obtained with an in-house finite-element analysis code developed in Matlab environment. The left graph of Figure 14 depicts the discretization of the domain in 16, 000 triangular T3 elements, obtained with software Gmsh [START_REF] Geuzaine | Gmsh: A 3-d finite element mesh generator with built-in pre-and post-processing facilities[END_REF].

The temperature field T (z) for two example realizations of the conductivity random field is shown in In Figure 16.

Because a single run of the considered "exact" model is computationally expensive (requires approximately 16 sec with an Intel(R) Xeon(R) CPU E3-1225 v3 processor), evaluation of the Sobol' indices using MCS is impractical in this problem. Therefore, the following variance-based sensitivity analysis relies solely on the meta-modeling approach. Figure 17 shows the ten highest total Sobol' indices, obtained by post-processing the coefficients of LRA and PCE meta-models built with an ED of size N = 2, 000. The shown LRA-and PCE-based indices are practically identical and are thus, considered to represent the reference solution. (In the preceding applications, we showed that the indices based on LRA and PCE meta-models converge to the reference solution obtained with a MCS approach.) For the considered ED, a similar agreement is observed between the LRA-and PCE-based first-order indices. Because contributions from interaction effects are insignificant, the differences between the total and the first-order indices are negligible; the latter are thus not shown.

A B D ∇T • n = 0 ∇T • n = 0 ∇T • n = 0 T = 0 ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡
In the following, we examine the convergence of LRA-and PCE-based indices to their reference values with increasing ED size. For the four most significant variables, i.e. ξ 2 , ξ 6 , ξ 1 and ξ 7 , Figure 18 shows the differences between the total indices and their reference values, while N varies between 100 and 2,000. With the LRA approach, the indices of ξ significant variables according to the absolute values of ρ S i , computed with the validation set. The ranking of the input variables in terms of the Spearman correlation based on the actual model is identical to their ranking in terms of the total Sobol' indices based on the meta-models.

This observation provides a validation of the sensitivity-analysis results obtained with the metamodels under the premise that the model under consideration behaves nearly monotonically.

Conclusions

In this paper, we introduce a novel approach for computing the Sobol' sensitivity indices by post-processing the coefficients of meta-models belonging to the class of low-rank approximations (LRA). The derived analytical expressions concern the particular case when LRA metamodels are built with polynomial functions due to the simplicity and versatility these offer. The proposed method can be particularly efficient in problems with high-dimensional input because:

(i) the number of unknowns in LRA grows only linearly with the input dimension and (ii) the construction of LRA relies on a series of least-square minimization problems of small size, independent of the input dimension. The novel approach is presented in parallel to revisiting Sobol' sensitivity analysis by means of polynomial chaos expansions (PCE), which also relies on post-processing coefficients of meta-models with polynomial bases.

Following the detailed description of the construction of LRA in a non-intrusive manner and the derivation of the associated sensitivity indices, the accuracy and efficiency of the novel approach is investigated in example applications involving analytical models with a rank-one structure and finite-element models deriving from structural mechanics and heat conduction. Despite the strong promise of the proposed approach for the global sensitivity analysis of computationally costly high-dimensional models, further investigations are needed in order to establish the capacities and limitations of LRA as a general meta-modeling tool. Because evaluation of the LRA coefficients relies a series of alternated minimizations along separate dimensions, issues of convergence may be encountered in certain cases. Investigations into these issues are currently underway.

  In the above equation, Var[Y] represents the empirical variance of the set comprising the model responses at the ED. Accordingly, the algorithm exits the r-th correction step if either I r reaches a maximum allowable value, denoted by I max , or ∆ err r becomes smaller than a prescribed threshold, denoted by ∆ err min .Konakli and Sudret (2015a) showed that the choice of I max and ∆ err min may have a significant effect on the accuracy of LRA; based on a number of example investigations, they proposed to use I max = 50 and ∆ err min ≤ 10 -6 . Updating step: After the completion of a correction step, the algorithm moves to an updating step, in which the set of coefficients b = {b 1 . . . b r } is obtained as the solution of the minimization problem: each updating step, the size of vector b is increased by one. In the r-th updating step, the value of the new element b r is determined for the first time, whereas the values of the existing elements {b 1 , . . . , b r-1 } are recomputed (updated). Construction of a rank-R approximation in the form of Eq. (45) requires repeating pairs of a correction and an updating step for r = 1, . . . , R. The algorithm is summarized below. Algorithm 1: Non-intrusive construction of a rank-R approximation of M with polynomial basesr = 1, . . . , R, repeat steps (a)-(d): (a) Initialize: v (i) r (x i ) = 1, i = 1, . . . , M ; I r = 0; ∆ err r = > ∆ err min . (b) While ∆ err r > ∆ err min and I r < I max , repeat steps i-iv: i. Set I r ← I r + 1.

  namely a common benchmark function in uncertainty quantification, of dimension M = 20, and a structural-mechanics model of dimension M = 5; the aforementioned analytical models have a rank-one structure. The subsequent applications involve finite-element models, namely a structural-mechanics model of dimension M = 10 and a heat-conduction model, with thermal conductivity described by a random field, of dimension M = 53. The LRA-based Sobol' indices are compared to respective PCE-based indices and reference indices based on the actual model.The latter are computed either analytically or by using a Monte Carlo Simulation (MCS) approach with large samples according to[START_REF] Janon | Asymptotic normality and efficiency of two Sobol index estimators[END_REF]. The computations of the PCE-and MCS-based indices are performed with the software UQLab[START_REF] Marelli | UQLab: a framework for uncertainty quantification in MATLAB[END_REF][START_REF] Marelli | UQLab user manual -sensitivity analysis[END_REF].

  LRA and PCE meta-models of the function using Legendre polynomials after an isoprobabilistic transformation into standard variables distributed uniformly over[-1, 1].We first assess the comparative accuracy of LRA and PCE in estimating the mean and standard deviation of the response, respectively denoted by µ Y and σ Y . The LRA-and PCEbased estimates of µ Y and σ Y are obtained in terms of the coefficients of the meta-models, as described in Sections 4.2 and 5.2. Table1lists the analytical solutions for µ Y and σ Y based on the actual model together with their LRA-and PCE-based estimates for two EDs of size N = 200 and N = 500 obtained with Sobol sequences (see the Appendix for the parameters of the respective meta-models). The relative errors ε of the estimates are given in parentheses. We observe that the mean is estimated with high accuracy in all cases, but LRA outperform PCE in the estimation of the standard deviation.

  Figures 1 and 2 respectively show the first-order and total Sobol' indices of the variables X 1 , . . . , X 5 considering the same EDs of size N = 200 and N = 500 as above. The indices of the remaining variables X 6 , . . . , X 20 are practically zero. The small differences between the first-order and total indices indicate only minor effects from interactions between the various variables. For N = 200, the sensitivity indices based on the meta-models are fairly close to the reference ones, whereas for N = 500, the agreement is nearly perfect. The values of the indices depicted in the two figures are listed in the Appendix.

Figure 1 :

 1 Figure 1: Sobol function: Comparison of LRA-and PCE-based first-order Sobol' indices to their analytical values for experimental designs obtained with Sobol sequences.

Figure 2 :

 2 Figure 2: Sobol function: Comparison of LRA-and PCE-based total Sobol' indices to their analytical values for experimental designs obtained with Sobol sequences.

Figure 3 :Figure 4 :

 34 Figure 3: Sobol function: Errors of LRA-and PCE-based first-order Sobol' indices for experimental designs obtained with Sobol sequences.

Figure 5 :

 5 Figure 5: Sobol function: Errors of LRA-and PCE-based total Sobol' indices for experimental designs obtained with maximin LHS (20 replications).

Figure 6 :Figure 7 :Figure 8 :

 678 Figure 6: Beam deflection: Comparison of LRA-and PCE-based first-order Sobol' indices to their reference values for experimental designs obtained with Sobol sequences.

Figure 9 :

 9 Figure 9: Truss structure.

For

  two EDs of size N = 50 and N = 200 obtained with Sobol sequences, Table 7 compares the estimates of the response mean and standard deviation obtained from the coefficients of LRA and PCE meta-models (see the Appendix for the meta-model parameters) with their respective values based on the actual model. The latter, which represent reference values, are computed with a MCS sample comprising n = 10 6 points. The relative errors ε of the LRA-and PCE-based estimates with respect to the reference values are given in parentheses. N = 50 is sufficient to obtain excellent estimates of both quantities with the LRA approach; for the same ED, the PCE approach is slightly inferior in the estimation of the standard deviation. Table 7: Truss deflection: Mean and standard deviation of response for experimental designs obtained with Sobol sequences. EDs of size N = 50 and N = 200, we next compare the Sobol' indices obtained by post-processing the LRA and PCE coefficients to respective reference values. The latter are

Figure 10 :Figure 11 :

 1011 Figure 10: Truss deflection: Comparison of LRA-and PCE-based total Sobol' indices to their reference values for experimental designs obtained with Sobol sequences.

Figure 12 :Figure 13 :

 1213 Figure 12: Truss deflection: Errors of LRA-and PCE-based total Sobol' indices for experimental designs obtained with maximin LHS (20 replications).

Figure 14 :

 14 Figure 14: Heat conduction: Domain and boundary conditions (left); finite-element mesh (right).

Figure 15 :

 15 Figure 15: Heat conduction: Shapes of the first 20 basis functions in the EOLE discretization (from left-top to bottom-right row-wise).

  6 and ξ 1 practically reach their reference values at N = 200, while the indices of ξ 2 and ξ 7 reach their reference values at N = 500. The PCE-based indices demonstrate a slightly slower convergence. Similar results are obtained for the first-order indices (not shown herein).

Figure 16 :Figure 17 :Figure 18 :

 161718 Figure 16: Heat conduction: Example realizations of the temperature field.

Figure 19 :

 19 Figure 19: Heat conduction: Sensitivity measures based on the Spearman rank correlation.

  The LRA-based indices are compared to PCE-based ones and to respective reference indices, representing the exact solution or Monte-Carlo estimates based on the actual model response. It is first shown that the indices obtained with the meta-models, using only a few tens or hundreds evaluations of the original model, provide excellent approximations of the exact indices or their Monte-Carlo estimates using a large number of evaluations of the original model. The heatconduction example then demonstrates that the meta-modeling approach can be used to estimate the Sobol' indices in a case when their computation with the original model is not affordable; such situations are often encountered by practitioners in real-life problems. In the examined applications, the LRA-based indices tend to outperform the PCE-based ones by converging faster to the reference solution.

  . . , x n val } is a sufficiently large set of realizations of the input vector, called validation set. The estimate of the relative generalization error, denoted by err G , is obtained by normalizing Err G with the empirical variance of Y val = {M(x 1 ), . . . , M(x n val )}.

	However, meta-models are typically used in cases when only a limited number of model
	evaluations is affordable. It is thus desirable to obtain an estimate of Err G by relying solely on
	the ED. One such error estimate is the empirical error Err E , given by:

Table 1 :

 1 Sobol function: Mean and standard deviation of response for experimental designs obtained with Sobol sequences.

			N = 200	N = 500
		Analytical LRA (ε%)	PCE (ε%)	LRA (ε%)	PCE (ε%)
	µ Y	1.000	1.005 (0.5)	0.998 (-0.2)	1.000 (0.0)	0.995 (-0.5)
	σ Y	0.3715	0.3820 (2.8) 0.3424 (-7.8) 0.3715 (0.0) 0.3536 (-4.8)

Table 2 :

 2 Sobol function: Relative generalization errors of meta-models based on Sobol sequences.

	N 100 8.08 • 10 -2 5.46 • 10 -2 err LRA G err PCE G 200 2.57 • 10 -2 3.64 • 10 -2 500 2.32 • 10 -3 1.45 • 10 -2 1,000 4.68 • 10 -4 6.34 • 10 -3 2,000 2.03 • 10 -4 2.48 • 10 -3

Table 3 :

 3 Beam deflection: Distributions of input variables.

	Variable Distribution mean CoV
	b (m)	Lognormal	0.15	0.05
	h (m)	Lognormal	0.3	0.05
	L (m)	Lognormal	5	0.01
	E (MPa) Lognormal 30,000 0.15
	P (KN)	Lognormal	10	0.20

Table 4 :

 4 Beam deflection: Mean and standard deviation of response for experimental designs obtained with Sobol sequences.

			N = 30	N = 50
		Analytical	LRA (ε%)	PCE (ε%)	LRA (ε%)	PCE (ε%)
	µ U (mm)	2.677	2.678 (0.0)	2.675 (-0.1)	2.677 (0.0)	2.673 (-0.

Table 6 :

 6 Truss deflection: Distributions of input random variables.

	Variable	Distribution	mean	CoV
	A 1 (m)	Lognormal	0.002	0.10
	A 2 (m)	Lognormal	0.001	0.10
	E 1 , E 2 (MPa)	Lognormal 210,000 0.10
	P 1 , . . . , P 6 (KN)	Gumbel	50	0.15

Table 8 :

 8 Truss deflection: Relative generalization errors of meta-models based on Sobol sequences..

	N 30 2.67 • 10 -2 3.60 • 10 -2 err LRA G err PCE G 50 2.83 • 10 -3 1.25 • 10 -2 100 2.07 • 10 -3 2.57 • 10 -3 200 1.22 • 10 -3 1.23 • 10 -4 500 1.17 • 10 -3 9.80 • 10 -6

  Table 9 shows the relative generalization

Table 9 :

 9 Heat conduction: Relative generalization errors of meta-models based on Sobol sequences.

	N 100 2.46 • 10 -2 3.89 • 10 -2 err LRA G err PCE G 200 1.38 • 10 -2 2.88 • 10 -2 500 9.68 • 10 -3 9.15 • 10 -3 1,000 8.19 • 10 -3 2.19 • 10 -3 2,000 7.72 • 10 -3 9.58 • 10 -4

Table 10 :

 10 Sobol function: Parameters and error estimates of LRA and PCE meta-models.

	LRA err CV3 N = 200 1 6 4.89 • 10 -2 0.5 9 1.50 • 10 -2 PCE R p q p t err * LOO N = 500 1 14 7.11 • 10 -3 0.5 13 5.49 • 10 -3

Table 11 :

 11 Sobol function: First-order Sobol' indices (see Figure1).

			N = 200	N = 500
	Variable Analytical LRA	PCE	LRA	PCE
	X 1	0.6037	0.5820 0.6150 0.6032 0.6124
	X 2	0.2683	0.2798 0.2707 0.2680 0.2638
	X 3	0.0671	0.0650 0.0682 0.0665 0.0636
	X 4	0.0200	0.0172 0.0165 0.0203 0.0197
	X 5	0.0055	0.0059 0.0028 0.0055 0.0055

Table 12 :

 12 Sobol function: Total Sobol' indices (see Figure2).

			N = 200	N = 500
	Variable Analytical LRA	PCE	LRA	PCE
	X 1	0.6342	0.6144 0.6279 0.6337 0.6419
	X 2	0.2944	0.3077 0.2844 0.2941 0.2877
	X 3	0.0756	0.0737 0.0812 0.0750 0.0716
	X 4	0.0226	0.0196 0.0191 0.0230 0.0220
	X 5	0.0062	0.0068 0.0051 0.0062 0.0063

Table 13 :

 13 Beam deflection: Parameters and error estimates of LRA and PCE meta-models.

	LRA err CV3 N = 30 1 2 1.21 • 10 -4 1 2 3.13 • 10 -2 PCE R p q p t err * LOO N = 50 1 3 3.14 • 10 -7 1 2 1.56 • 10 -3

Table 14 :

 14 Beam deflection: First-order Sobol' indices (see Figure6).

			N = 30	N = 50
	Variable Reference LRA	PCE	LRA	PCE
	P	0.4383	0.4351 0.4623 0.4381 0.4403
	h	0.2499	0.2501 0.2328 0.2490 0.2485
	E	0.2479	0.2486 0.2480 0.2467 0.2476
	b	0.0282	0.0274 0.0217 0.0274 0.0276
	L	0.0105	0.0100 0.0181 0.0099 0.0099

Table 15 :

 15 Beam deflection: Total Sobol' indices (see Figure7).

			N = 30	N = 50
	Variable Reference LRA	PCE	LRA	PCE
	P	0.4589	0.4565 0.4703 0.4597 0.4595
	h	0.2661	0.2666 0.2445 0.2657 0.2635
	E	0.2633	0.2651 0.2593 0.2632 0.2625
	b	0.0299	0.0298 0.0228 0.0298 0.0296
	L	0.0107	0.0109 0.0203 0.0108 0.0110

Table 16 :

 16 Truss deflection: Parameters and error estimates of LRA and PCE meta-models. LOO N = 50 1 2 6.26 • 10 -3 0.25 2 4.25 • 10 -2 N = 200 1 3 1.52 • 10 -3 0.75 4 3.99 • 10 -4

	R p	LRA err CV3	q	PCE p t err *

Table 17 :

 17 Truss deflection: First-order Sobol' indices.

			N = 50	N = 200
	Variable Reference LRA	PCE	LRA	PCE
	A 1	0.3664	0.3698 0.3913 0.3666 0.3662
	E 1	0.3662	0.3688 0.3626 0.3661 0.3664
	P 3	0.0777	0.0755 0.0776 0.0771 0.0766
	P 4	0.0770	0.0725 0.0644 0.0747 0.0756
	P 2	0.0383	0.0357 0.0344 0.0385 0.0368
	P 5	0.0380	0.0367 0.0370 0.0370 0.0368
	E 2	0.0137	0.0136 0.0118 0.0115 0.0126
	A 2	0.0138	0.0131 0.0141 0.0120 0.0124
	P 6	0.0059	0.0037 0.0051 0.0048 0.0048
	P 1	0.0060	0.0037 0.0017 0.0047 0.0046

Table 18 :

 18 Truss deflection: Total Sobol' indices (see Figure10).

			N = 50	N = 200
	Variable Reference LRA	PCE	LRA	PCE
	A 1	0.3712	0.3743 0.3913 0.3712 0.3713
	E 1	0.3696	0.3733 0.3626 0.3707 0.3715
	P 3	0.0776	0.0769 0.0776 0.0785 0.0777
	P 4	0.0773	0.0739 0.0644 0.0761 0.0767
	P 2	0.0374	0.0364 0.0344 0.0392 0.0374
	P 5	0.0374	0.0374 0.0370 0.0377 0.0373
	E 2	0.0126	0.0138 0.0118 0.0117 0.0128
	A 2	0.0126	0.0134 0.0141 0.0122 0.0126
	P 6	0.0048	0.0038 0.0051 0.0049 0.0049
	P 1	0.0047	0.0038 0.0017 0.0048 0.0047

Table 19 :

 19 Heat conduction: Parameters and error estimates of LRA and PCE meta-models. .84 • 10 -2 0.25 3 4.08 • 10 -2 N = 2, 000 1 2 7.76 • 10 -3 0.75 3 1.22 • 10 -3

		R p	LRA err CV3	q	PCE p t err * LOO
	N = 200	1 1 1

Table 20 :

 20 Heat conduction: First-order Sobol' indices.

		N = 200	N = 2, 000
	Variable LRA	PCE	LRA	PCE
	ξ 2	0.4034 0.4227 0.4138 0.4132
	ξ 6	0.1563 0.1548 0.1597 0.1593
	ξ 1	0.1432 0.1398 0.1445 0.1433
	ξ 7	0.0885 0.0900 0.0825 0.0813
	ξ 5	0.0643 0.0667 0.0594 0.0593
	ξ 9	0.0277 0.0278 0.0260 0.0265
	ξ 4	0.0214 0.0188 0.0243 0.0243
	ξ 3	0.0154 0.0123 0.0162 0.0160
	ξ 13	0.0165 0.0179 0.0149 0.0139
	ξ 10	0.0091 0.0072 0.0095 0.0094

Table 21 :

 21 Heat conduction: Total Sobol' indices (see Figure17).

		N = 200	N = 2, 000
	Variable LRA	PCE	LRA	PCE
	ξ 2	0.4103 0.4227 0.4209 0.4206
	ξ 6	0.1601 0.1548 0.1637 0.1636
	ξ 1	0.1468 0.1398 0.1481 0.1476
	ξ 7	0.0908 0.0900 0.0847 0.0842
	ξ 5	0.0660 0.0667 0.0610 0.0619
	ξ 9	0.0285 0.0278 0.0267 0.0275
	ξ 4	0.0220 0.0188 0.0250 0.0256
	ξ 3	0.0158 0.0123 0.0166 0.0173
	ξ 13	0.0170 0.0179 0.0153 0.0147
	ξ 10	0.0093 0.0072 0.0098 0.0101