One of t h e major drawbacks of Monte Carlo m e t h o d s is t h e i r insatiable demand for computer t i m e . Although t h e y a r e very efficient in t e r m s of t h e t i m e required by t h e analyst o r modeler to s e t u p a n appropriate s c h e m e for estimation a n d evaluation, this efficiency is t r a d e d against c o m p u t e r t i m e and, eventually, storage capacity.
There a r e a few basic rules t h a t c a n help t o make Monte Carlo tech- niques m o r e efficient in t e r m s of c o m p u t e r use.
a. M i n i m i z e the n u m b e r of trials
A reduction of t h e n u m b e r of trial ru.ns necessary t o identify a s e t of p a r a m e t e r vectors for a certain class of model response c a n be achieved in several ways. First, a given estimation problem c a n be split i n t o several cycles of trial r u n s in a n iterative way. Each cycle is analyzed before the next one is s t a r t e d . This eventually allows corrections t o be made, t h e ranges t h a t a r e t o be sampled t o be redefined, c o n s t r a i n t conditions t o be modified, e t c . After a relatively small n u m b e r of t r i a l r u n s (which certainly will depend o n t h e n u m b e r of unknowns estimated simultaneously) one m i g h t , for example, find a clear clustering of t h e "good" vectors in t h e p a r a m e t e r space already.
If, consequently, certain regions i n t h e p a r a m e t e r space s e e m "empty" ( i n
t e r m s of solutions), they c a n be discarded (by redefining t h e ranges sampled) to improve t h e efficiency of t h e sampling. Another example would be con- s t r a i n t conditions, which a r e always violated. This should lead t o t h e recon- sideration of t h e s e conditions and t h e p a r a m e t e r ranges sampled (here they might have t o be extended), o r a modification of t h e whole model s t r u c t u r e itself. Clearly, if after a first screening of t h e p a r a m e t e r space all model responses are off t h e i r t a r g e t in a systematic way (as in t h e example above), a n increase in t h e n u m b e r of t r i a l s will probably n o t be worth while.
Some intelligent check on t h e n u m b e r of r u n s c a n be made by defining complex stop rules for a cycle instead of simply using a fixed n u m b e r of tri- als. S u c h stop rules, for example, can monitor t h e m e a n s , standard devia- tions, and ranges of parameters of a c e r t a i n response class, and stop t h e esti- mation if new samples no longer change these values, i.e. when the e s t i m a t e s converge. Table 2 refers t o t h e example described above.
TABLE 2 Convergence of parameter estimates with increasing number of samples (independent cycles).
Number a b
of samples Mean Minimum Maximum Mean Minimum Maximum
b. * e e d u p the trial r u n s
Since a simulation program m a y r u n several thousand times in a Monte Carlo framework, streamlining t h e code will pay off. This includes, for exam- ple, t h e inactivation of all s t a t e m e n t s t h a t are n o t essential for t h e determi- nation of performance criteria. Examples might be auxiliary output variables, which a r e not used i n t h e testing procedure. Also, parts of t h e model t h a t a r e unchanged within a cycle of trial r u n s (for instance, setting up the geometry of t h e lake in t h e second application example, Section 3.2) should not be exe- cuted m o r e t h a n once in such a cycle. This, of course, requires more pro- gramming effort t h a n simply calling t h e entire model a s a subroutine of t h e Monte Carlo program - a compromise between programming effort and com- puter resource utilization has t o be found.
A somewhat simpler possibility is to abandon a r u n as soon as i t is obvi- ous (even during run-time) t h a t a given constraint condition will be violated.
Since this may happen within the first few time steps, savings in computer time can be considerable.
c. Reduce input/output
As even a small simulation program, when r u n several hundred or thousand times, can produce an absolutely incomprehensible mountain of output, the reduction of output is essential for more than one reason. First, t h e r e will rarely be enough space t o s t o r e i t all; second, nobody is going t o look a t i t all anyway; and third, 1/0 is time-consuming. Therefore, i t is essen- tial to reduce output t o a minimum and do whatever processing has to be done with the output (e.g. classification, and calculation of certain statistics) within the Monte Carlo program. Again, t h e r e is a trade-off between the size a program can have on a certain machine, setting an upper limit t o what can be done simultaneously, on-line, and storage capacity. Designing "intelligent"
programs for the automatic analysis of Monte Carlo runs is probably t h e most demanding - and most challenging
-
p a r t of t h e technique.Similarly, input should clearly also be reduced to the absolute minimum. The most obvious examples a r e time-variable inputs or forcings t o a dynamic simulation model, which should n o t be read a t each time step of each trial, but only once for a cycle of trials, and then stored in an appropri- a t e form within the program. Again, this calls for a compromise between time and core requirements.
d. Think f i r s t
As trivial as this last "rule" might seem, i t is probably t h e most impor- t a n t one. It is most tempting to just l e t t h e program r u n (specifically when computer t i m e is a free commodity) - and then to discover a little bug, somewhere, t h a t makes thousands of r u n s worthless. Time spent in carefully considering the estimation scheme will certainly pay off in the long r u n . For example, if the parameter ranges sampled are fairly large, most complex models a r e bound t o "crash" sooner or l a t e r - unless care is taken of zero divides, overflows, and underflows. Also, since operating systems tend to fail sometimes, provisions should be made t h a t , in case of t h e unavoidable crash, only a minimum amount of information is lost, and an estimation cycle can be restarted. The Morite Carlo approach is very forgiving and helpful in this respect, as sample runs can always be pooled.
3 APPLICATION EXA?dPLES
3.1 Hypothesis Testing: A Marine Pelagic Food-web Example*
The study of environmental systems as ecological a n d physicochemical as well a s socioeconomic entities requires a high degree of simplifying for- malism. However, a detailed understanding of a systems function and response t o various changes for t h e explicit purpose of systems management and planning still requires complex hypotheses, or models, which can hardly be subjected t o rigorous t e s t s without the aid of computers. Systems simula- tion is a powerful tool for subjecting complex hypotheses t o rigorous t e s t s of their logical s t r u c t u r e , as well as a possible means for rejecting or corrob- orating t h e underlying hypotheses.
The complexity and variability of environmental systems, the scarcity of appropriate observations and experiments, problems in t h e interpretation of empirical data, and the lack of a well established, comprehensive theoretical background make i t difficult to t e s t any possible conceptualization, o r hypothesis, describing a given system. A formal approach t o hypothesis test- ing, based on numerical simulation and Monte Carlo methods, which explic- itly considers t h e above constraints, is proposed in this section.
Based on a data s e t from the North Sea, a series of hypotheses on t h e s t r u c t u r a l relations and t h e dynamic function of t h e pelagic food web is for- mulated in t e r m s of numerical models. Hypotheses of various degrees of aggregation and abstraction a r e tested by comparing singular statements (predictions) deduced from t h e proposed hypotheses ( t h e models) with t h e observations. The basic processes of primary production, consumption, and remineralization, driven by light, heat, and advection/diffusion, a r e described in systems models ranging in complexity from two compartments to many compartments and species groups. Yearly cycles of systems behavior are simulated with each of t h e proposed models. A comparative analysis of t h e response of each of t h e models allows conclusions t o be drawn on t h e adequacy of t h e alternative hypotheses, including t h e i r "unknowns" or initial conditions (i.e. This analysis also allows one t o reject inadequate constructs, and provides some guidance on how t o improve a cer- tain hypothesis, even in t h e presence of a high degree of uncertainty.
Universal s t a t e m e n t s , describing those properties of a system t h a t are invariant in space and time, may be called models, whether they are of a n informal (e.g. verbal o r mental) or a formalized mathematical s t r u c t u r e . Such models, viewed a s scientific theories, have t o be t e s t a b l e . When one feeds or substitutes a s e t of specific singular s t a t e m e n t s into t h e model ( t h e initial conditions, which, in t h e case of a mathematical model, also include t h e model parameters in a general sense, as discussed in Section 2.2), i t m u s t be possible t o deduce or predict testable singular s t a t e m e n t s (i.e. possible observations or the outcome of possible experiments). Disagreement between t h e prediction deduced from t h e hypothesis or model a n d t h e available
*This section is largely based on Fedra (1981a, b).
observations would then require rejection of the hypothesis, modification and improvement, or the search for alternative hypotheses, which would then have to be subjected to the same procedure. This method, which would basi- cally represent the strategy of scientific research proposed by Popper (e.g.
1959), labeled falsificationism by critics such as Feyerabend (1975) and Laka- tos (1978), however, has a major drawback when applied to complex simula- tion models or dynamic hypotheses describing ecological systems, in t h a t t h e so-called initial conditions to be used with the basic s t r u c t u r e of t h e theory to deduce the testable predictions are not exactly known. In one simple example given by Popper (1959). where he refers to a mechanical experiment (breaking a piece of thread), the initial conditions to be specified are simple enough: a weight and the characteristics of t h e thread (e.g. material, diame- t e r etc.), which a r e measurable without considerable error (it is significant t h a t many examples used in epistemological analyses refer t o relatively sim- ple physical systems). Measurements "without" error, however, are not usu- ally possible when we a r e dealing with the complex aggregates conceptualized as "units" in large-scale systems thinking and models. This can certainly be seen as the result of two basic shortcomings, one in t h e measurement tech- niques available, another in the formulation of the models themselves: if t h e models require unknowns as inputs, they are not well formulated. The l a t t e r is certainly a generic shortcoming of environmental models and the underly- ing theoretical understanding.
The same line of argument can be followed with regard to the observa- tion used for comparison with model output in hypothesis testing. The break- ing of a thread, t h e singular prediction in Popper's example, is readily observ- able. It either happens, or does not. In most environmental applications, however, we have to compare predictions with measurements (as a rule, sam- ples) of t h e system, which always include some measurement e r r o r , t h a t i s t o say, these are ranges. Also, in environmental systems the degree of abstrac- tion and aggregation is quite different for measurements and for model con- ceptualization. Therefore, t h e observations and measurements can serve only a s samples of the properties or the state of the units conceptualized. As these units a r e generally heterogeneous (in terms of their measurable properties) and a r e generally characterized by a high degree of variability, further uncer- tainty h a s to be dealt with in the hypothesis-testing procedure.
Retaining t h e logical s t r u c t u r e of testing a proposed hypothesis, b u t including a t t h e same tirne t h e appropriate (or r a t h e r unavoidable) way of describing uncertain "initial conditions" as well as the expected outcome of t h e experiment, involves the following procedure. It is possible to describe regions in their respective spaces. Instead of the two vectors, we have to deal with s e t s of vectors with certain statistical properties and probability struc- tures.
To t e s t any specific hypothesis, we now examine whether, for a s e t of admissible initial conditions (i.e. the parameters), predictions (members of the s e t of allowable outcomes) can be made. The rejection of a hypothesis, whenever no allowable outcome can be generated, is based on a statistical argument, as the number of possible initial conditions forming the admissible s e t is infinite, and only samples can be examined. Also, t h e s e t of admissible rejection of inductive reasoning does not provide much help, but in practice hypotheses (and simulation models) a r e rarely generated randomly but a r e always based on empirical knowledge. However, the process of testing and rejecting a given hypothesis can also provide some diagnostic information about t h e causes of failure and about possible ways to improve the hypothesis.
One possibility is strict parsimony: t o s t a r t with the simplest possible conceptualization, or the least complex model one can formulate bona fide, which still may capture the relevant features of t h e system in view of t h e problem studied. Certainly, each hypothesis tested should be an honest candi- date for success: "What then is the point of setting up a [Poisson] model like a skittle, just to knock it down again?" (Finch 1981). If this simple version fails to give an acceptable behavior over t h e allowable parameter ranges, t h e model s t r u c t u r e is modified. Complexity is increased by adding elements and more complex process descriptions to the model (Figure 6), until a satisfac- tory behavior can be achieved. However, t h e r e is in any case more than one way t o increase t h e complexity of a model. A general formalization of this
"adding of complexity" seems to be most difficult, if not impossible. Some gui- dance for this process can be expected from t h e analysis of a series of errors, a s will be shown below. Also, as I a m only considering "conceptual" models (as opposed to purely statistical models, they are based on physical processes and only include t e r m s directly interpretable in t h e "real w o r l d ) , additional observations can be exploited in many cases. Knowledge accumulated from t h e study of similar systems may also be helpful in changing a given model s t r u c t u r e .
Building up complexity and iteratively subjecting each version or level of the model to extensive tests should allow one t o learn about t h e way struc- tural changes influence model response. At t h e same time, t h e intricate con- nection between s t r u c t u r e and t h e parameters has t o be emphasized, since model behavior is certainly responsive to both. As changes in the model s t r u c t u r e will, in almost every case, al.so necessitate changes in the parame- t e r s (their numbers, admissible ranges, and interpretation), comparisons of different versions a r e quite difficult. Although t h e approach described below is clearly far from being ideal, any a t t e m p t a t a formalization of t h e modeling
Model 1
Model 2
Model 3
FIGURE 6 Flow diagrams for t h e models compared: P, phosphate; A, phytoplankton; D.
detritus; 2 , zooplankton; Z1, herbivores; Z2, carnivores.
process seems preferable to a purely arbitrary and subjective procedure.
3 . 1 . 1 7he Empirical Background: &scribing the Environmental S y s t e m Considering t h e above constraints, t h e direct use of t h e raw data avail- able on any ecosystem seems to be r a t h e r inappropriate and difficult for t h e testing of complex and highly aggregated dynamic hypotheses. Consequently, we have to derive from t h e available data a description of t h e system and t h e processes we want to study a t a more appropriate level of abstraction and aggregation. This description, which already has to be formulated in t e r m s of t h e hypothesis to be tested, should take advantage of all t h e available infor- mation, and a t the same time provide a n estimate of t h e reliability of this information a t the required level of abstraction.
To illustrate t h e approach, a data s e t from t h e southern North Sea was used. Most of t h e information utilized s t e m s from t h e yearly reports of t h e Biological Station Helgoland, and describes physicochemical as well as biolog- ical variables a t t h e sampling station "Helgoland-Reede" for t h e period
1964-79.
Figure 7 summarizes t h e data used. The driving environmental variables, water t e m p e r a t u r e and radiation, were found suffiekntly smooth and well
behaved for a direct utilization of t h e long-term averages, approximated by simple sine waves. Data for n u t r i e n t s (PO4-P) and algae ( m e a s u r e d a's chloro- phyll a as well as in t e r m s of carbon, recalculated from counts) showed con- s i s t e n t yearly p a t t e r n s . However, when t h e year-to-year variations ( a s well a s t h e implicit sampling e r r o r s ) a r e included, t h e high variability of t h e obser- vations a s well as t h e difficulty in averaging over t i m e (several years) becomes obvious. Although t h e average phytoplankton dynamics show a sin- gle, but extended peak around July/August, t h e individual years exhibit a t l e a s t two peaks in t h e s u m m e r . As a r e s u l t of t h e i r variable timing, t h e peaks a r e averaged o u t when one looks a t t h e long-term m e a n . Also, t h e long-term m e a n is about one order of magnitude below t h e spiky peaks of t h e d a t a for t h e individual year. Little information was available on zooplankton biomass.
However, some additional information from independent experimentation, mainly on primary production, was also found. Also, t h e (time-variable) ratio of phytoplankton carbon to chlorophyll was used for t h e models described below, and approximated by a simple exponential curve.
Among t h e invariable, generalizable conditioris derived from t h e obser- vations a r e t h e following:
1. P r i m a r y producers a r e below a level of 4.0 mg m-3 chlorophyll d u r - ing t h e first t h r e e m o n t h s of the year.
2. Between Julian days 120 a n d 270 t h e r e is a t l e a s t a twofold i n c r e a s e in biomass.
3. There have t o be a t least two peaks within t h a t period, with a reduc-