Search

Reaction–diffusion models describing the movement, reproduction and death of individuals within a population are key mathematical modelling tools with widespread applications in mathematical biology. A diverse range of such continuum models have been applied in various biological contexts by choosing different flux and source terms in the reaction–diffusion framework. For example, to describe collective spreading of cell populations, the flux term may be chosen to reflect various movement mechanisms, such as random motion (diffusion), adhesion, haptotaxis, chemokinesis and chemotaxis. The choice of flux terms in specific applications, such as wound healing, is usually made heuristically, and rarely is it tested quantitatively against detailed cell density data. More generally, in mathematical biology, the questions of model validation and model selection have not received the same attention as the questions of model development and model analysis. Many studies do not consider model validation or model selection, and those that do often base the selection of the model on residual error criteria after model calibration is performed using nonlinear regression techniques. In this work, we present a model selection case study, in the context of cell invasion, with a very detailed experimental data set. Using Bayesian analysis and information criteria, we demonstrate that model selection and model validation should account for both residual errors and model complexity. These considerations are often overlooked in the mathematical biology literature. The results we present here provide a clear methodology that can be used to guide model selection across a range of applications. Furthermore, the case study we present provides a clear example where neglecting the role of model complexity can give rise to misleading outcomes.

Adaptive finite element models have allowed researchers to test hypothetical relationships between the local mechanical environment and the healing of bone fractures. However, their predictive power has not yet been demonstrated by testing hypotheses ahead of experimental testing. In this study, an established mechano-biological scheme was used in an iterative finite element simulation of sheep tibial osteotomy healing under a hypothetical fixation regime, “inverse dynamisation”. Tissue distributions, interfragmentary movement and stiffness across the fracture site were compared between stiff and flexible fixation conditions and scenarios in which fixation stiffness was increased at a discrete time-point. The modelling work was conducted blind to the experimental study to be published subsequently. The simulations predicted the fastest and most direct healing under constant stiff fixation, and the slowest healing under flexible fixation. Although low fixation stiffness promoted more callus formation prior to bridging, this conferred little additional stiffness to the fracture in the first 5 weeks. Thus, while switching to stiffer fixation facilitated rapid subsequent bridging of the fracture, no advantage of inverse dynamisation could be demonstrated. In vivo data remains necessary to conclusively test this treatment protocol and this will, in turn, provide an evaluation of the model’s performance. The publication of both hypotheses and their computational simulation, prior to experimental testing, offers an appealing means to test the predictive power of mechano-biological models.

Stochasticity is a key characteristic of intracellular processes such as gene regulation and chemical signalling. Therefore, characterizing stochastic effects in biochemical systems is essential to understand the complex dynamics of living things. Mathematical idealizations of biochemically reacting systems must be able to capture stochastic phenomena. While robust theory exists to describe such stochastic models, the computational challenges in exploring these models can be a significant burden in practice since realistic models are analytically intractable. Determining the expected behaviour and variability of a stochastic biochemical reaction network requires many probabilistic simulations of its evolution. Using a biochemical reaction network model to assist in the interpretation of time-course data from a biological experiment is an even greater challenge due to the intractability of the likelihood function for determining observation probabilities. These computational challenges have been subjects of active research for over four decades. In this review, we present an accessible discussion of the major historical developments and state-of-the-art computational techniques relevant to simulation and inference problems for stochastic biochemical reaction network models. Detailed algorithms for particularly important methods are described and complemented with Matlab®implementations. As a result, this review provides a practical and accessible introduction to computational methods for stochastic models within the life sciences community.

Numerical modelling has gained popularity in many science and engineering streams due to the economic feasibility and advanced analytical features compared to conventional experimental and theoretical models. Food drying is one of the areas where numerical modelling is increasingly applied to improve drying process performance and product quality. This investigation applies a three dimensional (3-D) Smoothed Particle Hydrodynamics (SPH) and Coarse-Grained (CG) numerical approach to predict the morphological changes of different categories of food-plant cells such as apple, grape, potato and carrot during drying. To validate the model predictions, experimental findings from in-house experimental procedures (for apple) and sources of literature (for grape, potato and carrot) have been utilised. The subsequent comaprison indicate that the model predictions demonstrate a reasonable agreement with the experimental findings, both qualitatively and quantitatively. In this numerical model, a higher computational accuracy has been maintained by limiting the consistency error below 1% for all four cell types. The proposed meshfree-based approach is well-equipped to predict the morphological changes of plant cellular structure over a wide range of moisture contents (10% to 100% dry basis). Compared to the previous 2-D meshfree-based models developed for plant cell drying, the proposed model can draw more useful insights on the morphological behaviour due to the 3-D nature of the model. In addition, the proposed computational modelling approach has a high potential to be used as a comprehensive tool in many other tissue morphology related investigations.

In this paper, we investigate an inverse differential game approach to modelling the mid-air collision avoidance behaviours of birds. We propose a general method for estimating the cost-functional parameters of a noncooperative differential game from partial-state measurements of an open-loop Nash equilibrium. We apply the method to data of birds performing mid-air collision avoidance. Our analysis suggests that a differential game model provides a close description of the observed bird behaviours, and could provide new insights for the design of collision avoidance strategies for unmanned aircraft.

Iterative computational models have been used to investigate the regulation of bone fracture healing by local mechanical conditions. Although their predictions replicate some mechanical responses and histological features, they do not typically reproduce the predominantly radial hard callus growth pattern observed in larger mammals. We hypothesised that this discrepancy results from an artefact of the models' initial geometry. Using axisymmetric finite element models, we demonstrated that pre-defining a field of soft tissue in which callus may develop introduces high deviatoric strains in the periosteal region adjacent to the fracture. These bone-inhibiting strains are not present when the initial soft tissue is confined to a thin periosteal layer. As observed in previous healing models, tissue differentiation algorithms regulated by deviatoric strain predicted hard callus forming remotely and growing towards the fracture. While dilatational strain regulation allowed early bone formation closer to the fracture, hard callus still formed initially over a broad area, rather than expanding over time. Modelling callus growth from a thin periosteal layer successfully predicted the initiation of hard callus growth close to the fracture site. However, these models were still susceptible to elevated deviatoric strains in the soft tissues at the edge of the hard callus. Our study highlights the importance of the initial soft tissue geometry used for finite element models of fracture healing. If this cannot be defined accurately, alternative mechanisms for the prediction of early callus development should be investigated.

Fibroblasts and their activated phenotype, myofibroblasts, are the primary cell types involved in the contraction associated with dermal wound healing. Recent experimental evidence indicates that the transformation from fibroblasts to myofibroblasts involves two distinct processes: the cells are stimulated to change phenotype by the combined actions of transforming growth factor β (TGFβ) and mechanical tension. This observation indicates a need for a detailed exploration of the effect of the strong interactions between the mechanical changes and growth factors in dermal wound healing. We review the experimental findings in detail and develop a model of dermal wound healing that incorporates these phenomena. Our model includes the interactions between TGFβ and collagenase, providing a more biologically realistic form for the growth factor kinetics than those included in previous mechanochemical descriptions. A comparison is made between the model predictions and experimental data on human dermal wound healing and all the essential features are well matched.

Biological networks are highly modular and contain a large number of clusters, which are often associated with a specific biological function or disease. Identifying these clusters, or modules, is therefore valuable, but it is not trivial. In this article we propose a recursive method based on the Louvain algorithm for community detection and the PageRank algorithm for authoritativeness weighting in networks. PageRank is used to initialise the weights of nodes in the biological network; the Louvain algorithm with the Newman-Girvan criterion for modularity is then applied to the network to identify modules. Any identified module with more than k nodes is further processed by recursively applying PageRank and Louvain, until no module contains more than k nodes (where k is a parameter of the method, no greater than 100). This method is evaluated on a heterogeneous set of six biological networks from the Disease Module Identification DREAM Challenge. Empirical findings suggest that the method is effective in identifying a large number of significant modules, although with substantial variability across restarts of the method.

Introduction & Aims Optimising fracture treatments requires a sound understanding of relationships between stability, callus development and healing outcomes. This has been the goal of computational modelling, but discrepancies remain between simulations and experimental results. We compared healing patterns vs fixation stiffness between a novel computational callus growth model and corresponding experimental data. Hypothesis We hypothesised that callus growth is stimulated by diffusible signals, whose production is in turn regulated by mechanical conditions at the fracture site. We proposed that introducing this scheme into computational models would better replicate the observed tissue patterns and the inverse relationship between callus size and fixation stiffness. Method Finite element models of bone healing under stiff and flexible fixation were constructed, based on the parameters of a parallel rat femoral osteotomy study. An iterative procedure was implemented, to simulate the development of callus and its mechanical regulation. Tissue changes were regulated according to published mechano-biological criteria. Predictions of healing patterns were compared between standard models, with a pre-defined domain for callus development, and a novel approach, in which periosteal callus growth is driven by a diffusible signal. Production of this signal was driven by local mechanical conditions. Finally, each model’s predictions were compared to the corresponding histological data. Results Models in which healing progressed within a prescribed callus domain predicted that greater interfragmentary movements would displace early periosteal bone formation further from the fracture. This results from artificially large distortional strains predicted near the fracture edge. While experiments showed increased hard callus size under flexible fixation, this was not reflected in the standard models. Allowing the callus to grow from a thin soft tissue layer, in response to a mechanically stimulated diffusible signal, results in a callus shape and tissue distribution closer to those observed histologically. Importantly, the callus volume increased with increasing interfragmentary movement. Conclusions A novel method to incorporate callus growth into computational models of fracture healing allowed us to successfully capture the relationship between callus size and fixation stability observed in our rat experiments. This approach expands our toolkit for understanding the influence of different fixation strategies on healing outcomes.

We present a mathematical model of glioma spread based on cellular movement by receptor-mediated haptotaxis, local proteolysis of healthy tissue components by glioma-derived proteinases, malignant proliferative enhancement and host up-regulation of specific key extracellular matrix (ECM) components in response to the invading glioma. We subsequently consider the nature of glioma-host interactions as predicted by our model in order to test the hypothesis given in (Knott et al. (1998) that production of adhesive ECM components by the brain in response to the invading glioma may have the counter-intuitive effect of enhancing glioma invasion by assisting haptotactic migration. We suggest that host production of certain adhesive ECM chemicals can have a profound effect on both glioma invasion speed and the character of the glioma-host interface. In particular, we conclude that up-regulation of host ECM production in the vicinity of the glioma may produce a less diffuse glioma, providing clearer demarcation between glioma and healthy tissue, and thus improving the possibility of surgical resection within reasonable bounds.

Background Directed cell migration is essential for normal development. In most of the migratory cell populations that have been analysed in detail to date, all of the cells migrate as a collective from one location to another. However, there are also migratory cell populations that must populate the areas through which they migrate, and thus some cells get left behind while others advance. Very little is known about how individual cells behave to achieve concomitant directional migration and population of the migratory route. We examined the behavior of enteric neural crest-derived cells (ENCCs), which must both advance caudally to reach the anal end and populate each gut region. Results The behaviour of individual ENCCs was examined using live imaging and mice in which ENCCs express a photoconvertible protein. We show that individual ENCCs exhibit very variable directionalities and speed; as the migratory wavefront of ENCCs advances caudally, each gut region is populated primarily by some ENCCs migrating non-directionally. After populating each region, ENCCs remain migratory for at least 24 hours. Endothelin receptor type B (EDNRB) signaling is known to be essential for the normal advance of the ENCC population. We now show that perturbation of EDNRB principally affects individual ENCC speed rather than directionality. The trajectories of solitary ENCCs, which occur transiently at the wavefront, were consistent with an unbiased random walk and so cell-cell contact is essential for directional migration. ENCCs migrate in close association with neurites. We showed that although ENCCs often use neurites as substrates, ENCCs lead the way, neurites are not required for chain formation and neurite growth is more directional than the migration of ENCCs as a whole. Conclusions Each gut region is initially populated by sub-populations of ENCCs migrating non-directionally, rather than stopping. This might provide a mechanism for ensuring a uniform density of ENCCs along the growing gut.

Hypertrophic scars arise when there is an overproduction of collagen during wound healing. These are often associated with poor regulation of the rate of programmed cell death(apoptosis) of the cells synthesizing the collagen or by an exuberant inflammatory response that prolongs collagen production and increases wound contraction. Severe contractures that occur, for example, after a deep burn can cause loss of function especially if the wound is over a joint such as the elbow or knee. Recently, we have developed a morphoelastic mathematical model for dermal repair that incorporates the chemical, cellular and mechanical aspects of dermal wound healing. Using this model, we examine pathological scarring in dermal repair by first assuming a smaller than usual apoptotic rate for myofibroblasts, and then considering a prolonged inflammatory response, in an attempt to determine a possible optimal intervention strategy to promote normal repair, or terminate the fibrotic scarring response. Our model predicts that in both cases it is best to apply the intervention strategy early in the wound healing response. Further, the earlier an intervention is made, the less aggressive the intervention required. Finally, if intervention is conducted at a late time during healing, a significant intervention is required; however, there is a threshold concentration of the drug or therapy applied, above which minimal further improvement to wound repair is obtained.

The accuracy of data derived from linked-segment models depends on how well the system has been represented. Previous investigations describing the gait of persons with partial foot amputation did not account for the unique anthropometry of the residuum or the inclusion of a prosthesis and footwear in the model and, as such, are likely to have underestimated the magnitude of the peak joint moments and powers. This investigation determined the effect of inaccuracies in the anthropometric input data on the kinetics of gait. Toward this end, a geometric model was developed and validated to estimate body segment parameters of various intact and partial feet. These data were then incorporated into customized linked-segment models, and the kinetic data were compared with that obtained from conventional models. Results indicate that accurate modeling increased the magnitude of the peak hip and knee joint moments and powers during terminal swing. Conventional inverse dynamic models are sufficiently accurate for research questions relating to stance phase. More accurate models that account for the anthropometry of the residuum, prosthesis, and footwear better reflect the work of the hip extensors and knee flexors to decelerate the limb during terminal swing phase.

The 3-dimensional coordinates of alpha-carbon atoms of proteins are used to distinguish the protein structural classes based on recurrence quantification analysis (RQA). We consider two independent variables from RQA of coordinates of alpha-carbon atoms, %determ1 and %determ2, which were defined by Webber et al. [C.L. Webber Jr., A. Giuliani, J.P. Zbilut, A. Colosimo, Proteins Struct. Funct. Genet. 44 (2001) 292]. The variable %determ2 is used to define two new variables, %determ21 and %determ22. Then three variables %determ1, %determ21 and %determ22 are used to construct a 3-dimensional variable space. Each protein is represented by a point in this variable space. The points corresponding to proteins from the α, β, α+β and α/β structural classes position into different areas in this variable space. In order to give a quantitative assessment of our clustering on the selected proteins, Fisher's discriminant algorithm is used. Numerical results indicate that the discriminant accuracies are very high and satisfactory.

Aims Fractures that are not rigidly fixed heal via the formation of a fibro-cartilaginous callus, which is progressively converted to bone. We examine the hypothesis that callus morphology is optimal to achieve a defined limiting stiffness of the bone-callus construct. We demonstrate the use of this principle in designing the starting geometry for iterative computational models of the fracture healing process. Methods For a given fracture geometry, optimal callus dimensions are determined by iteratively increasing the callus size until a defined stiffness is attained, as determined by finite element (FE) analysis of the bone-callus construct. This process is performed using various canonical shapes and material properties for the callus. The optimal morphologies thus predicted are compared to typical histological and radiological observations. The impact of callus morphology on the outcomes of an iterative FE model of fracture healing is investigated using experimental and "optimal" callus dimensions. This model is similar to those published previously, and plots the maturation of tissues in the callus until bony union. We analyse the predicted healing outcomes when optimal, sub-optimal and supra-optimal calluses are used as initial states, by evaluating against existing experimental data. Results Using an FE modelling scheme and comparing its results with typical histology, we establish that callus initially forms an optimally-determined morphology. Secondly, we illustrate the impact of sub-optimal and supra-optimal calluses on the outcomes of an evolutionary FE model, to establish the efficacy of employing optimal callus dimensions in computational models of fracture healing. Conclusions Defining principles for optimal callus morphology provides a feasible starting point for computational models of fracture healing. This may be useful in avoiding modelling artefacts when early-stage experimental data is insufficient. Such models may provide new insights into how different fracture stabilisation methods manipulate the healing process.

The mechanism by which mechanical stimulus of reparative tissues directs the pattern of healing of a bone fracture is not understood. Several hypotheses have been developed that predict the ossification pattern during healing for a given ambient mechanical environment. These have remained unproved because of the absence of data on stress fields in the reparative tissue of real fractures. The present study examines the predictive performance of the most recent hypothesis that was proposed by Claes and Heigele (J. Biomech. 32 (1999) 255), against measured and calculated data from the clinical fracture reported by Gardner et al. (J. Biomech. 33 (2000) 415). The hypothesis was used to predict ossification from the preceding stress and strain fields present in the FEM of Gardner et al. at four temporal stages during healing. Predictions were then compared with the observed differentiation and maturation. During early healing of the interfragmentary gap region, the hypothesis correctly predicted the formation of connective tissue and fibrocartilage, and during later healing it correctly predicted the beginning of endochondral ossification. At the periphery of the periosteal callus, the hypothesis correctly predicted intramembraneous ossification during early healing, and also its thickening by endochondral ossification during later healing. However, the hypothesis incorrectly predicted intramembraneous ossification during early healing of the main periosteal callus, although in later healing it correctly predicted endochondral ossification.

The exact mechanism by which mechanical stimulus regulates the healing process of a bone fracture is not understood. This has led to the development of several hypotheses that predict the pattern of differentiation of tissue during healing that may arise from characteristic fields of stress or strain at the fracture. These have so far remained unproved because data on stress fields in actual fracture tissue have been unavailable until recently. Thus the present study examines the predictive performance of the hypothesis proposed in J Orthop Res 6 (1988) 736, against measured and calculated data reported in J Biomech 33 (2000) 415, using a 2D FEM of a clinical fracture. The hypothesis was used to predict the influence of stress fields present in the Gardner et al. tissues at four temporal stages during healing. These predictions were then correlated with callus-size, rate of endochondral ossification and ossification pattern subsequently observed by Gardner et al. in the clinical fracture. Results corroborate the hypothesis that high octahedral shear stresses may increase the size of the callus during the initial phase of healing, and they also suggest that this may be true during the later stages of the fracture fixation period. However, compressive dilatational stresses were not found to inhibit endochondral ossification, as suggested by the hypothesis. Although high shear stresses were found in regions indicative of fibrous tissue as postulated by the hypothesis, this was not found to be the case for high tensile dilatational stresses. Also, contour diagrams of Osteogenic index (I) indicated only limited correlation with callus maturation and the pattern of healing. Therefore, the hypothesis was not wholly successful in predicting healing pattern.

Using six kinds of lattice types (4 x 4, 5 x 5, and 6 x 6 square lattices; 3 x 3 x 3 cubic lattice; and 2+3+4+3+2 and 4+5+6+5+4 triangular lattices), three different size alphabets (HP, HNUP, and 20 letters), and two energy functions, the designability of protein structures is calculated based on random samplings of structures and common biased sampling (CBS) of protein sequence space. Then three quantities stability (average energy gap), foldability, and partnum of the structure, which are defined to elucidate the designability, are calculated. The authors find that whatever the type of lattice, alphabet size, and energy function used, there will be an emergence of highly designable (preferred) structure. For all cases considered, the local interactions reduce degeneracy and make the designability higher. The designability is sensitive to the lattice type, alphabet size, energy function, and sampling method of the sequence space. Compared with the random sampling method, both the CBS and the Metropolis Monte Carlo sampling methods make the designability higher. The correlation coefficients between the designability, stability, and foldability are mostly larger than 0.5, which demonstrate that they have strong correlation relationship. But the correlation relationship between the designability and the partnum is not so strong because the partnum is independent of the energy. The results are useful in practical use of the designability principle, such as to predict the protein tertiary structure.