Size matters: An analytical study on the role of tissue size in spatiotemporal distribution of morphogens unveils a transition between different Reaction-Diffusion regimes

The reaction-diffusion model constitutes one of the most influential mathematical models to study distribution of morphogens in tissues. Despite its widespread use, the effect of finite tissue size on model-predicted spatiotemporal morphogen distributions has not been completely elucidated. In this study, we analytically investigated the spatiotemporal distributions of morphogens predicted by a reaction-diffusion model in a finite 1D domain, as a proxy for a biological tissue, and compared it with the solution of the infinite-domain model. We explored the reduced parameter, the tissue length in units of a characteristic reaction-diffusion length, and identified two reaction-diffusion regimes separated by a crossover tissue size estimated in ∼3.3 characteristic reaction-diffusion lengths. While above this crossover the infinite-domain model constitutes a good approximation, it breaks below this crossover, whereas the finite-domain model faithfully describes the entire parameter space. We evaluated whether the infinite-domain model renders accurate estimations of diffusion coefficients when fitted to finite spatial profiles, a procedure typically followed in Fluorescence Recovery After Photobleaching (FRAP) experiments. We found that the infinite-domain model overestimates diffusion coefficients when the domain is smaller than the crossover tissue size. Thus, the crossover tissue size may be instrumental in selecting the suitable reaction-diffusion model to study tissue morphogenesis.


Introduction 35
In their transition towards maturity, tissues are crucially regulated by molecules known as 36 morphogens, whose precise spatiotemporal distribution triggers the downstream changes in 37 protein expression responsible for the exact differentiation patterns. Nevertheless, tissues are not 38 only an inert scaffold upon which morphogens spread, but they are also fully responsible for the 39 morphogen uptake or their transformation by means of specific biochemical reactions. The 40 problem of how a morphogen propagates over a tissue while it is being eliminated was 41 mathematically encoded in the exquisite reaction-diffusion model by the great Alan Turing, who 42 coined the "morphogen" neologism to illustrate its character of "form generator" [1]. 43 The reaction-diffusion model constitutes one of the most influential quantitative approaches 44 within developmental biology. From the aforementioned Turing´s seminal article and the study 45 from Gierer and Meinhardt [2], a progressive wealth of reaction-diffusion models were developed, 46 paving the way to become an essential and pivotal concept to understand tissue morphogenesis 47 Here we briefly summarize the well-known reaction-diffusion model assuming an infinite domain 70 and its analytic solution [10,11]. Within this model, it is assumed that the dynamics of the 71 morphogen are faster than the proliferation rate of the tissue cells and, as a consequence, 72 advective effects can be neglected. Otherwise, an advective term could be included to the model 73 [17]. Since during developmental process tissues usually organize along a particular axis [18,19], 74 this model is studied in a one dimensional setting [10,11]. It is assumed that the morphogen 75 concentration 1 ( , ) can diffuse with a diffusion coefficient D and is linearly degraded with a 76 rate k. It is considered that there is no morphogen at the beginning, that is, the initial condition is: Where the morphogen source at the tissue origin, in nondimensional units, = 0, is: 97 And a morphogen sink at infinite in the nondimensionalized units is now: 99 lim →∞ ( , ) = 0 (Eq. 9) 100 Which leads to this solution: 101 (Eq. 10) 102

The reaction-diffusion model in finite domains: an analytical solution 104
The previous model variant entails an infinite domain (Eqs. 4 and 9). Since biological tissue sizes 105 require a finite domain, we decided to replace the condition imposed by Eq. 4 with: 106 Where L is the length of the tissue. To our knowledge, the general solution for any given L is yet 108 unknown. 109 We defined the quantity = √ , which is the only model parameter. This quantity represents the 110 tissue length L in units of the characteristic reaction-diffusion length , defined as = √ D [20,21]. 111 Thus, the second boundary condition for this model in nondimensionalized units is: 112 This equation replaces Eq. 9 in the section 2.1 assuming the finitude of the tissue. 114 We found the analytical solution of the general model for finite tissues (Eqs. 6-8 and 12) in the 115 nondimensionalized units to be as follows (see Supplementary information for the 116 demonstration): 117 Moreover, we also found the solution for different boundary conditions such as assuming a non-

128
We decided to compare the reaction-diffusion model assuming a finite tissue versus an infinite 129 domain. With the selected nondimensionalization, the latter does not have any free parameters. 130 In contrast, the finite model has only one free parameter, R, which represents the tissue size in 131 units of the characteristic length of the morphogen profile . By using our analytical solution for 132 the model of finite tissues (Eq. 13), we explored the predicted morphogen spatial profiles at 133 different tissue sizes (i.e., varying R) and compared them with those calculated from the previously 134 known solution assuming an infinite domain (Eq. 10), at three different time points (Fig. 1). We 135 observed that the morphogen concentrations predicted by the model assuming an infinite domain 136 are higher than those predicted by the model assuming a finite domain ( Fig. 1 A, B). For large 137 enough tissue lengths, morphogen profiles predicted by both models are indistinguishable at each 138 time point, as expected ( Fig. 1 C and Fig 1 D). Hence, the previously reported model assuming an 139 infinite domain is a reasonable description of the dynamics of morphogen profiles for larger 140 tissues. However, when addressing a tissue whose length is of the order of the characteristic 141 length  or smaller, the model introduced in the present work is a more accurate description.

142
Moreover, we observed that large tissues lead to morphogen spatial distributions temporarily 143 separated. In contrast, spatial distributions at different time points are indistinguishable in shorter 144 tissues, suggesting that they already approached the steady state (Fig. 1A). This result would 145 indicate that the larger the tissue, the longer the time necessary to reach the morphogen spatial 146 distribution at the steady state (see also sections 2.4 and 2.6).

Steady state morphogen spatial distributions 159
The morphogen spatial distribution assuming an infinite domain at the steady state ( (ε)) 160 is well known [10,11] and with our nondimensionalization it is the following exponential spatial 161 decay: 162 ( ) = − (Eq. 14) 163 We calculated the steady state solution for our model of finite tissues, Increasing the tissue size in this model modifies the steady state profile, augmenting the maximum 167 concentration at the origin and leading to a transition from a linear to an exponential curve ( Fig.  168 2), in agreement with the results observed at any time ( Fig. 1). Precisely, to estimate the limit 169 when the tissue size tends to zero, we calculated the Taylor series expansion of the steady state 170 solution (Eq. 15) on R to the first order. As is constrained by R, we subsequently obtained the 171 Taylor series expansion for the resulting expression on to the first order: 172 The limit when the tissue size tends to infinite was calculated: 174 Where the second term of Eq. 18 vanishes when the time  tends to infinity. Therefore, the  The analysis of the dependency of 10 with the tissue size can also be made before the morphogen 219 distribution achieves the steady state. Although we could not find an analytical expression for this 220 observable in the general case, we explored this dependency numerically (Fig. 4). We observed 221 that, for each tissue size, 10 increases in time until it reaches a plateau, which indicates that the 222 spatial profile stabilizes in the steady state. Moreover, the time needed to reach the plateau 223 monotonically increases with the tissue size until R ~ Rc. For larger tissue sizes, the time to reach 224 the plateau converges to the prediction of the model for infinite domains (Fig. 4). This result is 225 consistent with the fact that morphogen spatial distributions at different times are overlapped in 226 smaller tissues and separated in larger ones (Fig. 1). This is a consequence of the second term of 227 Eq. 18: the larger the tissue, the longer waiting times are required to vanish the exponential in the 228 second term. 229 , also experiences a transition near the crossover tissue size until 265 converging to √2 (Fig. 5 A). 266 For tissues smaller than the crossover size, the mean time to achieve the steady state and its error 267 in each position strongly depend on tissue size (Fig. 5 B). On the contrary, for tissue sizes higher 268 than the crossover tissue size, both magnitudes become independent of the size (Fig. 5 C). 269 Importantly, for tissues smaller than the crossover size, the steady state will be reached 270 significantly faster than the prediction from the model assuming an infinite domain. For tissues 271 larger than the crossover size, both models agree in the time to achieve steady state (Fig. 5 B and  272 C). 273 each of the datasets obtained using Eq. 15. We used Eq. 14 in the original coordinate as fitting 299 function and λ as the free parameter. We obtained the predicted value of λ as a function of R = L 300 (Fig. 6). For large values of R, the predicted λ is approximately 1, which is in agreement with the 301 value actually used to generate the data. In contrast, for values of R smaller than Rc, the predicted 302 value of λ deviate from 1, converging to 0 for small values of R. We concluded that both models 303 can be used to infer the kinetic parameters D and k from FRAP experiments, provided that tissue 304 sizes are higher than Rc. On the contrary, for tissues smaller than this crossover value, the model 305 assuming finite domains is the best alternative.  Our results showed that the morphogen spatial distributions predicted by our model assuming determining the spatial position along the tissue where the morphogen concentration is 10 % of 362 the source concentration (10), we geometrically characterized the steady state spatial distribution. crossover tissue size Rc (Fig. 3 and Fig. 7C and 7D). For tissues longer than Rc, the distributions are 365 exponential-like and cannot be distinguished from those predicted from the model assuming an 366 infinite domain (Fig. 2 and Fig. 7B). In this regime of the parameter space, the mean and standard 367 deviation of the time to reach the steady state (evaluated at the tissue origin) do not change much 368 with the tissue size and converged towards the corresponding values from the model assuming an 369 infinite domain (Fig. 5 and Fig. 7F). When comparing the morphogen concentrations predicted by 370 both models we found that the difference between them is mostly negligible (Fig. 8 A and B). 371 Hence, the model assuming an infinite domain can be considered a good approximation of the 372 model assuming finite domains for tissue sizes larger than Rc. 373  In contrast, for tissues smaller than Rc, the distributions tend to be linear and are clearly separated 382 from those predicted with the model assuming an infinite domain ( Fig. 2 and Fig. 7A). 383 Furthermore, the time to reach the steady state strongly depends on the tissue size in this regime 384 ( Fig. 5 and Fig. 7E). In particular, the error of using the model assuming an infinite domain 385 increases when  tends to R and the smaller the tissue the higher the error accumulated over the entire tissue (Fig. 8 A

398
The crossover tissue size provides a straightforward criterion to decide when to use any of the two 399 models presented here. As an example, the characteristic length of Wg was estimated in 6 m in 400 In conclusion, we found two reaction-diffusion regimes for large and small tissues, separated by a 418 crossover tissue size. While above this crossover the infinite-domain model constitutes a good 419 approximation, it breaks below this crossover, whereas the finite-domain model faithfully 420 describes the entire parameter space. Further studies will be needed to unveil the spatiotemporal 421 distribution of morphogens in tissues whose size is not fixed.