• Keine Ergebnisse gefunden

Analysis and Design of Simulation Experiments for the Approximation of Models

N/A
N/A
Protected

Academic year: 2022

Aktie "Analysis and Design of Simulation Experiments for the Approximation of Models"

Copied!
16
0
0

Wird geladen.... (Jetzt Volltext ansehen)

Volltext

(1)

NOT FOR QUOTATION WITHOUT PERMISSION OF THE AUTHOR

ANALYSIS AND DESIGN OF SIMULATION EXPERIMENTS FOR THE APPROXIMATION OF MODELS

V. Fedorov

July 1983 WP-83-71

Working Papers are interim reports on work of the International Insti- tute for Applied Systems Analysis and have received only limited review. Views or opinions expressed herein do not necessarily represent those of the Institute or of its National Member Organiza- tions.

INTERNATIONAL INSTITUTE FOR APPLIED SYSTEMS ANALYSIS 236 1 Laxenburg, Austria

(2)

PREFACE

Much of IIASA's work is concerned with modeling large a n d complex sys- t e m s . However, the resulting mathematical models t e n d to become very com- plicated a n d unwieldy, making it very difficult t o identify the k e y r e l a t i o n s h p s between t h e variables. To overcome this problem, i t is often necessary t o approximate these "primary models" by more t r a n s p a r e n t "secondary models".

In this p a p e r , Valeri Fedorov discusses these problems and presents a s t a - tistical procedure for secondary model construction. The use of t h e method is illustrated b y application to one of t h e IIASA energy models.

Andrzej Wierzbicki C h a i r m a n S y s t e m and D e c i s i o n S c i e n c e s

(3)

ANALYSIS AND DESIGN OF SIMULATION EXPERIEdlENTS FOR THE APPROXIMATION OF MODELS

V. Fedorov

1. INTRODUCTION

Modern computers make it possible to construct and run complicated mathematical models of complex systems (e.g., economic systems, ecological systems) which involve hundreds of inputs and equations. The links be tween the different variables (inputs and outputs) and equations in these models are usually very difficult to follow, and this is complicated by the fact that the models are continuously being updated and improved by the incorporation of new mathematical features. Sometimes models consist of modules (elements) prepared by different scholars, and this is one reason why mathematical models (or, more accurately, their computerized counterparts) occasionally become "mysterious" even to their authors. Analytical techniques prove to be useless in analyzing the properties of these models. Since it is not possible to obtain the required results in t h s way, it is natural to try another approach:

one possibility is to carry out experiments on the mathematical models them- selves. We shall call these s i m u l a t i o n e x p e r i m e n t s .

(4)

The question of the effectiveness of t h e experiments and whether the chosen model adequately describes t h e empirical data arises a t the very begin- ning of such investigations. To study t h s , "models" of the models are often constructed. In what follows the terms s e c o n d a r y m o d e l and p r i m a r y m o d e l will be used in an attempt to avoid confusion.

The construction of secondary models can also be stimulated by the fact that t h e primary models are frequently too detailed for the specific investiga- tions that the researcher wishes to perform. For instance, t o describe t h e behavior of a primary model over a relatively small range of input values it might be sufficient t o use a polynomial approximation of the model. One attractive feature of t h s approximation is that it then becomes possible to develop fast real-time interactive software. T h ~ s type of software can be extremely useful to decision makers because it allows them to scan a lot of variants in a relatively short time.

The secondary model should reflect the structure of both the primary model and the experiment. As in other experimental situations, it is possible to construct a number of primary models of one system which are all based on different principles and suited to a different type of experiment. Everybody can recall cases in which the same system has been described by either a sto- chastic model or a deterministic model, depending on the experiment planned.

For the sake of simplicity we shall restrict ourselves to deterministic pri- mary models. Assume that the primary model connects three sets of variables z, w , a n d y :

w = + ( ~ I Y ) (1)

where z E R ~ , w ER' and y€Rq. Usually vector z is composed of control vari- ables and ill-defined variables, while y comprises variables whose values are known with relatively h g h precision. The way in which these groups are defined

(5)

will, of course, depend on the researcher

We shall now explain some of the terms and notation used in the following sections. The result obtained by evaluating function (1) for given xi and yi ,

N .

will be called the (simulation) measurement. The set

CN

= tzi j wlll be taken to represent the design of an experiment, while the set

Z N

= tyi ,yi ,xi will be defined as a ( s i m d a t i o n ) ezperiment. (We shall generally omit the word "simu- lation" in what follows.)

In most cases, the dimension of the output or response vector w is small wlxle the dimensions of vectors z and y can be as large as several hundreds.

But (and this is one of the main assumptions) it is assumed that the responses y depend "strongly" upon only a few "significant" components of vector 3:.

The goal of simulation experiments is to identify these significant com- ponents (variables) and to construct some approximation r](z) of the response function $(z ,y). It should be emphasized that the function r](z) does not depend on variables y because the latter are assumed to be known relatively precisely.

2. STATISTICAL BACKGROUND

2.1. Screening Experiments

The aim of screening experiments is to detect the truly significant factors in a large collection of possibly significant factors (see, for instance, Li, 1962;

Meshalkin, 1970; Satterthwhaite, 1959). To get a n idea of t h e methods used, we shall consider one of the simplest approaches.

Suppose t h a t we have

(6)

where the yi are measurements, the are unknown parameters, b , I z , l c , , a=

-

l , m , and the E~ a r e independent random errors with zero

mean and variance u2. It is assumed that s parameters ( s l m ) are nonzero, where the value of s is known.

Consider a random design

tN

constructed in the following way. Each measurement is carried out under random conditions z such that xai = b , or

zai = C , with probability 0 . 5 . Assume that . . .

,a,,

is the solution of the

following extremal problem:

where p a can be zero or one. Let pN(u2)=1-bN be the probability that the nonzero parameters have been identified correctly. Then

For

N

1 m there are regular deterministic designs with the property lim 6(a2) = 1

.e-0

In other words, the use of such random designs makes it possible to reduce the number of measurements (which is usually essential) although the researcher pays for it with the resulting value of PN = 1-BN

<

1 (compare (3) and (4)).

Most of the efforts in the theory of screening experiments have been d r e c t e d towards minimizing the number of observations N necessary under given b N .

2.2. Design of Regression Experiments Consider the regression model

where ?P is the vector of unknown parameters, f ( x ) is the vector of basis

(7)

functions, the E ~ , i =

-

l , N are random errors with zero means and variances of E [ E ~ ] = h - l ( z i ) , and x i € X where X is the operability region. The precision of the best linear unbiased estimator

2

of the parameters 19 is determined by the information matrix:

where pi = .ri / N and ri is the number of measurements necessary under con- N

dition x i , such t h a t

C r i =

N . Recall that the covariance matrix ~ ( 2 ) of the

i = l

estimator

3

equals M - ' ( ( ~ ) . Then

is the optimal design for regression model (5) under criterion + ( M ) . + ( M ) is usually a monotonic convex function of a positive semidefinite matrix M , for example, + ( M ) = In det M-' or * ( M ) = t r M - ' .

2.3. Regression Analysis

This branch of statistics is composed of two main areas. The first is con- cerned with pure numerical problems, for instance, the extremal problem (2) or the following extremal problem:

where 1/, is usually motonically increasing. The second a r e a deals with the sta- tistical properties of the estimators obtained using methods similar to (2) or

In conclusion we should note t h a t the areas of statistics described above are well developed in terms of both theory and available software.

(8)

3. STATISTICAL METHOD FOR SECONDARY MODEX CONSTAUCTION

It has already been pointed out that, although the values of variables y are known more precisely than those of variables x , we still never know the exact values of the y. If the problem is approached deterministically, then our knowledge of y takes the form of a s e t of ranges b , I y, l c, , a =

F.

Under

the probabilistic approach, on the other hand, t h s information is given in the form of a distribution function F ( y ) , w h c h assigns a confidence level to each possible value of y.

Because the exact values of y, are not available the researcher should really calculate t h e function $(z,y) using different sets of values y l j , . .

.

, yqj to see how i t fluctuates. Obviously for high-dimensional y much effort can be wasted in trying to consider every possible $(z,yi); in practice it is generally sufficient to take t h e averaged behavior together with some confidence inter- val.

If the researcher is only interested in specific aspects of t h s averaged behavior (for instance, extreme points), stochastic optimization techniques can be used (see, for instance, Ermoliev, 1976). In cases where more detailed description is necessary, an approach based on the methods described briefly in Section 2 would be more appropriate.

The main steps in building a secondary model are summarized below.

( I ) The variables included in the primary model are divided into two vectors, y and z . The permissible sets for y and z are ascertained ( y U , z EX), and t h e possibility of calculating $(z,y) for points from X x r is considered.

(2) For each calculation (measurement), the values of yi are assigned using a random numberm generator with density function F ( y ) , y ~ r . The structure of F ( y ) is usually chosen in accordance with the Bayesian approach. In

(9)

the simplest case a variable y, can take values b , or c, with equal proba- bility. It should now be obvious that any measurement yij can be described by the following regression model:

where

( 3 ) It is, of course, unrealistic to hope to find ~ ( z ) analytically in practice, but i t may be possible to obtain a suitable approximation of ~ ( z ) . Very simple approximations are usually employed a t t h s stage, e.g.,

where m can be several hundred. I t is obvious that this approximation will be very rough. But then the goal is very modest: we only wish to identify the significant variables. In the simplest case the necessary calculations can be performed by two standard statistical programs: a procedure for generating random designs and another for stepwise regression (the latter should be present in any modern package of statistical programs). This yields the numbers a l , .

. .

,a, and the estimates $, of significant parame- ters. Here "significance" has its usual statistical meaning (broadly speak- ing, a parameter is significant if its estimated value is larger than its stan- dard deviation). As pointed out elsewhere (see, for instance, Devyatkina and Tereokhin, 1981), the classical statistical methods for testing the sig- nificance of parameters (for instance, the F-test) are not appropriate to stepwise procedures; simple permutation tests should be used instead.

(10)

(4) I t is now assumed that a comparatively small number (10-20) of significant variables are known from previous steps. Let

zT

= ( x , ~ . . . . ,

X,~)ERCR".

In region j? we can use the more sophisticated approximation (compare with (8)):

q ( z , * ) ~ - * ~ f ( z )

.

The basis functions f ( z ) a r e selected using some a priori information on t h e behavior of q ( z , $ ) ; a multidimensional second-degree polynomial f T ( z ) = ( l , Z 1 , . . . .ES, FIE2, . . . ,ZN,") often gives satisfactory results. If all variables a r e of equal interest then the simulation experiments should be carried out in accordance with D-optimal design,

+

= In d e t ~ - ' (see, for instance, Fedorov, 1972). If there is some need for interpolation or extrapolation, then t h e design should minimize

+ =

max w ( z ) d [q(ZN,s)]

z E X

where d [q(z",Z)]

=

f T ( E ) ~ - l f (ZN) is the variance of the estimator

ZTf

(Z),

and w (z) is a weight function reflecting the researcher's interests. If t h e systematic discrepancy between q ( E , d ) and gTf (E) is negligible t h e n t h e function (h-l(z)

+

f T ( E ) ~ - l f ( F ) ) " ~ represents the standard deviation of t h e forecast.

The problem of optimal design has received much attention and there a r e now catalogues of optimal designs for certain standard situations (see, for instance, Brodsky e t al., 1 982) as well as some r a t h e r advanced software.

( 5 ) The final step involves t h e use of secondary models t o analyze the system.

The type of problems t h a t can be solved are indicated in Figures 1-3. Fig- u r e 1 illustrates the possibility of testing or refining the primary model by

(11)

applying the secondary model t o initial d a t a . Figure 2 explains the possi- bility of combined analysis of two systems. Figure 3 illustrates how pri- m a r y models c a n be compared through approximation by t h e same secon- d a r y model.

Initial d a t a

Figure 1. Testing a primary model by applying a secondary model t o initial d a t a .

Figure 2. The combined analysis of two systems.

t

)

C

Primary model Secondary model

I I 11

>

Secondary model I

Primary model

I -

,

(12)

Figure 3. Comparison of primary models through approximation by the same secondary model.

+

Secondary model +

rl ( x

la)

+

System Comparison

To illustrate steps 1-4 in t h e procedure described above we shall consider a simple numerical example. The simplified version of the lIASA energy supply model MESSAGE-]] was chosen a s a primary model. T h s is a linear program- ming model:

Primary model I

min

cT

U

A u C b . u % O Primary model

+d2

Twenty elements of matrix A were chosen as components of vector x . These are the first twenty entries in Table 1. Analysis of a priori information showed t h a t all of these elements could vary in a 20% range. Ten other elements of matrix A were selected a s components of vector y and assumed to have a vari- ation of 1%. The response function was +(x ,y)=min

cTu.

I1 r] ( x lfl)

T

The experiment was conducted using a two-level random design for both z and y (see step 3 from the previous section). To improve the statistical proper- ties of this design the randomization was carried out under the constraint t h a t

+ Secondary model

(13)

Table 1. Components of z and y .

the correlation between vectors zj and zk should be less than 0.1 for all j and k .

The number of measurements was determined from (3) with m =20, s =6, and 6=0.03, and turned out to be 15. Of course, this number is only a rough estimate because the assumptions under which (3) holds are not completely fulfilled.

ue hh... a te ... .l a

ue hh... a te ... .2a

ue hh... a te

....

3a uehh ... a te .... 4a uehh ... a te ... .5a

ue hh... b te .... l b

ue hh... b te ... .2b

ue hh... b te .... 3b

I

Direct electric heating The matrix elements represent

ue hh... b te .... 4b the amount of electricity used

ue hh... b te .... 5b x

ue hh... c te .... l c ue hh... c te

....

2c ue hh... c te

....

3c ue hh... c te .... 4c

ue hh... c te ... .5c u0ii

...

a to

....

l a

u0ii.. .a t0 ... .2a

1

Production of process The matrix elements represent u0ii

...

a t0

....

3a heat in industry by an the amount of electricity used u0ii

...

a t0 .... 4a electric furnace

u0ii

...

a t0

....

5a uo hh... a uh

...

a

t

Oil heating system The matrix elements represent ug hh... a uh

....

.a Gas heating system the amount of heat produced ue nh... a uh.. . .a Electric night storage per unit of input

u2ii.. .a ui.. .. .a

u2ii

...

b ui

...

b Y Oil furnaces for industrial u;?ii...c ui

...

c

I

process heat u8ii.. .a ui.. .. .a

u8ii

...

b ui

...

b Gas furnaces for u8ii

...

c ui..

...

c industrial process heat yu.si..a mu.si..b Market penetration constraint

on industrial solar heat production systems

7

The results of the stepwise regression analysis are represented by the solid points in Figure 4, where

1 2 3 4 5 6 7 8 D 10 11 12 13 14 15 16 17 18 ID 20 21 22 23 24 25 26 27 28 28 30

(14)

Figure 4. Results of the stepwise regression analysis

is the simplest measure of the discrepancy between a response function q ( x ) and its approximation. The value

002

can be considered as a n estimate of the variability of t h e response function q ( x ) (see (7)) w i t h n the operability region

The smooth decrease in :a indicates that the contribution of every com- ponent of t h e vector z is comparable to t h e "noise" arising from the variation of 7. To avoid this, the stepwise regression analysis was repeated for both x and y simultaneously. Th.e open circles in Figure 4 illustrate the results of this analysis. The final estimates of the regression coefficients were ordered by

(15)

their absolute value: the largest are given in Table 2. These coefficients correspond to the scale -1 I var z , 1, -1 I var yp I 1.

Table 2. Estimates of regression coefficients, ordered by their absolute value.

It is clear that the most significant variables are y l , y2, and y3, despite the fact that they have a variation of only 1%. In other words, comparatively small changes in these variables provide the greatest contribution to the variability of the function +(x , y ) . This indicates that the main priority should be to evalu- ate these variables.

Variable Coefficient

I am very grateful to E. Nurminski, S. Scherbor and M. Strubegger for their help in obtaining the numerical results, and also to H. Gasking for editing the paper.

y2 -10.95

Y1 -4.11

Y3 -1.43

z l l

0.025

z 1 3

-0.022

2 7

0.020

(16)

REFERENCES

Brodsky, V.Z., L.I. Brodsky, T.I. Golikova, E.D. Nikitina, and L.A. Panchenko.

( 1982) Tables of m t i m a l Designs. Metallurgy, Moscow (in Russian).

Devyatkina, G.N., and A.T. Tereokhin (1981) Statistical hypothesis tests in stepwise regression analysis based on the permutation method. In V.

Fedorov and V. Nalirnov (Eds.), L i n e a r a n d N o n l i n e a r P a r a m e t r i z a t i o n in E z p e r i m e n t a l Design P r o b l e m s , pp. 1 1 1-121. Part of the series P r o b l e m s in Q b e r n e t i c s , Moscow (in Russian).

Ermoliev, Y. M. (1976) S t o c h a s t i c P r o g r a m m i n g Methods. Nauka, Moscow (in Russian).

Fedorov, V.V. (1972) W t i m a l Design of E z p e r i m e n t s . Academic Press, New York.

Li, C.H. (1962) Sequential method for screening experiments. J o u r n a l of t h e A m e r i c a n S t a t i s t i c a l Association, Vol. 57, pp. 455-477.

Meshalkin, L.D. (1970) Validity of the random balance method. I n d u s t r i a l Laboratory, Vol. 36, pp. 316-318.

Satterthwhaite, F.E. (1959) Random balance experimentation. T e c h n o m e t r i c s , Vol. 1, pp. 111-138.

Referenzen

ÄHNLICHE DOKUMENTE

The overall objective of the CATNETS project is to determine the applicability of a decentralized economic self-organized mechanism for service and resource allocation to be used

This assumption is supported by the analysis of the timescale of the dye-protein correlation, which found, that the motion of the dye, that is due to the protein flexibility, leads

That is the final step in the formulation of model (1); screening experiments can be carried out now. 1) Input variables can be separated into groups with the help of

Finally, MS2 fragmentation scans including noise are generated using the fragmentation function of pyteomics [21] for peptides and the fragments as defined in the Modomics

[ 1 ] The collision zone of the obliquely subducting Nazca Ridge and the erosive Peruvian forearc is characterized by enhanced tectonic erosion, normal faulting, forearc uplift, and

Working Papers are interim reports on work of the International Institute for Applied Systems Analysis and have received only limited review. Views or opinions expressed herein

The alternative replication models (or at least the conservative model; see footnote 10) could not be refuted then because there was an important untested auxiliary assumption in

In other words, b is a contravariant functor, mapping the category AB of abelian groups into the category CAB of compact abelian groups—and vice versa.. Program