Reconciling diverse mammalian pigmentation patterns with a fundamental mathematical model

A biased random walk approach for modeling the collective chemotaxis of neural crest cells

Collective cell migration is a multicellular phenomenon that arises in various biological contexts, including cancer and embryo development. 'Collectiveness' can be promoted by cell-cell interactions such as co-attraction and contact inhibition of locomotion. These mechanisms act on cell polarity, pivotal for directed cell motility, through influencing the intracellular dynamics of small GTPases such as Rac1. To model these dynamics we introduce a biased random walk model, where the bias depends on the internal state of Rac1, and the Rac1 state is influenced by cell-cell interactions and chemoattractive cues. In an extensive simulation study we demonstrate and explain the scope and applicability of the introduced model in various scenarios. The use of a biased random walk model allows for the derivation of a corresponding partial differential equation for the cell density while still maintaining a certain level of intracellular detail from the individual based setting.

Fluctuations in tissue growth portray homeostasis as a critical state and long-time non-Markovian cell proliferation as Markovian

Tissue growth is an emerging phenomenon that results from the cell-level interplay between proliferation and apoptosis, which is crucial during embryonic development, tissue regeneration, as well as in pathological conditions such as cancer. In this theoretical article, we address the problem of stochasticity in tissue growth by first considering a minimal Markovian model of tissue size, quantified as the number of cells in a simulated tissue, subjected to both proliferation and apoptosis. We find two dynamic phases, growth and decay, separated by a critical state representing a homeostatic tissue. Since the main limitation of the Markovian model is its neglect of the cell cycle, we incorporated a refractory period that temporarily prevents proliferation immediately following cell division, as a minimal proxy for the cell cycle, and studied the model in the growth phase. Importantly, we obtained from this last model an effective Markovian rate, which accurately describes general trends of tissue size. This study shows that the dynamics of tissue growth can be theoretically conceptualized as a Markovian process where homeostasis is a critical state flanked by decay and growth phases. Notably, in the growing non-Markovian model, a Markovian-like growth process emerges at large time scales.

Swapping in lattice-based cell migration models

Cell migration is frequently modeled using on-lattice agent-based models (ABMs) that employ the excluded volume interaction. However, cells are also capable of exhibiting more complex cell-cell interactions, such as adhesion, repulsion, pulling, pushing, and swapping. Although the first four of these have already been incorporated into mathematical models for cell migration, swapping has not been well studied in this context. In this paper, we develop an ABM for cell movement in which an active agent can "swap" its position with another agent in its neighborhood with a given swapping probability. We consider a two-species system for which we derive the corresponding macroscopic model and compare it with the average behavior of the ABM. We see good agreement between the ABM and the macroscopic density. We also analyze the movement of agents at an individual level in the single-species as well as two-species scenarios to quantify the effects of swapping on an agent's motility.

Genomic analysis reveals a KIT-related chromosomal translocation associated with the white coat phenotype in yak

White coat pigmentation is a striking phenotype of many domesticated species and has various genetic controls. The Tianzhu White yak, an indigenous breed with a complete white coat, has fascinated Tibetans for centuries. However, the genetic basis of this trait remains unknown. Here, we conducted population genomics analysis and genome-wide association study based on the whole-genome sequencing data of 38 white and 59 non-white-coated yak. The results revealed the presence of KIT-linked Cs alleles characterized by the translocations between chromosomes 6 and 29 in all-white yak. Furthermore, structural variations showed additional duplications of the Cs alleles in white yak compared with colour-sidedness cattle. Interestingly, the Cs alleles associated with the white coat phenotype in yak were found to have introgressed from taurine cattle. Our findings unveil the shared genetic control of the white coat phenotype and its evolution in closely related bovine species.

Genetically engineered multicistronic allele of Pmel yielding highly specific CreERT2-mediated recombination in the melanocyte lineage

Genetic approaches that allow lineage tracing are essential to our future understanding of melanocytes and melanoma. To date, the approaches used to label melanocytes in mice have relied on random integration of transgenes driven by the promoters of the Tyrosinase and Dopachrome tautomerase genes, knock-in to the Dopachrome tautomerase locus or knock-in to the Mlana locus in a bacterial artificial chromosome. These strategies result in expression in other tissues such as telencephalon and other cell types such as nerves. Here we used homologous recombination in mouse embryonic stem cells to generate a targeted multicistronic allele of the Pmel locus that drives melanocyte-specific expression of CreERT2, nuclear localised H2B-Cerulean and membrane localised marcks-mKate2 allowing live imaging of melanocytes and activation of other conditional alleles. We combined this allele with R26R-EYFP mice allowing induction of EYFP expression on administration of tamoxifen or its metabolite 4-OHT. The fluorescent proteins H2B-Cerulean and marcks-mKate2 label the cell nucleus and plasma membrane respectively allowing live imaging and FACS isolation of melanoblasts and melanocytes as well as serving to provide an internal control allowing estimation of recombination efficiency after administration of tamoxifen. We demonstrate the utility of the transgene in embryonic and adult tissues.

Generation of Pmel-dependent conditional and inducible Cre-driver mouse line for melanocytic-targeted gene manipulation

Conditional and inducible gene targeting using Cre/loxP-mediated recombination is a powerful reverse genetics approach used to study spatiotemporal gene functions in specified cell types. To enable temporal gene manipulation in the melanocyte lineage, we established a novel inducible Cre-driver mouse line by targeting an all-in-one tetracycline/doxycycline (Dox)-inducible Cre expression cassette into the Pmel locus (Pmel^{P2A-TetON3G-TRE3G-iCre} ), a gene locus preferentially expressed in pigment cells. By crossing these Cre-driver mice with a strong Cre-reporter mouse line, Gt(ROSA)26Sor^{tm9(CAG-tdTomato)Hze} , we show the effectiveness of the Pmel^{P2A-TetON3G-TRE3G-iCre} mouse line in facilitating Dox-inducible Cre/loxP recombination in a wide variety of pigment cell lineages including hair follicle melanocytes and their stem cells. Furthermore, to demonstrate proof of concept, we ablated Notch signaling postnatally in the Pmel^{P2A-TetON3G-TRE3G-iCre} mice. In agreement with the previously reported phenotype, induced ablation of Notch signaling in the melanocyte lineage resulted in premature hair graying, demonstrating the utility of the Pmel^{P2A-TetON3G-TRE3G-iCre} allele. Therefore, the Pmel^{P2A-TetON3G-TRE3G-iCre} mouse line is suitable for assessing gene functions in melanocytes using an in vivo inducible reverse genetics approach. Furthermore, we unexpectedly identified previously unrecognized PMEL-expressing cells in non-pigmentary organs in the mice, suggesting unanticipated functions of PMEL other than melanosome formation.

Reliable and efficient parameter estimation using approximate continuum limit descriptions of stochastic models

Stochastic individual-based mathematical models are attractive for modelling biological phenomena because they naturally capture the stochasticity and variability that is often evident in biological data. Such models also allow us to track the motion of individuals within the population of interest. Unfortunately, capturing this microscopic detail means that simulation and parameter inference can become computationally expensive. One approach for overcoming this computational limitation is to coarse-grain the stochastic model to provide an approximate continuum model that can be solved using far less computational effort. However, coarse-grained continuum models can be biased or inaccurate, particularly for certain parameter regimes. In this work, we combine stochastic and continuum mathematical models in the context of lattice-based models of two-dimensional cell biology experiments by demonstrating how to simulate two commonly used experiments: cell proliferation assays and barrier assays. Our approach involves building a simple statistical model of the discrepancy between the expensive stochastic model and the associated computationally inexpensive coarse-grained continuum model. We form this statistical model based on a limited number of expensive stochastic model evaluations at design points sampled from a user-chosen distribution, corresponding to a computer experiment design problem. With straightforward design point selection schemes, we show that using the statistical model of the discrepancy in tandem with the computationally inexpensive continuum model allows us to carry out prediction and inference while correcting for biases and inaccuracies due to the continuum approximation. We demonstrate this approach by simulating a proliferation assay, where the continuum limit model is the well-known logistic ordinary differential equation, as well as a barrier assay where the continuum limit model is closely related to the well-known Fisher-KPP partial differential equation. We construct an approximate likelihood function for parameter inference, both with and without discrepancy correction terms. Using maximum likelihood estimation, we provide point estimates of the unknown parameters, and use the profile likelihood to characterise the uncertainty in these estimates and form approximate confidence intervals. For the range of inference problems considered, working with the continuum limit model alone leads to biased parameter estimation and confidence intervals with poor coverage. In contrast, incorporating correction terms arising from the statistical model of the model discrepancy allows us to recover the parameters accurately with minimal computational overhead. The main tradeoff is that the associated confidence intervals are typically broader, reflecting the additional uncertainty introduced by the approximation process. All algorithms required to replicate the results in this work are written in the open source Julia language and are available at GitHub.

Incorporating domain growth into hybrid methods for reaction-diffusion systems

Reaction-diffusion mechanisms are a robust paradigm that can be used to represent many biological and physical phenomena over multiple spatial scales. Applications include intracellular dynamics, the migration of cells and the patterns formed by vegetation in semi-arid landscapes. Moreover, domain growth is an important process for embryonic growth and wound healing. There are many numerical modelling frameworks capable of simulating such systems on growing domains; however, each of these may be well suited to different spatial scales and particle numbers. Recently, spatially extended hybrid methods on static domains have been produced to bridge the gap between these different modelling paradigms in order to represent multi-scale phenomena. However, such methods have not been developed with domain growth in mind. In this paper, we develop three hybrid methods on growing domains, extending three of the prominent static-domain hybrid methods. We also provide detailed algorithms to allow others to employ them. We demonstrate that the methods are able to accurately model three representative reaction-diffusion systems accurately and without bias.

Synchronized oscillations in growing cell populations are explained by demographic noise

Understanding synchrony in growing populations is important for applications as diverse as epidemiology and cancer treatment. Recent experiments employing fluorescent reporters in melanoma cell lines have uncovered growing subpopulations exhibiting sustained oscillations, with nearby cells appearing to synchronize their cycles. In this study, we demonstrate that the behavior observed is consistent with long-lasting transient phenomenon initiated and amplified by the finite-sample effects and demographic noise. We present a novel mathematical analysis of a multistage model of cell growth, which accurately reproduces the synchronized oscillations. As part of the analysis, we elucidate the transient and asymptotic phases of the dynamics and derive an analytical formula to quantify the effect of demographic noise in the appearance of the oscillations. The implications of these findings are broad, such as providing insight into experimental protocols that are used to study the growth of asynchronous populations and, in particular, those investigations relating to anticancer drug discovery.

Apparent piebald variants in quolls (Dasyurus): examples of three recent cases in the spotted-tailed quoll Dasyurus maculatus

Pelage patterning plays an important role in animal behaviour. Variation in pelage patterns can change with pigment distribution and quantity in individuals. We present three cases of apparent piebaldism – a condition where the body is patchily unpigmented – in the spotted-tailed quoll Dasyurus maculatus. Using a comprehensive dataset of historical descriptions (from Dunlop et al. 2020), we conclude that these cases represent the first description of this phenotype in the genus Dasyurus, but acknowledge capture and testing of these individuals is required to be certain of the cause. Little is known about the implications of pelage patterning in quoll species and further investigation is required to understand the evolutionary and functional role of unaffected and unpigmented variants.

The blending region hybrid framework for the simulation of stochastic reaction-diffusion processes

The simulation of stochastic reaction–diffusion systems using fine-grained representations can become computationally prohibitive when particle numbers become large. If particle numbers are sufficiently high then it may be possible to ignore stochastic fluctuations and use a more efficient coarse-grained simulation approach. Nevertheless, for multiscale systems which exhibit significant spatial variation in concentration, a coarse-grained approach may not be appropriate throughout the simulation domain. Such scenarios suggest a hybrid paradigm in which a computationally cheap, coarse-grained model is coupled to a more expensive, but more detailed fine-grained model, enabling the accurate simulation of the fine-scale dynamics at a reasonable computational cost. In this paper, in order to couple two representations of reaction–diffusion at distinct spatial scales, we allow them to overlap in a ‘blending region’. Both modelling paradigms provide a valid representation of the particle density in this region. From one end of the blending region to the other, control of the implementation of diffusion is passed from one modelling paradigm to another through the use of complementary ‘blending functions’ which scale up or down the contribution of each model to the overall diffusion. We establish the reliability of our novel hybrid paradigm by demonstrating its simulation on four exemplar reaction–diffusion scenarios.

Modelling collective cell migration: neural crest as a model paradigm

A huge variety of mathematical models have been used to investigate collective cell migration. The aim of this brief review is twofold: to present a number of modelling approaches that incorporate the key factors affecting cell migration, including cell-cell and cell-tissue interactions, as well as domain growth, and to showcase their application to model the migration of neural crest cells. We discuss the complementary strengths of microscale and macroscale models, and identify why it can be important to understand how these modelling approaches are related. We consider neural crest cell migration as a model paradigm to illustrate how the application of different mathematical modelling techniques, combined with experimental results, can provide new biological insights. We conclude by highlighting a number of future challenges for the mathematical modelling of neural crest cell migration.

Melanoma mutations modify melanocyte dynamics in co-culture with keratinocytes or fibroblasts

Melanocytic cell interactions are integral to skin homeostasis, and affect the outcome of multiple diseases, including cutaneous pigmentation disorders and melanoma. By using automated-microscopy and machine-learning-assisted morphology analysis of primary human melanocytes in co-culture, we performed combinatorial interrogation of melanocyte genotypic variants and functional assessment of lentivirus-introduced mutations. Keratinocyte-induced melanocyte dendricity, an indicator of melanocyte differentiation, was reduced in the melanocortin 1 receptor (MC1R) R/R variant strain and by NRAS.Q61K and BRAF.V600E expression, while expression of CDK4.R24C and RAC1.P29S had no detectable effect. Time-lapse tracking of melanocytes in co-culture revealed dynamic interaction phenotypes and hyper-motile cell states that indicated that, in addition to the known role in activating mitogenic signalling, MEK-pathway-activating mutations may also allow melanocytes to escape keratinocyte control and increase their invasive potential. Expanding this combinatorial platform will identify other therapeutic target mutations and melanocyte genetic variants, as well as increase understanding of skin cell interactions.

Persistent exclusion processes: Inertia, drift, mixing, and correlation

In many biological systems, motile agents exhibit random motion with short-term directional persistence, together with crowding effects arising from spatial exclusion. We formulate and study a class of lattice-based models for multiple walkers with motion persistence and spatial exclusion in one and two dimensions, and use a mean-field approximation to investigate relevant population-level partial differential equations in the continuum limit. We show that this model of a persistent exclusion process is in general well described by a nonlinear diffusion equation. With reference to results presented in the current literature, our results reveal that the nonlinearity arises from the combination of motion persistence and volume exclusion, with linearity in terms of the canonical diffusion or heat equation being recovered in either the case of persistence without spatial exclusion, or spatial exclusion without persistence. We generalize our results to include systems of multiple species of interacting, motion-persistent walkers, as well as to incorporate a global drift in addition to persistence. These models are shown to be governed approximately by systems of nonlinear advection-diffusion equations. By comparing the prediction of the mean-field approximation to stochastic simulation results, we assess the performance of our results. Finally, we also address the problem of inferring the presence of persistence from simulation results, with a view to application to experimental cell-imaging data.

Modeling and analysis of melanoblast motion

Melanoblast migration is important for embryogenesis and is a key feature of melanoma metastasis. Many studies have characterized melanoblast movement, focusing on statistical properties and have highlighted basic mechanisms of melanoblast motility. We took a slightly different and complementary approach: we previously developed a mathematical model of melanoblast motion that enables the testing of biological assumptions about the displacement of melanoblasts and we created tests to analyze the geometric features of cell trajectories and the specific issue of trajectory interactions. Within this model, we performed simulations and compared the results with experimental data using geometric tests. In this paper, we developed the associated mathematical model and the main focus is to study the crossings between trajectories with new theoretical results about the variation of number of intersection points with respect to the crossing times. Using these results it is possible to study the random nature of displacements and the interactions between trajectories. This analysis has raised new questions, leading to the generation of strong arguments in favor of a trail left behind each moving melanoblast.

Modelling chase-and-run migration in heterogeneous populations

Cell migration is crucial for many physiological and pathological processes. During embryogenesis, neural crest cells undergo coordinated epithelial to mesenchymal transformations and migrate towards various forming organs. Here we develop a computational model to understand how mutual interactions between migrating neural crest cells (NCs) and the surrounding population of placode cells (PCs) generate coordinated migration. According to experimental findings, we implement a minimal set of hypotheses, based on a coupling between chemotactic movement of NCs in response to a placode-secreted chemoattractant (Sdf1) and repulsion induced from contact inhibition of locomotion (CIL), triggered by heterotypic NC–PC contacts. This basic set of assumptions is able to semi-quantitatively recapitulate experimental observations of the characteristic multispecies phenomenon of “chase-and-run”, where the colony of NCs chases an evasive PC aggregate. The model further reproduces a number of in vitro manipulations, including full or partial disruption of NC chemotactic migration and selected mechanisms coordinating the CIL phenomenon. Finally, we provide various predictions based on altering other key components of the model mechanisms.

Neural crest migration with continuous cell states

Models of cranial neural crest cell migration in cell-induced (or self-generated) gradients have included a division of labour into leader and follower migratory states, which undergo chemotaxis and contact guidance, respectively. Despite validated utility of these models through experimental perturbation of migration in the chick embryo and gene expression analysis showing relevant heterogeneity at the single cell level, an often raised question has been whether the discrete cell states are necessary, or if a continuum of cell behaviours offers a functionally equivalent description. Here we argue that this picture is supported by recent single-cell transcriptome data. Motivated by this, we implement two versions of a continuous-state model: (1) signal choice and (2) signal combination. We find that the cell population migrates further than in the discrete-state model and than in experimental observations. We further show that the signal combination model, but not the signal choice model, can be successfully adjusted to experimentally plausible regimes by reducing the chemoattractant consumption parameter. Thus we show an equivalently plausible, experimentally motivated, model of neural crest cell migration.

Extending the Mathematical Palette for Developmental Pattern Formation: Piebaldism

Here, we present a theoretical investigation with potential insights on developmental mechanisms. Three biological factors, consisting of two diffusing factors and a cell-autonomous immobile transcription factor are combined with different feedback mechanisms. This results in four different situations or fur patterns. Two of them reproduce classical Turing patterns: (1) regularly spaced spots, (2) labyrinth patterns or straight lines with an initial slope in the activation of the transcription factor. The third situation does not lead to patterns, but results in different homogeneous color tones. Finally, the fourth one sheds new light on the possible mechanisms leading to the formation of piebald patterns exemplified by the random patterns on the fur of some cows' strains and Dalmatian dogs. Piebaldism is usually manifested as white areas of fur, hair, or skin due to the absence of pigment-producing cells in those regions. The distribution of the white and colored zones does not reflect the classical Turing patterns. We demonstrate that these piebald patterns are of transient nature, developing from random initial conditions and relying on a system's bistability. We show numerically that the presence of a cell-autonomous factor not only expands the range of reaction diffusion parameters in which a pattern may arise, but also extends the pattern-forming abilities of the reaction-diffusion equations.

The invasion speed of cell migration models with realistic cell cycle time distributions

Cell proliferation is typically incorporated into stochastic mathematical models of cell migration by assuming that cell divisions occur after an exponentially distributed waiting time. Experimental observations, however, show that this assumption is often far from the real cell cycle time distribution (CCTD). Recent studies have suggested an alternative approach to modelling cell proliferation based on a multi-stage representation of the CCTD. In this paper we investigate the connection between the CCTD and the speed of the collective invasion. We first state a result for a general CCTD, which allows the computation of the invasion speed using the Laplace transform of the CCTD. We use this to deduce the range of speeds for the general case. We then focus on the more realistic case of multi-stage models, using both a stochastic agent-based model and a set of reaction-diffusion equations for the cells' average density. By studying the corresponding travelling wave solutions, we obtain an analytical expression for the speed of invasion for a general N-stage model with identical transition rates, in which case the resulting cell cycle times are Erlang distributed. We show that, for a general N-stage model, the Erlang distribution and the exponential distribution lead to the minimum and maximum invasion speed, respectively. This result allows us to determine the range of possible invasion speeds in terms of the average proliferation time for any multi-stage model.

The auxiliary region method: a hybrid method for coupling PDE- and Brownian-based dynamics for reaction-diffusion systems

Reaction-diffusion systems are used to represent many biological and physical phenomena. They model the random motion of particles (diffusion) and interactions between them (reactions). Such systems can be modelled at multiple scales with varying degrees of accuracy and computational efficiency. When representing genuinely multiscale phenomena, fine-scale models can be prohibitively expensive, whereas coarser models, although cheaper, often lack sufficient detail to accurately represent the phenomenon at hand. Spatial hybrid methods couple two or more of these representations in order to improve efficiency without compromising accuracy. In this paper, we present a novel spatial hybrid method, which we call the auxiliary region method (ARM), which couples PDE- and Brownian-based representations of reaction-diffusion systems. Numerical PDE solutions on one side of an interface are coupled to Brownian-based dynamics on the other side using compartment-based 'auxiliary regions'. We demonstrate that the hybrid method is able to simulate reaction-diffusion dynamics for a number of different test problems with high accuracy. Furthermore, we undertake error analysis on the ARM which demonstrates that it is robust to changes in the free parameters in the model, where previous coupling algorithms are not. In particular, we envisage that the method will be applicable for a wide range of spatial multi-scales problems including filopodial dynamics, intracellular signalling, embryogenesis and travelling wave phenomena.

Reactions, diffusion, and volume exclusion in a conserved system of interacting particles

Complex biological and physical transport processes are often described through systems of interacting particles. The effect of excluded volume on these transport processes has been well studied; however, the interplay between volume exclusion and reactions between heterogenous particles is less well studied. In this paper we develop a framework for modeling reaction-diffusion processes which directly incorporates volume exclusion. We consider simple reactions (unimolecular and bimolecular) that conserve the total number of particles. From an off-lattice microscopic individual-based model we use the Fokker-Planck equation and the method of matched asymptotic expansions to derive a low-dimensional macroscopic system of nonlinear partial differential equations describing the evolution of the particles. A biologically motivated, hybrid model of chemotaxis with volume exclusion is explored, where reactions occur at rates dependent upon the chemotactic environment. Further, we show that for reactions that require particle contact the appropriate reaction term in the macroscopic model is of lower order in the asymptotic expansion than the nonlinear diffusion term. However, we find that the next reaction term in the expansion is needed to ensure good agreement with simulations of the microscopic model. Our macroscopic model allows for more direct parametrization to experimental data than existing models.

Modeling persistence of motion in a crowded environment: The diffusive limit of excluding velocity-jump processes

Persistence of motion is the tendency of an object to maintain motion in a direction for short time scales without necessarily being biased in any direction in the long term. One of the most appropriate mathematical tools to study this behavior is an agent-based velocity-jump process. In the absence of agent-agent interaction, the mean-field continuum limit of the agent-based model (ABM) gives rise to the well known hyperbolic telegraph equation. When agent-agent interaction is included in the ABM, a strictly advective system of partial differential equations (PDEs) can be derived at the population level. However, no diffusive limit of the ABM has been obtained from such a model. Connecting the microscopic behavior of the ABM to a diffusive macroscopic description is desirable, since it allows the exploration of a wider range of scenarios and establishes a direct connection with commonly used statistical tools of movement analysis. In order to connect the ABM at the population level to a diffusive PDE at the population level, we consider a generalization of the agent-based velocity-jump process on a two-dimensional lattice with three forms of agent interaction. This generalization allows us to take a diffusive limit and obtain a faithful population-level description. We investigate the properties of the model at both the individual and population levels and we elucidate some of the models' key characteristic features. In particular, we show an intrinsic anisotropy inherent to the models and we find evidence of a spontaneous form of aggregation at both the micro- and macroscales.

Simulation of melanoblast displacements reveals new features of developmental migration

To distribute and establish the melanocyte lineage throughout the skin and other developing organs, melanoblasts undergo several rounds of proliferation, accompanied by migration through complex environments and differentiation. Melanoblast migration requires interaction with extracellular matrix of the epidermal basement membrane and with surrounding keratinocytes in the developing skin. Migration has been characterized by measuring speed, trajectory and directionality of movement, but there are many unanswered questions about what motivates and defines melanoblast migration. Here, we have established a general mathematical model to simulate the movement of melanoblasts in the epidermis based on biological data, assumptions and hypotheses. Comparisons between experimental data and computer simulations reinforce some biological assumptions, and suggest new ideas for how melanoblasts and keratinocytes might influence each other during development. For example, it appears that melanoblasts instruct each other to allow a homogeneous distribution in the tissue and that keratinocytes may attract melanoblasts until one is stably attached to them. Our model reveals new features of how melanoblasts move and, in particular, suggest that melanoblasts leave a repulsive trail behind them as they move through the skin.

Spatially extended hybrid methods: a review

Many biological and physical systems exhibit behaviour at multiple spatial, temporal or population scales. Multiscale processes provide challenges when they are to be simulated using numerical techniques. While coarser methods such as partial differential equations are typically fast to simulate, they lack the individual-level detail that may be required in regions of low concentration or small spatial scale. However, to simulate at such an individual level throughout a domain and in regions where concentrations are high can be computationally expensive. Spatially coupled hybrid methods provide a bridge, allowing for multiple representations of the same species in one spatial domain by partitioning space into distinct modelling subdomains. Over the past 20 years, such hybrid methods have risen to prominence, leading to what is now a very active research area across multiple disciplines including chemistry, physics and mathematics. There are three main motivations for undertaking this review. Firstly, we have collated a large number of spatially extended hybrid methods and presented them in a single coherent document, while comparing and contrasting them, so that anyone who requires a multiscale hybrid method will be able to find the most appropriate one for their need. Secondly, we have provided canonical examples with algorithms and accompanying code, serving to demonstrate how these types of methods work in practice. Finally, we have presented papers that employ these methods on real biological and physical problems, demonstrating their utility. We also consider some open research questions in the area of hybrid method development and the future directions for the field.

Understanding interactions between populations: Individual based modelling and quantification using pair correlation functions

Understanding the underlying mechanisms that produce the huge variety of swarming and aggregation patterns in animals and cells is fundamental in ecology, developmental biology, and regenerative medicine, to name but a few examples. Depending upon the nature of the interactions between individuals (cells or animals), a variety of different large-scale spatial patterns can be observed in their distribution; examples include cell aggregates, stripes of different coloured skin cells, etc. For the case where all individuals are of the same type (i.e., all interactions are alike), a considerable literature already exists on how the collective organisation depends on the inter-individual interactions. Here, we focus on the less studied case where there are two different types of individuals present. Whilst a number of continuum models of this scenario exist, it can be difficult to compare these models to experimental data, since real cells and animals are discrete. In order to overcome this problem, we develop an agent-based model to simulate some archetypal mechanisms involving attraction and repulsion. However, with this approach (as with experiments), each realisation of the model is different, due to stochastic effects. In order to make useful comparisons between simulations and experimental data, we need to identify the robust features of the spatial distributions of the two species which persist over many realisations of the model (for example, the size of aggregates, degree of segregation or intermixing of the two species). In some cases, it is possible to do this by simple visual inspection. In others, the features of the pattern are not so clear to the unaided eye. In this paper, we introduce a pair correlation function (PCF), which allows us to analyse multi-species spatial distributions quantitatively. We show how the differing strengths of inter-individual attraction and repulsion between species give rise to different spatial patterns, and how the PCF can be used to quantify these differences, even when it might be impossible to recognise them visually.

Multidisciplinary approaches to understanding collective cell migration in developmental biology

Mathematical models are becoming increasingly integrated with experimental efforts in the study of biological systems. Collective cell migration in developmental biology is a particularly fruitful application area for the development of theoretical models to predict the behaviour of complex multicellular systems with many interacting parts. In this context, mathematical models provide a tool to assess the consistency of experimental observations with testable mechanistic hypotheses. In this review, we showcase examples from recent years of multidisciplinary investigations of neural crest cell migration. The neural crest model system has been used to study how collective migration of cell populations is shaped by cell–cell interactions, cell–environmental interactions and heterogeneity between cells. The wide range of emergent behaviours exhibited by neural crest cells in different embryonal locations and in different organisms helps us chart out the spectrum of collective cell migration. At the same time, this diversity in migratory characteristics highlights the need to reconcile or unify the array of currently hypothesized mechanisms through the next generation of experimental data and generalized theoretical descriptions.

The patterns of birthmarks suggest a novel population of melanocyte precursors arising around the time of gastrulation

Systematic work in the mouse and chicken has mapped out two neural crest-derived pathways of melanocyte precursor migration. With these in mind, this study reappraises the patterns of congenital pigmentary disorders in humans and identifies three recurrent patterns consistent across genetically different diseases. Only two of these are seen in diseases known to be melanocyte cell-autonomous. The segmental pattern correlates well with the classical dorsolateral population from animal studies, demonstrating respect of the midline, cranio-caudal axial mixing, unilateral migration and involvement of key epidermally derived structures. Importantly however, the melanocyte precursors responsible for the non-segmental pattern, which demonstrates circular, bilateral migration centred on the midline, and not involving key epidermally derived structures, have not been identified previously. We propose that this population originates around the time of gastrulation, most likely within the mesoderm, and ultimately resides within the dermis. Whether it contributes to mature melanocytes in non-disease states is not known; however, parallels with the patterns of acquired vitiligo would suggest that it does. The third pattern, hypo- or hyperpigmented fine and whorled Blaschko's lines, is proposed to be non-cell-autonomous.

A Multi-stage Representation of Cell Proliferation as a Markov Process

The stochastic simulation algorithm commonly known as Gillespie’s algorithm (originally derived for modelling well-mixed systems of chemical reactions) is now used ubiquitously in the modelling of biological processes in which stochastic effects play an important role. In well-mixed scenarios at the sub-cellular level it is often reasonable to assume that times between successive reaction/interaction events are exponentially distributed and can be appropriately modelled as a Markov process and hence simulated by the Gillespie algorithm. However, Gillespie’s algorithm is routinely applied to model biological systems for which it was never intended. In particular, processes in which cell proliferation is important (e.g. embryonic development, cancer formation) should not be simulated naively using the Gillespie algorithm since the history-dependent nature of the cell cycle breaks the Markov process. The variance in experimentally measured cell cycle times is far less than in an exponential cell cycle time distribution with the same mean.

Here we suggest a method of modelling the cell cycle that restores the memoryless property to the system and is therefore consistent with simulation via the Gillespie algorithm. By breaking the cell cycle into a number of independent exponentially distributed stages, we can restore the Markov property at the same time as more accurately approximating the appropriate cell cycle time distributions. The consequences of our revised mathematical model are explored analytically as far as possible. We demonstrate the importance of employing the correct cell cycle time distribution by recapitulating the results from two models incorporating cellular proliferation (one spatial and one non-spatial) and demonstrating that changing the cell cycle time distribution makes quantitative and qualitative differences to the outcome of the models. Our adaptation will allow modellers and experimentalists alike to appropriately represent cellular proliferation—vital to the accurate modelling of many biological processes—whilst still being able to take advantage of the power and efficiency of the popular Gillespie algorithm.

Variable species densities are induced by volume exclusion interactions upon domain growth

In this work we study the effect of domain growth on spatial correlations in agent populations containing multiple species. This is important as heterogenous cell populations are ubiquitous during the embryonic development of many species. We have previously shown that the long-term behavior of an agent population depends on the way in which domain growth is implemented. We extend this work to show that, depending on the way in which domain growth is implemented, different species dominate in multispecies simulations. Continuum approximations of the lattice-based model that ignore spatial correlations cannot capture this behavior, while those that explicitly account for spatial correlations can. The results presented here show that the precise mechanism of domain growth can determine the long-term behavior of multispecies populations and, in certain circumstances, establish spatially varying species densities.

Using approximate Bayesian computation to quantify cell-cell adhesion parameters in a cell migratory process

In this work, we implement approximate Bayesian computational methods to improve the design of a wound-healing assay used to quantify cell–cell interactions. This is important as cell–cell interactions, such as adhesion and repulsion, have been shown to play a role in cell migration. Initially, we demonstrate with a model of an unrealistic experiment that we are able to identify model parameters that describe agent motility and adhesion, given we choose appropriate summary statistics for our model data. Following this, we replace our model of an unrealistic experiment with a model representative of a practically realisable experiment. We demonstrate that, given the current (and commonly used) experimental set-up, our model parameters cannot be accurately identified using approximate Bayesian computation methods. We compare new experimental designs through simulation, and show more accurate identification of model parameters is possible by expanding the size of the domain upon which the experiment is performed, as opposed to increasing the number of experimental replicates. The results presented in this work, therefore, describe time and cost-saving alterations for a commonly performed experiment for identifying cell motility parameters. Moreover, this work will be of interest to those concerned with performing experiments that allow for the accurate identification of parameters governing cell migratory processes, especially cell migratory processes in which cell–cell adhesion or repulsion are known to play a significant role.

Cell motility is a central process in wound healing and relies on complex cell-cell interactions. A team of mathematicians led by Ruth Baker and Kit Yates at the University of Oxford utilised computer simulations to re-design wound-healing assays that efficiently identify cell motility parameters. New experimental designs through computer simulation can more accurately identify cell motility parameters by expanding the size of the domain upon which the experiment is performed, as opposed to increasing the number of experimental replicates. The results describe time and cost-saving alterations for an experimental method for evaluate complex cell-cell interactions.

Coordination by Cdc42 of Actin, Contractility, and Adhesion for Melanoblast Movement in Mouse Skin

The individual molecular pathways downstream of Cdc42, Rac, and Rho GTPases are well documented, but we know surprisingly little about how these pathways are coordinated when cells move in a complex environment in vivo. In the developing embryo, melanoblasts originating from the neural crest must traverse the dermis to reach the epidermis of the skin and hair follicles. We previously established that Rac1 signals via Scar/WAVE and Arp2/3 to effect pseudopod extension and migration of melanoblasts in skin. Here we show that RhoA is redundant in the melanocyte lineage but that Cdc42 coordinates multiple motility systems independent of Rac1. Similar to Rac1 knockouts, Cdc42 null mice displayed a severe loss of pigmentation, and melanoblasts showed cell-cycle progression, migration, and cytokinesis defects. However, unlike Rac1 knockouts, Cdc42 null melanoblasts were elongated and displayed large, bulky pseudopods with dynamic actin bursts. Despite assuming an elongated shape usually associated with fast mesenchymal motility, Cdc42 knockout melanoblasts migrated slowly and inefficiently in the epidermis, with nearly static pseudopods. Although much of the basic actin machinery was intact, Cdc42 null cells lacked the ability to polarize their Golgi and coordinate motility systems for efficient movement. Loss of Cdc42 de-coupled three main systems: actin assembly via the formin FMNL2 and Arp2/3, active myosin-II localization, and integrin-based adhesion dynamics.

Loss of strumpellin in the melanocytic lineage impairs the WASH Complex but does not affect coat colour

The five-subunit WASH complex generates actin networks that participate in endocytic trafficking, migration and invasion in various cell types. Loss of one of the two subunits WASH or strumpellin in mice is lethal, but little is known about their role in mammals in vivo. We explored the role of strumpellin, which has previously been linked to hereditary spastic paraplegia, in the mouse melanocytic lineage. Strumpellin knockout in melanocytes revealed abnormal endocytic vesicle morphology but no impairment of migration in vitro or in vivo and no change in coat colour. Unexpectedly, WASH and filamentous actin could still localize to vesicles in the absence of strumpellin, although the shape and size of vesicles was altered. Blue native PAGE revealed the presence of two distinct WASH complexes, even in strumpellin knockout cells, revealing that the WASH complex can assemble and localize to endocytic compartments in cells in the absence of strumpellin.

How domain growth is implemented determines the long-term behavior of a cell population through its effect on spatial correlations

Domain growth plays an important role in many biological systems, and so the inclusion of domain growth in models of these biological systems is important to understanding how these systems function. In this work we present methods to include the effects of domain growth on the evolution of spatial correlations in a continuum approximation of a lattice-based model of cell motility and proliferation. We show that, depending on the way in which domain growth is implemented, different steady-state densities are predicted for an agent population. Furthermore, we demonstrate that the way in which domain growth is implemented can result in the evolution of the agent density depending on the size of the domain. Continuum approximations that ignore spatial correlations cannot capture these behaviors, while those that account for spatial correlations do. These results will be of interest to researchers in developmental biology, as they suggest that the nature of domain growth can determine the characteristics of cell populations.