PhoonKK2006_GeoCongress会议_Modeling and Simulation of Stochastic Data

Modeling and Simulation of Stochastic DataK.-K.Phoon1,M.ASCE1Department of Civil Engineering,National University of Singapore;Blk E1A,#07-03,1Engineering Drive2; Singapore117576;PH+65-65166783;FAX+65-67791635;email:cvepkk@.sgAbstractStochastic data appear as basic components in reliability analysis and geostatistics.It is rarely emphasized that the multivariate probability distributions underlying random vectors(or processes,fields)are very difficult to construct theoretically,to estimate empirically,and to simulate numerically.This paper discusses the modeling and simulation of non-Gaussian random vectors,highlighting some useful geotechnical applications,important limitations,and outstanding challenges.Discussion is restricted to translation vectors,although emerging techniques that can potentially simulate a wider class of non-Gaussian random vectors would be briefly introduced for completeness. The translation approach is quite natural and takes advantage of the practicality,theoretical generality,and simulation speed associated with the multivariate Gaussian distribution.Nonetheless,there are fundamental limitations that must be recognized.The present focus on probabilistic analysis should be balanced by more research on statistical inference.It is acknowledged that statistical inference for correlated data is a challenging problem,but its practical significance is obvious.IntroductionStochastic data appear as basic components in reliability analysis and geostatistics.Reliability analysis has a significant impact on the development of modern design codes.Much of its success could be attributed to the advent of the first-order reliability method(FORM)–which provides a practical scheme of computing small probabilities of failure at high dimensional space spanned by the random variables in the problem.The collection of random variables provides a formal probabilistic model for the uncertainties in the input parameters and is often called a random vector.Ang and Tang(1984)presented numerous practical applications of FORM in their well known book on Probability Concepts in Engineering Planning and Design.Geostatistics arose from the need to estimate spatial patterns in geology and mining.The physical premise that makes such estimation possible is that points close in space tend to assume similar values.The autocorrelation function(or variogram)is a fundamental tool describing similarities in a statistical sense as a function of distance.The works of G.Matheron,D.G.Krige,and F.P. Agterberg are notable.Parallel developments also took place in meteorology(L.S.Gandin)and forestry(B. Matérn).Geostatistics is mathematically founded on the theory of random processes/fields developed by A.Y. Khinchin,A.N.Kolmogorov,P.Lévy,N.Wiener,and A.M.Yaglom,among others.The interested reader can refer to books by Cressie(1993)and Chilès and Delfiner(1999)for details.Vanmarcke(1983)remarked that all measurements involve some degree of local averaging and that random field models do not need to consider variations below a finite scale because they are smeared by averaging.Baecher&Christian(2003)provided an updated review of geostatisical applications in geotechnical engineering.Fenton and Griffiths have been studying soil-interaction problems within the context of a random field since the early nineties(e.g.,Griffiths&Fenton,1993, 1997,2001;Fenton&Griffiths,1997,2002,2003).It is rarely emphasized that the multivariate probability distributions underlying random vectors(or processes, fields)are very difficult to construct theoretically,to estimate empirically,and to simulate numerically.A cursory review of basic probability texts would reveal that non-Gaussian distributions are usually presented in the uni-variate form.It is not possible to generalize these uni-variate formulae to higher dimensions unless the random variables are independent.This assumption is obviously invalid in geostatistics and overly-restrictive in reliability analysis where input parameters can be correlated by physics(e.g.,sliding resistance increases with normal load,undrained shear strength increases with loading rate,etc.)or by curve-fit(e.g.,cohesion and friction angle in a linear Mohr-Coulomb failure envelope).In principle,multivariate distributions can be constructed using their conditional distributions as illustrated below for a general random vector with three components:)x (F )x |x (F )x ,x |x (F )x (F )x (F )x ,x (F )x ,x (F )x ,x ,x (F )x ,x ,x (F 112213112121321321==(1)in which F is the n-dimensional cumulative distribution function (n is the number of arguments in F).The catch is that closed-form equations for the conditional probabilities in Eq.(1)are not available except for the Gaussian case,a class of elliptical distributions (Azzalini and Capitanio 2003;Arslan 2004),and a handful of others (Devroye 1997).It is evident from the statistics literature that simulation of random samples from multivariate distributions is not trivial,much less problems relating to statistical estimation and model identification.In practice,it is further recognized that a general random vector could only be characterized reliably up to the marginal distribution and a measure of dependency such as the popular product-moment (Pearson)correlation coefficient.With the exception of the Gaussian case and generalizations such as elliptical (Arslan 2004)and skew elliptical (Azzalini and Capitanio 2003)distributions,a general random vector could not be defined uniquely based on such limited information.An extensive body of literature on simulation of random processes exists (much of which are fairly recent)to attest to the theoretical and computational difficulties involved (e.g.,Grigoriu 1984;Grigoriu 1988;Yamazaki and Shinozuka 1988;Fenton &Vanmarcke,1990;Deodatis and Micaletti 2001;Phoon et al.2002a;Phoon et al.2002b;Puig et al.2002;Sakamoto and Ghanem 2002;Field and Grigoriu 2004;Phoon et al.2005a).At present,it is accurate to say that the most popular method for constructing non-Gaussian random vectors is to apply some suitable memoryless nonlinear transforms to the standard Gaussian vector.Such models are known as translation vectors.The central importance of the multivariate Gaussian distribution underlying the standard Gaussian vector is perhaps not surprising given its analytical tractability,complete definition by mean vector and covariance matrix,and availability of fast Gaussian simulation methods.This paper discusses the modeling and simulation of non-Gaussian random vectors,highlighting some useful geotechnical applications,important limitations,and outstanding challenges.Discussion is restricted to translation vectors,although emerging techniques that can potentially simulate a wider class of non-Gaussian random vectors would be briefly introduced for completeness.Issues raised in this paper are relevant to reliability analysis and geostatistics,because they deal with basic theoretical and computational aspects of stochastic data.Correlated Gaussian Random VectorsThe multivariate Gaussian probability density function is available analytically and can be defined uniquely by a mean vector and a covariance matrix:()()m X C 'm X 212n211-e )2(C )X (f =(2)in which X =(X 1,X 2,…X n ) is a Gaussian random vector with n components,m is the mean vector,and C is the covariance matrix.For the 2D (simplest)case,the mean vector and covariance matrix are given by:==2221212121C m m m (3)in which m i and i =mean and standard deviation of X i ,respectively,and =product-moment (Pearson)correlation between X 1and X 2.The Gaussian probability density function for the 2D (or bivariate)case is illustrated in Fig.1a.(a)(b)Figure1.Bivariate Gaussian probability density function for:(a)correlated case and(b)uncorrelated case The practical usefulness of Eq.(2)is not well appreciated.First,the full multivariate dependency structure of a Gaussian random vector only depends on a covariance matrix(C)containing bivariate information(correlations) between all possible pairs of components.The practical advantage of capturing multivariate dependencies in any dimension(i.e.,any number of random variables)using only bivariate dependency information is obvious.Note that coupling two random variables is the most basic form of dependency and also the simplest to evaluate from empirical data.In fact,there are usually insufficient data to compute reliable dependency information beyond correlations in actual engineering practice.For an n-dimensional random vector,the covariance matrix C is populated by n(n-1)/2correlation coefficients.The important practical question here is whether these correlation coefficients arising from distinct pairs of components can be selected independently.The answer is no and this important but somewhat neglected issue would be discussed in the first sub-section below.The second useful aspect is the elliptical form(X-m) C-1(X-m)appearing in the exponent of Eq.(2)(manifested as elliptical contours in Fig.1a).The simulation of a general n-D Gaussian random vector is fast and robust because of this elliptical form.The computational steps are also simple and easy to implement as discussed in the second sub-section below.Generalizations of this elliptical form are available,but the practical advantages of applying more complex generalizations are still being studied at present(Azzalini and Capitanio2003;Arslan2004).From a practical viewpoint[only bivariate information needed,correlations are easy to estimate,voluminous n-D data compiled into mere n(n-1)/2coefficients,etc.]and a computational viewpoint(fast,robust,simple to implement),it is accurate to say that the Gaussian probability model and non-Gaussian probability models generated from the translation approach are already very good and sufficiently general for many practical scenarios.For stochastic data that cannot be modeled using the translation approach,the hunt for probability models with comparable practicality, theoretical power,and simulation speed is still on-going.Active threads of research in this direction are copula theory(Schweizer1991)and non-Gaussian Karhunen-Loeve expansion(Phoon et al.2005a).Positive DefinitenessBy definition,correlation coefficients take on numerical values between-1and1.Unfortunately,not all covariance matrices containing such entries are valid.This observation is not completely unexpected given that a correlation contains pairwise dependency information and a covariance matrix contains many entries arising from overlapping pairs[e.g.,(1,2),(1,3),(2,3)].It is possible to demonstrate this mathematically using the following definition for partial correlation(Bedford and Cooke2002):)1)(1(2n ...4232n ...413n...423n ...413n ...412n ...312 = (4)All partial correlations can be computed from correlations by iterating the above equation.For example,the first-order partial correlation is determined from Eq.(4)as:)1)(1(223213231312312 = (5)Because partial correlations also lie between -1and 1,it is clear that a relationship exists between 12, 13and 23:2313223213122313223213)1)(1()1)(1( + + (6)Hence,if 13=0.8and 23=0.1,then -0.517 12 0.677.For a 3D random vector with 13=0.8and 23=0.1,the last entry in the covariance matrix 12cannot take values in the full range between -1and 1.In practice,it is possible to check the entire covariance matrix in one step (rather than entry by entry as shown above)using the property of positive definiteness.A symmetric matrix C is a valid covariance matrix if it is positive definite,i.e.:a T C a >0for all a 0(7)in which a =(a 1,a 2,…,a n )T is a vector of real numbers.The practical significance of Eq.(7)is most readily appreciated by defining a random variable Y as a weighted linear sum of components in a random vector X,i.e.:i n 1i iX a Y ==(8)If C is the covariance matrix for X,then it can be shown that the variance of Y is:Y 2=a T C a (9) Hence,if Eq.(7)is not satisfied,some combination of weights could produce a random variable with negative variance!This contradicts the definition of the variance.Note that it is not easy to check the validity of a covariance matrix by inspection –only one entry is different between the valid and invalid examples shown below.ValidInvalid 18.03.055.08.012.09.03.02.012.055.09.02.0118.03.02.08.012.09.03.02.012.02.09.02.01Y =1.81X 1+0.05X 2+1.66X 3–0.53X 4Y 2=0.57 Y 2=-0.10Fortunately,a common numerical procedure called Cholesky factorization can be used as a rigorous test for positive definiteness.A lower triangular matrix L is a Cholesky factor for C if:C =L L (10)Cholesky factorization can be roughly appreciated as taking the “square root”of matrices.The Cholesky factor can be computed in EXCEL using the array formula MAT_CHOLESKY,which is provided in a free matrix and linearalgebra add-in at http://digilander.libero.it/foxes/index.htm .The covariance matrix is valid (positive definite)if all diagonal entries in the Cholesky factor are positive as shown below.Valid Cholesky factorInvalid Cholesky factor 06.072.042.055.0044.002.09.00098.02.0000123.144.135.055.0044.002.09.00098.02.00001Simulation of Correlated Gaussian DataIt can be verified that the following transform would produce a correlated Gaussian vector X with the desired covariance matrix from a standard Gaussian vector U:X =LU +m (11)A standard Gaussian vector contains uncorrelated components with zero mean and unit variance.The probability density function is shown in Fig.1b.The “Random Number Generation”tool under “Data Analysis”in EXCEL can produce standard Gaussian numbers.For concreteness,consider a bivariate Gaussian vector with covariance given by Eq.(3).The lower triangular Cholesky factor L is:=222110L (12)Hence,the desired correlated Gaussian vector is:2221221111m )1U U (X m U X + + =+ =(13)In practice,p realizations of an n ×1correlated vector can be easily computed from Eq.(11)by substituting U with an n ×p matrix containing uncorrelated Gaussian numbers and m with p column duplicates of the mean vector (m 1,m 2,…m n ) .Matrix multiplication can be carried out in one step using the array formula MMULT in EXCEL.Eq.(11)can be modified to achieve the desired covariance matrix even for small p (Phoon 2004a).Translation Non-Gaussian Random VectorsIt is widely accepted that a practical procedure for modeling dependent non-Gaussian vectors should involve at most the marginal distributions and the covariance matrix.The challenge is to find a multivariate distribution function that is consistent with this second-order information.It should be emphasized that the solution is not unique .Given the computational ease of simulating correlated Gaussian vectors,it is natural to construct non-Gaussian vectors using the following translation approach:Y i =F i -1 (X i )(14)in which Y i =non-Gaussian random variable with cumulative distribution function F i and X i =standard Gaussian random variable with cumulative distribution function (but X i and X j may be correlated).Lognormal examplePhoon et al.(2006)proposed a simple but robust probabilistic model for load-displacement curves from augered cast-in-place (ACIP)piles.The basic idea is to:(a)reduce measured load-displacement curves into 2parameters via hyperbolic fit and (b)normalize resulting hyperbolic curves using an interpreted failure load to reduce the datascatter.The remaining scatter can be modeled using an appropriate random vector for the curve-fitting parameters.The normalized hyperbolic curve considered in their study is expressed as:by a y Q QSTC +=(15)in which Q =applied load,Q STC =failure load interpreted using slope tangent method,a and b =curve-fitting parameters,and y =pile butt displacement.Note that the reciprocals of a and b are equal to the initial slope and asymptotic value,respectively.The slope tangent method defines the failure load at the intersection of the load-displacement curve with the initial slope line offset by (0.15+B/120)inches,in which B =shaft diameter in inches.The normalized curves for piles with depth to diameter ratios (D/B)greater than 20are shown in Fig.2.The statistics of the a and b parameters (without outlier c30-1)are summarized in ing the Anderson-Darling goodness-of-fit test,one could formally state that there is no evidence to reject the null hypothesis of lognormality for both curve-fitting parameters at the customary 5%level of significance (p AD >0.05).In addition,the Pearson correlation between a and b is strongly negative ( a,b =-0.67),i.e.initial slope (1/a)tends to increase when the asymptotic value (1/b)decreases or vice-versa (Fig.3c).To construct a correlated lognormal random vector for the above example,observe that Eqs.(13)and (14)simplify to:])1U U (exp[b )U exp(a 22212111 + + = + =(16)The parameters 1and 1in Eq.(16)can be calculated from the mean (m a )and standard deviation (s a )of a as follows:()21a 12a 2a 15.0)m ln(m /s 1ln = += (17)The parameters 2and 2can be calculated from the mean and standard deviation of b in the same way.Note that is not the same as the a,b .It is given by:]1)][exp(1)[exp(1)exp(222121b ,a = (18)Fig.4shows 25simulated normalized hyperbolic curves for a,b =-0.67and a,b =0.It is evident that Fig.4a resembles Fig.2.However,Fig.4b does not compare well with Fig.2,even though the mean,standard deviation,and lognormal distribution for both a and b are followed.This example highlights the importance of verifying the presence of correlations whenever there are more than two variables in the problem and incorporating correlation information correctly in the simulation process.Note that there are many ways of coupling two lognormal random variables into a correlated random vector.Additional examples are given in Fig.11below.The above translation approach is considered reasonable based on similarities between actual (Fig.2)and simulated (Fig.4a)data.Model identification becomes an issue when a much broader spectrum of non-translation approaches are considered.051015202530Displacement (mm)0.01.02.03.0Q /Q S T CD/B > 20c30-1Figure 2.Normalized hyperboliccurves.05101520a parameter (mm)00.10.20.3R e l a t i v e F r e q u e n cy(a)0.00.5 1.0 1.5 2.0b parameter 00.10.20.3R e l a t i v e F r e q u e n cy(b)Figure 3.(a)Histogram for a parameter (b)Histogram for b parameter,and (c)correlation between a and b0510********Displacement (mm)0.01.02.03.0Q /Q S T C(a)051015202530Displacement (mm)0.01.02.03.0Q /Q S T C (b)Figure 4.Simulated normalized hyperbolic curves for a,b =(a)-0.67and (b)0.Hermite PolynomialsThere are two computational challenges in the simulation of translation vectors.First,the cumulative distribution function F in Eq.(14)is usually not available in closed-form,even for common classical distributions (e.g.,beta,gamma).It can be tedious to take the inverse of such a function.Second,the Gaussian correlation coefficient in Eq.(13)can only be computed from the empirical correlation coefficient in the non-Gaussian data using an integral equation.Eqs.(16)and (18)presented in the previous sub-section are only applicable for the special lognormal case.In a state-of-the-art report,Phoon (2004b)proposed the application of Hermite polynomials to approximate Eq.(14):==0k i k ik i )X (H aY (19)in which the Hermite polynomials H j (.)are given by:)X (H k )X (H X )X (H X3X )X (H 1X )X (H X)X (H 1)X (H 1k k 1k 332210 + = = ===(20)Note that the first coefficient a 0(k =0)is simply the mean of Y.The last row of Eq.(20)shows that Hermite polynomials of any degree (k)can be computed efficiently using a simple recurrence relation that depends only on two preceding Hermite polynomials.This recurrence relation can be implemented directly using EXCEL.The numerical values of the coefficients,a k ,depend on the distribution of Y.They can also be computed readily using EXCEL even if the empirical distribution of X cannot be fitted to any classical probability distribution functions:1.Let y be a p ×1vector containing measured/simulated data from a cumulative distribution function,F(y).2.Let x = -1F(y)be a p ×1vector containing p realizations of a standard Gaussian random variable.The function-1 can be invoked using NORMSINV in EXCEL.3.Let h 0be a p ×1vector containing ones,h 1=x,h 2=x *x –1,…h q-1 =x *h q-2 – (q-2)h q-3,and H be a p ×qmatrix containing h 0,h 1,h 2,…h q-1 in the columns.The operator “ *”means element-wise matrix multiplication (MATLAB convention),i.e.,for matrix A =B *C,(i,j)element in A,a ij =b ij ×c ij .4.Let a be a q ×1vector containing the unknown Hermite coefficients (a 0,a 1,a 2,…a q-1) .This vector iscomputed by solving the following system of linear equations:(H H)a =H y (21)Eq.(21)can be solved easily using array formulae and matrix functions in EXCEL (TRANSPOSE,MMULT,MINVERSE)or SYSLIN from matrix and linear algebra add-in.Phoon et al.(2005b)demonstrated that reasonably accurate Hermite coefficients can be obtained from a relatively small sample size (p =20).Figures 5and 6demonstrate that the Hermite expansion [Eq.(19)]is very efficient even for strongly non-Gaussian distributions -it is possible to match very small probabilities (say 10-4)using only 6to 8terms.Once the Hermite coefficients are computed,Eq.(19)can be applied for simulation with significantly less computational effort than executing Eq.(14)tens of thousands of times in a low probability simulation exercise (typical of civil engineering problems).Hermite polynomials can also improve the robustness and efficiency of first-order reliability analysis (Phoon et al.2005b;Phoon and Honjo 2005)Correlation StructureThe difficulty of calculating the Gaussian correlation coefficient from the empirical correlation coefficient in the non-Gaussian data is well known in the structural reliability literature.Recent papers (Puig et al.2002;Sakamoto and Ghanem 2002)apply the Hermite expansion for simulation of non-Gaussian translation processes,because it presents a relatively simpler power series solution to this difficulty (Sakamoto and Ghanem 2002):1010101010100x F (x )1010101010100x F (x )Figure 5.Probability tails of Hermite expansions for lognormal distribution with =0and =0.5(left); =1.0(right)1010101010100x F (x )1010101010100x F (x )Figure 6.Probability tails of Hermite expansions for gamma distributions with scale parameter =1and shape parameter =0.5(left);shape parameter =2.0(right)= =0k k 2k 1k 21a a!k Y EY (22)The above equation should be computed recursively,because factorials are notorious for growing very rapidly (see Phoon 2004b).Note that the unknown in Equation (22)is and analytical solutions are not available in the general case.However,it is possible to solve Eq.(22)using numerical methods such as the nonlinear optimization Solver function in EXCEL.Puig et al.(2002)noted that optimization is preferred rather than applying a nonlinear equation solver because it is possible to impose the added constraint of positive definiteness.Figure 7shows that a Hermite expansion containing less than or equal to 4terms is sufficient to reproduce the relationship between Y1Y2and X1X2(= )for lognormal random vectors with two significantly different components ( 1=0.3and 2=1)and with highly skewed identical components ( 1=1and 2=1).Although translation random vectors are general and easy to construct,there are two fundamental theoretical limitations that must be recognized.First,assuming that the empirical correlation matrix is positive definite,there is no guarantee that the correlation matrix of the equivalent Gaussian random variables is positive definite,even if the correlations are determined exactly from Eq.(22).This is a fundamental theoretical constraint related to the characterization of multivariate dependence based on pairwise measures of dependence [Eq.(22)is solved one pair of components at a time].Hence,performing Cholesky factorization on the correlation matrix of the underlying Gaussian vector X is still necessary.Second,the expedient assumption of Y1Y2! X1X2only works well when 1and 2are small (using the lognormal distribution for illustration).This is shown in Fig.8a.It can be seen in Fig.5that the lognormal probability density function becomes increasingly skewed as increases.Hence,it may be concluded (not surprising)that this assumption does not work for strongly non-Gaussian distributions.More significantly,the bounds for Y1Y2could be less than 1and/or greater than -1.The lower bound of -1becomes increasingly difficult to achieve when increases (Fig.8b).The upper bound of 1can be achieved as long as the components in the random vector are identically distributed (Grigoriu 1998).For the general case of non-identical components with possibly large differences in ,both theoretical upper and lower bounds are unachievable (Fig.8d).This situation may be rare in structural problems,because uncertainties in loads and structural materials typically fall within a narrow band (coefficients of variation between say 10%and 20%).For geotechnical uncertainties,the coefficients of variation fall within a much broader band (Phoon and Kulhawy 1999).A possible example is a random vector containing the effective stress friction angle and the Young’s modulus.The value ofX1X2 Y 1Y 2(a) 1=0.3and 2=1X1X2 Y 1Y2(b) 1=1and 2=1Figure 7.Relationship between Y1Y2and X1X2for lognormal randomvectors.X1X2 Y 1Y 2(a) 1=0.3; 2=0.3X1X2 Y 1Y 2(b) 1=1; 2=1 X1X2 Y 1Y 2(c) 1=1; 2=0.3 X1X2 Y 1Y 2(d) 1=2; 2=0.3Figure 8.Some relationships between observed correlation ( Y1Y2)and correlation of underlying equivalent Gaussian random variables ( X1X2)(approximately equal to the coefficient of variation of the physical variable)can be as low as about 0.1for the first component and as high as about 1for the second component (Phoon and Kulhawy 1999).The limitation highlighted above has been proven mathematically for the general case by Arwade (2004).The practical ramification is that some observed data cannot be modeled using translation vectors.It is tempting to conclude otherwise because the simple cumulative distribution transform in Eq.(14)looks “obvious”.However,Eq.(14)is only obvious in the univariate case,which demonstrates quite clearly that extension to the multivariate case is non-trivial.Li et al.(2006)demonstrated that the non-Gaussian Karhunen-Loeve expansion does not follow the standard translation model and the above theoretical limitation does not apply.This emerging research on non-translation models is significant,because a broader class of random vectors/processes can be simulated.Simulation-based Statistical Inference for Correlated DataA natural probabilistic model for correlated spatial data is the random field (Vanmarcke 1977).While the random field provides a concise description of spatial variation,it poses considerable practical difficulties for statistical inference in view of its complicated data structure.Previous studies are mostly probabilistic,wherein a random field is posited to describe soil variations and the probability of exceeding some limit state is computed.The statistical question of how this random field can be inferred from limited observations is no less important,but attempts to address this aspect are currently limited.Christian and Baecher (2006)commented on this problem as well and recommended the application of Bayesian methods to assess the probability that the model is true,given observed data.All classical statistical tests are invariably based on the important assumption that the data are independent (Cressie 1993).When they are applied to correlated data,large bias will appear in the evaluation of the test statistics (Phoon et al.2003).The application of standard statistical tests to correlated soil data is therefore potentially misleading.Independence is a very convenient assumption that makes a large part of mathematical statistics tractable.It is well recognized that real data are usually dependent.Statisticians can go to great length to remove this dependency (Fisher 1935)or be contented with less powerful tests that are robust to departures from the independence assumption.In recent years,a third approach involving the direct modeling of dependency relationships into very complicated test statistics through Monte Carlo simulation has become feasible because desktop computing machines have become very powerful.One of the first geostatistical studies that systematically exploits powerful digital simulation techniques for statistical inference in the correlated case was conducted recently by Phoon et al.(2003).The key features of this study are briefly summarized below to highlight the significant potential in this approach and to encourage further practical developments in this somewhat neglected area of statistical inference.Simulation of Random ProcessesThe basic step in simulation-based statistical inference is the generation of realizations that are consistent with a given target random process/field.The computational techniques discussed previously on the construction of translation random vectors capable of matching target marginal distributions and covariance matrices can be applied to random processes/fields.The only computational aspect that requires some elaboration is that simulation of a process/field is usually more efficient in the frequency domain.This is illustrated below using the Gaussian process.Realizations belonging to a zero-mean stationary Gaussian process X(t)can be generated using the well established spectral approach as follows:)t f 2cos V t f 2sin U ()t (X k k M 1k k k k += =(23)in which f )f (S 2k k "= ,"f is the interval over which the spectral density function S(f)is discretized,f k =(2k-1)"f/2,and U k and V k are standard uncorrelated Gaussian random variables.The single exponential autocorrelation function is commonly used in geostatistics:R(#)=exp(-|#|),in which #is the distance between data points.The corresponding spectral density function is:。