Mechanisms of Ageing and Development, 60 (1991) 99- I 12 Elsevier ScientificPublishers Ireland Ltd.
99
MATHEMATICAL MODELING FOR THE AGING PROCESS: NORMAL, ABNORMAL AND SELF-TERMINATING PHENOMENA IN SPATIO-TEMPORAL ORGANIZATION
MASATOSHI MURASE and MITSUYOSHI MATSUO Tokyo Metropolitan Institute of Gerontology, 35-2 Sakaecho, ltabashi-ku, Tokyo 173 (Japan)
(Received August 22nd, 1990) (Revision received February 7th, 1991)
SUMMARY An elucidation of the aging process is attempted using a simple one-dimensional multicellular system, a prototype of living organisms. This model analysis has the advantage of making us investigate the two types of modes of their dynamical behavior: (i) the local modes of behavior of individual cells and (ii) the global modes of behavior of the total system. At first, each cell is assumed to have biochemically excitable kinetics for local modes, and then what kind of the global modes results with change in intercellular interaction is examined. With a simple interaction as possibly occurs in the early stage of differentiation, the model displays wellcoordinated spatio-temporal patterns. This may be interpreted as a normal state. With a more complex interaction as possibly occurs in the late stage of growth, however, the model produces much more erratic patterns. This may refer to an abnormal state. Interestingly, these abnormal patterns can be transformed into normal patterns, when the activity of some parts of this model is turned 'off'; the system can survive at the sacrifice of its parts. This makes us imagine that programmed cell death plays an important role in development during morphogenesis. When individual cells become less sensitive to intercellular signals and still possess intrinsic excitability, then the normal patterns are developed for a short while before being replaced by rather irregular patterns. As time proceeds, however, all activity of them disappears. We call this a 'self-terminating' phenomenon, which may refer to aging. This strongly suggests that the loss of the total system function, leading to death, results from a global mode of system failure but not from a local mode of subsystem failure.
Key words: Aging; Theoretical model; Spatio-temporal organization
0047-6374/91/$03.50 Printed and Published in lreland
© 1991 ElsevierScientificPublishers Ireland Ltd.
100
INTRODUCTION
One of the most important problems in modern biology is why we and most other higher animals age. Although various theories of aging have been proposed (for reviews see Refs. 1---4), ese theories are fundamentally classified into four categories: error and damage theories [5,6], genetic program theories [7], evolutionary theories [81 and system theories [9--13]. The error and damage theories deal with unexpected chance (or stochastic events) occurring at each functional step (e.g. DNA duplication, cell division) due to intrinsic and/or extrinsic perturbations. Errors or damages are accumulated to cause malfunction in biological systems leading to aging and death [6]. The program theories, on the contrary, assume the presence of specific aging genes or longevity genes. According to these theories, aging may be a phenomenon similar to development and differentiation during morphogenesis. The features of the program theories are derived from deterministic events rather than stochastic events. Actually these two theories may account for some aspects of aging. They, however, are incomplete arguments, because the error and damage theories appear to be unable to explain species-specific lifespan, and because the program theories may not clarify the origin of such aging genes. The evolutionary theories suggest that aging is the result of natural selection through which animals have evolved to delay the time of the onset of harmful gene expressions. These theories are, however, incomplete in that what mechanism keeps time. If the time-keeping process is ascribed to a genetic clock, these theories are equal to the program theories; while if the time-keeping process is ascribed to the accumulation of errors, they are incorporated into the error theories [see 1]. The system theories are based on the concept that a biological system consists of a network of interacting subsystems. Thus it is possible for such a composite system to exhibit global modes of failure which are not associated with local subsystem failures. In other words, malfunction related to aging cannot be observed in any subsystem in isolation. It can only be defined with respect to the dynamical behavior of the system as a whole. Although these theories seem to be attractive, no simulation of aging phenomena based on spatio-temporal dynamics has been performed. The system theories prompted us to develop a simplified network model consisting of many cells which have excitable kinetics and are arranged in one-dimensional space. In the present paper, we describe that the model exhibits a rich dynamical behavior in 'terms' of cell biology. A MODEL FOR EXCITABLE KINETICS
Living organisms exhibit many types of excitability with alternate activation and inhibition such as electrical excitation in nerve membranes [14], biochemical excitation in cell membranes [ 15] and mechanochemical excitation in cilia and flagella [ 16].
101
These excitable dynamics are considered to be an important regulatory mechanism for them. Although the detailed underlying molecular kinetics must differ in different systems, their qualitative behaviors are described by a similar class of differential equations (cf. Refs. 17 and 18). We are interested in developing a formal model directly described by a particular differential equation. Let us consider an excitable kinetics model as a biochemical system. The excitable kinetics is generally described by a cubic source function f ( x ) for a substrate of dimensionless and normalized concentration x (0 _ x < 1), namely f(x)
= A(x
(1)
- x l ) ( x - xc)(1 - x )
where A is a positive constant and x I < x c • 1. A rate equation for this model is, for example, dx
--
=
A(x
-
xt)(x
-
Xc)(1
-
x)
+ b(xo
-
(2)
x)
dt
where b and x0 are positive constants. Without the cubic source function fix), this model describes that x is supplied at constant rate b x o and degrades linearly proportional to its concentration. W h e n f ( x ) is incorporated into the model as described by equation (2), however, it exhibits two types of excitable properties. The first type of excitability arises, when b is fixed, say b = 0. There are three positive steady states, two stable and one unstable. Figure 1 shows the rate of
dx dt
b=O •
0 ~
x
l
c
I
X
I b--b
c
b>b
c
Fig. I. Two types of excitable kinetics. For b = 0, there are three steady states (i.e., x = x l, x,. and 1) of equation (2). For x < x,., the system reaches x -- x l, while for x > x c it reaches x = 1. This system thus exhibits one type of excitability. If, in addition, b is allowed to vary, the two states move close together, and meet at b = b,.. If b continues to increase, the system reaches the steady state at high level of x. Thus this system exhibits another type of excitability. Parameters are: x, = 0 and x,. = 0.3.
102
kinetics d x / d t as a function of x. For all x < Xo x -- xl, while for all x > x , , x - 1. This model thus accounts for excitability in the sense that superthreshold x values trigger an abrupt increase in x. The concentration xc is a threshold concentration. The second type of excitability appears, when b is allowed to increase as a parameter. Initially the two steady states approach to each other. At a certain threshold bc these two steady states meet, and finally disappear when b > b,. Thus, for b > be, there remains a single stable steady state with a high level of concentration X.
In addition to the above-mentioned excitable kinetics, we will incorporate a c t i v a tion and i n h i b i t i o n processes into the model. In modelling such two opposed pro-
cesses, we consider two distinct subsystems I and II: one for activation and the other for inhibition. Figure 2A shows the two opposed functions denoted by FI and F H, (A) F
I
0
x~z
1
(B) n
0
I S 1
I
I
I S2
I 1
x
Fig. 2. T w o o p p o s e d excitable s u b s y s t e m s with alternate a c t i v a t i o n a n d inhibition. (A): The same cubic source functions are a s s u m e d in subsystem I a n d II, t h o u g h one is the m i r r o r image of the other. (B): Hysteresis switching function is a s s u m e d to a c c o u n t for alternate a c t i v a t i o n and inhibition processes. N o w consider w h a t h a p p e n s as x is increased from 0. Initially n is e q u a l to I. W h e n x > S 2, n becomes 0. I f x is n o w decreased, n r e m e m b e r s 0 unti! x < Sj and n becomes 1. P a r a m e t e r s are: x I = 0. x , = 0.3. x ' l = 1, x ' , . = 0.7, S I = 0.2 a n d S 2 = 0.8.
103 which corresponds to activation and inhibition. Analogous mechanisms have been considered in different situations: (i) two distinct isozymes catalyze the same substrate, but in two different ways [19]; (ii) there are two distinct genes, and each produces a repressor that turns off the other [20,21]; and (iii) two physiological systems influence one another [13]. Our model is thus viewed as a simple prototype of coupled subsystems. The complete system of equations is:
dx - - = Ft(x)n(x) + Fit(x)(1 - n(x)) + b(xo - x) dt
(3a)
Ft(x) = A(x - xt)(x - x,.)(l - x)
(3b)
Fil(x) = A(x - x'O(x
(3c)
- x',)(-
x)
where Fi and Fn are source functions representing excitable kinetics in subsystems I and II, respectively; x,. and x',. are threshold concentrations for kinetics in subsystem I and II, respectively. FI defined by equation (3b) is the same form as fix) in equation (1), by which x increases. By contrast, Fn is the mirror image of FI, by which x decreases. When subsystems I and II are both activated at the same time, no effective change occurs. To prevent this situation, we introduce an activation-inhibition switch, n, in equation (3a). The switching function, n, is represented by discrete 1 and 0 values. When n = 1 only subsystem I is activated (associated with the increase in x), while when n = 0 only subsystem II is activated (associated with the decrease in x). Furthermore, we assume that activation and inhibition depends not only on x but also on its history (i.e. hysteresis). This hysteresis switch thus accounts for ' m e m o r y ' which ensures alternate activation and inhibition. The hysteresis switch is described by the following binary function of x (see Fig. 2B). For x < 0, if initially n = 0 for x > S1, n switches from 0 to 1 at x = $1 as follows:
n(x)= [;
0S
99
MATHEMATICAL MODELING FOR THE AGING PROCESS: NORMAL, ABNORMAL AND SELF-TERMINATING PHENOMENA IN SPATIO-TEMPORAL ORGANIZATION
MASATOSHI MURASE and MITSUYOSHI MATSUO Tokyo Metropolitan Institute of Gerontology, 35-2 Sakaecho, ltabashi-ku, Tokyo 173 (Japan)
(Received August 22nd, 1990) (Revision received February 7th, 1991)
SUMMARY An elucidation of the aging process is attempted using a simple one-dimensional multicellular system, a prototype of living organisms. This model analysis has the advantage of making us investigate the two types of modes of their dynamical behavior: (i) the local modes of behavior of individual cells and (ii) the global modes of behavior of the total system. At first, each cell is assumed to have biochemically excitable kinetics for local modes, and then what kind of the global modes results with change in intercellular interaction is examined. With a simple interaction as possibly occurs in the early stage of differentiation, the model displays wellcoordinated spatio-temporal patterns. This may be interpreted as a normal state. With a more complex interaction as possibly occurs in the late stage of growth, however, the model produces much more erratic patterns. This may refer to an abnormal state. Interestingly, these abnormal patterns can be transformed into normal patterns, when the activity of some parts of this model is turned 'off'; the system can survive at the sacrifice of its parts. This makes us imagine that programmed cell death plays an important role in development during morphogenesis. When individual cells become less sensitive to intercellular signals and still possess intrinsic excitability, then the normal patterns are developed for a short while before being replaced by rather irregular patterns. As time proceeds, however, all activity of them disappears. We call this a 'self-terminating' phenomenon, which may refer to aging. This strongly suggests that the loss of the total system function, leading to death, results from a global mode of system failure but not from a local mode of subsystem failure.
Key words: Aging; Theoretical model; Spatio-temporal organization
0047-6374/91/$03.50 Printed and Published in lreland
© 1991 ElsevierScientificPublishers Ireland Ltd.
100
INTRODUCTION
One of the most important problems in modern biology is why we and most other higher animals age. Although various theories of aging have been proposed (for reviews see Refs. 1---4), ese theories are fundamentally classified into four categories: error and damage theories [5,6], genetic program theories [7], evolutionary theories [81 and system theories [9--13]. The error and damage theories deal with unexpected chance (or stochastic events) occurring at each functional step (e.g. DNA duplication, cell division) due to intrinsic and/or extrinsic perturbations. Errors or damages are accumulated to cause malfunction in biological systems leading to aging and death [6]. The program theories, on the contrary, assume the presence of specific aging genes or longevity genes. According to these theories, aging may be a phenomenon similar to development and differentiation during morphogenesis. The features of the program theories are derived from deterministic events rather than stochastic events. Actually these two theories may account for some aspects of aging. They, however, are incomplete arguments, because the error and damage theories appear to be unable to explain species-specific lifespan, and because the program theories may not clarify the origin of such aging genes. The evolutionary theories suggest that aging is the result of natural selection through which animals have evolved to delay the time of the onset of harmful gene expressions. These theories are, however, incomplete in that what mechanism keeps time. If the time-keeping process is ascribed to a genetic clock, these theories are equal to the program theories; while if the time-keeping process is ascribed to the accumulation of errors, they are incorporated into the error theories [see 1]. The system theories are based on the concept that a biological system consists of a network of interacting subsystems. Thus it is possible for such a composite system to exhibit global modes of failure which are not associated with local subsystem failures. In other words, malfunction related to aging cannot be observed in any subsystem in isolation. It can only be defined with respect to the dynamical behavior of the system as a whole. Although these theories seem to be attractive, no simulation of aging phenomena based on spatio-temporal dynamics has been performed. The system theories prompted us to develop a simplified network model consisting of many cells which have excitable kinetics and are arranged in one-dimensional space. In the present paper, we describe that the model exhibits a rich dynamical behavior in 'terms' of cell biology. A MODEL FOR EXCITABLE KINETICS
Living organisms exhibit many types of excitability with alternate activation and inhibition such as electrical excitation in nerve membranes [14], biochemical excitation in cell membranes [ 15] and mechanochemical excitation in cilia and flagella [ 16].
101
These excitable dynamics are considered to be an important regulatory mechanism for them. Although the detailed underlying molecular kinetics must differ in different systems, their qualitative behaviors are described by a similar class of differential equations (cf. Refs. 17 and 18). We are interested in developing a formal model directly described by a particular differential equation. Let us consider an excitable kinetics model as a biochemical system. The excitable kinetics is generally described by a cubic source function f ( x ) for a substrate of dimensionless and normalized concentration x (0 _ x < 1), namely f(x)
= A(x
(1)
- x l ) ( x - xc)(1 - x )
where A is a positive constant and x I < x c • 1. A rate equation for this model is, for example, dx
--
=
A(x
-
xt)(x
-
Xc)(1
-
x)
+ b(xo
-
(2)
x)
dt
where b and x0 are positive constants. Without the cubic source function fix), this model describes that x is supplied at constant rate b x o and degrades linearly proportional to its concentration. W h e n f ( x ) is incorporated into the model as described by equation (2), however, it exhibits two types of excitable properties. The first type of excitability arises, when b is fixed, say b = 0. There are three positive steady states, two stable and one unstable. Figure 1 shows the rate of
dx dt
b=O •
0 ~
x
l
c
I
X
I b--b
c
b>b
c
Fig. I. Two types of excitable kinetics. For b = 0, there are three steady states (i.e., x = x l, x,. and 1) of equation (2). For x < x,., the system reaches x -- x l, while for x > x c it reaches x = 1. This system thus exhibits one type of excitability. If, in addition, b is allowed to vary, the two states move close together, and meet at b = b,.. If b continues to increase, the system reaches the steady state at high level of x. Thus this system exhibits another type of excitability. Parameters are: x, = 0 and x,. = 0.3.
102
kinetics d x / d t as a function of x. For all x < Xo x -- xl, while for all x > x , , x - 1. This model thus accounts for excitability in the sense that superthreshold x values trigger an abrupt increase in x. The concentration xc is a threshold concentration. The second type of excitability appears, when b is allowed to increase as a parameter. Initially the two steady states approach to each other. At a certain threshold bc these two steady states meet, and finally disappear when b > b,. Thus, for b > be, there remains a single stable steady state with a high level of concentration X.
In addition to the above-mentioned excitable kinetics, we will incorporate a c t i v a tion and i n h i b i t i o n processes into the model. In modelling such two opposed pro-
cesses, we consider two distinct subsystems I and II: one for activation and the other for inhibition. Figure 2A shows the two opposed functions denoted by FI and F H, (A) F
I
0
x~z
1
(B) n
0
I S 1
I
I
I S2
I 1
x
Fig. 2. T w o o p p o s e d excitable s u b s y s t e m s with alternate a c t i v a t i o n a n d inhibition. (A): The same cubic source functions are a s s u m e d in subsystem I a n d II, t h o u g h one is the m i r r o r image of the other. (B): Hysteresis switching function is a s s u m e d to a c c o u n t for alternate a c t i v a t i o n and inhibition processes. N o w consider w h a t h a p p e n s as x is increased from 0. Initially n is e q u a l to I. W h e n x > S 2, n becomes 0. I f x is n o w decreased, n r e m e m b e r s 0 unti! x < Sj and n becomes 1. P a r a m e t e r s are: x I = 0. x , = 0.3. x ' l = 1, x ' , . = 0.7, S I = 0.2 a n d S 2 = 0.8.
103 which corresponds to activation and inhibition. Analogous mechanisms have been considered in different situations: (i) two distinct isozymes catalyze the same substrate, but in two different ways [19]; (ii) there are two distinct genes, and each produces a repressor that turns off the other [20,21]; and (iii) two physiological systems influence one another [13]. Our model is thus viewed as a simple prototype of coupled subsystems. The complete system of equations is:
dx - - = Ft(x)n(x) + Fit(x)(1 - n(x)) + b(xo - x) dt
(3a)
Ft(x) = A(x - xt)(x - x,.)(l - x)
(3b)
Fil(x) = A(x - x'O(x
(3c)
- x',)(-
x)
where Fi and Fn are source functions representing excitable kinetics in subsystems I and II, respectively; x,. and x',. are threshold concentrations for kinetics in subsystem I and II, respectively. FI defined by equation (3b) is the same form as fix) in equation (1), by which x increases. By contrast, Fn is the mirror image of FI, by which x decreases. When subsystems I and II are both activated at the same time, no effective change occurs. To prevent this situation, we introduce an activation-inhibition switch, n, in equation (3a). The switching function, n, is represented by discrete 1 and 0 values. When n = 1 only subsystem I is activated (associated with the increase in x), while when n = 0 only subsystem II is activated (associated with the decrease in x). Furthermore, we assume that activation and inhibition depends not only on x but also on its history (i.e. hysteresis). This hysteresis switch thus accounts for ' m e m o r y ' which ensures alternate activation and inhibition. The hysteresis switch is described by the following binary function of x (see Fig. 2B). For x < 0, if initially n = 0 for x > S1, n switches from 0 to 1 at x = $1 as follows:
n(x)= [;
0S
