Abstract
Variability in cell properties can be an important driving mechanism behind spatiotemporal patterns in biological systems, as the degree of celltocell differences determines the capacity of cells to locally synchronize and, consequently, form patterns on a larger spatial scale. In principle, certain features of spatial patterns emerging with time may be regulated by variability or, more specifically, by certain constellations of celltocell differences. Similarly, measuring variability in a system (i.e. the spatial distribution of cellcell differences) may help predict properties of laterstage patterns.
Here we apply and compare different statistical methods of extracting such systematic celltocell differences in the case of patterns generated with a simple model system of an excitable medium and of experimental data by the slime mold Dictyostelium discoideum. We demonstrate with the help of a correlation analysis that these methods produce systematic (i.e. stationary) results for cell properties. Furthermore, we discuss possible applications of our method, in particular how these cell properties may serve as predictors of certain laterstage patterns.
Background
In biological pattern formation a process of selforganization and a breaking of spatial symmetry are sometimes related. In physics symmetry breaking is often triggered by random fluctuations, enabling the system to select a particular stable steady state. In biology, however, differences between the constituents of the system may in a sense predetermine the outcome of symmetry breaking and can in principle allow predicting the layout of resulting patterns. This possibility to translate cellular variability into features of patterns requires new methods of analyzing spatiotemporal data. One set of methods in the core of this endeavor, the reconstruction of variability distributions from data, is the topic of the present paper.
In theoretical studies variability (sometimes refered to as disorder) is now appreciated as a source of randomness that, similarly to noise, can interact with the nonlinearities of the system and systematically influence patterns. For noise such influences are well known: Noiseinduced transitions and even noiseenhanced structure formation have been explored in theoretical model systems [1,2] as well as in nature (see, e.g., [3] for an overview of activities in this field related to biology). Remarkably, this discussion of stochastic contributions acting constructively can be found on all spatial and temporal scales, from the genetic level to the molecular and cellular levels and up to the level of ecological pattern formation. If we discuss a multicellular system, variability can be thought of as the spatial distribution of systematic celltocell differences. In contrast to dynamic noise, biological variability is a static system property. In order to study the effect of variability on patterns in a biological system, one has to reconstruct variability distributions from observable spatiotemporal data. For better understanding the details involved in reconstruction, it is convenient to complement the design of new analysis techniques by sample data generated by mathematical model systems and then putting similar restrictions on these sample data as in the case of an actual experiment. With the help of such sample data one can test, how well the analysis tools are capable of handling experimental data [4,5]. A model system, for which a regulation by variability could be a principal mechanism (see, e.g., [6,7]), is the slime mold Dictyostelium discoideum. In this example of biological pattern formation individual cells aggregate under the influence of the chemotactic signal cAMP and form a multicellular organism [8]. Under nutrient deprivation single cells start to emit cAMP into the environment. The molecules are detected by neighboring cells via highly specific surface receptors [9]. This initiates the intracellular autocatalytic synthesis of additional cAMP by the enzyme adenylat cyclase (ACA) and the subsequent segregation into the environment. Time delayed receptor desensitization and stopping of ACA activity are involved in the following refractory period. Extracellular cAMP is degraded by membranebound and segregated phosphodiesterase, which is on the other hand regulated by its inhibitor. The coupling of the underlying reaction kinetics with diffusion results in wave propagation. As long as the local cAMP concentration increases in time, the cells react with positive chemotaxis, resulting in the periodic movement perpendicular to the wave front, i.e. towards the origin of the chemical signal. The process of aggregation is accompanied by target patterns and spiral waves of cAMP [10,11] and, at a later stage, cell streams leading from the periphery to the center of the territory [1214]. At a later stage the cells have accumulated to a cell mound in the center of the aggregation territory. The advanced stages of the developmental cycle include transitions into successive types of multicellular aggregates and finally the formation of a differentiated fruiting body with spores capable of germination (see [15] for detailed information on D. discoideum life cycle).
We believe that in the case of D. discoideum variability is responsible for certain stages of symmetry breaking in the usual course of the developmental cycle (local pattern initiation, spatial distribution of cell streams, distribution and proportions of differentiated cell types). A strong support for our hypothesis that indeed cell properties can affect the collective patterns in D. discoideum has come from the recent observation that direction and magnitude of a cell's response to a signal pulse is indeed a cell property, which remains constant in time [16] and therefore falls into the general scheme outlined within the present paper. In that work the behavior of single cells under periodic cAMP signals has been analyzed and it is observed that the characteristics of the gradient sensing response of an individual cell at a certain time point strongly correlates with that of the same cell at a later time point.
The structure of our paper is as follows: First we will formulate spatiotemporal analysis tools capable of extracting distributions of cellular properties (Section 2). Next we will study for a simple model system of excitable media (Section 3) how well predefined distributions can be reconstructed from the simulated spatiotemporal patterns with the help of these analysis tools (Section 5.1) and lastly we will apply these tools to measured wave patterns for D. discoideum (Section 4) and show that the reconstructed distributions are indeed systematic: distributions found for a particular time window strongly correlate with those distributions of the same observable obtained for other time windows (Section 5.2). In the conclusions (Section 6) we place the process of variability reconstruction into the larger context of prediction schemes.
Spatiotemporal analysis tools
CA fluctuation number
Cellular automata (CA) are a useful mathematical approach for studying, by numerical simulation, the global patterns arising in a system on the basis of certain local interactions [17,18]. In a series of previous studies spatiotemporal filters have been formulated translating CAtype neighborhood constellations into quantitative estimates of a certain system property [19]. We have used these tools to study the phenomenon of spatiotemporal stochastic resonance [20], to quantify synchronization properties of biological patterns [21] and to develop an algorithm for evaluating independently the contributions of measurement noise and internal noise to a spatiotemporal data set [22]. All these applications were based upon the temporal change of spaceaveraged observables.
Here we will formulate and study a spatially explicit variant of this fluctuation analysis and apply it to patterns from excitable media. These tools can form a basis for identifying and quantifying variability in spatiotemporal data sets from biological systems.
Let denote a twodimensional spatial data set, i.e. a square matrix of size N with components a_{ij }∈ Σ, where Σ is the set of possible states. The restriction to a square matrix is only imposed for notational convenience. A time sequence of such matrices (or "images") is a set {(t) ; t = 1, 2, ..., N_{T}}, again with some normalized (dimensionless) time t and the number of images N_{T }in the sequence. We formulate a spatiotemporal filter, which approximates the contribution of noise to an observed dynamics by studying the relative movement of neighbors of a particular cell at a time t, i.e. by looking at changes of the quantities in
which serves as a means of separating directed and undirected (eventually stochastic) change of the state of a cell. If the discretization of the spatiotemporal data set in space (due to the finite cell size) and time (due to the finite number of images) is small enough, directed and stochastic changes will have very different scales in time and space. For means of separation, one has to assume that the scales for the stochastic part will be smaller than the scales present in the discretization of the data set and the time scales of deterministic dynamics themselves. This leads to a (suffcient) condition for a manifestation of noise in a specific change at
where the last two inequalities are subsidiary conditions introduced for convenience and the sign function Sig [x] gives 1, 1 or 0, if x is larger than, less than or equal to zero, respectively. Each transition fulfilling the condition (2) gives a contribution
Averaging with respect to k leads to the final expression for the spatially explicit version Ω_{ij}(t) of the CA fluctuation number Ω(t):
where the term in the second row is either 0 or 1 filtering the dynamics according to the fluctuation condition (2). This quantity Ω_{ij}(t) attributes an estimate of the fluctuation to each point of the spatiotemporal data set.
Variants of the mutual information
The mutual information I ([23]; see also [24]) has been shown to characterize the degree of complexity and information transport in complex biological systems [25,26]. The crucial idea is to analyze correlations over a spatial or temporal distance due to an interaction (or transfer of information) between the elements.
The general form of the mutual information I compares pair probabilities p_{ab }of observing the symbol combination (ab) in the output sequence of some stochastic model with the corresponding quantity p_{a}p_{b }of a reference model of independent symbols:
where p_{a }is the probability of observing the symbol a. The quantities a and b in Eq. (5) run over the whole set Σ of possible symbol states. When maximum likelihood estimates are used for the probabilities p_{a }and p_{ab}, it is convenient to project the data under consideration onto a binary state space Σ = {0, 1}.
Variants of applying the mutual information to the analysis of some spatiotemporal pattern differ in the precise definition of the joint event (ab), which leads to the pair probability p_{ab}, and in the specific projection of the spatiotemporal pattern onto a binary state space.
Figure 1 illustrates three different ways of defining the event (ab) entering the definition (5) using experimental data of D. discoideum pattern formation. In variant (a) the temporal neighborhood is used, i.e. two incidents a and b form a joint event (ab), when they are observed at two consecutive time points. Note that this variant does not involve spatial neighborhoods. In contrast, in variant (b) the event (ab) is defined via spatial neighborhood in xdirection. Variant (c) uses an arbitrarily chosen spatial reference point to define the joint event (ab). Here we focus on variant (a). It should be noted, however, that in very noisy data the variants (b) and (c) may show a better performance in extracting variability. As in all these cases the mutual information is computed for the time courses of individual spatial points the result of analyzing a full spatiotemporal data set is a spatial matrix I_{ij}
Figure 1. Binarization of the spatiotemporal data. The aim is to achieve a format suitable for computing the mutual information. (A) shows the time course x_{ij }(t) of a single spatial point, together with the corresponding time series of temporal differences x_{ij}(t + 1)  x_{ij}(t) for 300 images from the earlystage patterns of D. discoideum. The temporal differences are inserted into a binarization filter δ_{ij}(t) = Sig(x_{ij}(t + 1)  x_{ij}(t)). If δ_{ij}(t) is zero, the value at the previous time step is used instead. For convenience one then passes to values 1 and 0, rather than 1 and 1. The lower row in (A) shows the resulting binarizations as bar code diagram. As described in the main text, variant (a) of the mutual information can be computed from this sequence of zeroes and ones. For variant (b) the binary patterns of two adjacent spatial points are considered (B). Variant (c) compares each spatial point with an arbitrarily selected reference point (C).
A simple model of an excitable medium
In order to study the link between intrinsic cellular properties and the local signature of these properties in the spatiotemporal patterns it is convenient to consider a simple model of an excitable medium given by a cellular automaton. In cellular automata the spatial discretization coincides with the discrete nature of cells and the observed states within the spatiotemporal pattern are reduced to few essential elements. In its simplest form, an excitable medium is a spatial arrangement of identical elements, for which (at least) three states exist, namely "quiescent" (excitable) (Q), "excited" (E) and "refractory" (R). A typical time sequence of states for a single cell is characterized by a switch from the quiescent state Q to the excited state E when a certain condition is fulfilled; the falling into the refractory state R after one time step and the remaining in R for a fixed period of time (called the refractory time). This sequence immediately leads to the formation of propagating wave fronts and, when disturbed, to spiral waves. In a more formal way, the three rules thus read as follows: (1) An element a_{ij }in the state Q changes to E, when E appears in its neighborhood (a_{ij}). (2) An element in the state E always proceeds to the state R in the next time step. (3) When a_{ij }= R, the element will be in the state Q in the next time step, if the element has been in the state R for a time Δt(R) equal to or larger than the refractory time τ.
The letters Q, E and R introduced above can now be used almost in a telegram style to describe the dynamics of the slime mold D. discoideum, where one has an ensemble of cells aggregating under the influence of a chemotactic cAMP signal that each cell is capable of producing and detecting: Sensing cAMP (which is, e.g., produced by one of its neighbors) the amoeba changes from Q to E (i.e., it emits cAMP itself), enters a rest phase R and starts to move perpendicular to the wave front. After a short period (of the order of minutes) the amoeba is sensitive to cAMP again, i.e., it again enters the stage Q. The existence of a refractory period, together with the interaction with neighbors (detection, emission and degradation of the cAMP signal) results in the formation of characteristic excitablemedia patterns, namely concentric rings and spiral waves [27].
From the different cellular automaton models of excitable media (see, e.g., [28,29]) we select a variant, which combines simple update rules involving few parameters with a quasicontinuous state space helpful in achieving similarity with experimental data [30]. In this model we implement a simple interaction rule for epidemic spreading of excitations from cell to cell, which allows to study the global dynamics on a quasicontinuous state space (which in an epidemic scenario would correspond to the immunization state of the elements). The state space has the form
where healthy elements are represented by 0 (corresponding to the Q state of our general model introduced above), infected cells can be in states 1, ..., n  1 (the E state), and the diseased (sick) elements are given by n (the R state).
The update rules which determine the state of a cell in the next time step are as follows:
where [x] is the integer remainder of x (i.e. the remaining integer value after removing the decimal part), a denotes the number of infected cells in the neighborhood of a point (ij) and b represents the corresponding number of sick cells in . The quantity s in the update rule for infected cells is the sum over all elements in . The remaining quantities k_{1}, k_{2 }and g are model parameters, which regulate the impact of infected and sick cells on neighbors and the excitability of a cell, respectively. In a certain range of the parameter space the model exhibits a dynamics, which is highly comparable to that of real excitable media, particularly early patterns in D. discoideum or the spatially extended BelousovZhabotinsky (BZ) reaction (cf. Figure 2).
Figure 2. Snapshot of a pattern generated by the cellular automaton described in Section 3 after 2000 time steps (A) and of excitation waves in D. discoideum (B). In (A) the excitability was set to g = 28, the other parameter values are k_{1 }= k_{2 }= 2 and n = 100. The array size is 200 × 200 cells. Elements in the excitable state (0) are displayed in black, refractory elements (n) are white. The gray values correspond to the excited state (1, ..., n  1). The section in (B) corresponds to 7.4 × 7.4 mm. The black wave fronts are cells, which have just detected the cAMP molecules. In the bright area the cells are moving chemotactically towards the respective spiral center.
We introduce variability in this system as distributions of g and of k_{2}. In both cases we use two different integer values (a high value g* and , respectively, and a low background value g^{B }and , respectively). A certain percentage of spatial sites is assigned the high value. The percentage of highvalue sites constitutes the strength of variability, while the spatial distribution of parameter values is the matrix we aim at reconstructing with the help of our spatiotemporal observables defined in the previous sections. The spatial matrices K_{ij }and G_{ij }for k_{2 }and g distributions, respectively, constitute our model implementation of cellcell differences (in two different cellular properties). In the following, ν_{k }and ν_{g }denote the respective percentage of and g* values in the parameter matrices. These quantities allow tuning the strength of variability.
Indeed, the patterns produced by the model depend systematically on the two sources of variability. Figure 3 shows typical snapshots of simulations for different strengths of variability.
Figure 3. Cellular automaton model of an excitable medium for different variablities. Snapshots of a 132 × 132 lattice after 1000 time steps are shown with n = 100, k_{1 }= 2, g^{B }= 28, = 2, g* = 40 and = 3. The percentage of highvalue elements for g has been varied between 0.1 and 3.1 percent in steps of 0.5 percent (top to bottom), while the percentage of highvalue elements for k_{2 }has been varied between 5 and 55 percent in steps of 5 percent (left to right).
At low ν_{g }the effect of ν_{k }is very pronounced: with increasing ν_{k }one has larger domains and a bias towards target waves (compared to the spiraldominated patterns at low ν_{k}). At higher values of ν_{g }the effect of ν_{k }on domain size is less pronounced, but the bias towards target waves remains. Note that our scheme of implementing variability allows us to pass to a higher variablity strength simply by inserting further highvalue sites into the previous grid. When passing from one value of variability to a higher value of variability we keep the previous pacemaker positions fixed and add the correct amount of new pacemakers randomly. Qualitatively speaking, we implement variability here as a pacemaker density, e.g., as the density of highly excitable cells. Note, however, that in our setup such pacemakers can differ from the other cells in two properties: g (corresponding to the excitability) and k_{2 }(corresponding to a sensitivity). These two contributions to the overall excitability are well known in the case of D. discoideum [31]. The matrices G_{ij }and K_{ij }store the positions of these pacemaker elements. The percentages of highvalue entries ( and g*) in these matrices correspond to the variability strengths varied in Figure 3 and discussed quantitatively in the following analyses.
In this way we can understand more clearly, how the reconstruction of cellcell differences, e.g., depends on the number of highexcitability cells in the system (rather than a new arrangement of these elements).
Experimental methods
D. discoideum cells of the axenic strain AX2 were cultivated from frozen spores (18°C) according to the standard procedure [32] in HL5medium at 21°C to a final cell density of 6·10^{6} cells/ml. To initiate starvation and therefore pattern formation the cells have been harvested by low speed centrifugation and washed twice with phosphate buffer (15.7 mM KH_{2}PO_{4}/Na_{2}HPO_{4}, pH 6.14). The resuspended cells were spread homogeneously onto water agar plates (1% Difco Bacto agar, 2 mM caffeine in phophate buffer, pH 6.14) to a density of 6.2·10^{5 }cells/cm^{2}. The supernatant was removed and the plates were incubated in darkness. After five to six hours of starvation cAMP waves become indirectly visible in darkfield [33,34]. Individual cells change their light scattering properties in dependence of the local excitation state, resulting in macroscopic wave pattern. Our darkfield optics was constructed according to [35]. Successive images were taken in equidistant time intervals of 3 sec (Hamamatsu C 3077, DTOpen Layers DT 3155 Mach Series Frame Grabber).
Results
Results for the model simulation
From Figure 3 we see that the two forms of variability described in Section 3 have a systematic influence on properties of the resulting patterns. Note that here we are not interested in studying the effect of variability on the patterns in detail. In such a case the mean value of each parameter should be kept constant. The interesting property for our present aim of reconstructing cellcell differences is the mere presence of a parameter distribution. With the more general framework in mind, it is of course noteworthy that (for this deterministic system) the parameter distribution determines the exact layout of the patterns. We believe that by virtue of the process of selforganization even in a nondeterministic (stochastic) system this correspondence can exist and inscribes a certain amount of predictability into the system.
In order to see, how well the spatiotemporal observables from Section 2 reconstruct the parameter distributions we compute the correlation coefficient between the matrices Obs_{ij }and Par_{ij}, where Obs_{ij }is either Ω_{ij }or I_{ij }and Par_{ij }can be G_{ij }or K_{ij}, as introduced in Section 3. For the array of data from Figure 3 we thus obtain the four correlation arrays given in Figure 4, which are highly systematic: The correlation between Ω_{ij }and K_{ij }remains essentially constant with ν_{g}, as it should be, but changes systematically with ν_{k}. The opposite is true for the correlation between Ω_{ij} and G_{ij }. The other observable, the mutual information I_{ij}, reacts strongest to G_{ij}, but also the comparison with K_{ij }shows a small but consistently negative correlation coefficient varying systematically within the array. In both cases the observed correlation coefficients can be quite high (up to 0.2 or 0.3) showing that our attempt of reconstructing the variability distributions with these spatiotemporal observables is rather successful.
Figure 4. Correlation coefficients between the matrices Obs_{ij }(Ω_{ij }and I_{ij}, respectively) and Par_{ij }(K_{ij }and G_{ij}, respectively) for the array of patterns shown in Figure 3. Here, ν_{g }and ν_{k }have been varied in the same ranges as in Figure 3, but with half the step size.
An important issue for practical applications of these tools is their robustness with respect to noise. In our model, noise can be implemented, e.g., by adding a random integer number to each x_{ij }in each time step with the subsidiary condition of remaining in the model's state space. Here we are using n = 100, as before, and the random integers added as noise are equally distributed between 10 and 10. Figure 5 shows an array of typical spatial snapshots under the influence of this noise.
Figure 5. Same as Figure 3, but with the influence of noise, as described in the main text. Parameter values are the same as in Figure 3.
The general features of the patterns remain the same as in Figure 3. This is important for our analysis, as the noise does not destroy the patterns completely. Note also that the general trends of variability effects remain similar to Figure 3. Details of the wave propagation and spiral formation, however, are strongly affected by the noise. Surprisingly, the reconstruction of the underlying variability distributions is not impeded by noise. One rather observes in many cases an enhancement of the reconstruction. This is particularly true for the fluctuation number Ω_{ij}, which needs a certain amount of noise to properly detect systematic differences between spatial neighbors. Noise in this case helps sample the possibility space of neighbor differences entering the filter defined in equation 4. The mutual information does not show such an enhancement, but is still capable of filtering out the noise contribution. It is also clearly seen that reconstruction of the two sources of variability differs in their sensitivity to noise: While excitability (G_{ij}) is generally well reconstructed in the noisy case, for the sensitivity (K_{ij}) we observe no enhancement compared to the noisefree case.
Figure 6 shows the corresponding arrays of correlation coefficients. This enhancement due to noise (which is reflected in higher correlation coefficients and a more systematic parameter dependence) can be seen more clearly, when the correlation coefficients are plotted as a function of one variability strength, while keeping the other constant (i.e. by looking at cross sections of the correlation arrays from Figures 4 and 6). This is summarized in Figure 7.
Figure 6. Same as Figure 4, but for the spatiotemporal patterns under the influence of noise, as shown in Figure 5.
Figure 7. Correlation coefficients between the reconstructed matrices Obs_{ij }and the parameter distributions Par_{ij }as a function of the two variability strengths ν_{g }(at fixed ν_{k }= 20) and ν_{k }(at fixed ν_{g }= 0.6), both for the patterns without noise (lefthand side) and with noise (righthand side). The curves are cross sections of the correlation arrays from Figures 4 and 6.
Even in our minimal system (two sources of variability, two observables) we already see that the observables may depend differently on different sources of variability. In reallife systems such observations may help setting up a combinatorical scheme for estimating, which source of variability contains the strongest signal for predicting laterstage patterns.
An interesting numerical experiment is to introduce a single highly excitable element in the system and then follow the reconstruction of this element in the (Ω, I) plane. In this experiment all K_{ij }are set to = 2 and all but one G_{ij }to g^{B }= 25, while a single position is assigned a different value g*, which is changed between g* = 60 and g* = 10. Figure 8 displays the corresponding reconstruction from simulated data in the (Ω, I) plane. The cloud of small black dots denotes the (Ω, I) values for all other (nonpacemaker) elements. Starting now from a very high value of excitability g* = 60 for the single pacemaker element we gradually reduce the value of g* and observe the corresponding trajectory in the (Ω, I) plane. Over a wide range of g* (both at high and low values) this single pacemaker element stands out in terms of these reconstructed values. Only when g* is close to the values of the background elements, the pacemaker is within the cloud formed by the other elements.
Figure 8. Reconstruction of a single pacemaker element in the (Ω, I) plane as a function of its pacemaker strength g*, which was changed between g* = 60 and 10. The excitability of the background elements was g^{B }= 25. All other parameters are the same as in Figures 2 and 3 (at = = 2).
Application to D. discoideum pattern formation
The patterns of an early stage of the D. discoideum life cycle, where cellcell communication leads to propagating waves, are an ideal field of application for the reconstruction methods of cellcell variability described in the previous sections. As the darkfield images do not allow observing individual cells directly, we apply the analysis tools to each pixel and assume that the observables provide estimates valied as avarage for the cells residing at this spot. Cell movement is neglected in this analysis. This is not unrealistic, as at this stage of pattern formation directed movement is small. Figure 9 shows snapshots of corresponding experimental data sets. In order to see, if the observables from Section 2 indeed yield reproducible results, even though they essentially evaluate the systematics of neighborhood fluctuations behind the overall selforganized dynamics in a spatiotemporal data set, we compute Ω_{ij }and I_{ij }for three different time intervals and then see, if the reconstructed matrices correlate over time. Figure 10 summarizes the general scheme. Inspite of their respective focus on smallscale fluctuations (cf. the definitions of I_{ij }and Ω_{ij }in Section 2) the two observables show a very systematic result, which suggests that the individual pixels possess a specific dynamic response, even though the system as a whole displays a selforganized pattern with a high spatial order on a larger scale (Figure 11): For consecutive intervals (1, 2) and (2, 3) the correlation coefficients are almost identical, while they are (in most cases) systematically reduced for a larger time difference (1, 3). The result from Figure 11 complements nicely the singlecell observations from [16]. While these authors look at individual D. discoideum cells under welldefined stimuli, we analyze statistically a very large ensemble of cells in the process of pattern formation. In this way, our result is a cellpopulation variant of the findings in [16]. It is surprising that the individual cell properties contribute strongly and systematically enough to show up in this analysis.
Figure 9. Snapshots of experimental data sets analyzed on their spatial distribution in cellcell differences (bar size 2 mm). Time points are indicated above the array of snapshots. In addition the spatial size the experimental data differ in their resolution: (A) and (B) 22.6 pixels/mm, (C) 68.2 pixels/mm, (D) 68.0 pixels/mm, (E) 38.6 pixels/mm, (F) 53.3 pixels/mm. As a rule of thumb at a cell density of 6.172·10^{5 }cells/cm^{2 }one can expect that 1 pixel contains approximately 12 cells in (A) and (B), 1 cell in (C) and (D), 4 cells in (E) and 2 cells in (F).
Figure 10. Schematic view of the procedure of individual cell property extraction. Data sets were devided into intervals of 200 images corresponding to 10 minutes in the experiments considered here. For each interval the observables Ω_{ij }and I_{ij }were determined. For this particular data set (denoted (C) in Figure 9) the corresponding observables (for a particular image segment) are shown below the time axis: first row – Ω_{ij}, second row – I_{ij }. One sees from these distributions that the two observables focus on smallscale fluctuations rather than the largescale features of the original patterns. In the subsequent analysis the three correlation coefficients of the reconstructed matrices (first and second time interval, second and third, and first and third) are computed both for Ω_{ij }and I_{ij}.
Figure 11. Correlation coefficients of Ω_{ij }and I_{ij}, respectively, between the different intervals of experimental data. On the lefthand side of each image segment one can see the absolute values of correlation coefficients between intervals one and two (Corr(1,2), cf. Figure 10). The columns on the righthand image parts show correlation coefficients between the intervals normalized to Corr(1,2). The notation (A) to (F) corresponds to that of Figure 9.
Conclusion and outlook
The aim of the paper is twofold: First, we want to introduce the general idea that spatial distributions of cellular properties may serve as a basis for predicting later stages of a pattern formation process; second, we provide some tools for reconstructing such distributions from spatiotemporal data sets. The following results have been shown: The fluctuation number and the mutual information are adequate tools for extracting spatial variations in cell properties from experimental data; different cell properties are picked up differently by these tools (a point we discuss in detail for excitability and sensitivity); reconstruction of such cell properties by the fluctuation number is under certain conditions enhanced by noise; applying these tools to experimental data of D. discoideum pattern formation we find that, in agreement with the singlecell results from [16], individual loci in the spatial pattern possess systematic properties, which contribute to the pattern formation process.
While we are, at this stage of the investigation, far from explicitly formulating such a prediction scheme and applying it to specific biological scenarii, we nevertheless believe that this new view on pattern formation is of relevance to many biological systems. We do not imply that the techniques lined out here are the only means of reconstructing variability distributions from data. In fact, the choice of tools will strongly depend on the specific system at hand. Criteria for selecting the analysis tools could be: Is noise an important factor in the data? Which forms of variability are expected to be relevant? Which is the dominant form of dynamics; which types of patterns are formed?
For generic models of excitable media and, to a certain extent, for experimental data on D. discoideum pattern formation we have shown that simple spatiotemporal observables can be employed to reconstruct the spatial distribution of biological properties involved the pattern formation process. How does one pass from these methodological findings to a practical application of the prediction scheme lined out in the introduction? As long as there is no general theory connecting such distribution with the layout with laterstage patterns one has to rely on heuristic protocols. It should be possible for some systems where a large enough number of data sets is available to train an artificial neural network or a system of selforganizing maps to the task of linking distribution patterns with properties of the laterstage selforganized state. Which features of selforganized patterns can in principle be thought of as predictees of earlystage distributions in cellular properties? When one focuses on excitable media clearly the specific locations of phase singularities (i.e. the terminal points of spiral wave fronts) are a candidate for such predictive schemes. Particularly because recent work on D. discoideum relates the density of such phase singularities with the strength of a genetic feedback loop at the single cell level [31], it would be of huge interest to understand how a distribution of cell properties in this specific scenario translates into a distribution of phase singularities.
The new view on biological pattern formation as a consequence, beyond the general rules of the process, of distributed properties goes along with the potential of predicting features of the laterstage patterns. It should be noted, however, that such predictions will always be of statistical nature and never onetoone correspondences. Nevertheless, for the range of biological systems, for which such predictions could be of interest (ranging from excitation waves in heart tissue to calcium waves on cellular membranes and waves of epidemic spreading of diseases), even the possibility of assigning e.g. the probability of a phase singularity to a particular position has interesting implications for intervention schemes and risk assessment.
References

GarciaOjalvo J, Sancho JM: Noise in spatially extended systems. New York: SpringerVerlag; 1999.

Jung P, MayerKress G: Spatiotemporal stochastic resonance in excitable media.
Phys Rev Lett 1995, 74:21302133. PubMed Abstract  Publisher Full Text

Beck F, Hütt MTh, Lüttge U, Eds: Nonlinear dynamics and the spatiotemporal principles of biology. Nova Acta Leopold. Volume 332. Halle (Saale); 2003.

Hütt MTh: Datenanalyse in der Biologie. Berlin, Heidelberg: SpringerVerlag; 2001.

Schreiber T: Interdisciplinary application of nonlinear time series methods.
Phys Rep 1999, 308:164. Publisher Full Text

Levine H, Aranson I, Tsimring L, Truong TV: Positive genetic feedback governs cAMP spiral wave formation in Dictyostelium.
Proc Natl Acad Sci USA 1996, 93:63826386. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Lauzeral J, Halloy J, Goldbeter A: Desynchronization of cells on the developmental path triggers the formation of spiral waves of cAMP during Dictyostelium aggregation.
Proc Natl Acad Sci USA 1997, 94:91539158. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Weijer CJ: Morphogenetic cell movement in Dictyostelium.
Semin Cell Dev Biol 1999, 10:609619. PubMed Abstract  Publisher Full Text

Dormann D, Kim JY, Devreotes PN, Weijer CJ: cAMP receptor affinity controls wave dynamics, geometry and morphogenesis in Dictyostelium.
J Cell Sci 2001, 114:25132523. PubMed Abstract  Publisher Full Text

Höfer T, Sherratt JA, Maini PK: Dictyostelium discoideum: cellular selforganization in an excitable biological medium.
Proc Biol Sci 1995, 259:249257. PubMed Abstract  Publisher Full Text

Palsson E, Cox EC: Origin and evolution of circular waves and spirals in Dictyostelium discoideum territories.
Proc Natl Acad Sci USA 1996, 93:11511155. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Levine H, Reynolds W: Streaming instability of aggregating slime mold amoebae.
Phys Rev Lett 1991, 66:24002403. PubMed Abstract  Publisher Full Text

Kessler DA, Levine H: Pattern formation in Dictyostelium via the dynamics of cooperative biological entities.
Phys Rev E 1993, 48:48014804. Publisher Full Text

Höfer T, Maini PK: Streaming instability of slime mold amoebae: An analytical model.
Phys Rev E 1997, 56:20742080. Publisher Full Text

Kessin RH: Dictyostelium. Evolution, cell biology, and the development of multicellularity. Cambridge: University Press; 2001.

Samadani A, Mettetal J, van Oudenaarden A: Cellular asymmetry and individuality in directional sensing.
Proc Natl Acad Sci USA 2006, 103:1154911554. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Wolfram S: Cellular automata as models of complexity.
Nature 1984, 311:419424. Publisher Full Text

BarYam Y: Dynamics of complex systems. Reading: AddisonWesley; 1997.

Hütt MTh, Neff R: Quantification of spatiotemporal phenomena by means of cellular automata techniques.
Physica A 2001, 289:498516. Publisher Full Text

Hütt MTh, Neff R, Busch H, Kaiser F: A method for detecting the signature of spatiotemporal stochastic resonance.
Phys Rev E 2002, 66:026117. Publisher Full Text

Rascher U, Hütt MTh, Siebke K, Osmond CB, Beck F, Lüttge U: Spatiotemporal variation of metabolism in a plant circadian rhythm: the biological clock as an assembly of coupled individual oscillators.
Proc Natl Acad Sci USA 1996, 98:1180111805. Publisher Full Text

Busch H, Hütt MTh: Scaledependence of spatiotemporal filters inspired by cellular automata.
Int J Bif Chaos 2004, 14:19571974. Publisher Full Text

Shannon CE: A mathematical theory of communication.
Bell System Tech J 1948, 27:379423.
623–656

Ebeling W, Freund J, Schweitzer F: Komplexe Strukturen, Entropie und Information. Stuttgart: TeubnerVerlag; 1998.

Herzel H, Trifonov EN, Weiss O, Grosse I: Interpreting correlations in biosequences.
Physica A 1998, 249:449459. Publisher Full Text

Weiss O, Herzel H: Correlations in protein sequences and property codes.
J Theor Biol 1998, 190:341353. PubMed Abstract  Publisher Full Text

Winfree AT: The geometry of biological time. New York: SpringerVerlag; 1980.

Mikhailov AS, Calenbuhr V: From cells to societies: models of complex coherent action. Springer Series in Synergetics. Berlin, Heidelberg: SpringerVerlag; 2006.

Gaylord RJ, Wellin PR: Computer simulations with Mathematica. Explorations in complex physical and biological systems. Santa Clara: TELOS; 1995.

Dewdney AK: Computer recreations: the hodgepodge machine makes waves.

Sawai S, Thomason PA, Cox EC: An autoregulatory circuit for longrange selforganization in Dictyostelium cell populations.
Nature 2005, 433:323326. PubMed Abstract  Publisher Full Text

Sussman M: Cultivation and synchronous morphogenesis of Dictyostelium under controlled experimental conditions.
Methods Cell Biol 1987, 28:929. PubMed Abstract

Alcantara F, Monk M: Signal propagation during aggregation in the slime mould Dictyostelium discoideum.
J Gen Microbiol 1974, 85:321334. PubMed Abstract

Tomchik KJ, Devreotes PN: Adenosine 3',5'monophosphate waves in Dictyostelium discoideum: a demonstration by isotope dilutionfluorography.
Science 1981, 212:443446. PubMed Abstract  Publisher Full Text

Gross JD, Peacey MJ, Trevan DJ: Signal emission and signal propagation during early aggregation in Dictyostelium discoideum.
J Cell Sci 1976, 22:645656. PubMed Abstract  Publisher Full Text