A programmable reaction-diffusion system for spatiotemporal cell signaling circuit design.

Cells self-organize molecules in space and time to generate complex behaviors, but we lack synthetic strategies for engineering spatiotemporal signaling. We present a programmable reaction-diffusion platform for designing protein oscillations, patterns, and circuits in mammalian cells using two bacterial proteins, MinD and MinE (MinDE). MinDE circuits act like "single-cell radios," emitting frequency-barcoded fluorescence signals that can be spectrally isolated and analyzed using digital signal processing tools. We define how to genetically program these signals and connect their spatiotemporal dynamics to cell biology using engineerable protein-protein interactions. This enabled us to construct sensitive reporter circuits that broadcast endogenous cell signaling dynamics on a frequency-barcoded imaging channel and to build control signal circuits that synthetically pattern activities in the cell, such as protein condensate assembly and actin filamentation. Our work establishes a paradigm for visualizing, probing, and engineering cellular activities at length and timescales critical for biological function.

The Unreasonable Effectiveness of Reaction Diffusion in Vertebrate Skin Color Patterning.

In 1952, Alan Turing published the reaction-diffusion (RD) mathematical framework, laying the foundations of morphogenesis as a self-organized process emerging from physicochemical first principles. Regrettably, this approach has been widely doubted in the field of developmental biology. First, we summarize Turing's line of thoughts to alleviate the misconception that RD is an artificial mathematical construct. Second, we discuss why phenomenological RD models are particularly effective for understanding skin color patterning at the meso/macroscopic scales, without the need to parameterize the profusion of variables at lower scales. More specifically, we discuss how RD models (a) recapitulate the diversity of actual skin patterns, (b) capture the underlying dynamics of cellular interactions, (c) interact with tissue size and shape, (d) can lead to ordered sequential patterning, (e) generate cellular automaton dynamics in lizards and snakes, (f) predict actual patterns beyond their statistical features, and (g) are robust to model variations. Third, we discuss the utility of linear stability analysis and perform numerical simulations to demonstrate how deterministic RD emerges from the underlying chaotic microscopic agents.

The Clock and Wavefront Self-Organizing model recreates the dynamics of mouse somitogenesis in vivo and in vitro.

During mouse development, presomitic mesoderm cells synchronize Wnt and Notch oscillations, creating sequential phase waves that pattern somites. Traditional somitogenesis models attribute phase waves to a global modulation of the oscillation frequency. However, increasing evidence suggests that they could arise in a self-organizing manner. Here, we introduce the Sevilletor, a novel reaction-diffusion system that serves as a framework to compare different somitogenesis patterning hypotheses. Using this framework, we propose the Clock and Wavefront Self-Organizing model that considers an excitable self-organizing region where phase waves form independent of global frequency gradients. The model recapitulates the change in relative phase of Wnt and Notch observed during mouse somitogenesis and provides a theoretical basis for understanding the excitability of mouse presomitic mesoderm cells in vitro.

Collective chemomechanical oscillations in active hydrogels.

We report on the collective response of an assembly of chemomechanical Belousov-Zhabotinsky (BZ) hydrogel beads. We first demonstrate that a single isolated spherical BZ hydrogel bead with a radius below a critical value does not oscillate, whereas an assembly of the same BZ hydrogel beads presents chemical oscillation. A BZ chemical model with an additional flux of chemicals out of the BZ hydrogel captures the experimentally observed transition from oxidized nonoscillating to oscillating BZ hydrogels and shows this transition is due to a flux of inhibitors out of the BZ hydrogel. The model also captures the role of neighboring BZ hydrogel beads in decreasing the critical size for an assembly of BZ hydrogel beads to oscillate. We finally leverage the quorum sensing behavior of the collective to trigger their chemomechanical oscillation and discuss how this collective effect can be used to enhance the oscillatory strain of these active BZ hydrogels. These findings could help guide the eventual fabrication of a swarm of autonomous, communicating, and motile hydrogels.

Spiking at the edge: Excitability at interfaces in reaction-diffusion systems.

Excitable media, ranging from bioelectric tissues and chemical oscillators to forest fires and competing populations, are nonlinear, spatially extended systems capable of spiking. Most investigations of excitable media consider situations where the amplifying and suppressing forces necessary for spiking coexist at every point in space. In this case, spikes arise due to local bistabilities, which require a fine-tuned ratio between local amplification and suppression strengths. But, in nature and engineered systems, these forces can be segregated in space, forming structures like interfaces and boundaries. Here, we show how boundaries can generate and protect spiking when the reacting components can spread out: Even arbitrarily weak diffusion can cause spiking at the edge between two non-excitable media. This edge spiking arises due to a global bistability, which can occur even if amplification and suppression strengths do not allow spiking when mixed. We analytically derive a spiking phase diagram that depends on two parameters: i) the ratio between the system size and the characteristic diffusive length-scale and ii) the ratio between the amplification and suppression strengths. Our analysis explains recent experimental observations of action potentials at the interface between two non-excitable bioelectric tissues. Beyond electrophysiology, we highlight how edge spiking emerges in predator-prey dynamics and in oscillating chemical reactions. Our findings provide a theoretical blueprint for a class of interfacial excitations in reaction-diffusion systems, with potential implications for spatially controlled chemical reactions, nonlinear waveguides and neuromorphic computation, as well as spiking instabilities, such as cardiac arrhythmias, that naturally occur in heterogeneous biological media.

Differences in enantiomeric diffusion can lead to selective chiral amplification.

A fundamental question associated with chirality is how mixtures containing equal amounts of interconverting enantiomers can spontaneously convert to systems enriched in only one of them. Enantiomers typically have similar chemical properties, but can exhibit distinct reactivity under specific conditions, and these differences can be used to bias the system's composition in favor of one enantiomer. Transport properties are also expected to differ for enantiomers in chiral solvents, but the role of such differences in chiral symmetry breaking has not been clarified yet. In this work, we develop a theoretical framework to show that asymmetry in diffusion properties can trigger a spontaneous and selective symmetry breaking in mixtures of enantiomers. We derive a generic evolution equation for the enantiomeric excess in a chiral solvent. This equation shows that the relative stability of homochiral domains is dictated by the difference of diffusion coefficients of the two enantiomers. Consequently, deracemization toward a specific enantiomeric excess can be achieved when this difference is large enough. These results hold significant implications for our understanding of chiral symmetry breaking.

Harmonic resonance and entrainment of propagating chemical waves by external mechanical stimulation in BZ self-oscillating hydrogels.

Smart polymer materials that are nonliving yet exhibit complex "life-like" or biomimetic behaviors have been the focus of intensive research over the past decades, in the quest to broaden our understanding of how living systems function under nonequilibrium conditions. Identification of how chemical and mechanical coupling can generate resonance and entrainment with other cells or external environment is an important research question. We prepared Belousov-Zhabotinsky (BZ) self-oscillating hydrogels which convert chemical energy to mechanical oscillation. By cyclically applying external mechanical stimulation to the BZ hydrogels, we found that when the oscillation of a gel sample entered into harmonic resonance with the applied oscillation during stimulation, the system kept a "memory" of the resonant oscillation period and maintained it post stimulation, demonstrating an entrainment effect. More surprisingly, by systematically varying the cycle length of the external stimulation, we revealed the discrete nature of the stimulation-induced resonance and entrainment behaviors in chemical oscillations of BZ hydrogels, i.e., the hydrogels slow down their oscillation periods to the harmonics of the cycle length of the external mechanical stimulation. Our theoretical model calculations suggest the important roles of the delayed mechanical response caused by reactant diffusion and solvent migration in affecting the chemomechanical coupling in active hydrogels and consequently synchronizing their chemical oscillations with external mechanical oscillations.

RNA polymerase sliding on DNA can couple the transcription of nearby bacterial operons.

DNA transcription initiates after an RNA polymerase (RNAP) molecule binds to the promoter of a gene. In bacteria, the canonical picture is that RNAP comes from the cytoplasmic pool of freely diffusing RNAP molecules. Recent experiments suggest the possible existence of a separate pool of polymerases, competent for initiation, which freely slide on the DNA after having terminated one round of transcription. Promoter-dependent transcription reinitiation from this pool of posttermination RNAP may lead to coupled initiation at nearby operons, but it is unclear whether this can occur over the distance and timescales needed for it to function widely on a bacterial genome in vivo. Here, we mathematically model the hypothesized reinitiation mechanism as a diffusion-to-capture process and compute the distances over which significant interoperon coupling can occur and the time required. These quantities depend on molecular association and dissociation rate constants between DNA, RNAP, and the transcription initiation factor σ70; we measure these rate constants using single-molecule experiments in vitro. Our combined theory/experimental results demonstrate that efficient coupling can occur at physiologically relevant σ70 concentrations and on timescales appropriate for transcript synthesis. Coupling is efficient over terminator-promoter distances up to ∼1,000 bp, which includes the majority of terminator-promoter nearest neighbor pairs in the Escherichia coli genome. The results suggest a generalized mechanism that couples the transcription of nearby operons and breaks the paradigm that each binding of RNAP to DNA can produce at most one messenger RNA.

Stochastic Turing patterns of trichomes in Arabidopsis leaves.

Biological patterns that emerge during the morphogenesis of multicellular organisms can display high precision at large scales, while at cellular scales, cells exhibit large fluctuations stemming from cell-cell differences in molecular copy numbers also called demographic noise. We study the conflicting interplay between high precision and demographic noise in trichome patterns on the epidermis of wild-type Arabidopsis thaliana leaves, as a two-dimensional model system. We carry out a statistical characterization of these patterns and show that their power spectra display fat tails-a signature compatible with noise-driven stochastic Turing patterns-which are absent in power spectra of patterns driven by deterministic instabilities. We then present a theoretical model that includes demographic noise stemming from birth-death processes of genetic regulators which we study analytically and by stochastic simulations. The model captures the observed experimental features of trichome patterns.

The present and future of Turing models in developmental biology.

The Turing model (or reaction-diffusion model), first published in 1952, is a mathematical model that can account for autonomy in the morphogenesis of organisms. Although initially controversial, the model has gradually gained wider acceptance among experimental embryologists due to the accumulation of experimental data to support it. More recently, this model and others based on it have been used not only to explain biological phenomena conceptually but also as working hypotheses for molecular-level experiments and as internal components of more-complex 3D models. In this Spotlight, I will provide a personal perspective from an experimental biologist on some of the recent developments of the Turing model.

Ultraliser: a framework for creating multiscale, high-fidelity and geometrically realistic 3D models for in silico neuroscience.

Ultraliser is a neuroscience-specific software framework capable of creating accurate and biologically realistic 3D models of complex neuroscientific structures at intracellular (e.g. mitochondria and endoplasmic reticula), cellular (e.g. neurons and glia) and even multicellular scales of resolution (e.g. cerebral vasculature and minicolumns). Resulting models are exported as triangulated surface meshes and annotated volumes for multiple applications in in silico neuroscience, allowing scalable supercomputer simulations that can unravel intricate cellular structure-function relationships. Ultraliser implements a high-performance and unconditionally robust voxelization engine adapted to create optimized watertight surface meshes and annotated voxel grids from arbitrary non-watertight triangular soups, digitized morphological skeletons or binary volumetric masks. The framework represents a major leap forward in simulation-based neuroscience, making it possible to employ high-resolution 3D structural models for quantification of surface areas and volumes, which are of the utmost importance for cellular and system simulations. The power of Ultraliser is demonstrated with several use cases in which hundreds of models are created for potential application in diverse types of simulations. Ultraliser is publicly released under the GNU GPL3 license on GitHub (BlueBrain/Ultraliser). SIGNIFICANCE: There is crystal clear evidence on the impact of cell shape on its signaling mechanisms. Structural models can therefore be insightful to realize the function; the more realistic the structure can be, the further we get insights into the function. Creating realistic structural models from existing ones is challenging, particularly when needed for detailed subcellular simulations. We present Ultraliser, a neuroscience-dedicated framework capable of building these structural models with realistic and detailed cellular geometries that can be used for simulations.

Molecular Insights into Chemical Reactions at Aqueous Aerosol Interfaces.

Atmospheric aerosols facilitate reactions between ambient gases and dissolved species. Here, we review our efforts to interrogate the uptake of these gases and the mechanisms of their reactions both theoretically and experimentally. We highlight the fascinating behavior of N2O5 in solutions ranging from pure water to complex mixtures, chosen because its aerosol-mediated reactions significantly impact global ozone, hydroxyl, and methane concentrations. As a hydrophobic, weakly soluble, and highly reactive species, N2O5 is a sensitive probe of the chemical and physical properties of aerosol interfaces. We employ contemporary theory to disentangle the fate of N2O5 as it approaches pure and salty water, starting with adsorption and ending with hydrolysis to HNO3, chlorination to ClNO2, or evaporation. Flow reactor and gas-liquid scattering experiments probe even greater complexity as added ions, organic molecules, and surfactants alter the interfacial composition and reaction rates. Together, we reveal a new perspective on multiphase chemistry in the atmosphere.

Sensing the shape of a cell with reaction diffusion and energy minimization.

Some dividing cells sense their shape by becoming polarized along their long axis. Cell polarity is controlled in part by polarity proteins, like Rho GTPases, cycling between active membrane-bound forms and inactive cytosolic forms, modeled as a "wave-pinning" reaction-diffusion process. Does shape sensing emerge from wave pinning? We show that wave pinning senses the cell's long axis. Simulating wave pinning on a curved surface, we find that high-activity domains migrate to peaks and troughs of the surface. For smooth surfaces, a simple rule of minimizing the domain perimeter while keeping its area fixed predicts the final position of the domain and its shape. However, when we introduce roughness to our surfaces, shape sensing can be disrupted, and high-activity domains can become localized to locations other than the global peaks and valleys of the surface. On rough surfaces, the domains of the wave-pinning model are more robust in finding the peaks and troughs than the minimization rule, although both can become trapped in steady states away from the peaks and valleys. We can control the robustness of shape sensing by altering the Rho GTPase diffusivity and the domain size. We also find that the shape-sensing properties of cell polarity models can explain how domains localize to curved regions of deformed cells. Our results help to understand the factors that allow cells to sense their shape-and the limits that membrane roughness can place on this process.

Rapid formation of uniformly layered materials by coupling reaction-diffusion processes with mechanical responsiveness.

Straightforward manufacturing pathways toward large-scale, uniformly layered composites may enable the next generation of materials with advanced optical, thermal, and mechanical properties. Reaction-diffusion systems are attractive candidates to this aim, but while layered composites theoretically could spontaneously arise from reaction-diffusion, in practice randomly oriented patches separated by defects form, yielding nonuniformly patterned materials. A propagating reaction front can prevent such nonuniform patterning, as is the case for Liesegang processes, in which diffusion drives a reaction front to produce layered precipitation patterns. However, while diffusion is crucial to control patterning, it slows down transport of reactants to the front and results in a steady increase of the band spacing as the front advances. Here, we circumvent these diffusive limitations by embedding the Liesegang process in mechanically responsive hydrogels. The coupling between a moving reaction front and hydrogel contraction induces the formation of a self-regulated transport channel that ballistically carries reactants toward the area where patterning occurs. This ensures rapid and uniform patterning. Specifically, large-scale ([Formula: see text]5-cm) uniform banding patterns are produced with tunable band distance (d = 60 to 160 µm) of silver dichromate crystals inside responsive gelatin-alginate hydrogels. The generality and applicability of our mechanoreaction-diffusion strategy are demonstrated by forming patterns of precipitates in significantly smaller microscopic banding patterns (d = 10 to 30 µm) that act as self-organized diffraction gratings. By circumventing the inherent limitations of diffusion, our strategy unlocks the potential of reaction-diffusion processes for the manufacturing of uniformly layered materials.

Slow expanders invade by forming dented fronts in microbial colonies.

Most organisms grow in space, whether they are viruses spreading within a host tissue or invasive species colonizing a new continent. Evolution typically selects for higher expansion rates during spatial growth, but it has been suggested that slower expanders can take over under certain conditions. Here, we report an experimental observation of such population dynamics. We demonstrate that mutants that grow slower in isolation nevertheless win in competition, not only when the two types are intermixed, but also when they are spatially segregated into sectors. The latter was thought to be impossible because previous studies focused exclusively on the global competitions mediated by expansion velocities, but overlooked the local competitions at sector boundaries. Local competition, however, can enhance the velocity of either type at the sector boundary and thus alter expansion dynamics. We developed a theory that accounts for both local and global competitions and describes all possible sector shapes. In particular, the theory predicted that a slower on its own, but more competitive, mutant forms a dented V-shaped sector as it takes over the expansion front. Such sectors were indeed observed experimentally, and their shapes matched quantitatively with the theory. In simulations, we further explored several mechanisms that could provide slow expanders with a local competitive advantage and showed that they are all well-described by our theory. Taken together, our results shed light on previously unexplored outcomes of spatial competition and establish a universal framework to understand evolutionary and ecological dynamics in expanding populations.

Modeling protein dynamics in Caenorhabditis elegans embryos reveals that the PLK-1 gradient relies on weakly coupled reaction-diffusion mechanisms.

SignificanceIntracellular gradients have essential roles in cell and developmental biology, but their formation is not fully understood. We have developed a computational approach facilitating interpretation of protein dynamics and gradient formation. We have combined this computational approach with experiments to understand how Polo-Like Kinase 1 (PLK-1) forms a cytoplasmic gradient in Caenorhabditis elegans embryos. Although the PLK-1 gradient depends on the Muscle EXcess-5/6 (MEX-5/6) proteins, we reveal differences in PLK-1 and MEX-5 gradient formation that can be explained by a model with two components, PLK-1 bound to MEX-5 and unbound PLK-1. Our combined approach suggests that a weak coupling between PLK-1 and MEX-5 reaction-diffusion mechanisms dictates the dynamic exchange of PLK-1 with the cytoplasm, explaining PLK-1 high diffusivity and smooth gradient.

Bridging scales in a multiscale pattern-forming system.

Self-organized pattern formation is vital for many biological processes. Reaction-diffusion models have advanced our understanding of how biological systems develop spatial structures, starting from homogeneity. However, biological processes inherently involve multiple spatial and temporal scales and transition from one pattern to another over time, rather than progressing from homogeneity to a pattern. To deal with such multiscale systems, coarse-graining methods are needed that allow the dynamics to be reduced to the relevant degrees of freedom at large scales, but without losing information about the patterns at small scales. Here, we present a semiphenomenological approach which exploits mass conservation in pattern formation, and enables reconstruction of information about patterns from the large-scale dynamics. The basic idea is to partition the domain into distinct regions (coarse grain) and determine instantaneous dispersion relations in each region, which ultimately inform about local pattern-forming instabilities. We illustrate our approach by studying the Min system, a paradigmatic model for protein pattern formation. By performing simulations, we first show that the Min system produces multiscale patterns in a spatially heterogeneous geometry. This prediction is confirmed experimentally by in vitro reconstitution of the Min system. Using a recently developed theoretical framework for mass-conserving reaction-diffusion systems, we show that the spatiotemporal evolution of the total protein densities on large scales reliably predicts the pattern-forming dynamics. Our approach provides an alternative and versatile theoretical framework for complex systems where analytical coarse-graining methods are not applicable, and can, in principle, be applied to a wide range of systems with an underlying conservation law.

Performance characteristics allow for confinement of a CRISPR toxin-antidote gene drive for population suppression in a reaction-diffusion model.

Gene drive alleles that can bias their own inheritance could engineer populations for control of disease vectors, invasive species and agricultural pests. There are successful examples of suppression drives and confined modification drives, but developing confined suppression drives has proven more difficult. However, CRISPR-based toxin-antidote dominant embryo (TADE) suppression drive may fill this niche. It works by targeting and disrupting a haplolethal target gene in the germline with its gRNAs while rescuing this target. It also disrupts a female fertility gene by driving insertion or additional gRNAs. Here, we used a reaction-diffusion model to assess drive performance in continuous space, where outcomes can be substantially different from those in panmictic populations. We measured drive wave speed and found that moderate fitness costs or target gene disruption in the early embryo from maternally deposited nuclease can eliminate the drive's ability to form a wave of advance. We assessed the required release size, and finally we investigated migration corridor scenarios. It is often possible for the drive to suppress one population and then persist in the corridor without invading the second population, a potentially desirable outcome. Thus, even imperfect variants of TADE suppression drive may be excellent candidates for confined population suppression.

Modelling count data with partial differential equation models in biology.

Partial differential equation (PDE) models are often used to study biological phenomena involving movement-birth-death processes, including ecological population dynamics and the invasion of populations of biological cells. Count data, by definition, is non-negative, and count data relating to biological populations is often bounded above by some carrying capacity that arises through biological competition for space or nutrients. Parameter estimation, parameter identifiability, and making model predictions usually involves working with a measurement error model that explicitly relating experimental measurements with the solution of a mathematical model. In many biological applications, a typical approach is to assume the data are normally distributed about the solution of the mathematical model. Despite the widespread use of the standard additive Gaussian measurement error model, the assumptions inherent in this approach are rarely explicitly considered or compared with other options. Here, we interpret scratch assay data, involving migration, proliferation and delays in a population of cancer cells using a reaction-diffusion PDE model. We consider relating experimental measurements to the PDE solution using a standard additive Gaussian measurement error model alongside a comparison to a more biologically realistic binomial measurement error model. While estimates of model parameters are relatively insensitive to the choice of measurement error model, model predictions for data realisations are very sensitive. The standard additive Gaussian measurement error model leads to biologically inconsistent predictions, such as negative counts and counts that exceed the carrying capacity across a relatively large spatial region within the experiment. Furthermore, the standard additive Gaussian measurement error model requires estimating an additional parameter compared to the binomial measurement error model. In contrast, the binomial measurement error model leads to biologically plausible predictions and is simpler to implement. We provide open source Julia software on GitHub to replicate all calculations in this work, and we explain how to generalise our approach to deal with coupled PDE models with several dependent variables through a multinomial measurement error model, as well as pointing out other potential generalisations by linking our work with established practices in the field of generalised linear models.

A mathematical model of signalling molecule-mediated processes during regeneration of osteochondral defects after chondrocyte implantation.

Treating bone-cartilage defects is a fundamental clinical problem. The ability of damaged cartilage to self-repair is limited due to its avascularity. Left untreated, these defects can lead to osteoarthritis. Details of osteochondral defect repair are elusive, but animal models indicate healing occurs via an endochondral ossification-like process, similar to that in the growth plate. In the growth plate, the signalling molecules parathyroid hormone-related protein (PTHrP) and Indian Hedgehog (Ihh) form a feedback loop regulating chondrocyte hypertrophy, with Ihh inducing and PTHrP suppressing hypertrophy. To better understand this repair process and to explore the regulatory role of signalling molecules on the regeneration process, we formulate a reaction-diffusion mathematical model of osteochondral defect regeneration after chondrocyte implantation. The drivers of healing are assumed to be chondrocytes and osteoblasts, and their interaction via signalling molecules. We model cell proliferation, migration and chondrocyte hypertrophy, and matrix production and conversion, spatially and temporally. We further model nutrient and signalling molecule diffusion and their interaction with the cells. We consider the PTHrP-Ihh feedback loop as the backbone mechanisms but the model is flexible to incorporate extra signalling mechanisms if needed. Our mathematical model is able to represent repair of osteochondral defects, starting with cartilage formation throughout the defect. This is followed by chondrocyte hypertrophy, matrix calcification and bone formation deep inside the defect, while cartilage at the surface is maintained and eventually separated from the deeper bone by a thin layer of calcified cartilage. The complete process requires around 48 months. A key highlight of the model demonstrates that the PTHrP-Ihh loop alone is insufficient and an extra mechanism is required to initiate chondrocyte hypertrophy, represented by a critical cartilage density. A parameter sensitivity study reveals that the timing of the repair process crucially depends on parameters, such as the critical cartilage density, and those describing the actions of PTHrP to suppress hypertrophy, such as its diffusion coefficient, threshold concentration and degradation rate.