Ecosystem control theory

J. theor. Biol. (1986) 121, 417-437

Ecosystem Control Theory BRUCE HANNON

Geography Department, University of Illinois and Illinois Natural History Survey, Urbana, Illinois 61801, U.S.A. (Received 7 August 1985, and in revised form 18 March 1986) The flow accounting approach is reviewed for ecosystems and expanded to a linear dynamic description. The resulting differential equations are unstable. If they are to accurately model real ecosystems, the mechanisms of stabilization must be found. Two general stabilizing schemes are discussed and applied to an oyster reef ecosystem. The first is the anticipatory control process where the system is prepared (prepares) in advance of some known external future change. The second is the reactionary or feedback control where the ecosystem is seen as taking the appropriate stabilizing steps after an unexpected external change is made. Feedback controls are used in two specific ways to stabilize the oyster reef system. The means to test the applicability of the linear dynamic approach is discussed.

In Fig. 1 P is known as the production matrix. Its elements are the amount of commodity i which is used by process j in the given time period. For example, Po could be the daily amount of algal biomass (i) consumed by a particular herbivore (j). This is a multicommodity system since different, non-commensurable commodities are listed along each of the rows. The row sums may be calculated (since they are all the same commodity and, we assume, are regarded the same by all consumers). But, in general, the sums down the columns cannot be formed. Commodities of different qualities, even though measured in the same units (e.g. g-carbon) cannot be meaningfully added together. The inputs to omnivores and detritivores, for example, are of different qualities, both chemically and in nutritional meaning to the consumer. The diagonal elements in P are the self-use terms which are for example, own-waste consumption by rabbits and the consumption of decomposers by decomposers and cannibalism. The system in Fig. 1 is shown without joint products, that is, each process (column) is assumed to produce only one commodity. Joint output system theory is discussed in Hannon (1985a) and (Costanza & H a n n o n , 1986). The vector r in Fig. 1 is composed of two types of commodities: (1) the net export (export minus import) of each produced commodity, during the given period; and (2) the amount of stockreplaced during the given period (we use the basal metabolism as a surrogate for the commodity flows which are used in rebuilding the depreciated stock). The As is a vector of net changes in the commodity stocks during the time period (zero at steady state). 417 0022-5193/86/160417+21 $03.00/0

© 1986 Academic Press Inc. (London) Ltd

418

B. HANNON

Arnounfo.fi /'~

I

'I

FIG. 1. Ecosystemflow accounting diagram. The elements of the total output vector p are the row sum of p, r, and As. Excluded from this row sum is the amount of the commodity i which is given off as heat in non-basal respiration. These waste heat flows are the elements of the vector w. If the vectors w and p are added, the result is the total input vector pN. This vector should equal the mass or energy (enthalpy) measure of the total input. But pN otherwise has no meaning in flow analysis because o f the quality difference between the waste flow and the total output. The vector w, by definition are the flows which cannot be used by the system. In case some of the heat given off by the system's respiratory processes are used by, for example, the bacteria in the system, then that portion of the heat should be assigned the status of a commodity. The analysis becomes a multi input-multi output one (Costanza & H a n n o n , 1986). The vector e is the net input" vector, the amounts of commodities absorbed by the system but not produced by it, during the given time period. If the system boundary is carefully drawn, then all flows across this boundary are either produced by the system under investigation, or they are not. If they are produced-commodities, they are classed as part of r (with the exception of w); if they are not produced by the system, they are classed as a part of e. Inside the system, flows are either between compartments or not. If they are between compartments, they are classed with p. If they are not flows between compartments, they are either contributing to the rebuilding of the stocks which decomposed during the time period (assigned to r) or they are part of the waste flow, w or they are part of the net export (export minus import). The residual flow is the net gain (or loss) in stocks, As. The entire flow assignment for each compartment can be mathematically expressed in the following manner. Since s is a vector of stocks, st is the amount of commodity i held as a stock in the system. Letting the time period become infinitesimal, the commodity balance along the row becomes S, = P , - E P0 - r , - - p~ - W , - E J

Pu-r,

(1)

J

where s~ is ds~/dt. Note again that the non-basal waste heat vector is excluded from the balance in equation (1). This exclusion indicates that w has zero value to the ecosystem when compared to p and r. Let us assume that the flow between any two components in P divided by the stocks o f the reciving component is a constant matrix G'. That is: g'u = p j s ) . Equation

ECOSYSTEM

CONTROL

THEORY

419

(1) can be rewritten ~, = p;-Y~ ~usj - ri

(2)

J

or in vector-matrix form s=p-(~s-r.

(3)

If we assume further that the ratio of the stock to flow for each component is a constant diagonal matrix B, then the elements of this stock/flow matrix are b, = sdpi and equation (3) can be further rewritten as Bp = p - C ; B p - r = (I- GB)p-

r

and if we combine the two earlier assumptions, gu = Pu/Pj, and = B-X(I - G ) p -

B-~r.

(4)

This is the "standard" form for the dynamic equation of the flow analysis approach. It is sometimes referred to as the recipient-controlled dynamic flow equation. The stability of future p when the system is subjected to changes in r or e, depends entirely on the matrix which multiplies p in equation (4). If the real part of the eigenvalues of that matrix are negative, p is stable (Luenberger, 1979; p. 158). As shown below, the eigenvalues of the matrix are always positive and future p will always respond to changes in r or e in an unstable manner. The matrices B and G can be prescribed functions of time and still be classed as linear equations, hence a rich set of solutions still applies. The matrix B in economics is called the capital-output ratio and is sometimes assumed constant in economic calculations..In the framework used here, B also represents a "turnover" time, that time needed for the complete change of the stock of a particular commodity. The turnover time plays an important role in predicting changes in concentrations of toxic materials injected into an ecosystem (see Appendix). Equation (4) can be made to respond stably by modifying r to include an anticipatory control or a feedback control in addition to the desired change in r. The feedback control can be a function of p, a n d / o r the derivative or integral of p. These controls can be thought of in two ways: first, they give the ecosystem manager the appropriate co-ordinated steps to be taken in advance or following the desired changes in r or e; secondly, the ecosystem could be thought of as possessing certain types of these controls. In this case, control theory could be coupled with ecosystem experiment to discern the mathematical form of the natural controls. The details of these control theory approaches are the focus of the next section of this paper. Before developing the control theory approaches, one should examine the steady state version of the theory and the subsequent development of a set of important weighting factors. These factors are an aid in the use of control theory and in directing the form of any non-linear dynamic equations.

420

8. HANNON

If li = 0, equation (4) becomes p = (I - G ) - l r

(5)

which is independent of the stock level. Equation (5) is the static or steady state relationship between the net and the total output. In an ecosystem with no material flows across the boundary, equation (5) indicates the amount of total output required to make up for the structural depreciation (under these conditions, r is the basal metabolism). Allowing the addition of the same commodity along the row of p and r, actually requires that each consumer of that commodity value is the same. That is, a gram of algal biomass, for example, is equally useful to each of its consumers. In economics this condition is known as perfect competition: producer and consumer units are small and numerous and possess complete information about the variables which are important to them. Such a condition is the result of a long series of adjustments and it is assumed to be the condition at the steady state. In that case, there exists a set of weights, one for each commodity, ratios of which give the relative values of the commodities. If these weights are multiplied times the commodity flows, the result can be added down the columns of P. That is, the inputs to a particular compartment can be added, providing that they are first weighted correctly. This addition allows a solution for the weights. Let these weights be e. Then the following equationst hold e,p,j +ej = Ejpj i

or in vector-matrix form, ~P+e=~t' where ^ means a diagonalized vector, and --

- e P^

--1

(I-G)

--1

.

(6)

The weighted sum down the column of P plus the net input e is equal to the weighted total output, p. This leads us to the solution for the e as equation (6). These weights are the direct and indirect amount of the net input required for the production of, or embodied in, each of the commodities. In this role they are useful in determining what changes in r will produce zero changes in the net input requirements. For example, suppose that an estuarine manager wishes to increase the fishing rate but the algae are already absorbing nitrogen at maximum rate. Since nitrogen is the scarce net input, the manager must import other commodities (bacteria, herbivores, etc.) in such a way as to require no change in the present net input. In essence, these added quantities contain the needed embodied nitrogen which balances the increased fishing demand for it. The ~ values reveal how much added nitrogen is needed and how much would be effectively embodied in the balancing imports. t In a personal communique, M. Higashi of the University of Georgia, has shown how the ~ can be derived from the balance equation directly.

ECOSYSTEM

CONTROL

THEORY

421

From equation (6), it can be seen that if the net input per unit output, and the G, are both constant, then the e are constant. The first requirement seems consistent with the constancy of G, the combined requirement is that all inputs per unit output be constant. There is also the possibility that more than one type of net input exists for a given ecosystem (e.g. sunlight and calcium). It is even possible that both of them are simultaneously scarce. In Costanza & H a n n o n (1986), we argued that the relative values of such net inputs could be calculated by reference to the global ecosystem where only one net input (sunlight) exists. The non-solar net inputs to the ecosubsystem discussed above, can all be valued (equation (6)) in terms of sunlight since they are internal commodity flows in the global system. With this method, all net inputs to eco-subsystems can be given relative value in terms of the sunlight embodied in them by the largest of ecosystems. Equation (4) is classed as a linear dynamic set of differential equations because B and G are constant. Even if these matrices were functions of time, linear solution methods are available for equation (4). If these matrices were functions of the p, equations (4) becomes non-linear and analytic solutions are rare. Non-linear equations are usually solved discretely on the computer and an objective function is usually required to narrow the range of possible solutions (e.g. to maximize the weighted biomass stored in the ecosystem). The ~ would be used as weights in such solution seeking. Because linear differential equations in the form of equation (4) have not been experimentally proven inadequate for modeling ecosystems, their initial use is preferred over the more complicated non-linear versions. Linear models are always accurate over some range of the variables. Slightly different versions are often pieced together to form an approximate solution over the whole non-linear range of the variables. For example, the speed control in the m o d e m American automobile is a series of overlapping, linear microprocesser controls. The balance equation (6) can be recast in terms of a vector of concentrations of a tracer or perhaps of a toxic material. The details are given in the Appendix. Tracers can therefore be used directly to verify the linearity assumptions made in this paper.

Controlling the Ecosystem The instability of equation (4) leaves us in a quandary. Another linear modeling systems approach, generally referred to as compartment analysis (Costanza & Hannon, 1986), is stable in a predictive mode. But as we have shown, this form of analysis requires summing down the columns P and e, which is not a meaningful step, unless the P matrix is premultiplied by E. The only legitimate linear solution for p in terms of r seems to center on equation (4). One resolution of the stability problem is to impute some sort of control structure to the ecosystem. We observe that somehow the components of an ecosystem are capable o f stable response to changes in the net input or output. The stabilizing process is manifested as " s e l f " induced changes in the net inputs and net outputs for the linear approach. A non-linear stabilizing process would also include changes

422

e. HANNON

in B and G. The exact nature of how the components of an ecosystem sense when to take what action is unknown. However, using experimental ecosystems, it may be possible to quantify the coefficients in a mathematical control theory solution to the stability problem. Thus, the style and limits of the natural controls could be determined. Hopefully, through experiments on full scale natural systems, the same styles and limits can be estimated. Then, using the same theory, man-made augmenting control can be applied by an ecosystem manager when the desired changes exceed the natural controls. The situation is somewhat analogous to the engineer who can predict the deformation behavior of a complex metal part, based on the data obtained from a simple tension test on a small cylinder of the metal. The engineer does not understand the true microscopic nature of the resistance of the material to deformation. Such behavior is captured by the tension test. An early attempt to use control theory with ecosystem management is found in Olsen (1961) who used analog computers to model the flow of trace nuclides through ecosystems. Lowes & Blackwell (1975) and Boling & Van Sickle (1975) describe in general the application of engineering control theory to ecosystems. Mulholland & Sims (1976) used the stable compartmental model to demonstrate in a general way how control theory could be used to achieve desired flow levels with linear and non-linear models. Kercher (1983) applied control theory to stable, dynamic, linear models to find their frequency response to harmonic input. Goh (1979) demonstrates why non-linear control theory can be applied to only one or two variable systems.

Control Theory Applied to Ecosystems There are two basic types of ecosystem control: anticipatory and reactionary. Anticipatory control (Hannon, 1985b, c) refers to those changes in the net inputs and outputs which happen in advance of the anticipated major change in net input or output. The signals for the preparation are from outside the system. For example, the lengthening of the daily photoperiod may signal a terrestrial ecosystem to begin the release from storage of certain chemicals which initiate the leaf forming process. By mid-season, the vegetation is prepared for the major input of sunlight and the ecosystem is said to have reacted stably. As another example, suppose that an ecosystem manager wishes to allow an increase in the fishing rate. He begins by calculating the rate of imports which he must add to the system over a designated time in advance of the chosen day for the proposed fishing increase. He then adds these imports at the calculated rate over the time period to build up the stocks of each of the components in the system. He suspends the import program and begins" the fishing increase at the end of the preparatory period, and the system responds stably. The new stock level is exactly appropriate for the new production levels which are needed to produce the original outputs plus the new outputs needed to support the increased fishing. The reactionary controls are applied after the desired change in net input or output are enacted. These controls are almost infinite in their possible mathematical form but in general they could take one or a combination of three forms. The controls could be based on a composite history of the flows, or on the flow rates

ECOSYSTEM

CONTROL

THEORY

423

themselves (flow control) or on the rate of change of the flows (acceleration control). These reactionary controls can be thought of as embedded in the ecosystem or as a tool to guide the ecosystem manager toward stable response to change. It seems very likely that ecosystems do possess some levels of control. If the mathematical examples of controls can be verified by ecosystem experiment, then those particular ecosystems will have revealed the nature of their internal control processes. Hopefully, such revelations can be generalized to other systems, and the limits of the controls can be determined. Then ecosystem managers would have a much better idea of the natural limits to exploitation and when and how to intervene from outside the system to maintain stability. It seems possible that two levels of control could exist in ecosystem. The first of these would be the macro or system level control, related to the scarce net input and controlling mainly the rate and form of succession. The second level of control would be the micro level, related to the individual or perhaps to the species and controlling the rate and form of genetic change. The anticipatory controls seem to belong to the macro level and the reactionary controls appear to operate at the micro level of the ecosystem. ANTICIPATORY

CONTROLS

Equation (4) is exponentially unstable. If the ecosystem were at steady state when a unit increase in one of the net outputs began, the p would soon show exponential increases or decreases. By choosing some finite time T which immediately precedes the planned unit increase, the stocks of the system can (theoretically) be increased through a specific import program. If the imports are sized and proportioned correctly, they can be stopped at the end of time T (i.e. at t = 0). The system flows would continue steadily, having been increased to accommodate the planned unit export. A more thorough explanation of this type of control is given in Hannon (1985b, c). The calculation of these conditioning flows, r c, for the conditioning period T, begins with the general solution of equation (4) p=-I~

eM"-*'B-lr(~-) d~"

where M = B - I ( I - G ) , and where p ( - T), the total outputs at the conditioning period, are zero, and so equation (7) represents the We also know that an infinitesimal time after t = 0 p(0 ÷) = -

I_

e-MTB-~r c d~-

(7)

beginning of the change in p. (8)

T

where r c is a constant conditioning net output for - T <~ t <~0. Since we also want a steady state to prevail from t = 0 onward, equation (5) expresses the p v s r relationship during the conditioning time period. Integrating equation (8) with equation (5) gives r c = B ( I - eMr)-lB-lr(0+), -- T<~ t<~0 (9)

424

a. HANNON

the solution for the constant conditioning net input, with r(0 +) as the constant desired increase in net input, beginning at t = 0. Substituting equation (9) into equation (7) gives p = [I--eM(t+r)][I --eMZ]-l(I-- G) -1 r(0+),

- T ~< t<~0

(10)

as the solution for the change in the total outputs, in the conditioning period, for the transition between the initial and final steady states. Equations (9) and (10) have been applied the oyster reef ecosystem model (Dame & Patten, 1981). The diagram of the oyster reef system is shown in Fig. 2.

Phytoplorikton Et suspended orgonic moiler

~'°~l

I ~°°~'"

l°°~'"

(2000)

(69-21

~5.791 '~

o.~7

Deposiled detritus (1000)

(Resuspension)

:

0-64

Deposit feeders (16-27)

~'

817¢ " -4.~4-,..,><,1 ~.: f

5"56*I

L

* Resplrohon

l/o66

1 3"~8. J ~~"L '2"~' J

MicrO3biOtO T '~°'

l o.o.

~

t'-

MeiOfO4UnO

** Mortolify flows leovlngthe system

FIG. 2. The oyster reef ecosystem. Flow units are kcal/m2/day. Stock units: kcal/m 2. The organization of these data into the accounting framework is shown in Table 1. It was decided somewhat arbitrarily, to divide the respiratory heat into equal parts for each component, to represent the basal metabolism and waste heat. No suitable set of references could be found for typical data on these important flows. The resulting total requirements matrix is given in Table 2. The vector pN is included only to show balance. The desired net output change was a unit export of the predator, compartment number (6). the change is far too large for the flows in this syste, but it was chosen for convenience of presentation. Because we are dealing with the solution of a linear

ECOSYSTEM CONTROL THEORY

425

TABLE 1

Oyster reef input-output flow matrix (P), along with vectors for net export+stock replacement respiration (r), total output excluding waste heat (p), waste heat (w), and total output including waste heat (pN) P

Oysters (1) Detritus (2) Microbiota (3) Meiofauna (4) Deposit feeders (5) Predators (6) Net Input (e)

(1)

(2)

(3)

(4)

(5)

(6)

r

P

w

pn

0 0 0 0 0 0 41-47

15.79 0 0 4.24 1.91 0-33 0

0 8.17 0 0 0 0 0

0 7.27 1.21 0 0 0 0

0 0.64 1.21 0.66 0 0 0

0-51 0 0 0 0.17 0 0

10.44 6.19 0 0 0 0-05

26.74 22.27 2.42 4.90 2-08 0.38

14.73 0 5.75 3.58 0.43 0.30

41-47 22.27 8.17 8.4 2.51 0-68

model, the results can be scaled proportionately with a reduction in this desired change. The transition total output change and the conditioning imports are shown in Fig. 3. Because the r c is a set of co-ordinated net inputs, the total output changes rise smoothly from the original steady state to the final one. That final state includes the unit increase in predator export and the comensurate increase in all of the total outputs. Note that the increase in the total output of the predator is 1.07 kcal/m2/day: one unit for the increased export and 0.07 units of associated total output due to the increased " d e m a n d s " from component 2, the total output of which had to increase by 4-79 kcal/m2/day. This demand is a somewhat artificial one since the detritus is an abiotic component. But the links in the linear framework are rigid ones. On the other hand, the increased activity of the predator export is surely associated with an increased output of detritus material. So the flow from component 6 to 2 is not a mutually agreed upon exchange. Viewed in this way the resulting change seems appropriate in form if not in quantity. These new steady state total output changes are shown in column 6 of Table 2.

TABLE 2

Matrix ( I - G)-~ for the oyster reef ecosystem for the flow analysis approach: used in equation (5). Intensities, e, calculated from equation (6) (1)

(2)

(3)

(4)

(5)

(6)

Oysters (1) Detritus (2) Microbiota (3) Meiofauna (4) Deposit feeders (5) Predators (6)

1.000 0-000 0.000 0.000 0-000 0-000

2.60 3.56 0.385 0.782 0.330 0.053

8-76 12.02 2.30 2.64 1-11 0.178

6.01 8.25 1.140 2.811 0-765 0.124

7.80 10.7 1.81 2-67 1-99 0.158

4.83 4.79 0.813 1.19 0.890 1.07

(e)

1.55

4.03

13.60

9.34

12.11

7.50

426

B.

HANNON Output change, P (skCOl/mZ/day

4.83 Oyster (1) 4.79 Detritus (2) 1-19 Meiofauna (3) 1.07 Predator (6) I

0.89 Dep.Feeder (5) 0.81 Microbrota (3)

? ' = - 5 0 days

Conddlonmg period

Time,(doys)

i

End of original steady state

~

Import

rate:kcol/m2/day

Begin new steady state

4-5

0.02

0.12

i

with 1-0 kcol/m 2/day increase m predator export

Components

7.2

I,,

0.15

5.4

Conditioning import:equation (8)

FIG. 3. Anticipatory control for the oyster reef ecosystem. Prepared for a 1-0 kcal/m2/day increase in predator export by a 50-day conditioning period.

Next, suppose that we wished to increase predator export without changing the need for any increase in the net input. Suppose further that we will balance the predator export by a coincident import of deposit feeders (it could theoretically be done with any of the components or a cost-minimizing combination of them). The net input neutralizing ratio of the two components is the ratio of their input intensities (e). From Table 2, we see that this ratio is 7.50/12-11 or 0.62; that is, if one unit of increased predator export is desired, then 0.62 units of deposit feeders must be imported to neutralize the net input demand of the export. So r(0 +) = {0, 0, 0, 0, -0.62, +1.0}, where ( + ) means export. From equation (5) the new steady state total output changes would be p = { 0 , - 1 . 8 , - 0 . 3 1 , - 0 . 4 6 , - 0 - 3 4 , +0.97}, and the conditioning net inputs during a 50 day conditioning period would be: rc= {0, +1.12, -0.002, +0.018, -0-041, -3.06}, ( + ) being export. All units for these three vectors are kcal/m2/day. An interesting net effect has occurred with predator total output. While the increased net output was 1.0, the increased total output was less, 0.97. This reduction is generated by the last two numbers in columns 5 and 6 in Table 2. The net input of 0.62 units of deposit feeders reduces the direct and indirect total output of the predator by 0.62 times 0.158, while the increased unit export of the predator requires 1.07 units of increased predator output. The net effect is a 0.97 unit increase in the predator total output, needed to supply a unit export increase of the predator.

ECOSYSTEM

CONTROL

THEORY

427

Net input conditioning can also be used to control desired changes in the net input, and conversely, the net input can be used to control the system for changes in the net output. These applications are discussed in Hannon (1985c).

REACTIONARY

CONTROLS

The anticipatory controls are a special set of net inputs or net outputs which are applied in advance of a planned or expected change in the ecosystem flow balance. In one sense they are not controls since they are removed when the desired set of stock levels have been achieved. At that time, all the former instability problems return. Theoretically, if the anticipated changes do not occur as planned, the future flows will become unstable. The reactionary controls take place at or after the time of the changes in net input or output. These controls are in fact, activated and maintained by the input or output changes. Acceleration type control If T = 0 in equation (10) then the general solution to the dynamic equation (4) is p = [ I - - e M ' ] ( I - G ) -1 r(0+),

t>~0

(11)

where M - - - B - I ( I - G ) which is unstable. Note that ( I - G ) - t r ( 0 +) is the final steady state total output vector. For stability, the e M' must decay to zero. Here it does not because the real part of at least one o f t h e eigenvalues of M is positive. This condition is guaranteed since B is a positive diagonal matrix and the diagonals of G are positive and less than one: therefore, the trace of M is positive and, since this trace is the sum of the eigenvalues, at least one of them must be positive. Stabilizing equation 11 requires changes in B a n d / o r G which make all of the eigenvalues of M become negative. This is done by adding terms to the control variables r, which in effect, simulate the control action of the natural system or which inform the ecosystem managers of their needed control action. To stabilize equation (4) let r = r(0 ÷) - Kp

(12)

where k is called the acceleration control matrix. Substituting equation (12) into equation (4) gives ii = Np - (B - K)-~r(0 +)

(13)

where N = ( B - K ) - I ( I - G ) . The matrix K must be chosen to make the eigenvalues o f N negative. The absolute values of the eigenvalues control the rate at which the flows P approach the new steady state. Therefore, the choice of the elements of K is constrained from above and below: insufficiently large elements and the instabilities remain; too large and the final steady state is approached unrealistically

428

s. HANNON

fast. Within these constraints, however, usually lies a considerable range of choice. A theorem by Froebenius? is an aid in the determination of K. Solving equation (13) and (12) for the changes in p and r due to r(0+), the desired net output p = [I - eN']( I -- G)-]r(0 +)

(14)

r = [ I - KeN'(B- K)-']r(0+).

(15)

and

Using the data given above for the oyster reef ecosystem and choosing (diagonalized for convenience) K = 1.01*B and setting r(0 +) = {0, 0, 0, 0, 0, 1.0} as before, gives, with equation (14) a reasonable view of the controlled system (see Fig. 4). 5-0

~

-

5.5

: t :\ 2~

30 2.5

/

lo

°.

O.C

..~l[~f 0

o_.-.-.~-.-:-==t=t=~7-~-:=-. i

i 5

I

I

[

I

I 10

I Time

I

]

i

I t5

-. .: : = : : i

I

i

~

i 20

i

i

='~ I

I 25

(doys)

F]G. 4. The total output for the oyster reef ecosystem responding to a unit increase in predator export under an acceleration type control.

Figure 5 gives the variation in the desired net output plus the control r from equation (15). The peak values of the r are many times higher here than they were in the anticipatory control case. The physical problems which would likely be encountered in trying to add these control flows to the system are formidable. The addition of detritus for example, would no doubt hinder the photosynthesis process and reduce the visual contact of various predator and prey. The control rates would have to be reduced (and the approach time to the new steady state lengthened) or the extraction rate of the predator could be reduced. If the K were any larger, the total output flows would have approached the new steady state more quickly and the corresponding r would have been larger. In the oyster reef system, the choice of K is so restricted that the acceleration control form is probably not a useful concept; the accuracy demands placed on the stock per unit flow levels (B) are beyond experimental means. Other styles of control must be found. ? To guarantee that the eigenvalues of N are negative, each of its diagonal elements must be negative and less t h a n the s u m of the absolute values of all other elements in its row or column (Varga, 1962).

ECOSYSTEM CONTROL THEORY i) - 1 0 L~\

I\ -20 ii /

/

-4o -I

5 m.- . . . .

10

~7 II-

\.10 /

;-3

_a,,,-~'"

fl'A--A - ' ~ 1" ~ ' ~ ' . / i /

/t -'\.

Time (days) 15 a----ACt~,~'~

20 -¢--

. i 9 ~ . ' / o . -°'"

25 - -0.1

--0.2

./;7 0/ / .(~/

,/.i /.;ii

429

- -0"3

" - -o.5

"°-d

.-o6

FIG. 5. The net output for the oyster reef ecosystem responding to a unit increase in predator export under an acceleration type control. A n e n o r m o u s p e n a l t y m u s t b e p a i d , in the f o r m o f d r a m a t i c a l l y i n c r e a s e d , n e e d for net i n p u t , if the c o n t r o l s are started with a d e l a y , after the d e s i r e d net o u t p u t c h a n g e h a s b e g u n . I f the c o n t r o l s l a g g e d the p e r t u r b a t i o n , the c h a n g e in total a n d net o u t p u t are given b y p = [I - eM'+N']

( I -- G ) - ' r ( 0

+)

(17)

w h e r e ~- is t h e c o n t r o l lag time a n d , r = [I + KNeM~*N'(I - C ) - ' ] r ( 0 +)

(18)

these two e q u a t i o n s are p l o t t e d in Figs 6 a n d 7. F r o m Fig. 6, the i n s t a b i l i t i e s in the total o u t p u t s c a n be seen d u r i n g the u n c o n t r o l l e d first 10 d a y s . In Fig. 7, t h e c o n t r o l costs c a n b e seen to e x c e e d 1000 k c a l / m 2 / d a y

5 5 --

Uncontrolled

2--"~

t~ l ~ o ' l ~ a "

I

Oil.ll_n.a.ll,ll.lbll~~T 1 "0 5 ~,

I Controlled

~ ~1~'~

-'1 101

a" a ' a ' a - a ' a

""2

4\.. -

15

i I i i i t ~ L ~ ~ i i i i i i 20 25 30 35

FIG. 6. The total output for the oyster reef ecosystem responding to a unit increase in predator export under an acceleration type control which is lagging the export by 10 days.

430

B. HANNON

1200

~od ?

I

80 60E

; 0 ~Predetor(6)= 1.0

0

-20 -

Z

-40 -60 -

-80 -100 •

Uncontrolled

Controlled

FIG. 7. The net output for the oyster reef ecosystem responding to a unit increase in predator export under an acceleration type control which is lagging the export by 10 days.

for the Detritus sector. It seems as though the control costs rise very rapidly as delays are allowed in their application. Note that here the control is acting on current information on the system condition and so the situation is always the_oretically stable in the long run, regardless of the lag time. In more realistic situations, the control would be acting on lagged descriptions of the system. In this case, the controls may be destabilizing, depending on how quickly the reaction to the controls takes effect. A general analytical solution for this situation could not be found. Such solutions could be rather easily worked out for individual cases on digital computers.

Flow type control To stabilize the dynamic system description (equation (4)), rewrite equation (12) based on the flow rather than the change in the flow r = r(0 ÷) + r0 - Qp

(19)

where ro is a vector of control constants which allow one to choose the level o f the outputs in the final steady state; Q is the flow control matrix, here again, chosen as a diagonal matrix for convenience. The elements of Q, like those of K, must be chosen with care. I f they are too large the system response is too fast and the resulting r is too large. Equation (19) is a combination of the changes in the net output due to the desired export r(0 ÷) and the control net output r o - Q p . Equation (19) can be rewritten in terms of the desired change in net output, r(0÷), and the total output, as r = [I + Q ( I - G)-~]r(0 +) - Qp

(20)

ECOSYSTEM

CONTROL

431

THEORY

5"0

-~ 4.5-

~_. 4.o 3.5 -

...-*"

3' 0 -

/,/

I

o.o"*

/* ^"~/

2.0.-

5

/~

•

•

<>.o -<>- 0/<>10"~2

I -....~_~l/4,

~= 2.5o

..-*'*'*-*-*-*-*-*-*-'*'-'*"* <>..<>_<>._<>_<>..,o.<>-o-~

.-*

~_ 1.oL-/o,~"

_ ._£_^--~-A.-~-^--^--,~-^--~-,~--^-L~.~^-_~,

12~-~:|-~-~:-~;z-z-lf~i=ii=l=l=~:.-.-.-.-.-~.--.~.-.-.

O'Oz~'t"

I

I

I

0

I

I

I

I

5

I

I

I

I

i

I

I

10

I

I

i

I

15

i

I

I

2o

25

Time (doys) FIG. 8. The total output for the oyster reef ecosystemresponding to a unit increase in predator export under a flow type control.

and the solution for the total output (equations (4) and (20) is p = (I - e L ' ] ( / - - G ) - ' r ( 0 +)

(21)

where L = B - ] ( I - G + Q ) . The elements of Q are chosen to make the eigenvalues of L negative. In developing these elements for the oyster reef, the Frobenius theorem on negative eigenvalues was again invoked, the elements of Q were not nearly as sensitive as those of K. The diagonal elements used in the examples which follow are: { - 10.5, -5.2, - 4 . 2 , - 1 . 1 , - 2 . 4 , -30}. They were chosen to cause a reasonably uniform approach by all components to their respective steady states. The "best" choice of Q and K requires an optimality theory discussed below. Equations (20) and (21) are applied to the oyster reef ecosystem and plotted in Figs 8 and 9. Time(doys) 0 5 101 , ~ .3 ~ ¢ / ~ I 4

.

2

10 i

15 I

2O t

25

o~'-o~..,,-."

o"'/~.~ 6 / ° "

-~o~"

o/

i/•

:22" FZG. 9. The net output for the oyster reef ecosystem responding to a unit increase in predator export under a flow type control.

B. HANNON

432

It is difficult to meaningfully compare the results in Figs 8 and 9 (flow control) with those in Figs 4 and 5 (acceleration control) because the control levels are chosen rather arbitrarily--first to eliminate the instabilities and second, to make the total output levels somewhat similar. Finally, suppose that an ecosystem is limited by the net input and the increased extraction of a net output is desired. By coupling the desired export with an import of a less useful product, the need for increased net input can be averted. For example, suppose the Phytoplankton input rate was limiting the oyster reef ecosystem. Any increase in the predator export would require an increased input into the Phytoplankton, unless the export was offset by an import of say, detritus. The amount of detritus needed can be found using the ~. The ratio of these intensities (see Table 2) in 7.50/4-03 or 1.86: that is, if 1-86 units of detritus is imported while 1.0 units of predator is exported, the demand for the Phytoplankton remains unchanged. The system will still respond unstably and in this example, the flow type control is used. Actually, the e are also a function of. the control, Q, rather than the simpler form used here (equation (6)). While this procedure is an approximate one, the results presented in Figs I0 and II show that it is rather accurate in this case. It is approximately correct at the beginning and becomes increasingly accurate as the flows approach the new steady state. The dynamic E equation is presented in Hannon (1985a). I f the K or Q are selected to cause N or L to satisfy the special conditions for a continuous Markov process (Freedman, 1971), the /jth element of e N' or e Lt in equation (14) or (21) can be interpreted as the probability that a particle in j came from i by time t. This view may be of some value in estimating the amount of toxic materials which collects in the components, but the optimal version of the controls msst be determined first. Changes in the net input can act as controls for desired changes in the net output and conversely, the net outputs can act as controls for changes in the net input (Hannon, 1985b, c). For example, if the net output is sinusoidal (e.g. sunlight), one can easily show that frequency of the controlling net output (of the sunlight absorbing

"° I

,

.-~ 0.5 0

....

....

.....'~'" I

-0-5

./

~

o

~

A..-A ~A--dg,-" "-"~,-- A - - ~

/~

6

~5

--'--Ae-- &-- &-- & -- ~ . - A - - A - - " - " ~ ~ "

Time (days)

-1.0 -1-5 -2.0

v~o~

O~ O~ O~O_ 0__0_/._<>_O_ 0._0 -I-84

FIG. 10. The total output for the oyster reef ecosystem responding to a unit increase in predator export and a 1.86 unit increase in detritus import under a flow type control. No change in the net input.

ECOSYSTEM

CONTROL

433

THEORY

lo I

o~. o.." ,Comp.(2) / o / Ultimutely: 5 l - = ~°~o.~,.~'1-86 Import /

~

0

/ ,/ T

5

-6

Comp. (6) Ultimolely: 1"0 Export ' ~ / ~

P,

-

-10 -15 -20 -2 5

/

/

/

-o--~-o_~._o

5

E -5 -

O~CL... 0

/

10 a ~ = , ~ " " 15a/A ~

20

25

&/a/

/~/

Components (1,3,4,5) Ultimolely: Zero

/

-50 FIG. 11. The net output for the oyster reef ecosystem responding to a unit increase in predator export and a 1.86 unit increase in detritus import under a flow type control. No change in the net input.

component) must have the same frequency but it must lag its input change. For example, if each of the net inputs are ei = ai+c~ sin ¢o~t

(22)

where e is the net input vector, a and c are constant vectors and ¢o is the frequency vector. Equations (13) and (22) combine to give the lag angle, • for each of the net outputs sin ~ , = {1 + [(N

W),/W/wi]2}- 1 / 2

(23)

where N is (B - K)-m(I - G) for acceleration type control, and where w is the assumed constant ratio of the net input to the total output, for each component. Note that constants in equation (22) do not affect the lag angle. However, the control matrix appears in the N and the lag time of the response is therefore a function of the control, K. This means that the response time to a diurnal cycle might be several days, weeks, or possibly longer. Kercher (1983) reported a similar result. The cyclic behavior of the absorbing component is not passed on to the other components because of the assumptions of a diagonal B matrix (stock per unit output), and of a diagonal control matrix. If either or both of these matrices were non-diagonal, the input fluctuation to the producer would be seen (lagged) in other components. For example, a measured cyclic response in a herbivore to a fluctuation in the net input to the producer, may mean that the control matrix actually has off-diagonal terms. In this case, the natural control may be of the flow control type, acting through the off-diagonal terms in the G matrix. In a splendid paper, Dwyer & Perez (1983) showed that an input frequency is not faithfully replicated in the outputs of a series of controlled microcosm experiments. This elegant result could mean that their microcosm has no interesting range

434

a. HANNON

of the flows over which the B and G coefficients are constant (the theory in this paper) or it might mean that these coefficients are variable functions of time, but not of the stocks or flows (i.e. still linear). It is also possible to control the system based on the integral o f the flow rates. Accordingly, the r o f equation (4) would contain the integral of the flow vector as the control element. This approach means that the control mechanism contains the flow history of the system. While difficult to handle mathematically, the approach completes a trio of linear control theories involving the flow derivative, the flow and the flow history. Hopefully, these three control forms, taken singly or in combination, constitute a rich enough variety to allow modeling of the most complex ecosystems over a reasonably long time. The choice of the control parameters such as the K or the Q matrices is still open. So far, the only criterion has been that the eigenvalues of the critical matrix be transformed to negative ones, thus stabilizing the system. But if the parameters were too large, the convergence to the new steady state would be too rapid and the net inputs would be too large for feasible application. Singh (1979) demonstrates the use of the quadratic cost function to approximate an optimal choice of externally applied flow type control parameters. A generalization of this approach may be useful for classifying various optimality hypotheses regarding the behavior of ecosystems. It is difficult to overstate the importance of recognizing that the ecosystem has its own internal set of controls with certain limits on their ability to cope with disturbances. To blindly apply a human control criterion to optimize the system cropping without basing it on the natural controls is to invite ecological catastrophe. The natural controls, their limits and optimality criteria must be incorporated into such a process. The proper choice of the control parameters requires an accurate assessment of the collective behavior of the components of the system. The basic theoretical and experimental problem is to distinguish whether or not the criteria apply at the micro or macro level of the ecosystem: that is, do the criteria apply to individuals in the system (e.g. are they maximizing the probability of survival, or maximizing the value of their stored biomass?), or do the criteria apply only (also?) to the ecosystem as a whole (e.g. is the system striving to maximize the probability o f the survival o f at least some form o f life, or trying to maximize the value o f the collective biomass?). It is clear that some form o f optimizing behavior must be adopted. But the larger question is even more important: does the ecosystem possess its own controls? If so, what is the nature of these controls and what is their limit and how do we discover them? From the theory presented here, the control matrices must have units commensurate with the B or G matrices which they modify. But does the control act only through the existing flow pathways or on the existing stocks? While this seems possible, it is easy to imagine inter- and intra-species signals which flow along unique control pathways such as chemical, visual or audio signals which mark territories. How are these signals handled in the matrix modeling process? From even the most casual of ecosystem observation, we can see that some natural systems are often perturbed and they respond resiliently--clear indication that they do possess some kind of control mechanism. It seems equally clear that the first

ECOSYSTEM CONTROL THEORY

435

p r i o r i t y o f b o t h e c o s y s t e m theorists a n d e x p e r i m e n t a l i s t s c o u l d b e the d i s c o v e r y o f the c o n t r o l s a n d with t h e m , the r e v e l a t i o n o f the o p t i m i z i n g criterion. H a n n o n (1979, 1984, 1985) has s p e c u l a t e d r e p e a t e d l y on the n a t u r e o f such a criterion. U l a n o w i c z (1984) has e l e g a n t l y p o s e d m a x i m u m a s c e n d e n c y ( m a x i m u m i n f o r m a t i o n g a i n ) as an o p t i m a l i t y criterion. O d u m (1984) has p r o p o s e d the p r i n c i p l e o f m a x i m u m p o w e r as the d r i v i n g force b e h i n d the o r g a n i z a t i o n o f the e c o s y s t e m . These are e x a m p l e s o f m a c r o b e h a v i o r a l o p t i m a l i t y h y p o t h e s e s . It is s t r o n g l y suspected t h a t the o n l y w a y to u n r a v e l such c o m p l e x q u e s t i o n s is to s e e k a b l e n d i n g o f t h e o r y a n d m i c r o c o s m e x p e r i m e n t , p r e f e r a b l y in the s a m e l o c a t i o n , a i m e d at discovery o f t h e b a s i c p r i n c i p l e s o f e c o s y s t e m c o n t r o l a n d o p t i m i z a t i o n (see e s p e c i a l l y Hill & W i e g e r t , 1980; H e a t h , 1980). As o u r u n d e r s t a n d i n g o f such b e h a v i o r increases, we can e x p a n d the scale o f such e x p e r i m e n t s to m o r e c o m p l e x a n d n a t u r a l systems. E v e n t u a l l y , we s h o u l d d i s c o v e r the n a t u r e a n d the limits o f these n a t u r a l controls. T h e n if m a n desires m o r e t h a n the e c o s y s t e m c a n safely c o m p e n s a t e itself, for, the r e s p o n s i b i l i t y for the stability o f the system falls to him. O u r u n d e r s t a n d i n g o f the b e h a v i o r o f e c o s y s t e m s c o u l d a d v a n c e m o s t r a p i d l y if system t h e o r i s t s stated t h e i r p r e d i c t i o n s in e x p e r i m e n t a l l y falsifiable form a n d , e c o s y s t e m e x p e r i m e n t a l i s t s r e g u l a r l y c h a l l e n g e d - t h e s e p r e d i c t i o n s . S u c h an i d e a l s i t u a t i o n w o u l d r e q u i r e the theorists to l e a r n m o r e b i o l o g y a n d the e x p e r i m e n t a l i s t s to b e c o m e s o m e w h a t m o r e f a m i l i a r with m a t h e m a t i c s . REFERENCES

COSTANZA,R. &. HANNON,B. (1986). Multicommodity ecosystem analysis, to appear. BARBER, M. PA'VrEN, B. & FINN, J. (1979). Review and evaluation of 1-0 flow analysis for ecological. In: Compartmental Analysis of Ecosystem Models, Vol. 10 of Statistical Ecology. (Matis, J., Patten, B. & White, G. eds). Fairland, Md: International Cooperative Publishing House. BOEING, R. & VAN SICKLE,J. (1975). Control theory in ecosystem management. In: Ecological Analysis and Prediction. (Levin, S. ed.). pp. 219-229. Philadelphia: SIAM. DAME, R. & PATTEN,B. (1981). Analysis of energy flows in an intertidal oyster reef. Marine Ecology Progress Series 5, 115. DWYER, R. & PEREZ, K. (1983). An experimental examination of ecosystem linearization. Am. Nat. 121, 305. FINN, J. (1976). Measures of ecosystem structure and function derived from analysis of flows. J. theor. Biol. 56, 363. FREEDMAN~D. (1971). Markov chains. Ch. 5. San Francisco: Holden-Day. GOH, B. (1979). The usefulness of optimal control theory to ecological problems. In: Theoretical Systems Ecology. (Halfon, E. ed.). pp. 385-399. New York: Academic Press. HANNON, B. (1973). The structure of ecosystems. J. theor. Biol. 41, 535. HANNON, B. (1979). Total energy costs in ecosystems. J. theor. Biol. g0, 271. HANNON, B. (1984). Discounting in ecosystems. In: Integration of Economy and Ecology. (Jansson, A.-M. ed.). Proceedingsfrom the Wallenberg Symposium. pp. 73-84. Stockholm: University of Stockholm. HANNON, B. (1985a). Ecosystem flow analysis. Canadian Journal of Fisheries and Aquatic Sciences 213. Ecological Theory for Biologica. Oceanography. (Ulanowicz, R. & Platt, T. eds). pp. 97-118. HANNON, B. (1985b). Conditioning the ecosystem. Math. Biosc. 75, 23. HANNON, B. (1985c). Linear dynamic ecosystems. 3". theor. Biol. 116, 1: 89. HEATH,R. (1980). Are microcosms useful for ecosystem analysis? In: Microcosms in Ecological Research. (Giesy, J. ed.). Tech. Information Center, Conf-781101. pp. 138-163. Washington: US Dept of Energy. HERENDEEN, R. (1981). Energy intensities in ecological and ecological systems, J. theor. Biol. 91, 607. HILL, J. & WIEGERT, R. (1980). Microcosms in ecological modeling. In: Microcosms in Ecological Research. (Giesy, J. ed.). Tech. Information Center, Conf-781101. pp. 138-163. Washington: US Dept of Energy.

436

B. HANNON

KERCHER, J. (1983). Closed-form solutions to sensitivity equations in the frequency and time domains for linear models of ecosystems. Ecol. Mod. 18, 209. LUENEERGER, D. (1979). Introduction to Dynamic Systems: Theory, Models, and Application. New York: Wiley. LEVINE, S. (1977). Exploitation interactions and the structure of ecosystems. J. theor. Biol. 69, 345. LEVINE, S. (1980). Several measures of trophic structure applicable to complex food webs. J. theor. Biol. 83, 195. LOWES, A. R, BLACKWELL, C. (1975). Applications of modern control theory to ecological systems. In: Ecosystems Analysis and Prediction. (Levin, S. ed.). pp. 299-305. Philadelphia: SIAM. MULHOLLAND, R. & SIMS, C. (1976). Control theory and regulation of ecosystems. In: Systems Analysis and Simulation in Ecology. (Patten, B., ed.). pp. 373-390. New York: Academic Press. ODUM, H. (1984). Embodied energy, foreign trade and welfare of nations. In: Integration of Economy and Ecology. (Jansson, A.-M. ed). Proceedings from the Wallenberg Symposia, pp. 185-200. Stockholm: University of Stockholm. OLSEN, J. (1961). Analog computer models for movement of nuclides through ecosystems. In: Radioecology. Proc. First National Symp., Colo, State University, (Schultz & Klement, A. eds). pp. 121-125. New York: Rheinhold Publishing Co. O'NEILL, R. (1979). A review of linear compartmental analysis in ecosystem science. In: Compartmental Analysis of Ecosystem Models. (Matis, J., Patten, B. & White, G. eds). pp. 3-28, Fairland, Md: Internat. Cooperative Pub. House. PA'I-rEN, B., BOSSERMAN, R., FINN, J. & CALE, W. (1976). Systems Analysis and Simulation in Ecology 4, 457. S I N G H , M. (1979). Hierarchical methods in river pollution control. In: Theoretical Systems Ecology, (Halfon, E. ed.). pp. 419-451. New York: Academic Press. ULANOWICZ, R. (1984). Growth and Development: A phenomenological Perspective. Solomons, Md: Center for Environmental and Estuarine Studies, Chesapeake Biological Laboratory, University of Maryland. VAP.GA, R. (1962). Matrix lterative Analysis. New York: Prentice Hall.

APPENDIX

Experimental Verification o f the Linear M o d e l

I f a tracer is injected at a constant rate into a steady state ecosystem (O'Neill, 1979), the tracer concentrations will reach an equilibrium level in each of the components, assuming that the tracer does not interfere with system processes. I f the flow variables for each of the commodities are chosen to accurately represent the tracer flow, the equilibrium concentrations of the tracer can be calculated. Consider the following

I

~ c~P°

Processj

Here, cj is the concentration of the tracer q in the flow pj and in the stock of process j. The ¢b is the amount of the tracer injected directly into the stock o f j . The tracer balance across process j is rl

c~p~i + qj = cjpj i

ECOSYSTEM

CONTROL

THEORY

437

or

c = q ~ - l ( I - G ) -1 The vector c is the tracer concentration in the n process stocks. The connection between the concept of embodied net input (equation (6)) and the concentration of a tracer has been made. For example, if 41.47 ~g/m2/day of a tracer were injected into the oysters in the oyster reef ecosystem in Fig. 2, the equilibrium concentrations of the tracers would be the vector across the bottom of Table 2 (e) in ixg/kcal of stock. Of course, the multicommodity output form of the result should be used (Hannon, 1985a; Costanza & H a n n o n , 1986) for actual experimental verification. The equilibrium equation must also contain the control matrices to be an accurate predictor. The comparison of data from repeated tracer equilibrium experiments should yield the form of the control matrices. By separate introduction of a tracer into each of the n components of an ecosystem, n separate versions of the equilibrium concentration equation can be determined, given As/At, s and r and equation (3). With equation (23) (given • and e), there are n ( n + l ) equations and an equal number of unknowns (G and K or Q). The dynamic control equations can then be solved and the constancy of G and K or Q can be found by comparison with results using other tracers and other ecosystems. The results can also be compared with the predictions obtained from the use of optimality hypotheses.