Looking for the simplest models for adaptive branching, one has first to consider physical phenomena that can exhibit such a behaviour. Among them, diffusion-based models have been studied very early, for instance in relation to crystal growth. Diffusion can be coupled to reaction (after A.M. Turing) or to other phenomena like accretion, and involve a number of diffusible species. The biologically oriented models of Gierer and Meinhard [1] used two diffusible species. In [2,3] we described how to get adaptive branching with a single one, but the interesting range of parameters was limited and the phenomenon rather unstable. In the sequel, we will consider a model less directly relevant to physical epigenesis but very robust and perhaps better adapted to the study of branching in abstract biological spaces. We still view this model as of an epigenetic flavor in that we do not use any ‘preformation’ or genetic information.

Here we introduce the algorithm in a simple setting, which is not the most general one. The evolution takes place in the usual p-dimensional space ${ℝ}^p$. A fixed subset of it, finite or not, is called the target and noted C, with a probability defined on it. A second subset of ${ℝ}^p$, which we call the network, is going to be modified (it will ‘grow’) at each step of the algorithm. We shall note R_n the network after step n. Starting with an initial network R₀, each step will add a point to the network.

To start, one needs to define an initial network R₀, a target C with its probability distribution, and a real number $ϵ\in ]0,1[$. At each step, say i, we:

• – randomly draw a point, say a_i, from the target (using the given probability distribution);
• – look for b_i, the point of the network R_i−1 which is closest to a_i;
• – compute b′_i=ε·a_i+(1−ε)·b_i;
• – ‘graft’ b′_i on the set R_i−1 to get R_i=R_i−1∪{b′_i}.

The stopping criterion can be the number of steps, or some condition on the network. The case of a tie for the choice of b_i formally calls for a procedure, either stochastic or deterministic, to choose b_i, but those cases are very unprobable (we will not develop that issue here). The formula to compute b′_i uses vector (or affine) algebra; typically ε is taken small, making b′_i close to b_i but slightly displaced in the direction of a_i.

Since each point added to the network is related to a definite point of the former network, the network has the structure of a finite tree if R₀ was a single point (and a finite set of finite trees in general). Branching points can thus be defined as points with more than one point grafted, besides ‘continuation points’ (points with only one graft) and ‘tip points’ (points with no graft). One could include in the definition of the network the line segments joining points to their grafted sons and compute closest points on this object, but we will not do it.

Fig. 1 shows the result of the algorithm for a target being a line segment with uniform probability distribution, starting from a single point, and for ε=0.001. The network grows toward the target and branches several times in a way that looks adapted to the geometry of the target.

Starting from two points as the initial network can lead either to cooperation in a ‘market sharing’ way or to a ‘winner takes all’ competition, depending on the situation of the two seeds and the geometry of the target (Fig. 2).

Adaptation to more complex shapes of targets is also seen. In Fig. 3, the network adapts to two different targets (i.e. two components of the target) the one on top being a circle with uniform probability and the bottom one being a disk with uniform probability, with the disk getting a total of four times more probability than the circle. On that figure the target has not been drawn. Notice the different morphologies of the network in the two cases.

This algorithm can be used in higher dimensional spaces. Fig. 4 shows on the left a network growth from a point seed in ${ℝ}^6$ towards a 5-dimensional hypercube target. What is shown is the projection of the network onto a plane perpendicular to the target and which contains the seed. On the right of Fig. 4, the network is shown growing in ${ℝ}^{50}$ towards a 49-dimensional hypercube with the same conventions. Superimpositions due to the projection get of course more numerous in that case.

The former experimental results are supported by the following two theorems, which show the adaptive power of the algorithm. To state them, we slightly modify our definitions. Instead of defining first the target and then a probability on it, we start with a probability on ${ℝ}^p$ and define the target as the set of points of ${ℝ}^p$ for which any neighbourhood has positive (i.e. >0) probability (this set is known as the support of the probability measure). We shall write B(x,r) for the open ball $\{\mathrm{y}\in {ℝ}^p\mid \mathrm{d}(\mathrm{x},\mathrm{y})<\mathrm{r}\}$, where x is the centre and r is the radius, and d is the Euclidean distance.

The first theorem states that any neighbourhood of the target will be eventually approached by the network as close as we want with a probability as close to 1 as we want. To state it formally, we choose a ball B(x,r) intersecting the target, and a slightly larger concentric ball B(x,r′), together with a probability value 1−η as close to 1 as we wish. We then state that, under these conditions, there is always a definite number of steps after which the probability that some points of the network have entered the larger ball is at least the probability value chosen in advance.

The second theorem states that when the target is bounded, all its points can be simultaneously approached by points of the network as close as we want, with a probability as close to 1 as we want if we wait long enough. More formally, for any distance value r>0 and any probability value 1−ε, there is a definite number of steps after which the property that all the points of the target will be at a distance less than r from the network will hold with a probability greater than 1−ε.

Branching occurs when some grafting takes place at a point of the network, say A, which is not a tip. This implies that A is closer to some point of the target than any point on the branch already grafted on A.

Fig. 5 shows an instance of successive abortive bifurcations before a macroscopic branching appears. To see them, one should use small ε and magnify the network at a distance from the target where a macroscopic branching is likely to occur. Notice the constant qualitative pattern of the network at each abortive bifurcation: after one of the branches wins over the other one, i.e. ‘hides it’ from the target (the winner being closer than the other to any point of the target, the loser cannot grow anymore), the network bends back to the centre of the target. In some contexts like Evolution, such repetition of abortive attempts of differentiation might be looked for and perhaps used to detect the imminence of a major branching.

To understand why abortive branching happens repetitively before a major branching can remain, one can restrict the study to the simple 2-dimensional case with a linear target. The evolution of the two tips, say X and Y, that arise from a branching can be summarized by the evolution of the vector Y−X. Each of the two tips hides the other one and wins when the vector crosses one of two lines (we call them critical) related to the ends of the target. The evolution of the vector is stochastic, but the expectation of the dynamical system that governs it can be computed easily from the target geometry and the position of the tips. Fig. 6 shows the phase portraits of the dynamical system acting on the vector and the critical lines for two different distances to the target. When far from the target, all the trajectories eventually cross one of the winning lines. Getting closer to the target modifies the picture, and a growing number of trajectories are able to make their way without hitting either of the critical lines, meaning that the branching did not abort. The real vector field that acts on the couple of tips can be viewed as the sum of that deterministic (expected) field and a stochastic component that can be made small with ε small. The whole picture thus involves a bifurcation stochastically perturbed, the deterministic part of which depends on the position relative to the target. It fits with the idea of a ‘generalized catastrophe’ proposed for study by Thom in [4].

At each step of the algorithm, one has a to find a closest point within the network. But many points may become irreversibly inactive and it is computationally more efficient to periodically detect and mark inactive points to avoid testing them for proximity (one can also maintain a list of active points: several variants and data structures can be used). A network point x gets inactive when the active points y₁,…,y_k of the network collectively ‘hide it’ from all the target points, i.e. when for any target point c there is always a y_i such that d(y_i,c)<d(x,c). Sometimes a simple criterium can find some of the inactive points. For instance, when the target is a line segment with end points c₁ and c₂, for a point x in the network to be inactive, it is enough that there exists another point y in the network such that d(y,c₁)<d(x,c₁) and d(y,c₂)<d(x,c₂). This rule generalizes to polytopes (either for the target or a bounding box for it) and their extreme points in higher dimensions, but such information is not always easy to get in applications. More powerful rules have to adapt to complex target shapes and do not share the simplicity and versatility of the algorithm. Another reason for modifying the algorithm is the number of comparisons to be performed, even if active points have been selected.

We experimented with success a modified version of the algorithm, where a binary search tree is used to locate the closest points, and where no garbage collection is needed. The new rule is that after some branching occurs, one definitively assigns a part of the target to the future growth of each of the two branches. After each target point has been drawn, one thus first locates the branch to be grown by successively comparing the distances to couples of sons of the former branchings points, starting from the root (a classical binary search in a tree). Such algorithm remains very efficient even with large trees in high-dimension spaces. It also behaves quite like the former one, with good adaptive properties, but of course with no abortive bifurcations, since the subregions of the targets cannot move after each bifurcation.