Evaluating correlation coefficient for Nataf transformation

Probabilistic Engineering Mechanics 37 (2014) 1–6

Contents lists available at ScienceDirect

Probabilistic Engineering Mechanics journal homepage: www.elsevier.com/locate/probengmech

Evaluating correlation coefficient for Nataf transformation Qing Xiao School of Mechatronic Engineering and Automation, Shanghai University, Shanghai, China

art ic l e i nf o Article history: Received 11 December 2013 Accepted 27 March 2014 Available online 5 April 2014 Keywords: Nataf transformation Correlation coefficient Weierstrass approximation theorem Polynomial

a b s t r a c t In this paper, a novel approach is proposed to calculate the equivalent correlation coefficient ρz in the standard normal space for two correlated random variables with desired correlation coefficient ρx. According to Weierstrass approximation theorem, ρx is expressed as a polynomial function of ρz. For a given ρx, the associated ρz is evaluated by solving the polynomial equation. Especially, when one random variable is normal, ρx is proved to be a linear function of ρz. In order to check the proposed method, a Monte Carlo simulation method is put forward. Finally, three numerical examples are worked to demonstrate the proposed method. & 2014 Elsevier Ltd. All rights reserved.

1. Introduction For analysis of engineering systems under uncertainty, the input variables are preferable to be modeled as a random vector with each entry following respective probability distribution, that is, X ¼ ðx1 ; …; xi ; …; xm ÞT . To make the problem more tractable, a transformation of X into independent standard normal space is often invoked. If the joint probability density function (PDF) of X is known, Rosenblatt transformation [1] is viable. While, this method leads to m! different transformations according to the ordering of the inputs, and not all of them give rise to the same favorable numerical properties. Moreover, the joint PDF is seldom available in many practical applications. Thus, the Nataf transformation [2] makes a worthwhile alternative to normalize the inputs, which requires the marginal PDFs and the correlation matrix of the input random variables. The major obstacle for Nataf transformation is to evaluate the equivalent correlation matrix in the standard normal space. More specifically, it requires to evaluate the correlation coefficient ρz in the standard normal space for the correlation coefficient ρx of two correlated random variables. Heretofore, much research has been done to solve this problem, such as the root finding method [3,4], the empirical formulae given for commonly used distributions [2,5], the artificial neural method [6], the linear search method [7] and the polynomial normal transformation technique [8,9]. Among these approaches, the empirical formulae and the linear search method stand out for the generality and efficiency. Through numerical experiment, 49 empirical formulae for 10 kinds of distributions

E-mail addresses: [email protected], [email protected] http://dx.doi.org/10.1016/j.probengmech.2014.03.010 0266-8920/& 2014 Elsevier Ltd. All rights reserved.

have been given. Although accurate results are obtained by these formulae, not all probability distributions are covered. Li et al. propose a linear search method to find the value of ρz, of which the efficiency is slightly lower than the empirical formula [7]. This paper is devoted to the evaluation of ρz for random variables with arbitrary marginal distributions. According to Weierstrass approximation theorem, the original correlation coefficient ρx is expressed as a polynomial function of ρz. For two arbitrary probability distributions, the function relationship between ρx and ρz can be establish by performing the interpolation only once. The value of ρz is easily obtained by solving the polynomial equation. To verify the proposed method, a Monte Carlo simulation (MCS) for calculating ρz is also presented.

2. Nataf transformation The basic idea of Nataf transformation is to generate correlated random vector with a specified correlation matrix from independent standard normal deviates. This method is also known as the NORTA (NORmal To Anything) algorithm [10]. The transformation from standard normal deviates to random variables x with a specified cumulative distribution function (CDF) F(x) is as follows: FðxÞ ¼ ΦðzÞ

ð1Þ

x ¼ F  1 ½ΦðzÞ

ð2Þ

where ΦðÞ is the CDF of the standard normal variable z, F  1 ðÞ is the inverse CDF of x.

2

Q. Xiao / Probabilistic Engineering Mechanics 37 (2014) 1–6

Suppose X ¼ ðx1 ; …; xi ; …; xm ÞT is a random vector with correlation matrix 0 x 1 ρ1;1 ⋯ ρx1;j ⋯ ρx1;m B C ⋮ ⋮ ⋮ ⋮ C B ⋮ B x C x x B C RX ¼ B ρi;1 ⋯ ρi;j ⋯ ρi;m C B C ⋮ ⋮ ⋮ ⋮ C B ⋮ @ A ρxm;1 ⋯ ρxm;j ⋯ ρxm;m Based on the transformation in Eq. (2), X can be generated by a correlated standard normal vector Z ¼ ðz1 ; …zi ; …zm ÞT with correlation matrix 0 z 1 ρ1;1 ⋯ ρz1;j ⋯ ρz1;m B C ⋮ ⋮ ⋮ ⋮ C B ⋮ B z C z z B C RZ ¼ B ρi;1 ⋯ ρi;j ⋯ ρi;m C B C ⋮ ⋮ ⋮ ⋮ C B ⋮ @ A ρzm;1 ⋯ ρzm;j ⋯ ρzm;m As long as the correlation matrix RZ is known, Z can be easily generated by the following linear transformation: U ¼ L  1 Z2Z ¼ LU

ð3Þ

where U ¼ ðu1 ; …ui ; …um ÞT is an independent standard normal vector, L represents the lower triangular matrix obtained from Cholesky decomposition of RZ ¼ LLT

ð4Þ

The procedures of generating correlated random vector from standard normal deviates can be expressed as follows: 0 1 0 1 0 1 u1 z1 x1 B ⋮ C B ⋮ C B ⋮ C B C B C B C B C B Cx ¼ F i 1 ½Φðzi ÞB C B ui C ⟶ B zi C i ⟶ B xi C ð5Þ B C B C B C B ⋮ C B ⋮ C B ⋮ C @ A @ A @ A um zm xm RZ ¼ LLT ↓ Z ¼ LU

The principal problem for Nataf transformation is to determine a suitable correlation matrix RZ in the normal space, such that a desired correlation matrix RX of X is guaranteed. More specifically, it requires to calculate ρz ði; jÞ (i a j) of RZ for each entry ρx ði; jÞ of RX .

The lower triangular matrix from Cholesky decomposition of RZ is ! 1 0 p ffiffiffiffiffiffiffiffiffiffiffiffiffi L¼ ρz 1  ρ2z Via the procedures in Eq. (5), the two correlated random variable xi and xj are generated by xi ¼ F i 1 ½Φðui Þ xj ¼ F j 1 ½Φðρz ui þ

qffiffiffiffiffiffiffiffiffiffiffiffiffi 1  ρ2z uj Þ

ð7Þ

where ui and uj are mutually independent standard normal variables. The joint PDF of ui and uj is

ϕðui ; uj Þ ¼

1  ðu2i þ u2j Þ=2 e 2π

ð8Þ

Then, E½xi xj  is calculated as ZZ qffiffiffiffiffiffiffiffiffiffiffiffiffi 1  ðu2 þ u2 Þ=2 E½xi xj  ¼ dui duj F i 1 ½Φðui Þ  F j 1 ½Φðρz ui þ 1  ρ2z uj Þe i j 2π ð9Þ Eq. (6) is rewritten as ρx ¼ 

ZZ μi μj 1 F i 1 ½Φðui Þ þ si sj 2πsi sj

 F j 1 ½Φðρz ui þ

qffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 1  ρ2z uj Þe  ðui þ uj Þ=2 dui duj

ð10Þ The above equation defines the functional relationship between ρz and ρx:

ρx ¼ Gðρz Þ

ð11Þ

If the integral on the right-hand side of Eq. (10) can be carried out analytically, ρx ¼ Gðρz Þ can be easily obtained. However, in general, the integral is not tractable [3]. While, for three cases, the analytical expressions of GðÞ are obtainable, which are presented in Appendix A. However, some properties of the function GðÞ have been proved [2], which enable us to perform an efficient numerical method to find the value of ρz. Lemma 1.

ρx is a strictly increasing function of ρz.

Lemma 2. ρz ¼ 0 for ρx ¼ 0. Lemma 3. jρx j rjρz j.

3. The function relationship between

ρx and ρz

4. Monte Carlo simulation method

Suppose xi, xj are two correlated random variables with correlation coefficient ρx, which are generated by two correlated standard normal random variables zi and zj respectively. Let ρz denote the correlation coefficient between zi and zj. Using the first cross product moment of xi and xj, the functional relationship between ρx and ρz is established:

ρx si sj þ μi μj ¼ E½xi xj  ¼ ∬ F i 1 ½Φðzi ÞF j 1 ½Φðzj Þϕðzi ; zj ; ρz Þ dzi dzj

ð6Þ

where ϕðzi ; zj ; ρz Þ is the joint PDF of two correlated standard normal variables. μi, μj denote the means of xi, xj respectively, si, sj denote the standard deviations respectively. For the two dimensional standard normal vector Z ¼ ðzi ; zj ÞT , the correlation matrix is ! 1 ρz RZ ¼ ρz 1

This section presents a MCS method to evaluate specified ρx. According to Eq. (6), ρx is expressed as

ρx ¼

E½xi xj   μi μj

si sj

¼ Corrðxi ; xj Þ

ρz for a ð12Þ

The basic idea of the MCS method is straightforward. Lemma 1 indicates that GðÞ is a strictly increasing function, thus, there is a one to one correspondence between ρx and ρz. For a given value of ρx between xi and xj, ρz can be assessed by the following steps: (1) Select a set value of ρkz evenly spaced across the interval ½  1; 1 in steps of Δρz . (2) For each ρkz , generate N bivariate random vectors ðxi ; xj ÞT by the transformation in Eq. (7), evaluate the correlation coefficient of the samples, obtaining ρkx. (3) For a prescribed ρx, find the interval ½ρkx ; ρkx þ 1 , where ρx is located, i.e. ρkx r ρx r ρkx þ 1 . Then, ρz is contained in the interval

Q. Xiao / Probabilistic Engineering Mechanics 37 (2014) 1–6

½ρkz ; ρkz þ 1 , ρkz r ρz r ρkz þ 1 , and

ρz ¼ ρkz þ

ρ ρ

kþ1  z kþ1  x

ρ ρ ρ

k z ð  k x x

ρz can be calculated as follows:

ρkx Þ

ð13Þ

3

are the skewness and kurtosis of xj respectively, then, the function relationship GðÞ between x1 and x2 can be denoted as

ρx ¼ Gðγ i ; κ i ; γ j ; κ j ∣ρz Þ

ð16Þ

For normal distribution, it holds that the skewness kurtosis κ ¼ 3, and GðÞ for two normal variables is For two arbitrary probability distributions, the function relationship between ρx and ρz can be established using the pairs of points (ρkz ; ρkx ). However, in order to guarantee a precise value of ρz, the step size Δρz should be small, and the sample size N should be large, thus, the MCS method would not be efficient. In this paper, if the analytical expression of GðÞ is unobtainable, the MCS method is used to provided benchmark to check other approach. It is interesting to reinforce that this MCS method can accommodate the Nataf transformation involving discrete random variable. If the two random variables both follow continuous distributions, a much more efficient method in the next section can be employed.

Gð0; 3; 0; 3∣ρz Þ : ρx ¼ ρz

γ ¼ 0, the ð17Þ

Through numerical experiments, it is found that the nonlinearity of GðÞ would increase as the deviation of the skewness γi or γj from 0 and the kurtosis κi or κj from 3 increase, and a higher order polynomial should be employed to approximate the function GðÞ. Testing for various probability distributions (γ r 400, 1 r κ r9  107 ), it makes evident that GðÞ can be well approximated by a polynomial of degree less than 9, that is to say, only 9 pair values of (ρz ; ρx ) are required to establish the polynomial.

6. The case of one normal variable 5. Interpolation method based on Gauss–Hermite quadrature It can be seen from Eq. (10) that ρx is a continuous function with respect to ρz, which is located within the interval ½  1; 1. According to Weierstrass approximation theorem [11], the functional relationship between ρx and ρz can be approximated to any required accuracy by a polynomial in ρz:

ρx ¼ Gðρz Þ C am ρm z þ ⋯ þa1 ρz þa0

ð14Þ

If one random variable follow the normal distribution, GðÞ is a linear function. Suppose xj is a standard normal random variable with mean μj ¼ 0 and standard deviation sj ¼ 1, then, the inverse CDF of xj is 1 F j 1 ðÞ ¼ Φ ðÞ. From Eq. (7) it follows qffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffi 1 ð18Þ xj ¼ Φ ½Φðρz ui þ 1  ρ2z uj Þ ¼ ρz ui þ 1  ρ2z uj

Lemma 2 indicates that a0 ¼ 0. The polynomial approximation of the function GðÞ can be established as follows:

Since μj ¼ 0, sj ¼ 1, according to Eq. (10), ρx is calculated as ZZ qffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 1 ρx ¼ ð19Þ F i 1 ½Φðui Þ  ðρz ui þ 1  ρ2z uj Þe  ðui þ uj Þ=2 dui duj 2πsi

(1) Choose n values of ρkz over the interval ½ 1; 1. (2) Calculate the double integral in Eq. (10) by Gauss–Hermite quadrature via Eq. (38) in Appendix B, evaluate the value of ρkx. (3) Using n pairs of ðρkz ; ρkx Þ, establish a polynomial in ρz by Newton's interpolation method.

After several manipulations (the detailed mathematical derivation is presented in Appendix C), the following equation is obtained: Z 1 2 F i 1 ½Φðui Þui e  ui =2 dui  ρz ρx ¼ pffiffiffiffiffiffi ð20Þ 2 π si

For a specified ρx, the associated value of ρz is determined by solving the polynomial equation in Eq. (14) directly. The valid solution is restricted by the following conditions:

Denote

∣ρz ∣r 1

and

ρx ρz 4 0

Based on the polynomial in Eq. (14) and Lemma 3, be evaluated by a linear search method:

ð15Þ

ρz can also

(1) Determine the search interval ρz A ½a; b. According to Lemma 3, set a ¼ ρx , b¼ 1 for ρx Z 0, set a ¼ 1, b ¼ ρx for ρx r 0. Give the precision δ. (2) Denote ρnz ¼ ða þ bÞ=2, if ðb  aÞ o δ, stop and output ρz ¼ ρnz , otherwise evaluate ρnx at ρz ¼ ρnz by Eq. (14). (3) If ρnx ¼ ρx , stop and output ρz ¼ ρnz , otherwise set a¼ a, b ¼ ρnz for ρnx 4 ρx , set a ¼ ρnz , b ¼b for ρnx o ρx , and go to Step 2. It is worthy to note that the Pearson correlation coefficient ρx is invariant to the linear transformation. That is to say, for xni ¼ ai þbi xi and xnj ¼ aj þ bj xj (ai, bi, aj and bj are constants, bi bj 4 0), the function GðÞ for xni and xnj is identical to the one for xi and xj. Therefore, GðÞ is independent of the mean and standard deviation of xi and xj. Since each continuous distribution has its unique combination of mean, standard deviation, skewness and kurtosis, this implies that GðÞ is determined by the values of the skewness and kurtosis of the two distributions involved. Suppose γi and κi are the skewness and kurtosis of xi respectively, γj and κj

1 a ¼ pffiffiffiffiffiffi 2 π si

Z

F i 1 ½Φðui Þui e  ui =2 dui 2

ð21Þ

a can be evaluated by the Gauss–Hermite quadrature (Eq. (36) in Appendix B). Then

ρx ¼ aρz

ð22Þ

If xj is a normal random variable with mean μj and standard deviation snj , then, xj ¼ ðxnj  μnj Þ=sn is a standard normal variable, the function relationship GðÞ between xj and xi is a linear function as stated in Eqs. (20) and (22). Let ρx =ρnx denote the correlation coefficients between xj =xnj and xi n

ρnx ¼

n

E½xi xnj   μi μnj

si sj

n

¼

E½xi ðsnj xj þ μnj Þ  μi μnj

si snj

¼ ρx

ð23Þ

Therefore, ρnx ¼ ρx ¼ aρz , the function GðÞ remains linear for xi and xnj . GðÞ is a linear function if one random variable is normal.

7. Numerical example In this section, three examples are worked in Matlab to demonstrate the proposed method. Consider two random variables both following the Lognormal distribution, x1  lnNð0; 12 Þ, x2  lnNð0; 12 Þ, ρx is evaluated as

4

Q. Xiao / Probabilistic Engineering Mechanics 37 (2014) 1–6

The polynomial is obtained in the same way as the former case ð24Þ

e1

ρ z resulted from Eq. (24) is regarded as a benchmark. For the MCS method, 201 points of ρkz are selected over the interval ½ 1; 1 in steps of Δρz ¼ 0:01. 106 bivariate vectors ðx1 ; x2 ÞT are generated by the procedures in Section 4 to evaluate ρkx. For the interpolation method based on Gauss–Hermite quadrature, 9 points of ρkz are selected evenly over the interval ½  1; 1. For each ρkz , evaluate the double integral in Eq. (10) by an 11-point Gauss–Hermite quadrature as stated in Eq. (38), obtaining ρkx. Using 9 pairs of (ρkz , ρkx), an 8-order polynomial is established by Newton's interpolation formula:

ρx ¼ 1:4737  10  5 ρ8z þ 1:1851  10  4 ρ7z þ 8:0813 10  4 ρ6z þ 0:0048ρ5z þ 0:0242ρ4z þ 0:0970ρ3z þ 0:2910ρ2z þ 0:5820ρz þ 1:8763  10  14

ð25Þ

The function curves given in Eq. (24), MCS and Eq. (25) are plotted in Fig. 1. From Fig. 1, it is seen that curves obtained by MCS and the interpolation method exactly coincide with the true curve established by Eq. (24). Several values of ρx are selected. The polynomial equation in Eq. (25) is solved directly by the program in Matlab. As for MCS, ρz is calculated by Eq. (13). The results are summarized in Table 1. Inspection of Table 1 indicates that ρz is precisely evaluated. Here, another example is given in terms of the Lognormal distribution lnNð0; 12 Þ and the Exponential distribution Expð1Þ.

ρx ¼ 2:4001  10  6 ρ8z þ 1:2944  10  6 ρ7z þ 0:0001ρ6z  0:0003ρ5z  0:0018ρ4z þ 0:0255ρ3z þ 0:2272ρ2z þ 0:6890ρz ð26Þ Since the analytical expression of GðÞ is unobtainable. Results from MCS with Δρz ¼ 0:01 and N ¼ 106 are regarded as a benchmark. The function curves are depicted in Fig. 2. A linear search with the required precision δ ¼ 10  4 is employed to find the real roots of the polynomial equation in Eq. (26). Along with the results obtained by Eq. (13), the values of ρz are shown in Table 2. The results given by two methods are in close agreement. Finally, a numerical example related to the standard normal distribution Nð0; 12 Þ and Lognormal distribution lnNð0; 12 Þ is worked. According to Eq. (31) in Appendix A, the analytical expression of GðÞ is 1

ρx ¼ pffiffiffiffiffiffiffiffiffiffiρz

ð27Þ

e1

ρz resulted from Eq. (27) is regarded as a benchmark. 1 MCS Interpolation

0.5

x

ρx ¼

eρz  1

ρ

follows [12]:

0

1 0.8

True curve Interpolation MCS −0.5

0.6

−0.8

−1

−0.6

−0.4

−0.2

0

0.2

0.4

0.6

0.8

1

ρ

z

0.4 ρx

Fig. 2. The graphs of function Gðρz Þ (lnNð0; 12 Þ and Expð1Þ).

0.2

Table 2 ρz for lnNð0; 12 Þ and Expð1Þ.

0 −0.2 −0.4 −1

−0.8 −0.6 −0.4 −0.2

0 ρz

0.2

0.4

0.6

0.8

1

ρx

MCS

Linear search

 0.4  0.2 0.1 0.3 0.6 0.9

 0.749  0.324 0.139 0.385 0.698 0.967

 0.750  0.324 0.139 0.385 0.698 0.967

Fig. 1. The graphs of function Gðρz Þ (lnNð0; 12 Þ and lnNð0; 12 Þ). Table 3 Table 1

ρz for Nð0; 12 Þ and lnNð0; 12 Þ.

ρz for lnNð0; 12 Þ and lnNð0; 12 Þ. ρx

Exact value

MCS

Polynomial

 0.3  0.1 0.1 0.3 0.5 0.7 0.9

 0.725  0.189 0.159 0.416 0.620 0.790 0.935

 0.725  0.188 0.159 0.416 0.620 0.790 0.934

 0.725  0.189 0.159 0.416 0.620 0.790 0.935

ρx

Exact value

Polynomial

Li et al.'s method

 0.7  0.5  0.3  0.1 0.1 0.3 0.5 0.7

 0.918  0.655  0.393  0.131 0.131 0.393 0.655 0.918

 0.918  0.655  0.393  0.131 0.131 0.393 0.655 0.918

 0.918  0.655  0.393  0.131 0.131 0.393 0.655 0.918

Q. Xiao / Probabilistic Engineering Mechanics 37 (2014) 1–6

In this case, the function GðÞ is linear. The integral in Eq. (21) is evaluated by an 11-point Gauss–Hermite quadrature as stated in Eq. (36), yielding the answer a ¼0.7629

ρx ¼ 0:7629ρz

ð28Þ

The linear search method proposed by Li et al. [7] is also employed to evaluate ρz, the results are shown in Table 3. As shown in Table 3, both methods yield accurate values of ρz. The numerical experiment is performed on a PC with Intel DualCore 3.2 GHz operator on a 3 GB RAM. The linear search method [7] takes about 1.55 s (CPU time) to evaluate ρz for the 8 different values of ρx, while, using the interpolation method, it costs only 0.05 s, showing a much higher efficiency. For the linear search method [7], if the value of ρx changes, another linear search has to be performed to evaluate ρz. While, the proposed interpolation method can establish the function relationship of ρz and ρx, which can be used to calculate ρz for an arbitrary ρx repeatedly.

8. Conclusion By Nataf transformation, multivariate random vector can be generated from the standard normal deviates. To allow for the generation of random vector with arbitrary marginal distributions and specified correlation matrix, this paper develops two methodologies to evaluate the equivalent correlation coefficient ρz in the standard normal space: the MCS method and the interpolation method. The MCS method is introduced to provided benchmark in the case that the analytical expression between ρx and ρz is unobtainable, while the interpolation method employs a polynomial in ρz to approximate the function relationship between ρx and ρz, and ρz is evaluated for a desired ρx by solving the polynomial equation.

5

Appendix B R þ1 2 The univariate integral of the form  1 e  t gðtÞ dt can be calculated by a Gaussian-type formula [13]: Z þ1 n 2 e  t gðtÞ dt ¼ ∑ ωk gðt k Þ þ Rn 1

k¼1

Rn ¼

pffiffiffiffi n! π ð2nÞ g ðξÞ 1 o ξ o 1 n 2 ð2nÞ!

ð32Þ

where n is the number of the quadrature nodes. tk is the kth zero of the Hermite polynomial H n ðtÞ ¼ ð  1Þn et

2

n

d  t2 ne dt

The associated weights pffiffiffiffi 2n  1 n! π ωk ¼ 2 n ½H n  1 ðt k Þ2

ð33Þ

ωk are given by ð34Þ

Nodes and weights up to 20 are tabulated in Ref. [14]. As can be seen in Eq. (32), Rn is the residual error of the Gauss– R þ1 2 Hermite quadrature for calculating  1 e  t gðtÞ dt. As long as the function g(t) is not a discontinuous function or an oscillation function, and enough quadrature nodes are selected, the residual error would be very small. Using Eq. (32), the following formula can be easily derived: Z þ1 n ω pffiffiffi 1 2 k ffigð 2t k Þ gðtÞ  pffiffiffiffiffiffie  t =2 dt ¼ ∑ pffiffiffi ð35Þ π 2π 1 k¼1 For the integral in Eq. (21), the constant a is evaluated as Z 1 2 F i 1 ½Φðui Þui e  ui =2 dui a ¼ pffiffiffiffiffiffi 2 π si Z 1 1 2 ¼ F i 1 ½Φðui Þui  pffiffiffiffiffiffi e  ui =2 dui si 2π pffiffiffi pffiffiffi 1 n ωk  1 ¼ ∑ pffiffiffiffi F i ½Φð 2t k Þ 2t k ð36Þ

si k ¼ 1 π

Extension of Gauss–Hermite quadrature to double integral is as follows: Z

Acknowledgments

þ1

Z

1

The author thank the referees for the careful reading of the manuscript and helpful comments.

þ1 1

gðt 1 ; t 2 Þ 

pffiffiffi pffiffiffi 1  t2 =2  t 2 =2 1 n n e 1  e 2 dt 1 dt 2 ¼ ∑ ∑ ωk ωl gð 2t k ; 2t l Þ πk¼1l¼1 2π

ð37Þ

ρx in Eq. (10) is calculated by ρx ¼ 

μi μ j 1 þ si sj πsi sj

n

n

 ∑ ∑

k¼1l¼1

pffiffiffi

pffiffiffi

qffiffiffiffiffiffiffiffiffiffiffiffiffipffiffiffi

ωk ωl F i 1 ½Φð 2t k Þ  F j 1 ½Φðρz 2t k þ 1  ρ2z 2t l Þ

ð38Þ Appendix A

Appendix C

The function ρx ¼ Gðρz Þ can be obtained analytically in three cases: Uniform distribution Uðai ; bi Þ and Uniform distribution Uðaj ; bj Þ [3]:

ρz ¼ 2 sin

π 6

ρx



ð29Þ 2

Lognormal distribution lnNðai ; bi Þ and Lognormal distribution 2 lnNðaj ; bj Þ [12]: expðρz bi bj Þ 1 ρx ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ρz 2

expðbi Þ  1

Normal distribution 2 lnNðaj ; bj Þ: b

j ρx ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ρz 2

expðbj Þ  1

ð30Þ

2

expðbj Þ  1

2

Nðai ; bi Þ

and

Lognormal

distribution

ð31Þ

ρx ¼

1

ZZ

F i 1 ½Φðui Þ  ðρz ui þ

qffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 1  ρ2z uj Þe  ðui þ uj Þ=2 dui duj

2πsi ZZ 2 2 1 ¼ F i 1 ½Φðui Þ  ui  e  ðui þ uj Þ=2 dui duj  ρz 2πsi ZZ qffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 1 F i 1 ½Φðui Þ  uj  e  ðui þ uj Þ=2 dui duj  1  ρ2z þ 2πsi Z Z 2 1 2 F i 1 ½Φðui Þui  e  ui =2 dui  e  uj =2 duj  ρz ¼ 2πsi Z Z qffiffiffiffiffiffiffiffiffiffiffiffiffi 2 1 2 þ F i 1 ½Φðui Þ  e  ui =2 dui  uj e  uj =2 duj  1  ρ2z 2πsi Z pffiffiffiffiffiffi 1 2 F i 1 ½Φðui Þui  e  ui =2 dui  2π  ρz ¼ 2πsi Z qffiffiffiffiffiffiffiffiffiffiffiffiffi 1 2 F i 1 ½Φðui Þ  e  ui =2 dui  0  1  ρ2z þ 2πsi Z 1 2 F i 1 ½Φðui Þui  e  ui =2 dui  ρz ð39Þ ¼ pffiffiffiffiffiffi 2 π si

6

Q. Xiao / Probabilistic Engineering Mechanics 37 (2014) 1–6

References [1] Rosenblatt M. Remarks on a multivariate transformation. Ann Math Stat 1952;23(3):470–472. [2] Liu PL, Der Kiureghian A. Multivariate distribution models with prescribed marginals and covariances. Probab Eng Mech 1986;1(2):105–112. [3] Li ST, Hammond JL. Generation of pseudorandom numbers with specified univariate distributions and correlation coefficients. IEEE Trans Syst, Man Cybern 1975;SMC-5(5):557–561. [4] Chen HF. Initialization for norta: generation of random vectors with specified marginals and correlations. INFORMS J Comput 2001;13(4):312–331. [5] Der Kiureghian A, Liu PL. Structural reliability under incomplete probability information. J Eng Mech 1986;112(1):85–104. [6] Niaki STA, Abbasi B. Generating correlation matrices for normal random vectors in Norta algorithm using artificial neural networks. J Uncertain Syst 2008;2(3):192–201. [7] Li HS, Lü ZZ, Yuan XK. Nataf transformation based point estimate method. Chin Sci Bull 2008;53(17):2586–2592.

[8] Chen XY, Tung YK. Investigation of polynomial normal transform. Struct Saf 2003;25(4):423–445. [9] Headick TC. Fast fifth-order polynomial transforms for generating univariate and multivariate nonnormal distributions. Comput Stat Data Anal 2002;40 (4):685–711. [10] Cario MC, Nelson BL. Modeling and generating random vectors with arbitrary marginal distributions and correlation matrix. Department of Industrial Engineering and Management Sciences. Evanston: Illinois: Northwestern University; 1997. [11] Saxe K. Beginning functional analysis. New York, Berlin, Heidelberg: SpringerVerlag; 2002. [12] Mostafa M, Mahmoud M. On the problem of estimation for the bivariate lognormal distribution. Biometrika 1964;51(3/4):522–527. [13] Naylor JC, Smith AF. Applications of a method for the efficient computation of posterior distributions. J R Stat Soc Ser C (Appl Stat) 1982;31(3):214–225. [14] Salzer HE, Zucker R, Capuano R. Table of the zeros and weight factors of the first twenty Hermite polynomials. J Res Natl Bur Stand 1952;48(2):111–116.