Modiﬁed Dimension Reduction-Based Polynomial Chaos Expansion for Nonstandard Uncertainty Propagation and Its Application in Reliability Analysis

: This paper presents an algorithm for efﬁcient uncertainty quantiﬁcation (UQ) in the presence of many uncertainties that follow a nonstandard distribution (e


Introduction
Uncertainty, originating from the inherent randomness of a complex system, is common in first-principle models widely used in various engineering problems. Uncertainty quantification (UQ), which quantitatively studies the impact of uncertainty on a system's performance, is important since uncertainty in model parameters (e.g., external loads, geometry, and material properties in a complex system) can degrade the accuracy of model predictions, thus affecting prediction-based control design and analysis. Many UQ approaches have been developed and applied to engineering problems to improve a system's performance [1][2][3][4][5][6][7][8]. For example, reliability assessment has been explored for the controlled structures of viscoelastic dampers [3], and the stochastic responses of bridge structures under uncertainty have been studied [4]. Further, reliability-based design optimization (RBDO), used to consider uncertainty in process design, has been developed for structural analysis, including the sequential approximate optimization to deal with multimodal random variables [6] and the adaptive surrogate-based algorithm to replace expensive-toevaluate models [7]. A hybrid reliability-based optimization method was also proposed for thermal structure design [9]. However, several challenges remain unsolved. For example, a major challenge of existing UQ tools is how to quantify the effect of a large number of uncertainties on model responses in an accurate and computationally efficient manner.
There are several popular UQ methods, including Monte Carlo (MC) simulations [10,11], Taylor expansion-based methods [12,13], and polynomial chaos expansion (PCE) [14,15]. Among these, MC has been widely used since it is easy to implement. It only requires random generation of samples from the distribution of uncertainty and repetitive execution of the first principle models with each sample. Despite the apparent simplicity, MC and its associated approaches are insufficient because they often require a large number of samples to ensure UQ accuracy, thus increasing the computational cost. To reduce the computational cost, different variants of MC were developed, which include an active learning approach to combine Kriging and MC simulation (i.e., AK-MCS) [8]. Further, the Taylor expansion-based method is another way for UQ, which is effective for less complicated models [16,17]. This includes the first-order reliability method (FORM) and second-order reliability method (SORM) [12,18] in structural reliability. However, when models involve complex nonpolynomial terms, the UQ accuracy of the Taylor expansionbased method cannot be guaranteed [12,13,18].
Recently, PCE has been used in different areas including the communities of structural reliability and control theory [4,19]. This approach was proposed by Wiener [14] to quantify uncertainty that follows a normal distribution and later extended to uncertainty with several standard distributions (e.g., uniform) by Xiu [20,21]. The PCE-based methods approximate uncertainty and its effect on model responses with a spectral expansion as a function of PCE coefficients. The PCE coefficients of uncertainty in models are given by a prior probability density function (PDF), but coefficients of model outputs are unknown. Thus, the key of PCE is to solve the PCE coefficients of model responses. Once the coefficients of model responses (e.g., performance function) are obtained, they can be used to approximate the PDF of performance function [19,20] and to evaluate failure probability in reliability analysis [18,22,23]. Since each model response can be represented explicitly with random variables of input uncertainties, it is useful to obtain the statistics of the response. Depending on how the PCE coefficients of model responses are calculated, it is categorized into intrusive and nonintrusive methods [20,24,25].
For the intrusive method, the coefficients of model outputs are obtained by projecting the governing equations into PCE, which requires calculating the inner product between a basis polynomial and the equations described with a truncated PCE in a stochastic Galerkin approach. The inner product generates a family of nested deterministic models to represent the original stochastic models. Calculating the coefficients with the inner product can be computationally prohibitive when models have complex functions (e.g., nonpolynomial functions) and many uncertainties. Nonintrusive PCE is another way, which has been used in the reliability analysis of bridge, truss, and frame structures [4,19,26,27]. Similar to MC, the convergence and accuracy of nonintrusive PCE depend on samples used for UQ. These samples are often referred to as collocation points. While practical, the convergence and accuracy of the nonintrusive PCE-based UQ can be affected by the total number of collocation points and by how the collocation points are generated [28,29]. The nonintrusive methods can be further categorized as a regression-based approach and projection-based approach, depending on how the PCE expressions of output responses are calculated. The former estimates the PCE coefficients by minimizing the mean square error of the approximation of model responses, while each PCE coefficient for the latter is defined with a multidimensional integral. It should be noted that as the number of uncertainties increase, both the size of dimension in the integral and the total number of PCE terms increase. Thus, the projection-based approach can be computationally demanding.
To address these issues, our previous works [29,30] have presented algorithms to combine the PCE with the generalized dimension reduction method (gDRM) or M-variate DRM [31]. This approach approximates high-dimensional integrals in the spectral projection (SP) [28] with a few low-dimensional ones. These low-dimensional integrals involve at most M integration variables [29] that are much smaller than the number of uncertainties. To further reduce the computational cost, the lower-dimensional integrals resulting from the gDRM were calculated with a sampling-based method (i.e., Gaussian quadrature rules) [30]. However, our previous works only deal with uncertainty that follows a standard distribution as done in the classic PCE theory. Uncertainty may have nonstandard distributions. For example, lognormal distribution describes several natural phenomena, e.g., the bubble size distributions in the gas-liquid and gas-liquid-solid systems [32]. Because there is an explicit relationship between the lognormal and normal variables, we capitalize on the relationship following the work by Ghanem in [33] to efficiently deal with the lognormal uncertainty. Several tools were developed to approximate lognormal random variables [20,24,[34][35][36][37][38][39], which include the isoprobabilistic transformation, such as the Rosenblatt or Nataf transformation [37,38], that transforms the original random variable into independent normal variables. Besides, the lognormal distribution was described in [39] by constructing orthogonal polynomials using the Stieltjes-Wigert polynomials. Compared to existing techniques, our algorithm approximates a lognormal distribution with a normal random variable based on the explicit relationship [33] and couples the PCE-based approximation with the gDRM to quantify the effect of many uncertainties on model responses in a computationally efficient way.
In summary, this work expands our previous work to couple the gDRM with the PCE to quantify the effect of nonstandard uncertainties on model outputs. For algorithm validation, lognormal and Weibull distributions are chosen as the testbed, which are approximated with standard random variables [20,33]. In addition, four examples in structural reliability analysis were chosen to show the efficiency of the algorithm and to discuss the UQ accuracy by comparing the results with recent works in the literature.
This paper is organized as follows. A brief description of the PCE-based UQ is given in Section 2. Section 3 shows procedures to approximate nonstandard uncertainty, followed by the gDRM-based PCE to deal with high-dimensional UQ problems. For algorithm verification, four numerical examples in structural reliability analysis are presented in Section 4, and a brief conclusion is given in Section 5.

Background of Polynomial Chaos Expansion (PCE)
Suppose that a stochastic process can be defined as u = K(g), where g = (g 1 , . . . , g n ) ∈ R n is a Gaussian vector involving n uncertain parameters (n ≥ 1) and K is the function to describe the relationship between a model response u and g. It is assumed each parameter in g is independent, i.e., any correlation among uncertain parameters is not considered. In this work, uncertainty in model parameters, i.e., input random variable of the models, is hereafter referred to as parametric uncertainty.
To obtain a PCE expression of u, the first step is to rewrite each parametric uncertainty g i as a function of a Gaussian variable x i as in [20]: where g i,k are the PCE coefficients to estimate the i th parametric uncertainty, which is defined by the i th standard normal distribution, i.e., x i~N (0, 1). These coefficients g i,k are often assumed to be a given a priori or can be estimated with parameter estimation techniques [40]. In (1), H k is the one-dimensional Hermite polynomial basis function, and the polynomial order p defines the number of terms to approximate g i . Uncertainty in model parameters will introduce uncertainty in model responses. Thus, once the PCE coefficients of uncertainty in (1) are available, the PCE expansion of model response u can be subsequently defined with random variables as [20]: where x = (x 1 , . . . , x i , . . . , x n ) and H j,n (x) are the j th n-dimensional orthogonal polynomials defined by basis functions {H k (x i )} for each uncertainty. Explicit expressions of H j,n (x) can be found in [15]. In (2), u j are the deterministic PCE coefficients of model response u, which are unknown and have to be calculated. Further, m is the total number of terms to approximate u(x), which can be calculated with the number of parametric uncertainty n and the polynomial order p as in [20]: Unlike g i,k , PCE coefficients of model response, u j , are unknown and have to be calculated by projecting (2) onto each polynomial basis function H j,n and by using the spectral projection (SP), which can be defined as in [20]: where the angle brackets · represent the inner product operator for the n-dimensional random space R n defined by x. For example, the inner product of a function f (x) is defined as: where W (x) is the weighted function defined by the joint PDFs of random variables x. Note that the inner product of f (x) in (5) can be further specified by the expectation operator (4) can be further written as in [20]: where γ j = E H 2 j,n is the normalization factor. It is worth mentioning that the j th normalization factor of the Hermite polynomial basis can be easily calculated as γ j = j! [20]. In addition, the multidimensional integral in (6) is solved over the entire random domain R n of x, and W (x) here is the Gaussian weighted function. By substituting the PCE coefficients in (6) into (2), uncertainty in model responses can be rapidly estimated. For example, the statistical central moments µ u and σ 2 u , i.e., mean and variance of model response u, can be analytically calculated as in [20]: where H 0,n ≡ 1 is the constant polynomial. As seen in (7) and (8), the first PCE coefficient u 0 can be used to estimate the mean, while the rest of the coefficients (u j =0 ) and the expectation of the squared basis functions (E H 2 j =0,n ) can be used to estimate the uncertainty around the mean of u [20]. Similarly, other higher-order statistical central moments, such as skewness s u and kurtosis κ u , can be rapidly calculated with the high-order PCE coefficients as in [41].
Details to calculate these statistical central moments can be found in [41][42][43]. These central moments are used to evaluate the reliability index and failure probability for structural reliability analysis in moment methods [18,22]. Thus, they are chosen in this work to discuss the performance of different UQ algorithms in terms of UQ accuracy and computational efficiency in Section 4.

Approximation of Lognormal Uncertainties
Since lognormal uncertainty has an explicit relationship with Gaussian uncertainty, it will be transformed and estimated with the Gaussian variable as in [33]. For clarity, only a parametric uncertainty, which refers to an input random variable of the model, is considered, thus the subscript i in (1) will not be used.
Suppose a normally distributed uncertainty g(x) is approximated with a random variable x; and a lognormal uncertainty l(x) can be calculated by taking the exponential operator exp[·] of g(x), which results in a relationship as in [33]: In (11), g(x) will be referred to as the normalized Gaussian distribution, and the mean and the variance of g(x) are defined as µ g and σ 2 g , respectively. As done for the Gaussian uncertainty in (1), let one assume that l(x) can be approximated with a PCE expansion as in: where l k are the PCE coefficients, describing the original lognormal uncertainty, which will be derived based on normalized Gaussian variables. As defined in (1), H k is the Hermite polynomial basis function. Since the first coefficient of PCE in (12) represents the mean value as described in (7), l 0 can be easily given as µ l = e µ g +(σ 2 g /2) [33]. For other higher-order terms, we apply the SP to calculate PCE coefficients as a function of the normalized Gaussian random variables. Following the similar procedures to calculate the coefficients of model responses in (4), the resulting PCE expansion of the lognormal uncertainty in (12) can be represented as in [20,33]: Once the PCE coefficients of a lognormal uncertainty are available, the resulting uncertainty in model response as in (2) can be approximated following the similar procedures as in Section 2.
Compared to the PDF of the Gaussian random variable that has a symmetrical shape, the required polynomial order p in (13) for approximating the original lognormal uncertainty needs to be carefully selected to capture its asymmetric property, such as the long tail. As an example, Figure 1 shows the simulation results of a lognormal uncertainty, where three different polynomial orders (i.e., p = 2, 3, and 4) were used. Note that the PDFs shown in Figure 1 are estimated with the kernel density estimation function in MATLAB. The accuracy to approximate a lognormal distribution with different polynomial orders is compared to MC with 10 7 samples.
As seen in Figure 1, the accuracy to approximate a lognormal uncertainty can be affected by p and the magnitude of uncertainty. When uncertainty is small, the approximation of a lognormal uncertainty with different p values exhibits almost identical results as in Figure 1a. Whereas, when uncertainty is large, there is a noticeable difference between the PCE-based approximation and MC. However, when p was increased (>2), as seen in Figure 1b, the difference between the PCE-based approximation and MC is insignificant. Thus, we focus on two different values of p (2 vs. 3) in this work to show the performance of the UQ method, which will be discussed in Section 4. fected by and the magnitude of uncertainty. When uncertainty is small, the approximation of a lognormal uncertainty with different values exhibits almost identical results as in Figure 1a. Whereas, when uncertainty is large, there is a noticeable difference between the PCE-based approximation and MC. However, when was increased (>2), as seen in Figure 1b, the difference between the PCE-based approximation and MC is insignificant. Thus, we focus on two different values of (2 vs. 3) in this work to show the performance of the UQ method, which will be discussed in Section 4.

Approximation of Other Nonstandard Uncertainties
Compared to approximating a lognormal uncertainty with a PCE expression that has an explicit correlation to a Gaussian random variable, the approximation of other nonstandard uncertainties can be made with its probability distribution function [20]. To estimate a nonstandard uncertainty in PCE, suppose is an uncertainty with nonstandard distributions (e.g., a Weibull distribution), which has a PCE expression as: where are PCE coefficients, and are polynomial basis functions. Considering follows a nonstandard distribution, the calculation of is difficult since the explicit correlation between Y and x is unknown [20]. In this work, a technique in [20] is used to estimate an uncertainty that follows a nonstandard distribution. The resulting PCE coefficients in (14) can be calculated as: where is the cumulative distribution function (CDF) defined as ( ) = ( ) , and is the inverse of a cumulative distribution function of , i.e., . For example, let one suppose that a Gaussian random variable is used to approximate a Weibulldistributed uncertainty . Then, ( ) is the probability distribution function of ;

Approximation of Other Nonstandard Uncertainties
Compared to approximating a lognormal uncertainty with a PCE expression that has an explicit correlation to a Gaussian random variable, the approximation of other nonstandard uncertainties can be made with its probability distribution function [20]. To estimate a nonstandard uncertainty in PCE, suppose Y is an uncertainty with nonstandard distributions (e.g., a Weibull distribution), which has a PCE expression as: where Y k are PCE coefficients, and {Ψ k } are polynomial basis functions. Considering Y follows a nonstandard distribution, the calculation of Y k is difficult since the explicit correlation between Y and x is unknown [20]. In this work, a technique in [20] is used to estimate an uncertainty that follows a nonstandard distribution. The resulting PCE coefficients in (14) can be calculated as: (15) where F x is the cumulative distribution function (CDF) defined as F Y is the inverse of a cumulative distribution function of Y, i.e., F Y . For example, let one suppose that a Gaussian random variable x is used to approximate a Weibulldistributed uncertainty Y. Then, W (x) is the probability distribution function of x; Ψ k is the Hermite polynomial basis function; F x is the CDF of the Gaussian random variable x; and F −1 Y is the inverse of the CDF of Weibull-distributed uncertainty Y. Details about the derivation of (15) and its approximation are not given for brevity and can be found in [20].

Modified gDRM-Based PCE Using Quadrature Rules
When the PCE coefficients of uncertainty are available, such as (1) for the normal and (12) for lognormal uncertainty, SP is used to calculate the PCE coefficients of model response in (2), which requires calculating multivariate integrals as in (6). As in (3), the total number of terms (m) to accurately approximate uncertainty in model response in (2) is a function of the polynomial order (p) and the total number of uncertainties (n). When p is large (e.g., >2) and when the number of uncertainties increases, the number of terms in (2) can be greatly increased. This can increase the computational cost to calculate the PCE coefficients of model outputs. To reduce the computational cost, the gDRM will be used with quadrature rules to quickly calculate the PCE coefficients as in [30]. This approach is hereafter referred to as the modified gDRM-based PCE (or mgDRM-PCE).
Suppose that a continuous and differentiable function can be defined as f (x) = u(x)H j,n (x), which is a portion of the integrand in (6). Then, the multivariate integral in (6) is rewritten as: where W (x) represents the weighted function described by the joint PDFs of x. Note that the function f (x) is defined to easily calculate a specific PCE coefficient u j in (6), and this expression can also be reused for different j = 0, 1, . . . , m − 1.
Following the definition of gDRM [31] and the representation in (16), each of the unknown PCE coefficients of the model response in (6) can be calculated as in [29]: where f M−r in the expectation operator E[·] is a function of (M − r) random variables selected from x. In this way, several lower-dimensional functions, involving up to M random variables of x, can be used to approximate the original function f, which is written as in [31]: Since the gDRM considers combinations of the random variables in x, the total number of terms in (18) can be calculated as n M − r .
Notably, the last term in (18) As noted in [29,30], the accuracy to approximate a high-dimensional integral with several lower-dimensional ones increases as M increases. For example, compared to the BiDRM, a tri-variate DRM (TriDRM) approximates a high-dimensional integral in (16) with n 3 three-dimensional integrals, n 2 two-dimensional integrals, n 1 one-dimensional integrals, and a constant term f 0 . As such, extra effort is required to calculate n 3 additional three-dimensional integrals. Since the total number of lowerdimensional integrals increases, the computational time to approximate (16) may increase.
To address this issue, we will estimate the resulting lower-dimensional integrals with quadrature rules [44] to reduce the computational cost [30,31]. It should be noted that this approach reduces the computational difficulty to evaluate high-dimensional integrals using sampling methods, which appears to be methodologically similar to the anchored analysisof-variance (ANOVA) expansion in [45]. However, the anchored ANOVA expansion is dependent on the choice of the anchor point, which impacts the accuracy of the expansion and the truncation dimension [45], while the mgDRM does not require any strategy for the decomposition procedure in (18). A detailed description on the ANOVA expansion can be found in [45,46].
To calculate the expectation of one-variate function E[ f 1 ]-an integral only involves a single random variable, e.g., x d 1 of f 1 in (18), a one-dimensional quadrature rule can be defined as [30]: Processes 2021, 9, 1856 is a set of quadrature points (x q 1 d 1 ) and their corresponding weights (α q 1 d 1 ) used to compute the one-dimensional integral resulting from the gDRM. In this work, Gauss-Hermite quadrature rules [44] are used to estimate low-dimensional integrals in (17) in the gDRM step.
To generate quadrature points to approximate multidimensional integrals in (17), a full tensor product grid can be constructed based on the one-dimensional quadrature rules. Note that these lower-dimensional integrals resulting from the gDRM only involve a few integration variables. Thus, the approximation with the full tensor product grid can be finished in real time.
To better illustrate the quadrature points-based approximation, suppose the PCE coefficients of model outputs in (17) are calculated with the TriDRM. This implies that these lower-dimensional integrals to approximate the multidimensional integrals in SP involve at most three random variables in x (or M = 3). That is, a multidimensional integral two-dimensional, and n 3 three-dimensional integrals in total, which can be computed with the quadrature rules as in [30]: where θ i (i = 1, 2, 3) is the number of quadrature points, x is the weight of each quadrature point; and ⊗ means the tensor product. In summary, the total number of quadrature points for low-dimensional integrals can be defined as [44]. It is important to note that as the dimension of a random space increases, the number of integrals resulting from gDRM increases. In this case, the total number of evaluations of these integrals with quadrature points is also increased, which may increase the computational cost. This issue will be discussed with benchmark examples in Section 4.

Summary of the UQ Algorithm
Our algorithm quantifies the effect of nonstandard uncertainty on model responses by coupling PCE with the gDRM. A flowchart summarizing the procedures of the modified gDRM-based PCE to deal with nonstandard random variables is shown in Figure 2.
As in Figure 2, nonstandard uncertainties are firstly approximated with standard random variables. Once the PCE-based surrogate model of parametric uncertainty, i.e., an input random variable of the models is available, model responses can be approximated with coupled PCE coefficients as in (2). The PCE coefficients of model responses are calculated with (4). This produces high-dimensional integrals, which are not trivial to solve, especially when the number of uncertainties is large and when the model involves nonpolynomial terms. To address this, the gDRM, such as the BiDRM and TriDRM noted in Section 3.3, is used to convert the high-dimensional integral into few low-dimensional ones. To further decrease the computational cost, quadrature points generated with Gaussianquadrature rules are used to numerically solve these low-dimensional integrals. The approximation of lower-dimensional integrals is highlighted in green boxes in Figure 2.
Once the PCE coefficients of model responses are available, the statistical moments of Processes 2021, 9, 1856 9 of 21 model outputs (e.g., mean, variance, skewness, and kurtosis) can be quickly calculated with (7)~(10). culated with (4). This produces high-dimensional integrals, which are not trivial to solve, especially when the number of uncertainties is large and when the model involves nonpolynomial terms. To address this, the gDRM, such as the BiDRM and TriDRM noted in Section 3.3, is used to convert the high-dimensional integral into few low-dimensional ones. To further decrease the computational cost, quadrature points generated with Gaussian-quadrature rules are used to numerically solve these low-dimensional integrals. The approximation of lower-dimensional integrals is highlighted in green boxes in Figure  2. Once the PCE coefficients of model responses are available, the statistical moments of model outputs (e.g., mean, variance, skewness, and kurtosis) can be quickly calculated with (7)~(10).

Benchmark Examples in Structural Reliability Analysis
Four examples in structural reliability analysis are chosen to illustrate the performance of the UQ algorithm to deal with nonstandard uncertainties. The correlation among uncertainties is not considered, since our objective is to study the applicability of the algorithm to tackle many uncertainties. When uncertainties are correlated, the number of random variables to approximate uncertainties will be decreased, thus reducing the number of dimensions. To evaluate the UQ accuracy, statistical moments (mean, standard

Benchmark Examples in Structural Reliability Analysis
Four examples in structural reliability analysis are chosen to illustrate the performance of the UQ algorithm to deal with nonstandard uncertainties. The correlation among uncertainties is not considered, since our objective is to study the applicability of the algorithm to tackle many uncertainties. When uncertainties are correlated, the number of random variables to approximate uncertainties will be decreased, thus reducing the number of dimensions. To evaluate the UQ accuracy, statistical moments (mean, standard deviation, skewness, and kurtosis) are calculated and compared with other techniques, which include MC simulations and the least angle regression-based PCE.

Example 1: Linear Performance Function
As given in Figure 3, a linear performance function is used to study the plastic collapse mechanisms of a one-bay frame, which is mathematically described as in [47,48]: 6 (23) Following [47,48], X 1~X6 in (23) are assumed to be independent and lognormally distributed; the mean and standard deviation of each variable are listed in Table 1. In this case study, all parameters are considered as parametric uncertainties, resulting in a six-dimensional random space. Since uncertainty follows a lognormal distribution, Hermite polynomials were used as the basis functions to build the PCE-based surrogate model of the response Z in (23).
The first step to build a surrogate model of Z in (23) is to formulate the PCE expressions of lognormally distributed uncertainties, X 1~X6 . In this case study, (13) was used to build a PCE model for each uncertainty with normalized Gaussian variables that are related to the lognormal uncertainty as in Section 3.1. The second coefficients of the normalized Gaussian variables (i.e., σ g for X 1~X6 ) were calculated as 0.0998 for X 1~X4 and 0.2936 for X 5 and X 6 , respectively. These were also used to derive the explicit PCE expressions of lognormal uncertainties as in (13). Once the PCE expressions of uncertainties (X 1~X6 ) were available, the mgDRM-PCE in Section 3.3 was used to compute the PCE coefficients of Z in (23). Figure 3, a linear performance function is used to study the plastic collapse mechanisms of a one-bay frame, which is mathematically described as in [47,48]:

As given in
Following [47,48], ~ in (23) are assumed to be independent and lognormally distributed; the mean and standard deviation of each variable are listed in Table 1. In this case study, all parameters are considered as parametric uncertainties, resulting in a sixdimensional random space. Since uncertainty follows a lognormal distribution, Hermite polynomials were used as the basis functions to build the PCE-based surrogate model of the response in (23). Figure 3. Illustration of the one-bay frame for a plastic collapse mechanism [48]. Reprinted from [48], Copyright (2019), with permission from Elsevier. The first step to build a surrogate model of in (23) is to formulate the PCE expressions of lognormally distributed uncertainties, ~ . In this case study, (13) was used to build a PCE model for each uncertainty with normalized Gaussian variables that are related to the lognormal uncertainty as in Section 3.1. The second coefficients of the normalized Gaussian variables (i.e., for ~ ) were calculated as 0.0998 for ~ and 0.2936 for and , respectively. These were also used to derive the explicit PCE expressions of lognormal uncertainties as in (13). Once the PCE expressions of uncertainties ( ~ ) were available, the mgDRM-PCE in Section 3.3 was used to compute the PCE coefficients of in (23).
In this case study, different polynomial orders of ~ ( = 2 vs. 3) and different numbers of random variables in the gDRM ( = 2 vs. 3) were used. The UQ accuracy with different combinations of and was compared to the results of MC by calculating the relative error ( ), which is defined as in: Figure 3. Illustration of the one-bay frame for a plastic collapse mechanism [48]. Reprinted from [48], Copyright (2019), with permission from Elsevier. In this case study, different polynomial orders of X 1~X6 (p = 2 vs. 3) and different numbers of random variables in the gDRM (M = 2 vs. 3) were used. The UQ accuracy with different combinations of p and M was compared to the results of MC by calculating the relative error ( R ), which is defined as in: where e MC is a reference of a specific statistical moment calculated with MC, and e represents the corresponding statistical moment calculated with other UQ methods. The simulation results are summarized in Table 2 and compared with existing algorithms: a recent work in [48] that integrates the BiDRM with a high-order unscented transformation [49] and the least angle regression-based PCE (LAR-based PCE) [50,51]. For MC, 10 7 samples were used for each uncertainty to ensure UQ accuracy. In Table 2, mBiDRM-PCE means that the quadrature rules were used to estimate integrals from the BiDRM, which involve at most two integration variables, while mTriDRM-PCE means that the quadrature rules were applied to the TriDRM, which has at most three integration variables. The BiDRM converted a six-dimensional integral into 15 two-and 6 one-dimensional ones, while 20 three-, 15 two-, and 6 one-dimensional integrals were used in the TriDRM. In addition, five quadrature points of each dimension were generated following the Gauss-Hermite quadrature rules (i.e., θ i = 5), for which a full-tensor product grid was built to calculate the low-dimensional integrals from the BiDRM and the TriDRM.
In addition, for the LAR-based PCE, the LAR procedure was implemented to select a set of basis polynomials to build a sparse PCE and estimate the unknown PCE coefficients of the response Z in (23) with the least square regression. It is important to note the LARbased PCE in this work does not involve the adaptive solver as in [52] because our objective is to assess the accuracy of the surrogate model derived by the gDRM-based PCE. To obtain model evaluations of Z, 1000 pairs of random samples of parametric uncertainties were used. Thus, the total number of model runs is 1000. Further, the LAR method was used to identify the sensitive basis polynomials of output Z until the maximum correlation between the residual and the basis polynomials is less than a preselected effective value, which was set to 1.0 × 10 −6 . In this case study, two different polynomial orders (p = 2 and 3) were studied. Thus, the number of PCE terms of LAR-based PCE was 13 and 19, respectively. In contrast, the number of PCE terms for our algorithm was 28 and 84, respectively. As seen in Table 2, the difference in the mean value of Z (µ Z ) calculated with the mBiDRM-PCE and mTriDRM-PCE is insignificant, and their results are almost identical to MC. For the other higher-order statistical moments, it was found that the relative error is relatively larger when p was set to 2, as compared to the results when p was set to 3. This indicates that more polynomial terms to approximate a lognormal uncertainty can improve the UQ accuracy. Besides, compared to the LAR-based PCE, our algorithm shows identical results for the four statistical moments, which demonstrates the accuracy of the gDRM-based PCE. It was also found that, when p was set to 3, the mBiDRM-PCE outperforms the results in a recent work [48], for which only the BiDRM was used. This shows the advantage of combining the PCE with the gDRM, because PCE can not only provide an analytical expression for model responses as a function of random variables, but also explicitly quantify the impact of approximation on UQ accuracy.
For comparison, the simulation results of mBiDRM-PCE and mTriDRM-PCE with different polynomial orders and MC are shown in Figure 4, where the cumulative density function (CDF) of each UQ method was approximated with the built-in kernel density estimation function in MATLAB. As seen, the CDF of Z in (23) is almost identical to the results obtained with MC, when p was 3. This shows our algorithm can accurately quantify the effect of lognormal uncertainties on model responses.

Example 2: Roof Structure
The performance function of a roof truss structure as illustrated in Figure 5 was used to study the UQ accuracy of our algorithm to deal with different nonstandard uncertainties, when the system involves nonlinear polynomial functions. The performance function is described as in [19,53]: where q and l represent a vertical load and roof span that are used to calculate the nodal load P = ql/4, and A and E are the cross-sectional area and elastic modulus, respectively. In this case study, all parameters in the performance function in (25) are uncertain with predefined PDFs. This yields a six-dimensional random space (n = 6). Details about the statistical description of uncertainties are given in Table 3.
forms the results in a recent work [48], for which only the BiDRM was used. This shows the advantage of combining the PCE with the gDRM, because PCE can not only provide an analytical expression for model responses as a function of random variables, but also explicitly quantify the impact of approximation on UQ accuracy. For comparison, the simulation results of mBiDRM-PCE and mTriDRM-PCE with different polynomial orders and MC are shown in Figure 4, where the cumulative density function (CDF) of each UQ method was approximated with the built-in kernel density estimation function in MATLAB. As seen, the CDF of in (23) is almost identical to the results obtained with MC, when was 3. This shows our algorithm can accurately quantify the effect of lognormal uncertainties on model responses.

Example 2: Roof Structure
The performance function of a roof truss structure as illustrated in Figure 5 was used to study the UQ accuracy of our algorithm to deal with different nonstandard uncertainties, when the system involves nonlinear polynomial functions. The performance function is described as in [19,53]: where and represent a vertical load and roof span that are used to calculate the nodal load = /4, and A and E are the cross-sectional area and elastic modulus, respectively.
In this case study, all parameters in the performance function in (25) are uncertain with predefined PDFs. This yields a six-dimensional random space ( = 6). Details about the statistical description of uncertainties are given in Table 3.  [19,53]. Reprinted from [19], Copyright (2018), with permission from Elsevier. Like Example 1, Hermite polynomials were used as the basis functions, and the PCE expressions of lognormal uncertainties were constructed with normalized Gaussian variables. In addition, it is assumed that follows a Weibull distribution, for which the PCE expression was approximated with standard random variables using its PDF as described in Section 3.2. Additionally, to compute the unknown PCE coefficients of the performance function in (25), the mBiDRM and mTriDRM were used to reduce the computational cost by converting the high-dimensional integrals into low-dimensional ones, which were further approximated with quadrature points constructed by the Gaussian-Hermite quadrature rules. Note that the number of quadrature points for each dimension was set to 5 as done in Example 1. To quantify the UQ accuracy, the relative errors defined in (24) were calculated by referring to the results of MC. For MC, 10 7 samples were used to ensure UQ accuracy.
The results of Example 2 are summarized in Table 4. It is worth mentioning that the results from previous works (e.g., [19,53]) were not given. This is because the uncertainty in [19,53] follows a normal distribution and we intentionally introduced lognormal and Weibull distributions to test the performance of our method to deal with different types Figure 5. Illustration of a roof truss structure [19,53]. Reprinted from [19], Copyright (2018), with permission from Elsevier. Like Example 1, Hermite polynomials were used as the basis functions, and the PCE expressions of lognormal uncertainties were constructed with normalized Gaussian variables. In addition, it is assumed that l follows a Weibull distribution, for which the PCE expression was approximated with standard random variables using its PDF as described in Section 3.2. Additionally, to compute the unknown PCE coefficients of the performance function Z in (25), the mBiDRM and mTriDRM were used to reduce the computational cost by converting the high-dimensional integrals into low-dimensional ones, which were further approximated with quadrature points constructed by the Gaussian-Hermite quadrature rules. Note that the number of quadrature points for each dimension θ i was set to 5 as done in Example 1. To quantify the UQ accuracy, the relative errors R defined in (24) were calculated by referring to the results of MC. For MC, 10 7 samples were used to ensure UQ accuracy.
The results of Example 2 are summarized in Table 4. It is worth mentioning that the results from previous works (e.g., [19,53]) were not given. This is because the uncertainty in [19,53] follows a normal distribution and we intentionally introduced lognormal and Weibull distributions to test the performance of our method to deal with different types of nonstandard uncertainties. Following Example 1 in the previous section, however, the LAR-based PCE was chosen to compare the UQ accuracy with two different polynomial orders (p = 2 and 3). The simulation results are given in Table 4. The total number of PCE terms for our algorithm is 28 and 84, when the polynomial order for each uncertain parameter was set to 2 and 3, respectively, while the number of PCE terms of the LAR-based PCE is 28 and 77 for each polynomial order, respectively. Similar to Example 1, it was found that the UQ accuracy can be affected by the polynomial order (p) and the number of random variables in the gDRM (M). It is important to note that the mean values (µ Z ) of different PCE-based gDRM methods in Table 4 reached the same result (0.0064), since only two significant figures were shown. Specifically, as the number of p and M increases, the UQ accuracy can be improved. For example, when M was set to 3 for the mTriDRM and when p was set to 3, the relative error R of κ Z is one order of magnitude smaller than other methods (~0.0202%). It was also found that when the polynomial order p was set to 2, the UQ accuracy of LAR-based PCE is close to our method (mBiDRM and mTriDRM) for all four statistical moments. When p was set to 3, the UQ accuracy of mTriDRM-PCE has the same order of magnitude as the LAR-based PCE. This clearly shows the algorithm in this work can deal with nonlinear problems with nonstandard uncertainties and provide accurate results as one of the most popular techniques in the literature.
For comparison, Figure 6 shows the CDFs of the performance function Z in (25), with different combinations of p and M. As seen, the UQ results of our algorithm converge to MC, as the number of p and M increases. It is important to note that the accuracy of MC is related to the total number of samples used for UQ. When the model is highly nonlinear and when the number of uncertainties is large, a large number of samples is required, thereby increasing the computational cost for MC. This is further discussed with two more complicated examples below.
For comparison, Figure 6 shows the CDFs of the performance function in (25), with different combinations of and . As seen, the UQ results of our algorithm converge to MC, as the number of and increases. It is important to note that the accuracy of MC is related to the total number of samples used for UQ. When the model is highly nonlinear and when the number of uncertainties is large, a large number of samples is required, thereby increasing the computational cost for MC. This is further discussed with two more complicated examples below.

Example 3: Truss Structure with 13 Members
A 13-bar truss structure [48] was chosen to study the performance of the UQ algorithm for dealing with a combination of different types of uncertainties. A schematic of the truss structure is shown in Figure 7, which involves 8 nodes and 13 members.

Example 3: Truss Structure with 13 Members
A 13-bar truss structure [48] was chosen to study the performance of the UQ algorithm for dealing with a combination of different types of uncertainties. A schematic of the truss structure is shown in Figure 7, which involves 8 nodes and 13 members. Illustration of a 13-bar truss structure [48]. Reprinted from [48], Copyright (2019), with permission from Elsevier.
In this example, it is assumed that there are three external loads, i.e., , , and , imposed on nodes 6, 7, and 8, respectively, which follow a normal distribution. Additionally, the elastance and the sectional area for each member are assumed to be independent and lognormally distributed. The statistical description of each parameter is listed in Table 5. Finite element analysis was used to quantify uncertainty in the performance function of the 13-bar truss structure, which is mathematically defined as in [48]: where ℎ is the threshold describing the maximum allowable deflection, i.e., displacement, and ℎ is the vertical displacement on node 3. In reliability analysis, when the performance function exceeds a limit, e.g., zero in (26), it would be considered as a failure event. Additionally, the probability of failure events is referred to as the failure probability, which is often used for structural reliability analysis [54]. Thus, it is essential to calculate the statistics of the prediction h in a precise and computationally efficient way. Here we focus on evaluating the displacement h with different UQ methods. The proposed methods in this work were integrated with finite element analysis to assess the precise prediction h, and the performance was validated with other UQ techniques, such as MC, LAR-based PCE, and the method in [48].  . Illustration of a 13-bar truss structure [48]. Reprinted from [48], Copyright (2019), with permission from Elsevier.
In this example, it is assumed that there are three external loads, i.e., P 1 , P 2 , and P 3 , imposed on nodes 6, 7, and 8, respectively, which follow a normal distribution. Additionally, the elastance E and the sectional area A for each member are assumed to be independent and lognormally distributed. The statistical description of each parameter is listed in Table 5. Finite element analysis was used to quantify uncertainty in the performance function Z of the 13-bar truss structure, which is mathematically defined as in [48]: where h max is the threshold describing the maximum allowable deflection, i.e., displacement, and h is the vertical displacement on node 3. In reliability analysis, when the performance function Z exceeds a limit, e.g., zero in (26), it would be considered as a failure event. Additionally, the probability of failure events is referred to as the failure probability, which is often used for structural reliability analysis [54]. Thus, it is essential to calculate the statistics of the prediction h in a precise and computationally efficient way. Here we focus on evaluating the displacement h with different UQ methods. The proposed methods in this work were integrated with finite element analysis to assess the precise prediction h, and the performance was validated with other UQ techniques, such as MC, LAR-based PCE, and the method in [48]. As seen in Table 5, there are two lognormal uncertainties and three normal uncertainties. For the lognormal uncertainty, the normalized Gaussian variables and Hermite basis functions were used to build the PCE surrogate models in (13), while for normal uncertainty, the formulation of PCE models is defined in (1). Once these PCE models of uncertainties are available, the resulting surrogate model of h in (26) is described with unknown PCE coefficients that can be solved with the mgDRM involving five quadrature points in each dimension (i.e., θ i = 5). Two different polynomial orders (2 vs. 3) and two different values of M (2 vs. 3) were considered. Table 6 summarizes the simulation results. It is important to note that the total number of samples for each parametric uncertainty in MC was set to 10 7 to ensure UQ accuracy. For the implementation of LAR-based PCE, 1000 model runs were used to evaluate the truss structure model as in Figure 7. Further, the sparse PCE terms were determined with the LAR algorithm to find the sensitive basis polynomials with a predefined value (i.e., 1.0 × 10 −6 ) as done in previous examples. Note that the value was used to terminate the LAR algorithm when the maximum correlation between the residual and the basis polynomials was below the value. Details of the LAR procedures can be found in [50,51]. In addition, the resulting number of PCE terms for the LAR-based PCE was 21 and 40, when the polynomial order p was 2 and 3, respectively. In contrast, the total numbers of PCE coefficients of our method were 21 and 56, which were determined using (3).
As seen in Table 6, it was found that UQ accuracy can be affected by the polynomial order (p) and the total number of integration variables (M) used in the gDRM. For example, when the mTriDRM-PCE was used (p = 3 and M = 3), the relative error R of κ h is 0.1628% in Table 6. This is smaller than the results (~0.4273%) in [48]. As compared to the LAR-based PCE, our algorithm provides comparable results. This shows the potential of the gDRM-based PCE to deal with more complicated problems, especially when estimating the higher-order statistical moments since it was previously considered challenging [55].
To graphically study the UQ accuracy with different combinations of p and M, Figure 8 shows the estimated CDFs of the displacement h. As seen, the mTriDRM-PCE, when p was 3, can provide almost identical results, as compared to MC.
coefficients of our method were 21 and 56, which were determined using (3). As seen in Table 6, it was found that UQ accuracy can be affected by the polynomial order ( ) and the total number of integration variables ( ) used in the gDRM. For example, when the mTriDRM-PCE was used ( = 3 and = 3), the relative error of is ~0.1628% in Table 6. This is smaller than the results (~0.4273%) in [48]. As compared to the LAR-based PCE, our algorithm provides comparable results. This shows the potential of the gDRM-based PCE to deal with more complicated problems, especially when estimating the higher-order statistical moments since it was previously considered challenging [55]. To graphically study the UQ accuracy with different combinations of and , Figure 8 shows the estimated CDFs of the displacement ℎ. As seen, the mTriDRM-PCE, when was 3, can provide almost identical results, as compared to MC. As compared to the first two examples, UQ in this case study requires integrating the finite element analysis with the PCE, which may increase the computational burden for UQ. We further studied the computational time to calculate the displacement h with an office desktop (Core i5-8400 central processing unit (CPU) at 2.80 GHz). Using MC with 10 7 samples for each uncertainty, it was found that~31.03 min were required to simulate all nodes in Figure 7. For the mBiDRM-PCE, it took~21.36 and 48.06 s to simulate all nodes, when the polynomial order (p) was set to 2 and 3, respectively. In addition, for the mTriDRM-PCE, it was found that~26.28 and 67.31 s were required, when p was 2 and 3, respectively. The mTriDRM-PCE took longer than mBiDRM-PCE due to the larger number of integrals that need to be solved. Further, depending on the polynomial order p, the requisite number of terms in PCE to approximate the displacement h differs as shown in (3), which leads to the difference in computational times of each method (e.g., mBiDRMand mTriDRM-PCE).

Example 4: Planar Truss Structure with 23 Members
A planar truss structure with 23 members as in Figure 9 is considered in this case study. It is assumed that all parameters cannot be known with certainty, resulting in a ninedimensional random space. This allows us to validate the UQ accuracy of the proposed method for dealing with many uncertainties. Details of model parameters are given in Table 7. office desktop (Core i5-8400 central processing unit (CPU) at 2.80 GHz). Using MC with 10 7 samples for each uncertainty, it was found that ~31.03 min were required to simulate all nodes in Figure 7. For the mBiDRM-PCE, it took ~21.36 and 48.06 s to simulate all nodes, when the polynomial order ( ) was set to 2 and 3, respectively. In addition, for the mTriDRM-PCE, it was found that ~26.28 and 67.31 s were required, when was 2 and 3, respectively. The mTriDRM-PCE took longer than mBiDRM-PCE due to the larger number of integrals that need to be solved. Further, depending on the polynomial order p, the requisite number of terms in PCE to approximate the displacement ℎ differs as shown in (3), which leads to the difference in computational times of each method (e.g., mBiDRMand mTriDRM-PCE).

Example 4: Planar Truss Structure with 23 Members
A planar truss structure with 23 members as in Figure 9 is considered in this case study. It is assumed that all parameters cannot be known with certainty, resulting in a nine-dimensional random space. This allows us to validate the UQ accuracy of the proposed method for dealing with many uncertainties. Details of model parameters are given in Table 7. Figure 9. Illustration of a planar truss structure with 23 members [19,56]. Reprinted from [56], Copyright (2006), with permission from Elsevier.  Figure 9. Illustration of a planar truss structure with 23 members [19,56]. Reprinted from [56], Copyright (2006), with permission from Elsevier. In this case study, E is the elastic modulus; A 1 and A 2 are the cross-sectional area corresponding to the specific member; and P 1 to P 6 represent external loads imposed on nodes from 8 to 13. As done in Example 3, the finite-element analysis was used to calculate the performance function Z as in [19,56]: where v max represents a specific threshold, i.e., the maximum deflection, and v is the vertical displacement on node 4 that will be approximated with different UQ methods. Our objective in this example is to accurately approximate uncertainty in v such that the failure probability can be quantified for structural reliability analysis. As seen in Table 7, three parameters (E, A 1 , and A 2 ) follow a lognormal distribution, while the rest of the parameters follow a Weibull distribution. Thus, the resulting random space is nine-dimensional, i.e., n = 9. To derive a PCE expression of v and to compute the statistical moments, the finite element analysis was coupled with the proposed methods, i.e., mBiDRM-and mTriDRM-PCE, and two different polynomial orders, i.e., p = 2 and 3, were investigated. Note that in these methods, the BiDRM approximates a nine-dimensional integral with 36 two-and 9 one-dimensional integrals, while the TriDRM requires the evaluations of 84 three-, 36 two-, and 9 one-dimensional ones. These integrals were calculated with the Gauss-Hermite quadrature rule [44], where the number of quadrature points for each random variable was set to 5. In addition, the number of samples for MC was set to 10 7 in order to compute the relative errors R in (24) for comparison purposes. The UQ results and relative errors ( R ) are shown in Table 8. Since the distribution of uncertainties is different as compared to [19,56], the results in these works were not given. However, LAR-based PCE was simulated for two different polynomial orders (p = 2 and 3) for comparison. The simulation results are summarized in Table 8. For the LAR-based PCE, the total number of PCE terms was 55 and 212, when p was set to 2 and 3, respectively. In constrast, 55 and 220 PCE coefficients were used in our algorithm for the output response, when p was set to 2 and 3, respectively. Similar to previous examples, our method provides accurate results, as the number of polynomial order (p) and the integration variables (M) in gDRM increase. For example, the mTriDRM outperforms others when p was 3, since the relative errors are smaller as in Table 8. In addition, it was found that the mTriDRM-PCE can provide comparable results to the LAR-based PCE for all statistical moments, thus confirming the accuracy of the mTriDRM-PCE-based algorithm in this work. For comparison, the CDF of the displacement v on node 4 was approximated with the proposed method and is shown in Figure 10. As compared to MC, it was found that mTriDRM-PCE provided the most accurate result, when p was set to 3. This shows the proposed algorithm can deal with complicated problems involving many uncertainties.

Conclusions
Using the polynomial chaos expansion (PCE), an uncertainty quantification (UQ) algorithm is presented to deal with many uncertainties that follow nonstandard distributions. To build PCE-based surrogate models, standard random variables are used to identify the relationship between nonstandard and standard distributions. To reduce the computational cost for calculating the PCE coefficients of model outputs, a generalized dimension reduction method is used to transform a high-dimensional integral resulting from the spectral projection (SP) into a few lower-dimensional integrals, which can be rapidly solved with quadrature rules in real-time. To show the UQ accuracy of our algorithms, four examples of structural reliability analysis were used. As compared to Monte Carlo simulations and other works in the literature, our results show the superior performance of the algorithms in terms of UQ accuracy and computational time. This shows the potential of the algorithm to tackle UQ in more complicated engineering problems that require consideration of many uncertainties. However, when uncertainty has a nonstandard distribution with a large variance, more PCE coefficients of uncertainty might be needed in the proposed approach to ensure UQ accuracy. In this case, it can lead to intensive computational burden in the presence of many uncertainties. Future work will improve the proposed algorithm to address this issue by identifying orthogonal polynomial functions only for a given uncertainty. We further studied the computational time of the proposed method and compared the results to MC. Similar to the previous case study, it was found that the computational time for MC is larger, as compared to the gDRM-based method. For example,~55.54 min were required in this case study, when 10 7 samples were used for each parametric uncertainty for MC. However, the computational time for the mTriDRM-PCE was found to be~11.52 min, when p was set to 3, which is about 80% lower than MC. It is also important to note that the gDRM-PCE based UQ method, as compared to sampling-based MC, can provide an analytical expression of the displacement v as a function of uncertainties. This provides mathematical explanations to gain a deep understanding of the problem for improved reliability analysis. The analytical expression, for instance, can be combined with sensitivity analysis techniques to find the most sensitive random variable that can significantly affect the system's performance.

Conclusions
Using the polynomial chaos expansion (PCE), an uncertainty quantification (UQ) algorithm is presented to deal with many uncertainties that follow nonstandard distributions. To build PCE-based surrogate models, standard random variables are used to identify the relationship between nonstandard and standard distributions. To reduce the computational cost for calculating the PCE coefficients of model outputs, a generalized dimension reduction method is used to transform a high-dimensional integral resulting from the spectral projection (SP) into a few lower-dimensional integrals, which can be rapidly solved with quadrature rules in real-time. To show the UQ accuracy of our algorithms, four examples of structural reliability analysis were used. As compared to Monte Carlo simulations and other works in the literature, our results show the superior performance of the algorithms in terms of UQ accuracy and computational time. This shows the potential of the algorithm to tackle UQ in more complicated engineering problems that require consideration of many uncertainties. However, when uncertainty has a nonstandard distribution with a large variance, more PCE coefficients of uncertainty might be needed in the proposed approach to ensure UQ accuracy. In this case, it can lead to intensive computational burden in the pres-ence of many uncertainties. Future work will improve the proposed algorithm to address this issue by identifying orthogonal polynomial functions only for a given uncertainty.

Conflicts of Interest:
The authors declare no conflict of interest.