Ecological Modelling, 47 (1989) 291-302 Elsevier Science Publishers B.V., Amsterdam - Printed in The Netherlands
ANALYSIS OF MATHEMATICAL
291
PRINCIPLES
IN C R O P G R O W T H S I M U L A T I O N M O D E L S
P. RACSKO Computing Center of the Eotvos Lorand University, 1117 Bogdanffv u. l O/b, Budapest ( Hungao')
and M. SEMENOV Department of Mathematical Modelling in Ecology and Medicine, Computing Center, USSR Academy of Sciences, 40 Vavilova str., Moscow 117333 (U.S.S.R.)
(Accepted 9 March 1989)
ABSTRACT Racsko, P. and Semenov, M., 1989. Analysis of mathematical principles in crop growth simulation models. Ecol. Modelling, 47: 291-302. Computer simulation is one of the most effective methods in the analysis and estimation of agricultural crop production. However, very little is known about the growth mechanism of the plants above the level of physiological processes. The unknown mechanisms are substituted by statistical models and experimental data. This fact raises the question of portability of the models to an environment different from the one used for identification, Many analysts formulate heuristic principles and build models on them. A relatively simple adaptation principle which can be used effectively in building crop growth models and carry out simulation experiments is presented.
INTRODU CTION I n the p a s t few years, several t e a m s h a v e d e v e l o p e d m o d e l s d e s c r i b i n g the g r o w t h of the b i o m a s s of a g r i c u l t u r a l p l a n t s ( D e W i t et al., 1970, 1978: D u n c a n , 1971; M o l d a u , 1975; W a n g et al., 1977; W e i r et al., 1984: J o n e s a n d Kiniry, 1986; S e m e n o v et al., 1988). T h e m o d e l s h a v e b e e n d e s i g n e d to help b o t h the r e s e a r c h t e a m s to a n a l y s e the r e l e v a n t biological, ecological a n d t e c h n o l o g i c a l p r o c e s s e s a n d the f a r m e r s to select c r o p varieties a n d a g r i c u l t u r a l t e c h n o l o g y , a n d to e s t i m a t e the a n n u a l a n d l o n g - t e r m costs. T w o b a s i c classes of s i m u l a t i o n m o d e l s are - the t e r m ' s i m u l a t i o n ' is used as m a t h e m a t i c a l r e p r e s e n t a t i o n of the m o d e l static a n d d y n a m i c . G e n e r ally s p e a k i n g , static m o d e l s ( M o l d a u , 1975; H u e y et al., 1981; W e s t a n d 0304-3800/89/$03.50
© 1989 Elsevier Science Publishers B.V.
292
r . R A C S K O A N D M. S E M E N O V
C O z concenlTation
Radiation --9----
Temperature
~
t . t Precipitatlon
Soil water content
~ / ~ - ~
- ~ ~ ~ ~ ~ ---.~
Plant gro~rth
Available N,P,K
Fig. 1. A generalized plant growth model.
Lauenroth, 1982; Palmer, 1983) consist of a few (sometimes one) equation(s) which help to predict the yield or loss of yield due to weather conditions or technology applied. The input variables of the equations are genetic parameters, and soil and weather conditions. Most crop growth simulation models, however, represent a system which changes in time. These models are called dynamic. Dynamic simulation models describe the growth of biomass of the roots, stem, leaves and grain or the whole plant. Some crop growth models cover a relatively large number of physiological, physical and agrotechnical processes, as C E R E S Maize (Jones and Kiniry, 1986), T A M W (Maas and Arkin, 1985) and CROP (Semenov et al., 1988), while other processes are connected with the models by surrogate parameters ('stress factors') estimated by means of statistical methods. The general structure of growth models is presented in Fig. 1. It is clear that the basic process to be modelled is the plant growth. The model of plant growth has to govern and coordinate the submodels which describe the surrounding environment. Thus we have arrived at the basic problem of plant growth simulation how to design the model of the main process? In plant growth modelling as in a large number of simulation models from other areas - a natural approach is to apply statistical analyses and build regression models. The independent variables are the environmental and technological parameters, while the dependent variables are yield and plant biomass. The time series used for model identification are usually real field data or data measured in controlled experiments, where the effect of individual, manageable factors (such as water or nitrogen supply) are observed, while other factors are maintained constant.
CROP GROWTH SIMULATIONMODELS
293
Formally, most plant growth simulation models can be written as a systems of differential equations:
;~,(t) =f,(x,, )~2(t) =~(x,
. . . , Xn, ~, . . . . . ~,,, X,, . . . , Xk) .....
x . , ~1 . . . .
n(t) = L ( x , . . . . . x . .
. . . . .
, btm, X, . . . . .
Xk)
x, . . . . .
xk)
(1)
or the relevant difference equations with initial data: x,(0) = xO
(i= 1.....
n)
where x~(t) represent a quantitative characteristic of an organ of a plant (root, leaf, stem or grain biomass, or leaf, root, etc., surface); /~ deterministic data, genetic parameters, technological data (sowing data, quantity of nutrients applied, etc.); and k i uncontrollable or partly controllable factors such as weather parameters and quantity of nutrients in the soil. Let us call the functions fl, -.-, fn growth functions. Unfortunately, our present knowledge of plant physiology does not allow to find explicit growth functions, and perhaps we shall never know the real deterministic growth functions which transform microscopic physiological processes to macroscopic data. Researchers try to avoid this problem in different ways. The most popular and in some sense most effective method is to define fl . . . . . . /, on a basis of a large number of empirical data (CERES-Maize, TAMW, CREAMS, Knisel, 1980). Supposing a high level of professional knowledge and enough computer time, the model will provide the required level of precision of approximation of the empirical data not considering the instability of too large models due to computing errors. But without simulating the internal plant growth mechanisms, why should we believe that the model will work at strongly non-similar conditions? As a matter of fact, most crop simulation models fail to work at external conditions, which are largely different from the original environment. Another approach to the definition of function fi is the use of general biological principles. Further on we present an overview of some principles as they are used in plant growth models and give qualitative results of simulation experiments of a model with 'principles' in the description of the growth dynamics. We start with the discussion of some elementary principles, which are widely used in plant growth simulation. LIEBIG'S P R I N C I P L E A N D STRESS FACTORS
Eiebig, in the 19th century, introduced the principle of limiting factors (Liebig, 1885). It was constructed as an analogy to stechiometric relations in
294
P. RACSKO AND M. SEMENOV
chemistry. Liebig's principle says that the biological processes are determined by those factors in which quantity is at the minimum level. Let us give a formal definition. Let H be a biological process, formally a mapping of a subset of space of input variables to the space of output variables: n ( ~ , , ..., ~ , ) ~ (~,, ..., ~m) AS an example, let us consider the process of photosynthesis:
(7,)
.....
where ~, is the quantity of available water; ~2 the quantity of available carbon-dioxide; ~3. . . . . ~t are status variables describing the status of the plant and the environment; and ~7, is the quantity of new photosynthate produced during a time unit. Let FIp~h be the projection of Hph on the k t h coordinate axis. Liebig's principle says that: 7,=Ilph(f,,
-- . , f , ) = m i n {
rtph('
,),
...,
II~h(. ) and II2ph(-) are defined from the simple stechiometric equation of photosynthesis: 6 HzO + 6 CO 2 ~ C6Hl206 + 6 0 z Unfortunately (or fortunately ) , macroscopic biological processes very often are not governed by direct stechiometric laws. For example, it is a well known fact that nitroge n and water cannot be considered as independent limiting factors of the plant's growth. In the general case: . . . .
,
Now, let us analyse the concept of stress factors. Most crop models include a central growth block which simulates the crop growth assuming an optimal environment, and water, temperature, and nutrient blocks for computing independent water, temperature, and nutrient stress factors, which are zero to unity values. Then, the optimal production is multiplied by the minimum of stress factors. Thus, system (1) is simplified as follows: ~i=f~.(Xl, ..., Xn, ~ll,, ..., /lm)'Si(~k , . . . . . x,(0)
=x °
(i=1
.....
~kk)
(2)
n)
where
Si(~k' . . . . .
~k) = min{ S ~ ( ~ k l )
i , $2(~k2)
.....
i Sk(kk) }
and sk(kk) i are the stress factors for the ith state variable, generated by the uncontrollable parameter ~k.
CROP GROWTH
SIMULATION
MODELS
295
Now it is obvious, that the use of stress factors is nothing else but the formalization of Liebig's principle. It has been pointed out that macroscopic processes do not always follow Liebig's principle for those input factors taking part in microscopic processes. Consequently, in our opinion the concept of stress factors is one of the most vulnerable points in simulation models. Perhaps models with stress factors can be identified on time series if the uncontrollable factors are close to optimum, however they will have little power of explanation or extrapolation when they are far from optimum. For further facts see Odum (1971) and Hanson et al. (1985). STATISTICAL MODELS AND METHODS We know very little about the biological laws and mechanisms of the macroscopic development of living organisms; they may not even exist in the sense of physical and chemical laws. A quite natural way to find these 'laws' is to use statistical methods to fit the functions of system (1) to experimental data. The simplest way is the application of regression methods. Multiple regression analysis, however, works only in case of a very high level of aggregation of status variables and parameters. Empirical data suggest that the growth dynamic itself changes in time. Thus, if the model is supposed to give a detailed description of the growth dynamics, the vegetation period is divided into several time intervals in such a way that within an interval the crop has a 'homogeneous' law of growth. Empirical data and - first of all - common sense suggest that the growth function must be homogeneous within one phenophase. In fact, all crop growth simulation models use a phenophase-partitioning of the vegetation period and apply a different growth function to each period. Most growth functions are based on simple principles of distribution of new photosynthate, e.g., during the phase of vegetative growth, the biomass of the roots, stem and leaves must grow proportionally or, during the phase of the growth of reproductive organs, they have an absolute priority. In each phenophase the parameters of the growth functions are determined by statistical methods from experimental data. It is obvious that the mechanism (growth function), statistically fitted to experimental data, expresses the will and professional level of the analyst rather than the real growth mechanism of the plant - if there is any. Why should one suppose that these 'pseudo'-mechanisms are more or less invariant to the enviromental conditions, which are different from those used in the statistical analysis? One does not know the weight and significance of the terms in the statistical growth functions. Moreover, it is impossible to check all coefficients by experiments, as sometimes they are not natural or measurable data.
296
P. R A C S K O A N D M. S E M E N O V
HEURISTIC MODELS, HEURISTIC PRINCIPLES
One of the basic problems in modelling biological objects is a lack of knowledge about their behavior on the macroscopic level. The only positive information is the existence of a genetic program, which is nothing else but a set of specific reactions to the infinite number of different constellations of environmental conditions. However, we possess information - empirical and experimental data - on only a relatively small number of possibilities. Now, how to fill the gap between the small - at least finite - number of observed cases and the infinite number of possibilities? Scientists in similar cases apply general 'laws' or 'principles'. We all know several physical and chemical laws and a few biological ones. In biology, starting from Darwin, laws are based on the "presumption of existence of reason" (Gorban and Chlebopros, 1988), i.e. biological objects behave in such a way as if they followed a purpose. This assumption is very flexible. One can always check the assumed behavior against empirical data and modify and refine the concrete theory. This approach has been practiced by several authors when building plant growth models. Let us see some characteristic representatives - without pretending to have a complete collection of 'principles' used in plant growth models. Perhaps no author thinks that plants have 'reason', memory, objectives or a built-in computer to carry out sophisticated optimization programs with a large number of variables. Principles are only used as effective modelling tools. Let us give some examples: Vernadsky promoted the idea that the surface of green leaves of each plant tends to the maximum and the leaves are distributed in space in such a way that no light ray can avoid the microscopic mechanism of energy transformation (Vernadsky, 1967). Paltridge recommended the following principle (Paltridge, 1972): in each time step the maximum quantity of new assimilates is distributed to that layer of the above surface part of the plants or to that layer of roots which provides the maximum quantity of new assimilates, assuming that the enviromental conditions of the previous time step (day) remain unchanged. The same assumption was used by Sirotenko (1981). Moldau (1975) formulated another possible principle of optimal distribution of new assimilates to the leaves, roots and reproductive organs. According to Moldau the new photosynthate is distributed in such a way as to provide a maximum reproductive mass by the end of the vegetation period. The distribution of the new assimilates must provide the optimal rate Sr/S~, where Sr is the surface of the roots, and S~ the surface of the leaves. The optimal rate depends on the moisture content in the soil: S r / S t increases when the soil moisture content decreases.
('ROP GROWTtt SIMULATION MODELS
297
Cohen (1971) developed a mathematical model of growth and seed production in which growth is limited by time or limiting resources. The model is based on an optimal strategy of distribution of new assimilates. It has been proved that optimal strategy is switching over from total vegetative growth to total reproductive growth, assuming that the growing period is known and deterministic. For a random growing period the optimal transition strategy is gradual. Further on, we describe our crop growth simulation model, which was developed on the basis of a heuristic optimality principle. We postulated an optimality criterion which expresses the plant's ability to adapt to the changing environment, rather than its 'intention' to produce an optimal final result. CROP MODEL The model was developed for the simulation of the growth of spring barley. The input data are the daily average temperature, water content of the soil, and quantity of photosynthetically active radiation. The output data set consists of data describing the growth dynamics, i.e. dynamics of roots, stem, leaves and grains. The simulation model was identified for spring barley in the area of Kursk. This is a typical chernoziom area; all nutrients are present in sufficient quantity. This fact makes the nutrient block unnecessary. The identification was carried out on two independent data sets, observed in two different areas for several years.
Conceptual scheme The green biomass of the plant, mainly the leaves, produces new biomass under the impact of photosynthetically active radiation, provided other enviromental conditions such as available CO2 and water, necessary temperature etc. are present. The new assimilates are then distributed between the different organs of the plant and increase their biomass. Energy, necessary for maintaining life processes, such as transport of materials, biochemical processes etc., is generated from the same pool of assimilates by reduction. Thus, the total increase of biomass of the plant's organs is a result of two processes: assimilation and respiration. The intensity of these two processes depends on the biomass of the organs of the plant, their geometry, age of tissues, and external factors. The conceptual growth scheme is illustrated in Fig. 2. The unknown mechanism of the distribution of the new assimilates, however~ does not allow to 'close' the model. Now a heuristic adaptation
298
P. R A C S K O A N D M. S E M E N O V
)
I Photosynthesisand[ respiration
1
I
i I Vi
new
assimilates
[
environmental parameters
o(
$
$
$
$
organs
Fig. 2. Conceptual scheme of the plant growth.
principle will be introduced which describes a possible low of distribution. Of course, we cannot formally prove the truthfulness of the following principle, but as far as it was succesfully identified in several models (tree growth: Racsko, 1978; cotton growth: Sadulloev and Tarko, 1983; wheat and barley: Semenov et al., 1988) it can be accepted as a working modelling tool. Let us formulate the adaptation principle: The new assimilates, produced during a time step by the plant, are distributed in a way which provides the m a x i m u m growth during the next time step. We shall use the following notation: t - discrete time step (day); y(t) - intensity of new biomass production, i.e. the newly produced assimilates, minus respiration loss; x l ( t ) , X z ( t ) , X3(t), X4(t ) - leaves, stem, root and grain biomass on the unity of area. External parameters: vl(t ) - intensity or F A R ( c a l m - 2 d a y - l ) ; Vz(t ) - parameters, describing the water regime of the plant (day m 1); v3(t) - CO 2 concentration in the atmosphere (g m - 3); v4(t) - air temperature ( ° C).
FAR, photosynthetically active radiation.
299
CROP G R O W T H S I M U L A T I O N MODELS
N o w we formulate the balance equations of the plant growth: x i ( t + 1) = { x i ( t
) + ei(t ) y[X(t), V(t)]}(1 - wi)
with initial data: xi(0) = x °
(i=1 . . . .
,4)
where
X ( t ) = ( x a ( t ) , x 2 ( t ) , x3(t), x4(t)) V(t)=(v](t), v2(t), v3(t), "4(t)) and w i is the loss rate of leaves, stem, roots and grains. The coefficients e(t) give the distribution rate of new biomass, if we make the following restriction: 4
Y'. eift) = 1
(O ~ eift) <.%l, i= l , . . . , 4)
i--1
In order to solve the system of equations (2) we have to know the coefficients ei(t ) for every time step t. In most crop models e~(t)'s are given explicitly (Giliamin, 1981; Sirotenko, 1981). We divide the growth period into two phases. The first phase lasts until the appearance of reproductive organs, while the second period ends with the harvest. The beginning of the second phenophase is calculated from temperature data in temperature days. For the first phenophase, coefficients ei(t) are calculated using the following optimality principle: The new biomass in each time step is distributed between the leaves, stem and roots in such a way, that the plant could produce the m a x i m u m quantity of new assimilates during the next time period provided the enviromental conditions do not change. After the appearance of grains the m a x i m u m possible quantity of new assimilates is distributed to the grains, while the remaining biomass is added to the leaves, stem and roots, in accordance with the optimization principle. In formal language, the vector (el(t), . . . , e a ( t ) ) is the solution of the following problem of mathematical programming:
y [ X ( t + 1), V( t )] --, max el(t) +
e2(t ) + e3(t ) = 1 - e4(t)
ei(t)>~O e4(t)=~rrfin(1, L x 4 ( t ) / y ( t ) } , ,~
if if
t
or
t>t*
300
P. RACSKOAND M. SEMENOV
TABLE 1 Dependence of the yield on environmental factors: Simulation results for spring barley Phenophases
Main factors
Optimal conditions
Unfavorable conditions
double ridge
temperature soil water
temperate cold and wet weather
cold or hot and dry weather
tillering
radiation soil water
sunny and wet weather
drought or low light intensity
start grain filling
radiation temperature
temperate cold and sunny weather
hot weather or cold or drought
maturity
weak dependence on external conditions
Vegetatioe period
Reproductive period
where L is a genetically determined constant, t * the time of grain appearance, e is a small constant (say 0.01) providing positive biomass of grains at time t. The growth function y ( t ) is taken from (Moldau, 1971): avis t -
av 1(casl + ... +
Y = 1 + a v , + slb(1 + v)v2/(
C4S 4 ) 3s3) -(RIS
+ ...
+R4s4)
where a is the derivative of the curve of photosynthesis at very low light intensity; s~ surface of leaves, stem, roots and grains; c, are coefficients of photorespiration of leaves, stem, root and grains; b is the rate of molecular diffusion of CO2; y the coefficient connected with the stomatal resistance; and R i are coefficients of 'dark' respiration. Then y ( t ) is multiplied by a bell-shape function of the temperature, expressing the effect of non-optimal temperature on the photosynthesis. The plant's geometry competition and soil parameters are taken into consideration at the definition of environmental parameters and thus they are built into a formula (4) rather then being stress factors (for more details, see Semenov et al., 1988). Table I illustrates the qualitative results of the computer simulation experiments with different weather scenarios. CONCLUSION
In our opinion, the 'heuristic' approach to crop growth simulation is as powerful as statistical models in quantitative predictions and can be successfully used for comparative quantitative analysis even for extreme scenarios, where statistical and empirical models fail.
CROP GROWTH SIMULATION MODELS
301
REFERENCES Cohen, D., 1971. Maximizing final yield when growth is limited by time or by limiting resources. J. Theor. Biol., 33: 299-307. De Wit, C.T., Brouwer, R. and Penning de Vries, F.W.T., 1970. The simulation of photosynthetic systems., In: I. Setlik (Editor), Production and Measurement of Photosynthetic Productivity. Pudoc, Wageningen, The Netherlands, pp. 47-70. De Wit, C.T., Goudriaan, J., Van Laar, H.H., Penning de Vries, F.W.T., Rabl~inge, R.R., Van Keulen, H., Louwerse, W., Simba, L. and Jonge, C., 1978. Simulation of assimilation, respiration and transpiration of crops. Pudoc, Wageningen, The Netherlands. Duncan, W.G., 1971. SIMCOT, a simulator of cotton growth and yield. In: Proc. Workshop Tree Dynamics and Modeling, Duke University, Durham, NC, pp. 115 118. Giliamin, E.P., 1981. Optimizacia Operativnogo Raspredelenia Vodnih Resursov v Orosenii (Optimization of Distribution of Water Resources for Irrigation). Gidrometeoizdat, Leningrad, 276 pp. (in Russian with English summary). Gorban, A.N. and Chlebopros, R.G., 1988. Demon Darvina. Ideia Optimalnosti i Estestvennij Otbor (Darwin's Demon. The Idea of Optimality and the Natural Selection). Nauka, Moscow, 207 pp. (In Russian with English summary). Hanson, J.D., Parton, W.J. and Innis, G.S., 1985. Plant growth and production of grassland ecosystems: a comparison of modelling approaches. Ecol. Modelling, 29: 131-144. Huey, B.A., Wells, B.R., Chapman, S.L. and Boston, N.R., 1981. DD50 system for predicting management practices for rice. Rice Inf. AGR-6, 48 (Rev. 6). Arkansas Cooperative Extension Service, 150 pp. Jones, C.A. and Kiniry, J.R. (Editors), 1986. CERES-Maize. A simulation model of maize growth and development. Texas A&M University Press, College Station, TX, 194 pp. Knisel, W.G. (Editor), 1980. CREAMS - a field-scale model for chemicals, runoff and erosion from agricultural management systems. Conserv. Res. Rep. 26, U,S. Department of Agriculture, Washington, DC, 640 pp. Liebig, J., 1885. Die Grundsatze der Agricultur-Chemie. Brunnswick. Maas, S.J. and Arkin, G.F., 1980. TAMW: a wheat growth and development simulation model. Texas Agricultural Experiment Station, College Station, TX, 60 pp. Moldau, H., 1971. Model of plant production at limited water supply considering adaptation. Photosynthetica, 5: 16-21. Moldau, H., 1975. Optimalnoe raspredelenie assimilatov pri deficite vodi (Optimal distribution of assimilates at water shortage). Izv. Acad. Sci. Estonia Ser. Biol., 2 4 : 3 - 9 (in Russian with English summary). Odum, E.P., 1971. Fundamentals of Ecology (3rd Edition). Saunders, Philadelphia, PA, 574 PP. Palmer, J.H., 1983. Microcomputer program to aid soybean growers and others in variety selection. Agron. Abstr., p. 23. Paltridge, G.W., 1972. Experiments on a mathematical model of a pasture. Agric. Meteorol., 10: 39-44. Racsko, P., 1978. Imitacionnaja model dereva kak elementa lesnogo biogenocenoza (Simulation model of the tree growth, as the element of the forest biogeocenosis). Vopr. Kibern. 52 (1979), USSR Academy of Sciences, Moscow, pp. 73-110 (In Russian with English summary). Sadulloev, R.1. and Tarko, A.M., 1984. Matematiceskaja model rosta i razvitija chlopchatnika c ucjotom azothogo pitanija (Mathematical model of the growth and development of a cotton field including nitrogen uptake). Comput Centre, USSR Academy of Sciences, Moscow, 39 pp. (in Russian).
302
r'. R A C S K O A N D M. S E M E N O V
Semenov, M.A., Polenok, S.P. and Denise~ako, 1988. Model agrocenosa iarovoi kulturi (Model of the agrocenosis of spring crop cultures). Working paper, Computing Centre of the USSR Academy of Sciences, Moscow, 29 pp. Sirotenko, O.D., 1981. Matematiceskaja Modelirovanije Vodnoteplogogo Rezhima i Produktivnosti Agroecosystem (Mathematical Modelling of Water and Heat Regime and the Productivity of Agroecosystems). Gidrometeoizdat, Lengingrad, 167 pp. (in Russian with English summary). Vernadsky, V.I., 1967. Biosphera - Izbrannije Trudi po Biogeochimii (The Biosphere Selected Papers on Biogeochemistry). Misl, Moscow, 376 pp. (in Russian). Wang, Y., Gutierrez, A.P., Oster G. and Daxl, R., 1977. A population model for plant growth and development: coupling cotton-herbivore interaction. Can. Entomol., 109: 1359-1374. Weir, A.H., Bragg, R.L., Porter, J.R. and Rayner, J.H., 1984. A winter wheat crop simulation model without water or nutrient limitations. J. Agric. Sci. Cambridge, 102: 371-382. West, W.H. and Lauenroth, W.K., 1982. Estimating long-term forage production in rangelands. In: W.K. Lauenroth, G:V. Skogerboe and M. Flug (Editors), State of the Art in Ecological Modelling. Developments in Environmental Modelling, 5. Elsevier, Amsterdam, pp. 443-450.