Every cell contains several millions of diffusing and reacting biological molecules. The interactions between these molecules ultimately manifest themselves in all aspects of life, from the smallest bacterium to the largest whale. One of the greatest open scientific challenges is to understand how the microscopic chemistry determines the macroscopic biology. Key to this challenge is the development of mathematical and computational models of biochemistry with molecule-level detail, but which are sufficiently coarse to enable the study of large systems at the cell or organism scale. Two such models are in common usage: the reaction-diffusion master equation, and Brownian dynamics. These models are utterly different in both their history and in their approaches to chemical reactions and diffusion, but they both seek to address the same reaction-diffusion question. Here we make an in-depth study into the physical validity of these models under various biological conditions, determining when they can reliably be used. Taking each model in turn, we propose modifications to the models to better model the realities of the cellular environment, and to enable more efficient computational implementations. We use the models to make predictions about how and why cells behave the way they do, from mechanisms of self-organisation to noise reduction. We conclude that both models are extremely powerful tools for clarifying the details of the mysterious relationship between chemistry and biology.

Quantitative descriptions of reaction kinetics formulated at the stochastic mesoscopic level are frequently used to study various aspects of regulation and control in models of cellular control systems. For this type of systems, numerical simulation offers a variety of challenges caused by the high dimensionality of the problem and the multiscale properties often displayed by the biochemical model. In this thesis I have studied several aspects of stochastic simulation of both well-stirred and spatially heterogenous systems. In the well-stirred case, a hybrid method is proposed that reduces the dimension and stiffness of a model. We also demonstrate how both a high performance implementation and a variance reduction technique based on quasi-Monte Carlo can reduce the computational cost to estimate the probability density of the system. In the spatially dependent case, the use of unstructured, tetrahedral meshes to sample realizations of the stochastic process is proposed. Using such meshes, we then extend the reaction-diffusion framework to incorporate active transport of cellular cargo in a seamless manner. Finally, two multilevel methods for spatial stochastic simulation are considered. One of them is a space-time adaptive method combining exact stochastic, approximate stochastic and macroscopic modeling levels to reduce the simualation cost. The other method blends together mesoscale and microscale simulation methods to locally increase modeling resolution. / eSSENCE

Mathematical models are important tools in systems biology, since the regulatory networks in biological cells are too complicated to understand by biological experiments alone. Analytical solutions can be derived only for the simplest models and numerical simulations are necessary in most cases to evaluate the models and their properties and to compare them with measured data. This thesis focuses on the mesoscopic simulation level, which captures both, space dependent behavior by diffusion and the inherent stochasticity of cellular systems. Space is partitioned into compartments by a mesh and the number of molecules of each species in each compartment gives the state of the system. We first examine how to compute the jump coefficients for a discrete stochastic jump process on unstructured meshes from a first exit time approach guaranteeing the correct speed of diffusion. Furthermore, we analyze different methods leading to non-negative coefficients by backward analysis and derive a new method, minimizing both the error in the diffusion coefficient and in the particle distribution. The second part of this thesis investigates macromolecular crowding effects. A high percentage of the cytosol and membranes of cells are occupied by molecules. This impedes the diffusive motion and also affects the reaction rates. Most algorithms for cell simulations are either derived for a dilute medium or become computationally very expensive when applied to a crowded environment. Therefore, we develop a multiscale approach, which takes the microscopic positions of the molecules into account, while still allowing for efficient stochastic simulations on the mesoscopic level. Finally, we compare on- and off-lattice models on the microscopic level when applied to a crowded environment.

Self-organizing patterns arise in a variety of ways in nature, the complex patterning observed on animal coats is such an example. It is already known that the mechanisms responsible for pattern formation starts at the developmental stage of an embryo. However, the actual process determining cell fate has been, and still is, unknown. The mathematical interest for pattern formation emerged from the theories formulated by the mathematician and computer scientist Alan Turing in 1952. He attempted to explain the mechanisms behind morphogenesis and how the process of spatial cell differentiation from homogeneous cells lead to organisms with different complexities and shapes. Turing formulated a mathematical theory and proposed a reaction-diffusion system where morphogens, a postulated chemically active substance, moderated the whole mechanism. He concluded that this process was stable as long as diffusion was neglected; otherwise this would lead to a diffusion-driven instability, which is the fundamental part of pattern formation. The mathematical theory describing this process consists of solving partial differential equations and Turing considered deterministic reaction-diffusion systems.   This thesis will start with introducing the reader to the problem and then gradually build up the mathematical theory needed to get an understanding of the stochastic reaction-diffusion systems that is the focus of the thesis. This study will to a large extent simulate stochastic systems using numerical computations and in order to be computationally feasible a compartment-based model will be used. Noise is an inherent part of such systems, so the study will also discuss the effects of noise and morphogen kinetics on different geometries with boundaries of different complexities from one-dimensional cases up to three-dimensions.

Complex spatial structures in biology arise from random interactions of molecules. These molecular interactions can be studied using spatial stochastic models, such as Reaction Diffusion Master Equation (RDME), a mesoscopic model that subdivides the spatial domain into smaller, well mixed grid cells, in which the macroscopic diffusion-controlled reactions take place. While RDME has been widely used to study how fluctuations in number of molecules affect spatial patterns, simulations are computationally expensive and it requires a lower bound for grid cell size to avoid an apparent unphysical loss of bimolecular reactions. In this thesis, we propose Gaussian Reaction Diffusion Master Equation (GRDME), a novel model in the RDME framework, based on the discretization of the Laplace operator with Particle Strength Exchange (PSE) method with a Gaussian kernel. We show that GRDME is a computationally efficient model compared to RDME. We further resolve the controversy regarding the loss of bimolecular reactions and argue that GRDME can flexibly bridge the diffusion-controlled and ballistic regimes in mesoscopic simulations involving multiple species. To efficiently simulate GRDME, we develop Gaussian Next Subvolume Method (GNSM). GRDME simulated with GNSM up to six-times lower computational cost for a three-dimensional simulation, providing a significant computational advantage for modeling three-dimensional systems. The computational cost can be further lowered by increasing the so-called smoothing length of the Gassian jumps. We develop a guideline to estimate the grid resolution below which RDME and GRDME exhibit loss of bimolecular reactions. This loss of reactions has been considered unphysical by others. Here we show that this loss of bimolecular reactions is consistent with the well-established theory on diffusion-controlled reaction rates by Collins and Kimball, provided that the rate of bimolecular propensity is interpreted as the rate of the ballistic step, rather than the macroscopic reaction rate. We show that the reaction radius is set by the grid resolution. Unlike RDME, GRDME enables us to explicitly model various sizes of the molecules. Using this insight, we explore the diffusion-limited regime of reaction dynamics and discover that diffusion-controlled systems resemble small, discrete systems. Others have shown that a reaction system can have discreteness-induced state inversion, a phenomenon where the order of the concentrations differs when the system size is small. We show that the same reaction system also has diffusion-controlled state inversion, where the order of concentrations changes, when the diffusion is slow. In summary, we show that GRDME is a computationally efficient model, which enables us to include the information of the molecular sizes into the model.:1 Modeling Mesoscopic Biology 1.1 RDME Models Mesoscopic Stochastic Spatial Phenomena 1.2 A New Diffusion Model Presents an Opportunity For A More Efficient RDME 1.3 Can A New Diffusion Model Provide Insights Into The Loss Of Reactions? 1.4 Overview 2 Preliminaries 2.1 Reaction Diffusion Master Equation 2.1.1 Chemical Master Equation 2.1.2 Diffusion-controlled Bimolecular Reaction Rate 2.1.3 RDME is an Extention of CME to Spatial Problems 2.2 Next Subvolume Method 2.2.1 First Reaction Method 2.2.2 NSM is an Efficient Spatial Stochastic Algorithm for RDME 2.3 Discretization of the Laplace Operator Using Particle Strength Exchange 2.4 Summary 3 Gaussian Reaction Diffusion Master Equation 3.1 Design Constraints for the Diffusion Model in the RDME Framework 3.2 Gaussian-jump-based Model for RDME 3.3 Summary 4 Gaussian Next Subvolume Method 4.1 Constructing the neighborhood N 4.2 Finding the Diffusion Event 4.3 Comparing GNSM to NSM 4.4 Summary 5 Limits of Validity for (G)RDME with Macroscopic Bimolecular Propensity Rate 5.1 Previous Works 5.2 hmin Based on the Kuramoto length of a Grid Cell 5.3 hmin of the Two Limiting Regimes 5.4 hmin of Bimolecular Reactions for the Three Cases of Dimensionality 5.5 hmin of GRDME in Comparison to hmin of RDME 5.6 Summary 6 Numerical Experiments To Verify Accuracy, Efficiency and Validity of GRDME 6.1 Accuracy of the Diffusion Model 6.2 Computational Cost 6.3 hmin and Reaction Loss for (G)RDME With Macroscopic Bimolecular Propensity Rate kCK 6.3.1 Homobiomlecular Reaction With kCK at the Ballistic Limit 6.3.2 Homobiomlecular Reaction With kCK at the Diffusional Limit 6.3.3 Heterobiomlecular Reaction With kCK at the Ballistic Limit 6.4 Summary 7 (G)RDME as a Spatial Model of Collins-Kimball Diffusion-controlled Reaction Dynamics 7.1 Loss of Reactions in Diffusion-controlled Reaction Systems 7.2 The Loss of Reactions in (G)RDME Can Be Explained by Collins Kimball Theory 7.3 Cell Width h Sets the Reaction Radius σ∗ 7.4 Smoothing Length ε′ Sets the Size of the Molecules in the System 7.5 Heterobimolecular Reactions Can Only Be Modeled With GRDME 7.6 Zeroth Order Reactions Impose a Lower Limit on Diffusivity Dmin 7.6.1 Consistency of (G)RDME Could Be Improved by Redesigning Zeroth Order Reactions 7.7 Summary 8 Difussion-Controlled State Inversion 8.1 Diffusion-controlled Systems Resemble Small Systems 8.2 Slow Diffusion Leads to an Inversion of Steady States 8.3 Summary 9 Conclusion and Outlook 9.1 Two Physical Interpretations of (G)RDME 9.2 Advantages of GRDME 9.3 Towards Numerically Consistent (G)RDME 9.4 Exploring Mesoscopic Biology With GRDME Bibliography