View source for Design and Modeling for Computer Experiments

{{Cite book
|last1  = Fang
|first1 = Kai-Tai
|last2  = Li
|first2 = Runze
|last3  = Sudjianto
|first3 = Agus
|title  = [[Design and Modeling for Computer Experiments]]
|year   = 2006
|publisher = Chapman and Hall/CRC
}}


=Chapter 1: Introduction=

==Concepts/Definitions==

* '''Factor''' - controllable variable that is of interest in the experiment
** quantitative vs. qualitative
** '''quantitative''' - can be measured on numerical scale (e.g. T, P, ratio, rxn rate, etc.)
** '''qualitative''' - values are categories (e.g. operators, material type, etc.); also called categorical factor or indicator factor
** computer experiments: factor = input variable

* '''Experimental domain''' - hypercube of all possible factor values (also called input variable space)
** DC approach: this is the initial hypercube <math>H = \bigcap \left[ \alpha_i \leq x_i \leq \beta_i \right] </math>

* '''Run/trial''' - implementation of level-combination in experimental environment
** Computer experiments: no random error, trials are deterministic

* '''Response''' - result of a run/trial based on purpose of experiment
** Can be a function: functional response
** Chapter 7: computer experiments with functional responses
** Responses also called outputs

* '''Factorial design''' - set of level-combinations with purpose of estimating main effects, interaction effects among factors
** '''symmetric''' - all factors have same number of levels
** '''asymmetric''' - factors have diff. numbers of levels
** '''full factorial design''' - all level combinations appear
** '''fractional factorial desgin''' - subset of all level combinations

* '''ANOVA models''' - factorial designs are based on statistical model
** e.g. 1 factor experiment, q levels, expressed as:

<math>
y_{ij} = \mu + \alpha_j + \epsilon_{ij} = \mu_{j} + \epsilon_{ij} , j = 1 \dots q, i = 1 \dots n_j
</math>

*** <math>\mu</math> - overall mean of y
*** <math>\mu_j</math> - true value of response y at <math>x_j</math>
*** <math>\epsilon_{ij}</math> - random error in ith replication of jth level of x (<math>x_j</math>)
*** all errors <math>\epsilon_{ij}</math> assumed independently/identically distributed according to <math>N(0,\sigma^2)</math>
*** mean decomposed into <math>\mu_j = \mu + \alpha_j</math>, where <math>\alpha_j</math> is main effect of y at x_j, and satisfies <math>\sum_{j=1}^{q} \alpha_j = 0</math>

** e.g. 2 factor experiment (factor A, factor B):

<math>
\begin{align}
y_{ijk} &=& \mu + \alpha_i + \beta_j + ( \alpha \beta )_{ij} + \epsilon_{ijk} \\
& & i=1,\dots,p \\
& & j=1,\dots,q \\
& & k=1,\dots,K 
\end{align}
</math>

*** <math>\mu</math> = overall mean
*** <math>\alpha_i, \beta_j</math> = main effect of factor A/factor B at level i/level j
*** <math>\epsilon_{ijk}</math> = random error in kth trial at level combination <math>\alpha_i \beta_j</math>
*** <math>( \alpha \beta )_{ij}</math> = interaction between A and B at level combination <math>\alpha_i \beta_j</math>, under restrictions:
**** <math>\sum_{i=1}^{p} \alpha_i = 0</math>
**** <math>\sum_{j=1}^{q} \beta_j = 0</math>
**** <math>\sum_{i=1}^{p} ( \alpha \beta )_{ij} = \sum_{j=1}^{q} ( \alpha \beta )_{ij} = 0</math>

* Factorial design cost - for an s-factor experiment with <math>q_1,\dots,q_s</math> levels, the number of main effects plus interactions is <math>\prod_{j=1}^{s} q_j - 1</math> and this exponentially increases as s increases
** '''Sparsity principle''' - number of relatively important effects and interactions in factorial design is small
** '''Hierarchical ordering principle''' - lower order effects are more likely to be important than higher order effects; main effects more likely to be important than interactions; effects of same order are equally likely to be important

* '''Optimal design''' - given an underlying relationship, different optimization approaches may be taken
** General regression model: 
** <math>y_{k} = \sum_{j=1}^{m} \beta_j g_j ( x_{k1}, \dots, x_{ks} ) + \epsilon_{k} = \sum_{j=1}^{m} \beta_j g_j (\boldsymbol{x}_k) + \epsilon_k, k=1,\dots,n</math>
** <math>g_j(\cdot)</math> are specified or known functions
** <math>\epsilon</math> is random error, <math>E(\epsilon)=0</math> and <math>Var(\epsilon)=\sigma^2</math>
** This can be applied to several specific cases, e.g. linear model, quadratic model, etc., but g can also be nonlinear functions of x
** Rewriting in matrix form:
** <math>\mathbf{y} = \mathbf{G} \boldsymbol{\beta} + \epsilon</math>
** matrix G: design matrix
*** <nowiki>row 1 = [ g_1(x_1) \dots g_m(x_1) ]</nowiki>
*** <nowiki>row n = [ g_1(x_n) \dots g_m(x_n) ]</nowiki>
** matrix M: information matrix
*** <math>\mathbf{M} = \frac{1}{n} \mathbf{G^{\prime} G}</math>
** covariance matrix of least squares estimator:
*** <math>Cov(\hat{\beta}) = \frac{ \sigma^2 }{ n } \mathbf{M}^{-1}</math>
** Optimization techniques:
*** Want covariance matrix to be small as possible, which suggests maximizing <math>\mathbf{M}</math> with respect to dsign <math>D_n = \left[ \mathbf{x}k = (x_1,\dots,x_s) \right]</math>
** Many different proposed optimization criteria
** '''D-optimality''': maximize determinant of <math>\mathbf{M}</math>
*** Equivalent to minimizing generalized variance (determinant of covariance matrix)
** '''A-optimality''': minimize trace of <math>\mathbf{M}^{-1}</math>
*** equivalent to minimizing the sum of variances <math>\sum_{i=1}^{m} \hat{\beta}_i</math>
** '''E-optimality''': minimize largest eigenvalue of <math>\mathbf{M}^{-1}</math>

==Motivation==

* Computer model: <math>y = f(\mathbf{x})</math>
* physical experiments to understand relationship between response y and inputs <math>x_j</math> are too expensive or time consuming
* computer models important for investigating complicated physical phenomena
* one goal of computer experiments is to find an approximate model that is much simpler than the true but complicated model
* True model: <math>y = f(\mathbf{x})</math>
* Metamodel: <math>\hat{y} = f( \mathbf{x} )</math>

Comprehensive review papers:
* Sacks Welch Mitchell Wynn 1989
* Bates Buck Riccomango Wynn 1996
* Koehler Owen 1996

Computer vs. physical experiments
* involve larger numbers of variables compared to typical physical experiments
* larger experiment domain or design space employed to explore nonlinear functions
* computer experimens are deterministic

Uses for metamodels:
* preliminary study and visualization
* prediction and optimization
* sensitivitiy analysis
* probabilistic analysis (effect of input uncertainty on variability of output variable; reliability and risk assessment appliations)
* robust design and reliability-based design

Statistical approach for computer experiments:
* Design - find set of points in input space so that a model can be constructed; e.g. space filling designs
* Modeling - fitting highly adaptive models using various techniques; these are more complex, straightforward interpretations not available, use of sophisticated ANOVA-like global sensitivity analysis needed to interpret metamodel

General discussion of computer models and their use in industry (internal combustion engine application)
* Robust design and prbabilistic-based design optimization approaches have been proposed:
** Wu Wang 1998
** Du Chen 2004
** Kalagnanam Diwekar 1997
** Du Sudjianto Chen 2004
** Hoffman Sudjianto Du Stout 2003
** Yang et al 2001
** Simpson Booker Ghosh Giunta Koch Yang 2000
** Du et al 2004

* Factorial design widely used in industrial designs
** Montgomery 2001
** Wu Hamada 2000

Space-filling designs
* trying to minimize deviation between expensive/full model and metamodel
* stochastic approaches - e.g. latin hypercube sampling (LHS)
* deterministic approaches - e.g. uniform design

Chapter 2/3: LHS and UD designs

Koehler Owen 1996: different way to classify approaches to computer experiments
* "There are two main statistical approaches to computer experiments, one based on Bayesian statistics and a frequentist one based on sampling techniques."
** LHS, UD = frequentist experimental designs
** optimal LHS designs = Bayesian designs

==Modeling Techniques==

Metamodels: can be represented using linear combination of set of specific basis functions

===Univariate Functions===

Polynomial models:
* popular for computer experiments
* 2nd order polynomial models most popular
* <math>g(\mathbf{x}) = \beta_0 + \sum_{i=1}^{s} \beta_i x_i + \sum_{i=1}^{s} \sum_{j=1}^{s} \beta_{ij} x_i x_j</math>
* "response surfaces" refers to 2nd order polynomial models
** Myers Montgomery 1995
** Morris Mitchell 1995
* problems:
** unstable computations... bypassed by centering variables, e.g. replace <math>x_i</math> with <math>x_i - \overline{x_i}</math>
** collinearilty problem with high-order polynomials... bypassed by using orthogonal polynomial models
* splines - variation of polynomial models ddesigned to work in high collinearity/high order case, better than polynomials alone

Fourier basis models:
* True model is approximated using Fourier regression, set of periodic funcitons
* Number of terms increases exponentially with dimension
* In practice, one particular Fourier metamodel used
** Riccomango, Schwabe, Wynn 1997
* Wavelets:
** used to improve Fourier basis
** work esp. well when function being approximated is not smooth
** Chui 1992
** Daubechies 1992
** Antoniadis Oppenheim 1995

polynomials, splines, fourier bases, and wavelets are powerful for univariate functions, but lose effectiveness and applicability for multivariate functions

===Multivariate Functions===

Kriging model
* assumes that <math>y(\mathbf{x}) = \mu + z(\mathbf{x})</math
** <math>\mu</math> = overall mean of <math>y(\mathbf{x})</math>
** <math>z(\mathbf{x})</math> = Gaussian process with mean 0 and covariance function <math>Cov(z(\mathbf{x}_i),z(\mathbf{x}_j) = \sigma^2 R(\mathbf{x}_i,\mathbf{x}_j)</math>
** <math>\sigma^2</math> = unknown variance of <math>z(\mathbf{x})</math>
** R = correlation function with pre-specified functional form, some unknown parameters
** Typical correlation function: <math>r(\mathbf{x}_1, \mathbf{x}_2) = exp \left[ - \sum_{k=1}^{s} \Theta_i (x_{k1} -x_{k2} )^2 \right]</math>
** <math>\Theta_i</math> are unknowns

* Resulting metamodel can be written <math>g(\mathbf{x}) = \sum_{i=1}^{n} \beta_i r(\mathbf{x},\mathbf{x}_i)</math>
** which is of the general form of the linear combination of basis functions
** advantage of Krigging approach: constructs the basis directly using the correlation function
** under certain conditions, it can be shown that resulting metamodel from Kriging approach is the best linear unbiased predictor (see Ch. 5.4.1)
** Gaussian Kriging approach admits a Bayesian interpretation

Bayesian interpolation
* Proposed by:
** Currin Mitchell Morris Ylvisaker 1991
** Morris Mitchell Ylvisaker 1993
* advantage: can easily incorporate auxiliary information
** e.g.: Bayesian Kriging method
** Morris et al 1993

Neural networks (multilayer perceptron network)
* mathematical definition...
* nonlinear optimization for a least squares objective function
* training algorithms
* etc.

Radial basis function methods
* used for neural network modeling, closely relate to Kriging approach
* generally, for design inputs <math>\mathbf{x}_1, \dots, \mathbf{x}_n</math> and associated outputs <math>y_1, \dots, y_n</math>:
* K \left( \left\Vert \mathbf{x} - \mathbf{x}_i \right\Vert / \theta \right), i=1,\dots,n_i</math>
** <math>K(\cdot)</math> = kernel function
** <math>\theta</math> = smoothing parameter
* Resulting metamodel:
** <math>g(\mathbf{x}) = \sum_{i=1}^{n} \beta_i K \left( \left\Vert \mathbf{x} - \mathbf{x}_i \right\Vert / \theta \right)</math>
* if kernel function is taken to be the Gaussian kernel function (density function of normal distribution, the resulting metamodel has the same form as the Kriging metamodel

Local polynomial models
* Fan 1992
* Fan Gijbels 1996
* Concept: data point closer to <math>\mathbf{x}</math> carries more information about the value of <math>f(\mathbf{x})</math> than one that is further away
* regression function estimator: running local average
* improved version of local average: locally weighted average, e.g.
* <math>g(\mathbf{x}) = \sum_{i=1}^{n} w_i ( \mathbf{x} ) y_i</math>

'''Chapter 5''' - more explanations of these various modeling approaches, more modeling techniques, etc.

==Book Map==

Part II: design of computer experiments

Chapter 2: 
* Latin hypercube sampling
* its modifications 
** randomized orthogonal array
** symmetric Latin hypercube sampling
** optimal Latin hypercube designs

Chapter 3: 
* uniform design
* measures of uniformity
* modified L2-discrepancies
* algebraic approaches for constructing several classes of uniform design

Chapter 4
* stochastic optimization techniques for constructing optimal space-filling designs
* heuristic optimization algorithms
* high-quality space filling designs under variuos optimality criteria
* popular algorithms

Chapter 5: 
* introduction to various modeling techniques
* fundamental concepts
* logical progression from simple to increasingly complex models
** polynomial models
** splines
** kriging
** bayesian approaches
** neural networks
** local polynomials
* unified view of all models is provided
* Kriging is a central concept to this chapter

Chapter 6
* special techniques for model interpretatino
* generalizations of the traditional ANOVA (analysis of variance) for linear models
* ''highly recommended that readers study this chapter'', especially those interested in understanding sensitivity of input variables to the output
* beginning: traditional sum-of-squares decomposition (linear models)
* sequential sum of squares decomposition for general models
* Sobol functional decomposition (generalization of ANOVA decomposition)
* analytic functional decomposition (tensor product metamodel structures)
* computational technique: FAST (fourier amplitude sensitivity test)

Chapter 7
* computer experiments with functional responses
* response is in the form of a curve, where response data are collected over a range of time interval, space interval, or operation interval (e.g. experiment measuring engine noise at range of speeds)
* analysis of functional response in context of design of experiments is new area

=Chapter 2: Latin Hypercube Sampling and its Modifications=

=Chapter 3: Uniform Experimental Data=

=Chapter 4: Optimization in Construction of Designs for Computer Experiments=

=Chapter 5: Metamodeling=

==Review of Modeling Concepts==

===Mean square error and prediction error===

Mean square error (MSE):

<math>MSE(g) = \int_{T} \left[ f(\mathbf{x}) - g(\mathbf{x}) \right]^2 d\mathbf{x}</math>

No random error in computer experiments, so MSE = PE (prediction error)

Weighted mean square error (WMSE):

<math>WMSE(g) = \int_{T} \left[ f(\mathbf{x} - g(\mathbf{x}) \right]^2 w(\mathbf{x}) d\mathbf{x}</math>

<math>w(\mathbf{x}) \geq 0</math> is a weighted function, with <math>\int_{T} w(\mathbf{x}) d\mathbf{x} = 1</math>

The weighting function allows you to incorporate prior information about the distribution of <math>\mathbf{x}</math> over the domain

When no information about distribution of x over domain is available, it is assumed to be uniformly distributed

In cases where computer experiments are expensive (order of hours), impractical to evaluate prediction error directly

General strategy: ''estimate'' prediction error of metamodel g by using cross-validation procedure

for <math>i=1,\dots,n</math> let <math>g_{-i}</math> denote the metamodel based on the sample ''excluding'' <math>(\mathbf{x}_i, y_i)</math>.

===Cross Validation for Prediction Error of Expensive Models===

Cross validation score:

<math>CV_{n} = \frac{1}{2} \sum_{i=1}^{n} \left[ f(x_i) - g_{-i}(\mathbf{x}_i) \right]^2</math>

This gives a good estimate for the prediction error of g

Procedure  also called "leave-one-out cross validation"

If sample size n is large, and process of building nonlinear metamodel is time-consuming (Kriging and/or neural network models), using <math>CV_1</math> scores becomes computationally too demanding (b/c need to build n metamodels)

To reduce computational burden even further, can modify procedure

For pre-specified K, divide sample into K groups (equal sample size)

Let <math>g_{(-k)}</math> be metamodel built on sampel excluding observations in the kth group

Let <math>\mathbf{y}_{(k)},\mathbf{g}_{(k)}</math> be vector consisting of observed values and predicted values for the kth group using the <math>g_{(-k)}</math> respectively

This yields a K-fold cross validation score:

<math>
CV_{k} = \frac{1}{n} \sum_{k=1}^{K} \left( \mathbf{y}_{(k)} - \mathbf{g}_{(k)} \right)^{\prime} \left( \mathbf{y}_{(k)} - \mathbf{g}_{(k)} \right)
</math>

===Regularization Parameter===

For most modeling procedures: metamodel g depends on a regularization parameter, say <math>\lambda</math>

Cross validation score depends on regularization parameter, denoted by <math>CV(\lambda)</math> (either <math>CV_1</math> or <math>CV_K</math>)

Goal is to minimize cross validation score with respect to <math>\lambda</math>

<math>
\hat{\lambda} = min CV(\lambda)
</math>

Minimization done by searching over grid values of <math>\lambda</math>

Theoretical properties of cross validation:
* Li 1987

==General Metamodel Form==

Metamodels can be generally expressed in the form:

<math>
g(\mathbf{x}) = \sum_{j=0}^{L} B_j (\mathbf{x}) \beta_j
</math>

where <math>B_j</math> are a set of basis functions defined over experimental domain (hypercube)

This may be a polynomial basis function, covariance function (Kriging), radial basis function (neural networks), etc.

Outputs of computer experiments are deterministic: therefore metamodel construction is interpolation problem

Matrix notation:

<math>
\begin{array}{rcl}
\mathbf{y} &=& (y_1,\dots,y_n)^{\prime} \\
\boldsymbol{\beta} &=& (\beta_0,\dots,\beta_L)^{\prime} \\
\mathbf{b(x)} &=& ( B_0(\mathbf{x}),\dots,B_L(\mathbf{x}))^{\prime} 
\end{array}
</math>

and

<math>
\mathbf{B} = 
\left( 
\begin{array}{ccc}
B_0 (\mathbf{x}_1) & \dots & B_L (\mathbf{x}_1) \\
B_0 (\mathbf{x}_2) & \dots & B_L (\mathbf{x}_2) \\
\dots & \dots & \dots \\
B_0 (\mathbf{x}_n) & \dots & B_L (\mathbf{x}_n) 
\end{array} 
\right)
</math>














[[Category:Reading]]
[[Category:Textbooks]]